diff --git a/dust3d.pro b/dust3d.pro
index f6d53077..b220abc7 100644
--- a/dust3d.pro
+++ b/dust3d.pro
@@ -477,6 +477,9 @@ HEADERS += src/triangulatefaces.h
 SOURCES += src/booleanmesh.cpp
 HEADERS += src/booleanmesh.h
 
+SOURCES += src/imageskeletonextractor.cpp
+HEADERS += src/imageskeletonextractor.h
+
 SOURCES += src/main.cpp
 
 HEADERS += src/version.h
diff --git a/src/imageskeletonextractor.cpp b/src/imageskeletonextractor.cpp
new file mode 100644
index 00000000..f00a234a
--- /dev/null
+++ b/src/imageskeletonextractor.cpp
@@ -0,0 +1,104 @@
+#include "imageskeletonextractor.h"
+
+// This is an implementation of the following paper:
+// <A Fast Parallel Algorithm for Thinning Digital Patterns>
+// T. Y. ZHANG and C. Y. SUEN
+
+const int ImageSkeletonExtractor::m_targetHeight = 256;
+
+ImageSkeletonExtractor::~ImageSkeletonExtractor()
+{
+    delete m_image;
+    delete m_grayscaleImage;
+}
+
+void ImageSkeletonExtractor::setImage(QImage *image)
+{
+    delete m_image;
+    m_image = image;
+}
+
+QImage *ImageSkeletonExtractor::takeResultGrayscaleImage()
+{
+    QImage *resultImage = m_grayscaleImage;
+    m_grayscaleImage = nullptr;
+    return resultImage;
+}
+
+bool ImageSkeletonExtractor::firstSubiterationSatisfied(int i, int j)
+{
+    if (!isBlack(i, j))
+        return false;
+    auto blackNeighbors = countBlackNeighbors(i, j);
+    if (blackNeighbors < 2 || blackNeighbors > 6)
+        return false;
+    auto neighborTransitions = countNeighborTransitions(i, j);
+    if (1 != neighborTransitions)
+        return false;
+    if (isBlack(i + neighborOffsets[P2].first, j + neighborOffsets[P2].second) &&
+            isBlack(i + neighborOffsets[P4].first, j + neighborOffsets[P4].second) &&
+            isBlack(i + neighborOffsets[P6].first, j + neighborOffsets[P6].second)) {
+        return false;
+    }
+    if (isBlack(i + neighborOffsets[P4].first, j + neighborOffsets[P4].second) &&
+            isBlack(i + neighborOffsets[P6].first, j + neighborOffsets[P6].second) &&
+            isBlack(i + neighborOffsets[P8].first, j + neighborOffsets[P8].second)) {
+        return false;
+    }
+    return true;
+}
+
+bool ImageSkeletonExtractor::secondSubiterationSatisfied(int i, int j)
+{
+    if (!isBlack(i, j))
+        return false;
+    auto blackNeighbors = countBlackNeighbors(i, j);
+    if (blackNeighbors < 2 || blackNeighbors > 6)
+        return false;
+    auto neighborTransitions = countNeighborTransitions(i, j);
+    if (1 != neighborTransitions)
+        return false;
+    if (isBlack(i + neighborOffsets[P2].first, j + neighborOffsets[P2].second) &&
+            isBlack(i + neighborOffsets[P4].first, j + neighborOffsets[P4].second) &&
+            isBlack(i + neighborOffsets[P8].first, j + neighborOffsets[P8].second)) {
+        return false;
+    }
+    if (isBlack(i + neighborOffsets[P2].first, j + neighborOffsets[P2].second) &&
+            isBlack(i + neighborOffsets[P6].first, j + neighborOffsets[P6].second) &&
+            isBlack(i + neighborOffsets[P8].first, j + neighborOffsets[P8].second)) {
+        return false;
+    }
+    return true;
+}
+
+void ImageSkeletonExtractor::extract()
+{
+    m_grayscaleImage = new QImage(m_image->convertToFormat(
+        QImage::Format_Grayscale8).scaled(
+            QSize(m_targetHeight, m_targetHeight), Qt::KeepAspectRatio));
+    
+    while (true) {
+        std::vector<std::pair<int, int>> firstSatisfied;
+        for (int i = 1; i < (int)m_grayscaleImage->width() - 1; ++i) {
+            for (int j = 1; j < (int)m_grayscaleImage->height() - 1; ++j) {
+                if (firstSubiterationSatisfied(i, j))
+                    firstSatisfied.push_back(std::make_pair(i, j));
+            }
+        }
+        for (const auto &it: firstSatisfied)
+            setWhite(it.first, it.second);
+        std::vector<std::pair<int, int>> secondSatisfied;
+        for (int i = 1; i < (int)m_grayscaleImage->width() - 1; ++i) {
+            for (int j = 1; j < (int)m_grayscaleImage->height() - 1; ++j) {
+                if (secondSubiterationSatisfied(i, j))
+                    secondSatisfied.push_back(std::make_pair(i, j));
+            }
+        }
+        for (const auto &it: secondSatisfied)
+            setWhite(it.first, it.second);
+        if (firstSatisfied.empty() && secondSatisfied.empty())
+            break;
+        printf("firstSatisfied:%d\r\n", firstSatisfied.size());
+        printf("secondSatisfied:%d\r\n", secondSatisfied.size());
+    }
+}
diff --git a/src/imageskeletonextractor.h b/src/imageskeletonextractor.h
new file mode 100644
index 00000000..8792375a
--- /dev/null
+++ b/src/imageskeletonextractor.h
@@ -0,0 +1,83 @@
+#ifndef DUST3D_IMAGE_SKELETON_EXTRACTOR_H
+#define DUST3D_IMAGE_SKELETON_EXTRACTOR_H
+#include <QImage>
+#include <QObject>
+#include <vector>
+
+class ImageSkeletonExtractor : QObject
+{
+    Q_OBJECT
+public:
+    const std::vector<std::pair<int, int>> neighborOffsets = {
+        { 0, -1},
+        { 1, -1},
+        { 1,  0},
+        { 1,  1},
+        { 0,  1},
+        {-1,  1},
+        {-1,  0},
+        {-1, -1},
+    };
+    enum {
+        P2 = 0,
+        P3,
+        P4,
+        P5,
+        P6,
+        P7,
+        P8,
+        P9
+    };
+
+    ~ImageSkeletonExtractor();
+    void setImage(QImage *image);
+    void extract();
+    QImage *takeResultGrayscaleImage();
+private:
+    QImage *m_image = nullptr;
+    QImage *m_grayscaleImage = nullptr;
+    static const int m_targetHeight;
+    
+    bool isBlack(int i, int j)
+    {
+        return QColor(m_grayscaleImage->pixel(i, j)).black() > 0;
+    }
+    
+    bool isWhite(int i, int j)
+    {
+        return !isBlack(i, j);
+    }
+    
+    void setWhite(int i, int j)
+    {
+        m_grayscaleImage->setPixel(i, j, qRgb(255, 255, 255));
+    }
+    
+    int countNeighborTransitions(int i, int j)
+    {
+        int num = 0;
+        for (size_t m = 0; m < neighborOffsets.size(); ++m) {
+            size_t n = (m + 1) % neighborOffsets.size();
+            if (isWhite(i + neighborOffsets[m].first, j + neighborOffsets[m].second) &&
+                    isBlack(i + neighborOffsets[n].first, j + neighborOffsets[n].second)) {
+                ++num;
+            }
+        }
+        return num;
+    }
+    
+    int countBlackNeighbors(int i, int j)
+    {
+        int num = 0;
+        for (const auto &it: neighborOffsets) {
+            if (isBlack(i + it.first, j + it.second))
+                ++num;
+        }
+        return num;
+    }
+    
+    bool firstSubiterationSatisfied(int i, int j);
+    bool secondSubiterationSatisfied(int i, int j);
+};
+
+#endif