[gegl/npd-squashed] libs, operations: Add n-point deformation



commit 7bd297888509d38889424393caaf0a31e2e2ca15
Author: Marek Dvoroznak <dvoromar gmail com>
Date:   Wed Jun 26 23:32:22 2013 +0200

    libs, operations: Add n-point deformation

 .gitignore                      |    1 +
 Makefile.am                     |    1 +
 configure.ac                    |    9 +
 gegl.pc.in                      |    2 +-
 libs/npd/.gitignore             |    6 +
 libs/npd/Makefile.am            |   40 +++
 libs/npd/deformation.c          |  207 +++++++++++++
 libs/npd/deformation.h          |   48 +++
 libs/npd/graphics.c             |  556 +++++++++++++++++++++++++++++++++++
 libs/npd/graphics.h             |  128 ++++++++
 libs/npd/lattice_cut.c          |  278 ++++++++++++++++++
 libs/npd/lattice_cut.h          |   40 +++
 libs/npd/npd.h                  |   32 ++
 libs/npd/npd_common.c           |  619 +++++++++++++++++++++++++++++++++++++++
 libs/npd/npd_common.h           |  157 ++++++++++
 libs/npd/npd_gegl.c             |   68 +++++
 libs/npd/npd_gegl.h             |   33 ++
 libs/npd/npd_math.c             |   47 +++
 libs/npd/npd_math.h             |   58 ++++
 operations/external/Makefile.am |    8 +
 operations/external/npd.c       |  329 +++++++++++++++++++++
 po/POTFILES.in                  |    1 +
 22 files changed, 2667 insertions(+), 1 deletions(-)
---
diff --git a/.gitignore b/.gitignore
index 69c9757..d411137 100644
--- a/.gitignore
+++ b/.gitignore
@@ -47,3 +47,4 @@
 cscope.files
 cscope.out
 /*.project
+/nbproject
diff --git a/Makefile.am b/Makefile.am
index 31f43be..588a324 100644
--- a/Makefile.am
+++ b/Makefile.am
@@ -7,6 +7,7 @@ SUBDIRS=\
        libs \
        opencl \
        gegl \
+       libs/npd \
        seamless-clone \
        operations \
        bin \
diff --git a/configure.ac b/configure.ac
index f567090..a6f780b 100644
--- a/configure.ac
+++ b/configure.ac
@@ -1074,6 +1074,14 @@ if test "x$have_p2tc" != "xyes"; then
   AC_SUBST(P2TC_LIBS, "-L\$(top_builddir)/libs/poly2tri-c/poly2tri-c -lpoly2tri-c")
 fi
 
+##################
+# Check for libnpd
+##################
+have_libnpd="yes (internal)"
+
+AC_SUBST(NPD_CFLAGS, "-I\$(top_srcdir)/libs")
+AC_SUBST(NPD_LIBS, "-L\$(top_builddir)/libs/npd -lgegl-npd-\$(GEGL_API_VERSION)")
+
 #######################
 # Check for other items
 #######################
@@ -1139,6 +1147,7 @@ gegl/property-types/Makefile
 gegl/opencl/Makefile
 libs/Makefile
 libs/rgbe/Makefile
+libs/npd/Makefile
 libs/poly2tri-c/Makefile
 libs/poly2tri-c/poly2tri-c/Makefile
 libs/poly2tri-c/poly2tri-c/render/Makefile
diff --git a/gegl.pc.in b/gegl.pc.in
index 614bb78..f4b70fc 100644
--- a/gegl.pc.in
+++ b/gegl.pc.in
@@ -8,5 +8,5 @@ Name: GEGL
 Description: Generic Graphics Library 
 Version: @GEGL_REAL_VERSION@
 Requires: @GLIB_PACKAGES@ babl
-Libs: -L${libdir} -l PACKAGE_NAME@- GEGL_API_VERSION@
+Libs: -L${libdir} -l PACKAGE_NAME@- GEGL_API_VERSION@ -l PACKAGE_NAME@-npd- GEGL_API_VERSION@
 Cflags: -I${includedir}/@PACKAGE_NAME -@GEGL_API_VERSION@
diff --git a/libs/npd/.gitignore b/libs/npd/.gitignore
new file mode 100644
index 0000000..4e0b01c
--- /dev/null
+++ b/libs/npd/.gitignore
@@ -0,0 +1,6 @@
+/Makefile
+/Makefile.in
+/.deps/
+/.libs/
+/*.lo
+/*.la
\ No newline at end of file
diff --git a/libs/npd/Makefile.am b/libs/npd/Makefile.am
new file mode 100644
index 0000000..080bb9d
--- /dev/null
+++ b/libs/npd/Makefile.am
@@ -0,0 +1,40 @@
+include $(top_srcdir)/operations/Makefile-common.am
+
+######################################################
+# A shared library for n-point image deformation API #
+######################################################
+
+GEGL_NPD_publicdir = $(includedir)/gegl-$(GEGL_API_VERSION)/npd
+
+GEGL_NPD_public_HEADERS = \
+       npd_common.h     \
+       deformation.h    \
+       npd_math.h       \
+       graphics.h       \
+       lattice_cut.h    \
+       npd_gegl.h       \
+       npd.h
+
+GEGL_NPD_SOURCES = \
+       $(GEGL_NPD_public_HEADERS)      \
+       npd_common.c                    \
+       deformation.c                   \
+       npd_math.c                      \
+       graphics.c                      \
+       lattice_cut.c                   \
+       npd_gegl.c
+
+libgegl_npd_ GEGL_API_VERSION@_la_SOURCES = \
+       $(GEGL_NPD_public_HEADERS)      \
+       $(GEGL_NPD_SOURCES)
+
+libgegl_npd_ GEGL_API_VERSION@_la_LDFLAGS = \
+       -avoid-version -export-dynamic $(no_undefined)
+
+libgegl_npd_ GEGL_API_VERSION@_la_CFLAGS = \
+  $(AM_CFLAGS)
+
+lib_LTLIBRARIES = libgegl-npd- GEGL_API_VERSION@.la
+
+libgegl_npd_ GEGL_API_VERSION@_la_LIBADD = \
+       $(libgegl)
diff --git a/libs/npd/deformation.c b/libs/npd/deformation.c
new file mode 100644
index 0000000..87f4362
--- /dev/null
+++ b/libs/npd/deformation.c
@@ -0,0 +1,207 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#include "deformation.h"
+#include <math.h>
+
+#define NPD_COMPUTE_CENTROID(suffix, arg1, arg2, WEIGHTS, accessor)     \
+void                                                                    \
+npd_compute_centroid_##suffix (gint      num_of_points,                 \
+                               arg1,                                    \
+                               arg2,                                    \
+                               NPDPoint *centroid)                      \
+{                                                                       \
+  gfloat x_sum = 0, y_sum = 0, weights_sum = 0;                         \
+  gint i;                                                               \
+                                                                        \
+  /* first compute sum of all values for each coordinate */             \
+  for (i = 0; i < num_of_points; ++i)                                   \
+    {                                                                   \
+      x_sum       += (WEIGHTS) * points[i] accessor x;                  \
+      y_sum       += (WEIGHTS) * points[i] accessor y;                  \
+      weights_sum += (WEIGHTS);                                         \
+    }                                                                   \
+                                                                        \
+  /* then compute mean */                                               \
+  centroid->x = x_sum / weights_sum;                                    \
+  centroid->y = y_sum / weights_sum;                                    \
+}
+
+NPD_COMPUTE_CENTROID (of_overlapping_points,
+                      NPDPoint *points[],
+                      gfloat weights[],
+                      1,
+                      ->)
+NPD_COMPUTE_CENTROID (from_weighted_points,
+                      NPDPoint points[],
+                      gfloat weights[],
+                      weights[i],
+                      .)
+#undef NPD_COMPUTE_CENTROID
+
+void
+npd_compute_ARSAP_transformation (gint     num_of_points,
+                                  NPDPoint reference_points[],
+                                  NPDPoint current_points[],
+                                  gfloat   weights[],
+                                  gboolean ASAP)
+{
+  NPDPoint pc = {0, 0}, qc = {0, 0};
+  gfloat a = 0, b = 0, mu_part = 0, mu, r1, r2, x0, y0;
+  gint i;
+
+  /* p - points of reference pose */
+  npd_compute_centroid_from_weighted_points (num_of_points,
+                                             reference_points,
+                                             weights,
+                                            &pc);
+  /* q - points of current pose */
+  npd_compute_centroid_from_weighted_points (num_of_points,
+                                             current_points,
+                                             weights,
+                                            &qc);
+
+  /* get rotation */
+  for (i = 0; i < num_of_points; ++i)
+    {
+      gfloat px_minus_pcx = reference_points[i].x - pc.x;
+      gfloat py_minus_pcy = reference_points[i].y - pc.y;
+      gfloat qx_minus_qcx =   current_points[i].x - qc.x;
+      gfloat qy_minus_qcy =   current_points[i].y - qc.y;
+
+      a += weights[i]
+              * ((px_minus_pcx)
+              *  (qx_minus_qcx)
+              +  (py_minus_pcy)
+              *  (qy_minus_qcy));
+      b += weights[i]
+              * ((px_minus_pcx)
+              *  (qy_minus_qcy)
+              -  (py_minus_pcy)
+              *  (qx_minus_qcx));
+
+      mu_part += weights[i]
+              * ((px_minus_pcx)
+              *  (px_minus_pcx)
+              +  (py_minus_pcy)
+              *  (py_minus_pcy));
+    }
+
+  mu = 1;
+  if (ASAP) mu = mu_part;
+  else      mu = sqrt(a * a + b * b);
+
+  r1 =  a / mu;
+  r2 = -b / mu;
+
+  /* get translation */
+  x0 = qc.x - ( r1 * pc.x + r2 * pc.y);
+  y0 = qc.y - (-r2 * pc.x + r1 * pc.y);
+
+  /* transform points */
+  for (i = 0; i < num_of_points; ++i)
+    {
+      if (!current_points[i].fixed)
+        {
+          current_points[i].x =  r1 * reference_points[i].x
+                  + r2 * reference_points[i].y + x0;
+          current_points[i].y = -r2 * reference_points[i].x
+                  + r1 * reference_points[i].y + y0;
+        }
+    }
+}
+
+void
+npd_compute_ARSAP_transformations (NPDHiddenModel *hidden_model)
+{
+  gint     i;
+  
+  for (i = 0; i < hidden_model->num_of_bones; ++i)
+    {
+      NPDBone *reference_bones = &hidden_model->reference_bones[i];
+      NPDBone *current_bones   = &hidden_model->current_bones[i];
+      npd_compute_ARSAP_transformation (reference_bones->num_of_points,
+                                        reference_bones->points,
+                                        current_bones->points,
+                                        current_bones->weights,
+                                        hidden_model->ASAP);
+    }
+}
+
+void
+npd_deform_model (NPDModel *model,
+                  gint      rigidity)
+{
+  gint i;
+  for (i = 0; i < rigidity; ++i)
+    {
+      npd_deform_model_once (model);
+    }
+}
+
+void
+npd_deform_model_once (NPDModel *model)
+{
+  gint i, j;
+  
+  /* updates associated overlapping points according to this control point */
+  for (i = 0; i < model->control_points->len; ++i)
+    {
+      NPDControlPoint *cp = &g_array_index (model->control_points,
+                                            NPDControlPoint,
+                                            i);
+
+      for (j = 0; j < cp->overlapping_points->num_of_points; ++j)
+        {
+          npd_set_point_coordinates (cp->overlapping_points->points[j],
+                                     &cp->point);
+        }
+    }
+
+  npd_deform_hidden_model_once (model->hidden_model);
+}
+
+void
+npd_deform_hidden_model_once (NPDHiddenModel *hidden_model)
+{
+  gint i, j;
+  npd_compute_ARSAP_transformations (hidden_model);
+
+  /* overlapping points are not overlapping after the deformation,
+     so we have to move them to their centroid */
+  for (i = 0; i < hidden_model->num_of_overlapping_points; ++i)
+    {
+      NPDOverlappingPoints *list_of_ops
+              = &hidden_model->list_of_overlapping_points[i];
+      NPDPoint centroid;
+      
+      npd_compute_centroid_of_overlapping_points (list_of_ops->num_of_points,
+                                                  list_of_ops->points,
+                                                  NULL,
+                                                 &centroid);
+
+      for (j = 0; j < list_of_ops->num_of_points; ++j)
+        {
+          list_of_ops->points[j]->x = centroid.x;
+          list_of_ops->points[j]->y = centroid.y;
+        }
+    }
+}
diff --git a/libs/npd/deformation.h b/libs/npd/deformation.h
new file mode 100644
index 0000000..c363e1d
--- /dev/null
+++ b/libs/npd/deformation.h
@@ -0,0 +1,48 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#ifndef __NPD_DEFORMATION_H__
+#define        __NPD_DEFORMATION_H__
+
+#include "npd_common.h"
+
+void  npd_compute_centroid_from_weighted_points
+                                        (gint              num_of_points,
+                                         NPDPoint          points[],
+                                         gfloat            weights[],
+                                         NPDPoint         *centroid);
+void  npd_compute_centroid_of_overlapping_points
+                                        (gint              num_of_points,
+                                         NPDPoint         *points[],
+                                         gfloat            weights[],
+                                         NPDPoint         *centroid);
+void  npd_compute_ARSAP_transformation  (gint              num_of_points,
+                                         NPDPoint          reference_shape[],
+                                         NPDPoint          current_shape[],
+                                         gfloat            weights[],
+                                         gboolean          ASAP);
+void  npd_compute_ARSAP_transformations (NPDHiddenModel   *model);
+void  npd_deform_model                  (NPDModel         *model,
+                                         gint              rigidity);
+void  npd_deform_model_once             (NPDModel         *model);
+void  npd_deform_hidden_model_once      (NPDHiddenModel   *hidden_model);
+
+#endif /* __NPD_DEFORMATION_H__ */
diff --git a/libs/npd/graphics.c b/libs/npd/graphics.c
new file mode 100644
index 0000000..e5c5e5f
--- /dev/null
+++ b/libs/npd/graphics.c
@@ -0,0 +1,556 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#include "graphics.h"
+#include "glib.h"
+#include <math.h>
+#include "npd_math.h"
+#include "lattice_cut.h"
+#include <stdio.h>
+#include <string.h>
+
+void
+npd_create_mesh_from_image (NPDModel *model,
+                            gint      width,
+                            gint      height,
+                            gint      position_x,
+                            gint      position_y)
+{
+  gint      square_size = model->mesh_square_size;
+  NPDImage *image = model->reference_image;
+  gint      i, cy, cx, y, x, r, c, ow, oh;
+  NPDColor  pixel_color = { 0, 0, 0, 0 };
+  GArray   *squares;
+  gint     *sq2id;
+  gboolean *empty_squares;
+  GList    *tmp_ops = NULL;
+  GArray   *ops = g_array_new (FALSE, FALSE, sizeof (NPDOverlappingPoints));
+  GList   **edges;
+  NPDHiddenModel *hm = model->hidden_model;
+
+  /* create squares above the image */
+  width  = ceil ((gfloat) width  / square_size);
+  height = ceil ((gfloat) height / square_size);
+
+  squares = g_array_new (FALSE, FALSE, sizeof (NPDBone));
+  empty_squares = g_new0 (gboolean, width * height);
+  sq2id = g_new0 (gint, width * height);
+
+  i = 0;
+  for (cy = 0; cy < height; cy++)
+  for (cx = 0; cx < width; cx++)
+    {
+      gboolean is_empty = TRUE;
+
+      for (y = cy * square_size; y < (cy + 1) * square_size; y++)
+      for (x = cx * square_size; x < (cx + 1) * square_size; x++)
+        {
+          /* test of emptiness */
+          npd_get_pixel_color (image, x, y, &pixel_color);
+          if (!npd_is_color_transparent (&pixel_color))
+            {
+              is_empty = FALSE;
+              goto not_empty;
+            }
+        }
+
+#define NPD_SQ_ID cy * width + cx
+
+      not_empty:
+      if (!is_empty)
+        {
+          /* create a square */
+          NPDBone square;
+          npd_create_square (&square,
+                             position_x + cx * square_size,
+                             position_y + cy * square_size,
+                             square_size, square_size);
+          g_array_append_val (squares, square);
+          sq2id[NPD_SQ_ID] = i;
+          i++;
+        }
+      else
+        empty_squares[NPD_SQ_ID] = TRUE;
+    }
+
+  /* find empty edges */
+  edges = npd_find_edges (model->reference_image, width, height, square_size);
+
+  /* create provisional overlapping points */
+#define NPD_ADD_P(op,r,c,point)                                              \
+  if ((r) > -1 && (r) < (oh - 1) && (c) > -1 && (c) < (ow - 1) &&            \
+      edges[op] == NULL && !empty_squares[(r) * width + (c)])                \
+    {                                                                        \
+      tmp_ops = g_list_append (tmp_ops, GINT_TO_POINTER ((r) * width + (c)));\
+      tmp_ops = g_list_append (tmp_ops, GINT_TO_POINTER (point));            \
+      num++;                                                                 \
+    }
+
+  ow = width + 1; oh = height + 1;
+
+  for (r = 0; r < oh; r++)
+  for (c = 0; c < ow; c++)
+    {
+      gint index = r * ow + c, num = 0;
+      NPD_ADD_P (index, r - 1, c - 1, 2);
+      NPD_ADD_P (index, r - 1, c,     3);
+      NPD_ADD_P (index, r,     c,     0);
+      NPD_ADD_P (index, r,     c - 1, 1);
+      if (num > 0)
+        tmp_ops = g_list_insert_before (tmp_ops,
+                                        g_list_nth_prev (g_list_last (tmp_ops), 2 * num - 1),
+                                        GINT_TO_POINTER (num));
+    }
+#undef NPD_ADD_P
+
+  /* cut lattice's edges and continue with creating of provisional overlapping points */
+  tmp_ops = g_list_concat (tmp_ops, npd_cut_edges (edges, ow, oh));
+  for (i = 0; i < ow * oh; i++) g_list_free (edges[i]);
+  g_free (edges);
+
+  /* create model's bones */
+  hm->num_of_bones = squares->len;
+  hm->current_bones = (NPDBone*) (gpointer) squares->data;
+  g_array_free (squares, FALSE);
+
+  /* create model's overlapping points */
+  while (g_list_next (tmp_ops))
+    {
+      GPtrArray *ppts = g_ptr_array_new ();
+      gint count = GPOINTER_TO_INT (tmp_ops->data);
+      gint j = 0;
+
+      for (i = 0; i < count; i++)
+        {
+          gint sq, p;
+          tmp_ops = g_list_next (tmp_ops);
+          sq = GPOINTER_TO_INT (tmp_ops->data);
+          tmp_ops = g_list_next (tmp_ops);
+          p  = GPOINTER_TO_INT (tmp_ops->data);
+
+          if (!empty_squares[sq])
+            {
+              g_ptr_array_add (ppts, &hm->current_bones[sq2id[sq]].points[p]);
+              j++;
+            }
+        }
+
+      if (j > 0)
+        {
+          NPDOverlappingPoints op;
+          op.num_of_points = j;
+          op.points = (NPDPoint**) ppts->pdata;
+          g_ptr_array_free (ppts, FALSE);
+          op.representative = op.points[0];
+          g_array_append_val (ops, op);
+        }
+      tmp_ops = g_list_next (tmp_ops);
+    }
+
+  g_list_free (tmp_ops);
+  g_free (empty_squares);
+  g_free (sq2id);
+
+  hm->num_of_overlapping_points = ops->len;
+  hm->list_of_overlapping_points = (NPDOverlappingPoints*) (gpointer) ops->data;
+  g_array_free (ops, FALSE);
+
+  /* create reference bones according to current bones */
+  hm->reference_bones = g_new (NPDBone, hm->num_of_bones);
+  for (i = 0; i < hm->num_of_bones; i++)
+    {
+      NPDBone *current_bone = &hm->current_bones[i];
+      NPDBone *reference_bone = &hm->reference_bones[i];
+      NPDPoint *current_points, *ref_points;
+      gint j, n = current_bone->num_of_points;
+
+      reference_bone->num_of_points = n;
+      reference_bone->points = g_new (NPDPoint, n);
+      memcpy (reference_bone->points, current_bone->points, n * sizeof (NPDPoint));
+      reference_bone->weights = current_bone->weights;
+
+      current_points = current_bone->points;
+      ref_points = reference_bone->points;
+
+      for (j = 0; j < n; j++)
+        {
+          current_points[j].current_bone = current_bone;
+          current_points[j].reference_bone = reference_bone;
+          ref_points[j].current_bone = current_bone;
+          ref_points[j].reference_bone = reference_bone;
+          ref_points[j].x = current_points[j].x - position_x;
+          ref_points[j].y = current_points[j].y - position_y;
+          current_points[j].counterpart = &ref_points[j];
+          ref_points[j].counterpart = &current_points[j];
+        }
+    }
+
+/*
+// could be useful later
+  gint j;
+  for (i = 0; i < hm->num_of_overlapping_points; i++)
+  for (j = 0; j < hm->list_of_overlapping_points[i].num_of_points; j++)
+    {
+      NPDPoint *p = hm->list_of_overlapping_points[i].points[j];
+      p->overlapping_points = &hm->list_of_overlapping_points[i];
+      p->counterpart->overlapping_points = &hm->list_of_overlapping_points[i];
+    }
+*/
+}
+
+void
+npd_draw_mesh (NPDModel   *model,
+               NPDDisplay *display)
+{
+  NPDHiddenModel *hm = model->hidden_model;
+  gint i, j;
+
+  for (i = 0; i < hm->num_of_bones; i++)
+    {
+      NPDBone  *bone  = &hm->current_bones[i];
+      NPDPoint *first = &bone->points[0];
+      NPDPoint *p0, *p1 = NULL;
+
+      for (j = 1; j < bone->num_of_points; j++)
+        {
+          p0 = &bone->points[j - 1];
+          p1 = &bone->points[j];
+          npd_draw_line (display, p0->x, p0->y, p1->x, p1->y);
+        }
+      npd_draw_line (display, p1->x, p1->y, first->x, first->y);
+    }
+}
+
+gfloat
+npd_bilinear_interpolation (gfloat I0,
+                            gfloat I1,
+                            gfloat I2,
+                            gfloat I3,
+                            gfloat dx,
+                            gfloat dy)
+{
+  return (I0 * (1 - dx) + I1 * dx) * (1 - dy) +
+         (I2 * (1 - dx) + I3 * dx) *      dy;
+}
+
+void
+npd_bilinear_color_interpolation (NPDColor *I0,
+                                  NPDColor *I1,
+                                  NPDColor *I2,
+                                  NPDColor *I3,
+                                  gfloat    dx,
+                                  gfloat    dy,
+                                  NPDColor *out)
+{
+  out->r = npd_bilinear_interpolation (I0->r, I1->r, I2->r, I3->r, dx, dy);
+  out->g = npd_bilinear_interpolation (I0->g, I1->g, I2->g, I3->g, dx, dy);
+  out->b = npd_bilinear_interpolation (I0->b, I1->b, I2->b, I3->b, dx, dy);
+  out->a = npd_bilinear_interpolation (I0->a, I1->a, I2->a, I3->a, dx, dy);
+}
+
+gfloat
+npd_blend_band (gfloat src,
+                gfloat dst,
+                gfloat src_alpha,
+                gfloat dst_alpha,
+                gfloat out_alpha_recip)
+{
+  return (src * src_alpha +
+          dst * dst_alpha * (1 - src_alpha)) * out_alpha_recip;
+}
+
+void
+npd_blend_colors (NPDColor *src,
+                  NPDColor *dst,
+                  NPDColor *out_color)
+{
+#ifdef NPD_RGBA_FLOAT
+  gfloat src_A = src->a,
+         dst_A = dst->a;
+#else
+  gfloat src_A = src->a / 255.0,
+         dst_A = dst->a / 255.0;
+#endif
+  gfloat out_alpha = src_A + dst_A * (1 - src_A);
+  gfloat out_alpha_recip = 1 / out_alpha;
+
+  out_color->r = npd_blend_band (src->r, dst->r, src_A, dst_A, out_alpha_recip);
+  out_color->g = npd_blend_band (src->g, dst->g, src_A, dst_A, out_alpha_recip);
+  out_color->b = npd_blend_band (src->b, dst->b, src_A, dst_A, out_alpha_recip);
+#ifdef NPD_RGBA_FLOAT
+  out_color->a = out_alpha;
+#else
+  out_color->a = out_alpha * 255;
+#endif
+}
+
+void
+npd_texture_fill_triangle (gint       x1,
+                           gint       y1,
+                           gint       x2,
+                           gint       y2,
+                           gint       x3,
+                           gint       y3,
+                           NPDMatrix *A,
+                           NPDImage  *input_image,
+                           NPDImage  *output_image,
+                           NPDSettings settings)
+{
+  gint yA, yB, yC, xA, xB, xC;
+  gint tmp, y;
+  gint deltaXAB, deltaYAB;
+  gint deltaXBC, deltaYBC;
+  gint deltaXAC, deltaYAC;
+
+  gfloat slopeAB;
+  gfloat slopeBC;
+  gfloat slopeAC;
+
+  gfloat k, l;
+  gfloat slope1, slope2, slope3, slope4;
+
+  if (y1 == y2 && x1 > x2)
+    {
+      tmp = y1; y1 = y2; y2 = tmp;
+      tmp = x1; x1 = x2; x2 = tmp;
+    }
+  if (y1 == y3 && x1 > x3)
+    {
+      tmp = y1; y1 = y3; y3 = tmp;
+      tmp = x1; x1 = x3; x3 = tmp;
+    }
+  if (y2 == y3 && x2 > x3)
+    {
+      tmp = y2; y2 = y3; y3 = tmp;
+      tmp = x2; x2 = x3; x3 = tmp;
+    }
+
+  if (y1 <= y2)
+    {
+      if (y2 <= y3)
+        {
+          /* y1 <= y2 <= y3 */
+          yA = y1; yB = y2; yC = y3;
+          xA = x1; xB = x2; xC = x3;
+        }
+      else if (y1 <= y3)
+        {
+          /* y1 <= y3 < y2 */
+          yA = y1; yB = y3; yC = y2;
+          xA = x1; xB = x3; xC = x2;
+        }
+      else
+        {
+          /* y3 < y1 < y2 */
+          yA = y3; yB = y1; yC = y2;
+          xA = x3; xB = x1; xC = x2;
+        }
+    }
+  else
+    {
+      if (y1 <= y3)
+        {
+          /* y2 < y1 <= y3 */
+          yA = y2; yB = y1; yC = y3;
+          xA = x2; xB = x1; xC = x3;
+        }
+      else if (y2 <= y3)
+        {
+          /* y2 <= y3 < y1 */
+          yA = y2; yB = y3; yC = y1;
+          xA = x2; xB = x3; xC = x1;
+        }
+      else
+        {
+          /* y3 < y2 < y1 */
+          yA = y3; yB = y2; yC = y1;
+          xA = x3; xB = x2; xC = x1;
+        }
+    }
+
+  deltaXAB = xB - xA, deltaYAB = yB - yA;
+  deltaXBC = xC - xB, deltaYBC = yC - yB;
+  deltaXAC = xC - xA, deltaYAC = yC - yA;
+
+  slopeBC = (deltaYBC == 0 ? 0 : (gfloat) deltaXBC / deltaYBC);
+  slopeAC = (deltaYAC == 0 ? 0 : (gfloat) deltaXAC / deltaYAC);
+
+  if (deltaYAB == 0)
+    {
+      slopeAB = 0;
+      k = xA;
+      l = xB;
+
+      slope1 = slopeAB;
+      slope2 = slopeAC;
+      slope3 = slopeAC;
+      slope4 = slopeBC;
+    }
+  else
+    {
+      slopeAB = (gfloat) deltaXAB / deltaYAB;
+      k = xA;
+      l = xA;
+
+      if (slopeAB > slopeAC)
+        {
+          slope1 = slopeAC;
+          slope2 = slopeAB;
+          slope3 = slopeAC;
+          slope4 = slopeBC;
+        }
+      else
+        {
+          slope1 = slopeAB;
+          slope2 = slopeAC;
+          slope3 = slopeBC;
+          slope4 = slopeAC;
+        }
+    }
+
+  for (y = yA; y < yB; y++)
+    {
+      npd_draw_texture_line ((gint) round (k), (gint) round (l),
+                             y, A,
+                             input_image, output_image,
+                             settings);
+      k += slope1;
+      l += slope2;
+    }
+
+  for (y = yB; y <= yC; y++)
+    {
+      npd_draw_texture_line ((gint) round (k), (gint) round (l),
+                             y, A,
+                             input_image, output_image,
+                             settings);
+      k += slope3;
+      l += slope4;
+    }
+}
+
+void
+npd_texture_quadrilateral (NPDBone    *reference_bone,
+                           NPDBone    *current_bone,
+                           NPDImage   *input_image,
+                           NPDImage   *output_image,
+                           NPDSettings settings)
+{
+  /* p1 are points of domain, p2 are points of codomain */
+  NPDPoint *p1 = current_bone->points;
+  NPDPoint *p2 = reference_bone->points;
+
+  NPDMatrix *A = NULL;
+  npd_new_matrix (&A);
+
+  npd_compute_affinity (&p1[0], &p1[1], &p1[2],
+                        &p2[0], &p2[1], &p2[2], A);
+  npd_texture_fill_triangle ((gint) p1[0].x, (gint) p1[0].y,
+                             (gint) p1[1].x, (gint) p1[1].y,
+                             (gint) p1[2].x, (gint) p1[2].y,
+                             A, input_image, output_image,
+                             settings);
+
+  npd_compute_affinity (&p1[0], &p1[2], &p1[3],
+                        &p2[0], &p2[2], &p2[3], A);
+  npd_texture_fill_triangle ((gint) p1[0].x, (gint) p1[0].y,
+                             (gint) p1[2].x, (gint) p1[2].y,
+                             (gint) p1[3].x, (gint) p1[3].y,
+                             A, input_image, output_image,
+                             settings);
+
+  npd_destroy_matrix (&A);
+}
+
+void
+npd_draw_texture_line (gint        x1,
+                       gint        x2,
+                       gint        y,
+                       NPDMatrix  *A,
+                       NPDImage   *input_image,
+                       NPDImage   *output_image,
+                       NPDSettings settings)
+{
+  gint x, fx, fy;
+  gfloat dx, dy;
+  
+  for (x = x1; x <= x2; x++)
+    {
+      NPDPoint p, q;
+      NPDColor I0, interpolated, *final;
+
+      q.x = x; q.y = y;
+      npd_apply_transformation (A, &q, &p);
+
+      fx = floor (p.x);
+      fy = floor (p.y);
+
+      npd_get_pixel_color (input_image, fx, fy, &I0);
+      final = &I0;
+
+      /* bilinear interpolation */
+      if (settings & NPD_BILINEAR_INTERPOLATION)
+        {
+          NPDColor I1, I2, I3;
+
+          dx =  p.x - fx;
+          dy =  p.y - fy;
+
+          npd_get_pixel_color (input_image, fx + 1, fy,     &I1);
+          npd_get_pixel_color (input_image, fx,     fy + 1, &I2);
+          npd_get_pixel_color (input_image, fx + 1, fy + 1, &I3);
+          npd_bilinear_color_interpolation (&I0, &I1, &I2, &I3, dx, dy, &interpolated);
+          final = &interpolated;
+        }
+
+      /* alpha blending */
+      if (settings & NPD_ALPHA_BLENDING)
+        {
+          NPDColor dest;
+          npd_get_pixel_color (output_image, x, y, &dest);
+          npd_blend_colors (final, &dest, final);
+        }
+
+      npd_set_pixel_color (output_image, x, y, final);
+    }
+}
+
+gboolean
+npd_compare_colors (NPDColor *c1,
+                    NPDColor *c2)
+{
+  if (npd_equal_floats (c1->r, c2->r) &&
+      npd_equal_floats (c1->g, c2->g) &&
+      npd_equal_floats (c1->b, c2->b) &&
+      npd_equal_floats (c1->a, c2->a))
+    return TRUE;
+
+  return FALSE;
+}
+
+gboolean
+npd_is_color_transparent (NPDColor *color)
+{
+  if (npd_equal_floats (color->a, 0.0))
+    return TRUE;
+
+  return FALSE;
+}
diff --git a/libs/npd/graphics.h b/libs/npd/graphics.h
new file mode 100644
index 0000000..db174fc
--- /dev/null
+++ b/libs/npd/graphics.h
@@ -0,0 +1,128 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#ifndef __NPD_GRAPHICS_H__
+#define        __NPD_GRAPHICS_H__
+
+#include "npd_common.h"
+
+//#define NPD_RGBA_FLOAT
+
+struct _NPDColor {
+#ifdef NPD_RGBA_FLOAT
+  gfloat r;
+  gfloat g;
+  gfloat b;
+  gfloat a;
+#else
+  guint8 r;
+  guint8 g;
+  guint8 b;
+  guint8 a;
+#endif
+};
+
+typedef enum
+{
+  NPD_BILINEAR_INTERPOLATION = 1,
+  NPD_ALPHA_BLENDING         = 1 << 1
+} NPDSettings;
+
+void        npd_create_model_from_image       (NPDModel   *model,
+                                               NPDImage   *image,
+                                               gint        square_size);
+void        npd_create_mesh_from_image        (NPDModel   *model,
+                                               gint        width,
+                                               gint        height,
+                                               gint        position_x,
+                                               gint        position_y);
+void        npd_draw_model                    (NPDModel   *model,
+                                               NPDDisplay *display);
+void        npd_draw_mesh                     (NPDModel   *model,
+                                               NPDDisplay *display);
+
+gboolean    npd_load_image                    (NPDImage   *image,
+                                               const char *path);
+void        npd_destroy_image                 (NPDImage   *image);
+//void        npd_draw_image                    (NPDImage *image);
+void        npd_texture_fill_triangle         (gint        x1,
+                                               gint        y1,
+                                               gint        x2,
+                                               gint        y2,
+                                               gint        x3,
+                                               gint        y3,
+                                               NPDMatrix  *A,
+                                               NPDImage   *input_image,
+                                               NPDImage   *output_image,
+                                               NPDSettings settings);
+void        npd_texture_quadrilateral         (NPDBone    *reference_bone,
+                                               NPDBone    *current_bone,
+                                               NPDImage   *input_image,
+                                               NPDImage   *output_image,
+                                               NPDSettings settings);
+void        npd_draw_texture_line             (gint        x1,
+                                               gint        x2,
+                                               gint        y,
+                                               NPDMatrix  *A,
+                                               NPDImage   *input_image,
+                                               NPDImage   *output_image,
+                                               NPDSettings settings);
+void      (*npd_draw_line)                    (NPDDisplay *display,
+                                               gfloat      x0,
+                                               gfloat      y0,
+                                               gfloat      x1,
+                                               gfloat      y1);
+gfloat      npd_bilinear_interpolation        (gfloat      I0,
+                                               gfloat      I1,
+                                               gfloat      I2,
+                                               gfloat      I3,
+                                               gfloat      dx,
+                                               gfloat      dy);
+void        npd_bilinear_color_interpolation  (NPDColor   *I0,
+                                               NPDColor   *I1,
+                                               NPDColor   *I2,
+                                               NPDColor   *I3,
+                                               gfloat      dx,
+                                               gfloat      dy,
+                                               NPDColor   *out);
+gfloat      npd_blend_band                    (gfloat      src,
+                                               gfloat      dst,
+                                               gfloat      src_alpha,
+                                               gfloat      dest_alpha,
+                                               gfloat      out_alpha);
+void        npd_blend_colors                 (NPDColor   *src,
+                                               NPDColor   *dst,
+                                               NPDColor   *out_color);
+void      (*npd_get_pixel_color)              (NPDImage   *image,
+                                               gint        x,
+                                               gint        y,
+                                               NPDColor   *color);
+void      (*npd_set_pixel_color)              (NPDImage   *image,
+                                               gint        x,
+                                               gint        y,
+                                               NPDColor   *color);
+gboolean    npd_compare_colors                (NPDColor   *c1,
+                                               NPDColor   *c2);
+gboolean    npd_is_color_transparent          (NPDColor   *color);
+gboolean    npd_init_display                  (NPDDisplay *display);
+void        npd_destroy_display               (NPDDisplay *display);
+
+#endif /*__NPD_GRAPHICS_H__ */
diff --git a/libs/npd/lattice_cut.c b/libs/npd/lattice_cut.c
new file mode 100644
index 0000000..d8e7f72
--- /dev/null
+++ b/libs/npd/lattice_cut.c
@@ -0,0 +1,278 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#include "lattice_cut.h"
+#include "npd_common.h"
+#include "graphics.h"
+#include <glib.h>
+#include "npd_math.h"
+
+#define NPD_SWAP_INTS(i,j) { gint tmp = i; i = j; j = tmp; }
+
+/* only works for straight lines */
+gboolean
+npd_is_edge_empty (NPDImage *image,
+                   gint      X1,
+                   gint      Y1,
+                   gint      X2,
+                   gint      Y2)
+{
+  gint x, y;
+  NPDColor color;
+
+  if (Y1 > Y2) NPD_SWAP_INTS (Y1, Y2);
+  if (X1 > X2) NPD_SWAP_INTS (X1, X2);
+
+  for (y = Y1; y <= Y2; y++)
+  for (x = X1; x <= X2; x++)
+    {
+      npd_get_pixel_color (image, x, y, &color);
+      if (!npd_is_color_transparent (&color))
+        return FALSE;
+    }
+
+  return TRUE;
+}
+
+GList**
+npd_find_edges (NPDImage *image,
+                gint      count_x,
+                gint      count_y,
+                gint      square_size)
+{
+  gint i, j;
+  gint ow = count_x + 1,
+       oh = count_y + 1;
+  GList **empty_edges = g_new0 (GList*, ow * oh);
+
+  for (j = 1; j <= count_y; j++)
+  for (i = 1; i <= count_x; i++)
+    {
+#define NPD_TEST_EMPTY(from_op_x,from_op_y,to_op_x,to_op_y)                         \
+      if (npd_is_edge_empty (image,                                                 \
+                             (from_op_x) * square_size, (from_op_y) * square_size,  \
+                             (to_op_x)   * square_size, (to_op_y)   * square_size)) \
+        {                                                                           \
+          gint from_op_id = (from_op_y) * ow + (from_op_x),                         \
+               to_op_id   = (to_op_y)   * ow + (to_op_x);                           \
+          empty_edges[from_op_id] = g_list_append (empty_edges[from_op_id],         \
+                                                   GINT_TO_POINTER (to_op_id));     \
+          empty_edges[to_op_id]   = g_list_append (empty_edges[to_op_id],           \
+                                                   GINT_TO_POINTER (from_op_id));   \
+        }
+
+      if (j != count_y) NPD_TEST_EMPTY (i, j, i - 1, j);
+      if (i != count_x) NPD_TEST_EMPTY (i, j, i,     j - 1);
+#undef NPD_TEST_EMPTY
+    }
+
+  return empty_edges;
+}
+
+GList*
+npd_cut_edges (GList **edges,
+               gint    ow,
+               gint    oh)
+{
+  gint i, r, col, width = ow - 1;
+  GList *ops = NULL;
+  gint neighbors[4];
+
+  for (r = 0; r < oh; r++)
+  for (col = 0; col < ow; col++)
+    {
+      gint index = r * ow + col;
+      GList *op = edges[index];
+      gint num_of_neighbors = g_list_length (op);
+
+      if (num_of_neighbors == 0) continue;
+
+#define NPD_ADD_COUNT(count) ops = g_list_append (ops, GINT_TO_POINTER (count))
+#define NPD_ADD_P(r,col,point)                                                 \
+      if ((r) > -1 && (r) < (oh - 1) && (col) > -1 && (col) < (ow - 1))        \
+        {                                                                      \
+          ops = g_list_append (ops, GINT_TO_POINTER ((r) * width + (col)));    \
+          ops = g_list_append (ops, GINT_TO_POINTER (point));                  \
+        }
+
+      for (i = 0; i < num_of_neighbors; i++)
+        neighbors[i] = GPOINTER_TO_INT (g_list_nth_data (op, i));
+
+      if (num_of_neighbors == 1)
+        {
+          gboolean border = FALSE;
+          if (r == 0 || col == 0 || r == (oh - 1) || col == (ow - 1))
+            border = TRUE;
+
+          if (border) NPD_ADD_COUNT (1);
+          else        NPD_ADD_COUNT (4);
+          NPD_ADD_P (r - 1, col - 1, 2);
+          NPD_ADD_P (r - 1, col,     3);
+          NPD_ADD_P (r,     col,     0);
+          NPD_ADD_P (r,     col - 1, 1);
+          if (border)
+            ops = g_list_insert_before (ops,
+                                        g_list_nth_prev (g_list_last (ops), 1),
+                                        GINT_TO_POINTER (1));
+#undef NPD_ADD_P
+        }
+      else
+      if (num_of_neighbors == 2)
+        {
+          gboolean x_differs = FALSE, y_differs = FALSE;
+          gint a, b, c, d;
+
+#define NPD_OP_X(op) ((op) % ow)
+#define NPD_OP_Y(op) ((op) / ow)
+
+          if (NPD_OP_X (neighbors[0]) != NPD_OP_X (neighbors[1])) x_differs = TRUE;
+          if (NPD_OP_Y (neighbors[0]) != NPD_OP_Y (neighbors[1])) y_differs = TRUE;
+
+          if (x_differs && y_differs)
+            {
+              /* corner */
+              gint B = neighbors[0], C = neighbors[1];
+
+              if (NPD_OP_X (index) == NPD_OP_X (neighbors[0]))
+                { B = neighbors[1]; C = neighbors[0]; }
+
+              if (NPD_OP_Y (index) < NPD_OP_Y (C))
+                {
+                  if (NPD_OP_X (index) < NPD_OP_X (B))
+                    { /* IV. quadrant */  a = 2; b = 3; c = 1; d = 0; }
+                  else
+                    { /* III. quadrant */ a = 2; b = 3; c = 0; d = 1; }
+                }
+              else
+                {
+                  if (NPD_OP_X (index) < NPD_OP_X (B))
+                    { /* I. quadrant */   a = 2; b = 0; c = 1; d = 3; }
+                  else
+                    { /* II. quadrant */  a = 0; b = 3; c = 1; d = 2; }
+                }
+
+#define NPD_OP2SQ(op) (op == 0 ? ( r      * width + col) :                     \
+                      (op == 1 ? ( r      * width + col - 1) :                 \
+                      (op == 2 ? ((r - 1) * width + col - 1) :                 \
+                                 ((r - 1) * width + col))))
+#define NPD_ADD_P(square,point)                                                \
+              ops = g_list_append (ops, GINT_TO_POINTER (square));             \
+              ops = g_list_append (ops, GINT_TO_POINTER (point));
+
+              NPD_ADD_COUNT (3);
+              NPD_ADD_P (NPD_OP2SQ (a), a);
+              NPD_ADD_P (NPD_OP2SQ (b), b);
+              NPD_ADD_P (NPD_OP2SQ (c), c);
+              NPD_ADD_COUNT (1);
+              NPD_ADD_P (NPD_OP2SQ (d), d);
+            }
+          else
+            {
+              /* segment */
+              a = 0; b = 1; c = 2; d = 3;
+              if (y_differs) { a = 1; b = 2; c = 0; d = 3; }
+
+              NPD_ADD_COUNT (2);
+              NPD_ADD_P (NPD_OP2SQ (a), a);
+              NPD_ADD_P (NPD_OP2SQ (b), b);
+              NPD_ADD_COUNT (2);
+              NPD_ADD_P (NPD_OP2SQ (c), c);
+              NPD_ADD_P (NPD_OP2SQ (d), d);
+            }
+        }
+      else
+      if (num_of_neighbors == 3)
+        {
+          gint B = neighbors[0], C = neighbors[1], D = neighbors[2];
+          gint a = 2, b = 1, c = 3, d = 0;
+
+          if ((NPD_OP_X (B) != NPD_OP_X (C)) && (NPD_OP_Y (B) != NPD_OP_Y (C)))
+            {
+              /* B and C form corner */
+              B = neighbors[2]; D = neighbors[0]; /* swap B and D */
+
+              if ((NPD_OP_X (B) != NPD_OP_X (C)) && (NPD_OP_Y (B) != NPD_OP_Y (C)))
+                {
+                  /* (new) B and C form corner */
+                  C = neighbors[0]; D = neighbors[1]; /* swap C and D */
+                }
+            }
+
+          /* B and C form segment */
+          if (NPD_OP_X (B) == NPD_OP_X (C))
+            {
+              if (NPD_OP_X (B) < NPD_OP_X (D))
+                {
+                  /* |_
+                     |  */
+                  a = 2; b = 1; c = 3; d = 0;
+                }
+              else
+                {
+                  /* _|
+                      | */
+                  a = 3; b = 0; c = 2; d = 1;
+                }
+            }
+          else if (NPD_OP_Y (B) == NPD_OP_Y (C))
+            {
+              if (NPD_OP_Y (B) < NPD_OP_Y (D))
+                {
+                  /* _ _
+                      |  */
+                  a = 2; b = 3; c = 1; d = 0;
+                }
+              else
+                {
+                  /* _|_ */
+                  a = 1; b = 0; c = 2; d = 3;
+                }
+            }
+
+          NPD_ADD_COUNT (2);
+          NPD_ADD_P (NPD_OP2SQ (a), a);
+          NPD_ADD_P (NPD_OP2SQ (b), b);
+          NPD_ADD_COUNT (1);
+          NPD_ADD_P (NPD_OP2SQ (c), c);
+          NPD_ADD_COUNT (1);
+          NPD_ADD_P (NPD_OP2SQ (d), d);
+        }
+      else
+      if (num_of_neighbors == 4)
+        {
+          NPD_ADD_COUNT (1);
+          NPD_ADD_P (NPD_OP2SQ (0), 0);
+          NPD_ADD_COUNT (1);
+          NPD_ADD_P (NPD_OP2SQ (1), 1);
+          NPD_ADD_COUNT (1);
+          NPD_ADD_P (NPD_OP2SQ (2), 2);
+          NPD_ADD_COUNT (1);
+          NPD_ADD_P (NPD_OP2SQ (3), 3);
+        }
+    }
+#undef NPD_ADD_P
+#undef NPD_OP2SQ
+#undef NPD_OP_X
+#undef NPD_OP_Y
+#undef NPD_ADD_COUNT
+
+  return ops;
+}
diff --git a/libs/npd/lattice_cut.h b/libs/npd/lattice_cut.h
new file mode 100644
index 0000000..b331ba3
--- /dev/null
+++ b/libs/npd/lattice_cut.h
@@ -0,0 +1,40 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#ifndef __REFINE_H__
+#define        __REFINE_H__
+
+#include "npd_common.h"
+
+gboolean     npd_is_edge_empty    (NPDImage *image,
+                                   gint      X1,
+                                   gint      Y1,
+                                   gint      X2,
+                                   gint      Y2);
+GList**      npd_find_edges       (NPDImage *image,
+                                   gint      count_x,
+                                   gint      count_y,
+                                   gint      square_size);
+GList*       npd_cut_edges        (GList   **edges,
+                                   gint      ow,
+                                   gint      oh);
+
+#endif /* __REFINE_H__ */
diff --git a/libs/npd/npd.h b/libs/npd/npd.h
new file mode 100644
index 0000000..654c1ce
--- /dev/null
+++ b/libs/npd/npd.h
@@ -0,0 +1,32 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#ifndef __NPD_H__
+#define        __NPD_H__
+
+#include "npd_common.h"
+#include "graphics.h"
+#include "deformation.h"
+#include "npd_math.h"
+#include "lattice_cut.h"
+
+#endif /* __NPD_H__ */
+
diff --git a/libs/npd/npd_common.c b/libs/npd/npd_common.c
new file mode 100644
index 0000000..1b427b9
--- /dev/null
+++ b/libs/npd/npd_common.c
@@ -0,0 +1,619 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#include "npd_common.h"
+#include "npd_math.h"
+#include <math.h>
+#include <glib.h>
+#include <glib/gprintf.h>
+
+gint  npd_int_sort_function_descending (gconstpointer a,
+                                        gconstpointer b);
+
+void
+npd_init_model (NPDModel *model)
+{
+  NPDHiddenModel *hidden_model;
+  GArray         *control_points;
+  
+  /* init hidden model */
+  hidden_model        = g_new (NPDHiddenModel, 1);
+  model->hidden_model = hidden_model;
+  hidden_model->ASAP                      = FALSE;
+  hidden_model->MLS_weights               = FALSE;
+  hidden_model->MLS_weights_alpha         = 1;
+  hidden_model->num_of_bones              = 0;
+  hidden_model->num_of_overlapping_points = 0;
+
+  /* init control points */
+  control_points        = g_array_new (FALSE, FALSE, sizeof (NPDControlPoint));
+  model->control_points = control_points;
+  model->control_point_radius             = 6;
+  model->control_points_visible           = TRUE;
+
+  /* init other */
+  model->mesh_visible                     = TRUE;
+  model->texture_visible                  = TRUE;
+}
+
+void
+npd_destroy_hidden_model (NPDHiddenModel *hm)
+{
+  gint i;
+  for (i = 0; i < hm->num_of_overlapping_points; i++)
+    {
+      g_free (hm->list_of_overlapping_points[i].points);
+    }
+
+  g_free (hm->list_of_overlapping_points);
+
+  for (i = 0; i < hm->num_of_bones; i++)
+    {
+      g_free (hm->current_bones[i].weights);
+      g_free (hm->current_bones[i].points);
+      g_free (hm->reference_bones[i].points);
+    }
+  
+  g_free (hm->current_bones);
+  g_free (hm->reference_bones);
+}
+
+void
+npd_destroy_model (NPDModel *model)
+{
+  /* destroy control points */
+  g_array_free (model->control_points, TRUE);
+
+  /* destroy hidden model */
+  npd_destroy_hidden_model (model->hidden_model);
+  g_free (model->hidden_model);
+}
+
+/**
+ * Finds nearest (to specified position) overlapping points, creates a new
+ * control point at the position of overlapping points and assigns them to the
+ * control point.
+ * 
+ * @param  model
+ * @param  coord     specified position
+ * @return pointer to a newly created control point or NULL when there already
+ * is a control point at the position of nearest overlapping points
+ */
+NPDControlPoint*
+npd_add_control_point (NPDModel *model,
+                       NPDPoint *coord)
+{
+  gint                  num_of_ops, i, closest;
+  gfloat                min, current;
+  NPDOverlappingPoints *list_of_ops;
+  NPDControlPoint       cp;
+  NPDPoint             *closest_point;
+
+  list_of_ops = model->hidden_model->list_of_overlapping_points;
+  num_of_ops  = model->hidden_model->num_of_overlapping_points;
+
+  /* find closest overlapping points */
+  closest = 0;
+  min     = npd_SED (coord, list_of_ops[0].representative);
+
+  for (i = 1; i < num_of_ops; i++)
+    {
+      NPDPoint *op = list_of_ops[i].representative;
+      current      = npd_SED(coord, op);
+
+      if (min > current)
+        {
+          closest = i;
+          min     = current;
+        }
+    }
+
+  closest_point = list_of_ops[closest].representative;
+
+  /* we want to create a new control point only when there isn't any
+   * control point associated to the closest overlapping points - i.e. we
+   * don't want to have two (or more) different control points manipulating
+   * one overlapping points */
+  if (!npd_get_control_point_at (model, closest_point))
+    {
+      cp.point.weight       = closest_point->weight;
+      cp.overlapping_points = &list_of_ops[closest];
+
+      npd_set_point_coordinates (&cp.point, closest_point);
+      g_array_append_val (model->control_points, cp);
+
+      if (model->hidden_model->MLS_weights)
+        npd_compute_MLS_weights (model);
+
+      return &g_array_index (model->control_points,
+                             NPDControlPoint,
+                             model->control_points->len - 1);
+    }
+  else
+    return NULL;
+}
+
+/**
+ * Beware, when you use this function on previously stored pointers to control
+ * points it needn't to work properly, because g_array_remove_index can
+ * preserve pointers but change (move) data.
+ * In this situation use npd_remove_control_points instead.
+ */
+void
+npd_remove_control_point (NPDModel        *model,
+                          NPDControlPoint *control_point)
+{
+  gint i;
+  NPDControlPoint *cp;
+
+  for (i = 0; i < model->control_points->len; i++)
+    {
+      cp = &g_array_index (model->control_points, NPDControlPoint, i);
+
+      if (cp == control_point)
+        {
+          npd_set_control_point_weight (cp, 1.0);
+          g_array_remove_index (model->control_points, i);
+
+          if (model->hidden_model->MLS_weights)
+            npd_compute_MLS_weights (model);
+
+          return;
+        }
+    }
+}
+
+gint
+npd_int_sort_function_descending (gconstpointer a,
+                                  gconstpointer b)
+{
+  return GPOINTER_TO_INT (b) - GPOINTER_TO_INT (a);
+}
+
+void
+npd_remove_control_points (NPDModel *model,
+                           GList    *control_points)
+{
+  gint i;
+  NPDControlPoint *cp;
+  GList *indices = NULL;
+
+  /* first we find indices of control points we want to remove */
+  while (control_points != NULL)
+    {
+      for (i = 0; i < model->control_points->len; i++)
+        {
+          cp = &g_array_index (model->control_points, NPDControlPoint, i);
+          if (cp == control_points->data)
+            {
+              npd_set_control_point_weight (cp, 1.0);
+              indices = g_list_insert_sorted (indices,
+                                              GINT_TO_POINTER (i),
+                                              npd_int_sort_function_descending);
+            }
+        }
+
+      control_points = g_list_next (control_points);
+    }
+
+  /* indices are sorted in descending order, so we can simply iterate over them
+   * and remove corresponding control points one by one */
+  while (indices != NULL)
+    {
+      g_array_remove_index (model->control_points, GPOINTER_TO_INT (indices->data));
+      indices = g_list_next (indices);
+    }
+
+  if (model->hidden_model->MLS_weights)
+    npd_compute_MLS_weights (model);
+
+  g_list_free (indices);
+}
+
+void
+npd_remove_all_control_points (NPDModel *model)
+{
+  g_array_remove_range (model->control_points,
+                        0,
+                        model->control_points->len);
+}
+
+void
+npd_set_control_point_weight (NPDControlPoint *cp,
+                              gfloat           weight)
+{
+  npd_set_overlapping_points_weight(cp->overlapping_points, weight);
+}
+
+gboolean
+npd_equal_coordinates (NPDPoint *p1,
+                       NPDPoint *p2)
+{
+  return npd_equal_coordinates_epsilon(p1, p2, NPD_EPSILON);
+}
+
+gboolean
+npd_equal_coordinates_epsilon (NPDPoint *p1,
+                               NPDPoint *p2,
+                               gfloat    epsilon)
+{
+  if (npd_equal_floats_epsilon (p1->x, p2->x, epsilon) &&
+      npd_equal_floats_epsilon (p1->y, p2->y, epsilon))
+    {
+      return TRUE;
+    }
+  
+  return FALSE;
+}
+
+NPDControlPoint*
+npd_get_control_point_with_radius_at (NPDModel        *model,
+                                      NPDPoint        *coord,
+                                      gfloat           radius)
+{
+  gint i;
+  for (i = 0; i < model->control_points->len; i++)
+    {
+      NPDControlPoint *cp = &g_array_index (model->control_points,
+                                            NPDControlPoint,
+                                            i);
+      if (npd_equal_coordinates_epsilon (&cp->point,
+                                          coord,
+                                          radius))
+        {
+          return cp;
+        }
+    }
+
+  return NULL;
+}
+
+NPDControlPoint*
+npd_get_control_point_at (NPDModel *model,
+                          NPDPoint *coord)
+{
+  return npd_get_control_point_with_radius_at (model,
+                                               coord,
+                                               model->control_point_radius);
+}
+
+void
+npd_create_square (NPDBone *square,
+                   gint     x,
+                   gint     y,
+                   gint     width,
+                   gint     height)
+{
+  gint i;
+  square->num_of_points = 4;
+  square->points  = g_new (NPDPoint, 4);
+  square->weights = g_new (gfloat,   4);
+
+  square->points[0].x = x;         square->points[0].y = y;
+  square->points[1].x = x + width; square->points[1].y = y;
+  square->points[2].x = x + width; square->points[2].y = y + height;
+  square->points[3].x = x;         square->points[3].y = y + height;
+
+  for (i = 0; i < 4; i++)
+    {
+      square->weights[i] = 1.0;
+      square->points[i].weight = &square->weights[i];
+      square->points[i].fixed = FALSE;
+      square->points[i].index = i;
+    }
+}
+
+void
+npd_create_list_of_overlapping_points (NPDHiddenModel *hm)
+{
+  gint        i, j, num_of_bones;
+  NPDBone    *bone;
+  NPDPoint   *point;
+  GPtrArray  *list_of_ops;
+  GHashTable *coords_to_cluster;
+  
+  list_of_ops       = g_ptr_array_new ();
+  num_of_bones      = hm->num_of_bones;
+  coords_to_cluster = g_hash_table_new_full
+                          (g_str_hash, g_str_equal,
+                           g_free,     (GDestroyNotify) g_hash_table_destroy);
+
+  for (i = 0; i < num_of_bones; i++)
+    {
+      bone = &hm->current_bones[i];
+
+      for (j = 0; j < bone->num_of_points; j++)
+        {
+          point =  &bone->points[j];
+          add_point_to_suitable_cluster (coords_to_cluster,
+                                         point,
+                                         list_of_ops);
+        }
+    }
+
+  hm->list_of_overlapping_points = g_new (NPDOverlappingPoints,
+                                          list_of_ops->len);
+  hm->num_of_overlapping_points  = list_of_ops->len;
+
+  for (i = 0; i < list_of_ops->len; i++)
+    {
+      GPtrArray *op = g_ptr_array_index (list_of_ops, i);
+      hm->list_of_overlapping_points[i].points = (NPDPoint**) op->pdata;
+      hm->list_of_overlapping_points[i].num_of_points = op->len;
+      hm->list_of_overlapping_points[i].representative =
+              hm->list_of_overlapping_points[i].points[0];
+      
+      for (j = 0; j < op->len; j++)
+        {
+          NPDPoint *p = hm->list_of_overlapping_points[i].points[j];
+          p->overlapping_points = &hm->list_of_overlapping_points[i];
+          p->counterpart->overlapping_points = &hm->list_of_overlapping_points[i];
+        }
+
+      g_ptr_array_free (op, FALSE); /* we want to preserve the underlying
+                                       array */
+    }
+
+  /* free allocated memory */
+  g_hash_table_destroy (coords_to_cluster);
+  g_ptr_array_free (list_of_ops, TRUE);
+}
+
+#define NPD_FLOAT_TO_STRING(name_of_string, value)                             \
+/* must be freed */                                                            \
+name_of_string = g_new (gchar, 10);                                            \
+g_ascii_dtostr (name_of_string, 10, value);
+
+void
+add_point_to_suitable_cluster (GHashTable *coords_to_cluster,
+                               NPDPoint   *point,
+                               GPtrArray  *list_of_overlapping_points)
+{
+  gchar      *str_coord_x, *str_coord_y;
+  GHashTable *coord_y;
+  GPtrArray  *op;
+
+  NPD_FLOAT_TO_STRING (str_coord_x, point->x);
+  NPD_FLOAT_TO_STRING (str_coord_y, point->y);
+  
+  coord_y = g_hash_table_lookup (coords_to_cluster, str_coord_x);
+
+  if (coord_y == NULL)
+    {
+      /* coordinate doesn't exist */
+      coord_y = g_hash_table_new_full (g_str_hash,  /* is freed during   */
+                                       g_str_equal, /* destroying        */
+                                       g_free,      /* coords_to_cluster */
+                                       NULL);       /* hash table        */
+      g_hash_table_insert (coords_to_cluster, str_coord_x, coord_y);
+    }
+  
+  op = g_hash_table_lookup (coord_y, str_coord_y);
+  if (op == NULL)
+    {
+      op = g_ptr_array_new ();
+      g_hash_table_insert (coord_y, str_coord_y, op);
+      g_ptr_array_add (list_of_overlapping_points, op);
+    }
+  
+  g_ptr_array_add (op, point);
+}
+
+void
+npd_set_overlapping_points_weight (NPDOverlappingPoints *op,
+                                   gfloat                weight)
+{
+  gint i;
+  for (i = 0; i < op->num_of_points; i++)
+    {
+      (*op->points[i]->weight) = weight;
+    }
+}
+
+void
+npd_set_point_coordinates (NPDPoint *target,
+                           NPDPoint *source)
+{
+  target->x = source->x;
+  target->y = source->y;
+}
+
+/**
+ * Sets type of deformation. The function doesn't perform anything if supplied
+ * deformation type doesn't differ from currently set one.
+ *
+ * @param model
+ * @param ASAP          TRUE = ASAP deformation, FALSE = ARAP deformation
+ * @param MLS_weights   use weights from Moving Least Squares deformation method
+ */
+void
+npd_set_deformation_type (NPDModel *model,
+                          gboolean ASAP,
+                          gboolean MLS_weights)
+{
+  NPDHiddenModel *hm = model->hidden_model;
+
+  if (hm->ASAP == ASAP && hm->MLS_weights == MLS_weights) return;
+
+  if (MLS_weights)
+    npd_compute_MLS_weights (model);
+  else if (hm->MLS_weights)
+    npd_reset_weights (hm);
+
+  hm->ASAP = ASAP;
+  hm->MLS_weights = MLS_weights;
+}
+
+void
+npd_compute_MLS_weights (NPDModel *model)
+{
+  NPDHiddenModel       *hm = model->hidden_model;
+  NPDControlPoint      *cp;
+  NPDOverlappingPoints *op;
+  NPDPoint             *cp_reference, *op_reference;
+  gfloat                min, SED, MLS_weight;
+  gint                  i, j;
+
+  if (model->control_points->len == 0)
+    {
+      npd_reset_weights (hm);
+      return;
+    }
+
+  for (i = 0; i < hm->num_of_overlapping_points; i++)
+    {
+      op           = &hm->list_of_overlapping_points[i];
+      op_reference = op->representative->counterpart;
+      min          = INFINITY;
+
+      for (j = 0; j < model->control_points->len; j++)
+        {
+          cp = &g_array_index (model->control_points,
+                               NPDControlPoint,
+                               j);
+          cp_reference = cp->overlapping_points->representative->counterpart;
+
+          /* TODO - use geodetic distance */
+          SED = npd_SED (cp_reference,
+                         op_reference);
+          if (SED < min) min = SED;
+        }
+
+      if (npd_equal_floats (min, 0.0)) min = 0.0000001;
+      MLS_weight = 1 / pow (min, hm->MLS_weights_alpha);
+      npd_set_overlapping_points_weight (op, MLS_weight);
+    }
+}
+
+void
+npd_reset_weights (NPDHiddenModel *hm)
+{
+  NPDOverlappingPoints *op;
+  gint                  i;
+
+  for (i = 0; i < hm->num_of_overlapping_points; i++)
+    {
+      op  = &hm->list_of_overlapping_points[i];
+      npd_set_overlapping_points_weight (op, 1.0);
+    }
+}
+
+void
+npd_print_model (NPDModel        *model,
+                 gboolean         print_control_points)
+{
+  gint i;
+  g_printf ("NPDModel:\n");
+  g_printf ("control point radius: %f\n", model->control_point_radius);
+  g_printf ("control points visible: %d\n", model->control_points_visible);
+  g_printf ("mesh visible: %d\n", model->mesh_visible);
+  g_printf ("texture visible: %d\n", model->texture_visible);
+  g_printf ("mesh square size: %d\n", model->mesh_square_size);
+
+  npd_print_hidden_model (model->hidden_model, FALSE, FALSE);
+
+  if (print_control_points)
+    {
+      g_printf ("%d control points:\n", model->control_points->len);
+      for (i = 0; i < model->control_points->len; i++)
+        {
+          NPDControlPoint *cp = &g_array_index (model->control_points,
+                                                NPDControlPoint,
+                                                i);
+          npd_print_point (&cp->point, TRUE);
+        }
+    }
+}
+
+void
+npd_print_hidden_model (NPDHiddenModel *hm,
+                        gboolean        print_bones,
+                        gboolean        print_overlapping_points)
+{
+  gint i;
+  g_printf ("NPDHiddenModel:\n");
+  g_printf ("number of bones: %d\n", hm->num_of_bones);
+  g_printf ("ASAP: %d\n", hm->ASAP);
+  g_printf ("MLS weights: %d\n", hm->MLS_weights);
+  g_printf ("number of overlapping points: %d\n", hm->num_of_overlapping_points);
+  
+  if (print_bones)
+    {
+      g_printf ("bones:\n");
+      for (i = 0; i < hm->num_of_bones; i++)
+        {
+          npd_print_bone (&hm->current_bones[i]);
+        }
+    }
+  
+  if (print_overlapping_points)
+    {
+      g_printf ("overlapping points:\n");
+      for (i = 0; i < hm->num_of_overlapping_points; i++)
+        {
+          npd_print_overlapping_points (&hm->list_of_overlapping_points[i]);
+        }
+    }
+}
+
+void
+npd_print_bone (NPDBone *bone)
+{
+  gint i;
+  g_printf ("NPDBone:\n");
+  g_printf ("number of points: %d\n", bone->num_of_points);
+  g_printf ("points:\n");
+  for (i = 0; i < bone->num_of_points; i++)
+    {
+      npd_print_point (&bone->points[i], TRUE);
+    }
+}
+
+void
+npd_print_point (NPDPoint *point,
+                 gboolean  print_details)
+{
+  if (print_details)
+    {
+      g_printf ("(NPDPoint: x: %f, y: %f, weight: %f, fixed: %d)\n",
+              point->x, point->y, *point->weight, point->fixed);
+    }
+  else
+    {
+      g_printf ("(NPDPoint: x: %f, y: %f)\n",
+              point->x, point->y);
+    }
+}
+
+void
+npd_print_overlapping_points (NPDOverlappingPoints *op)
+{
+  gint i;
+  g_printf ("NPDOverlappingPoints:\n");
+  g_printf ("number of points: %d\n", op->num_of_points);
+  g_printf ("representative: ");
+  npd_print_point (op->representative, TRUE);
+  g_printf ("points:\n");
+  for (i = 0; i < op->num_of_points; i++)
+    {
+      npd_print_point (op->points[i], TRUE);
+    }
+}
diff --git a/libs/npd/npd_common.h b/libs/npd/npd_common.h
new file mode 100644
index 0000000..b3be82e
--- /dev/null
+++ b/libs/npd/npd_common.h
@@ -0,0 +1,157 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#ifndef __NPD_COMMON_H__
+#define        __NPD_COMMON_H__
+
+#include <glib.h>
+
+/* opaque types for independency on used display library */
+typedef struct _NPDImage    NPDImage;
+typedef struct _NPDColor    NPDColor;
+typedef struct _NPDDisplay  NPDDisplay;
+typedef struct _NPDMatrix   NPDMatrix;
+
+typedef struct _NPDPoint    NPDPoint;
+typedef struct _NPDBone     NPDBone;
+typedef struct _NPDOverlappingPoints NPDOverlappingPoints;
+
+struct _NPDPoint
+{
+  gfloat                x;
+  gfloat                y;
+  gboolean              fixed;
+  gfloat               *weight;            /* pointer to weight in array of weights */
+  gint                  index;
+  NPDBone              *current_bone;
+  NPDBone              *reference_bone;
+  NPDPoint             *counterpart;
+  NPDOverlappingPoints *overlapping_points;
+};
+
+struct _NPDBone
+{
+  gint                  num_of_points;
+  NPDPoint             *points;            /* array of points */
+  gfloat               *weights;           /* array of weights */
+};
+
+struct _NPDOverlappingPoints
+{
+  gint                  num_of_points;
+  NPDPoint             *representative;    /* pointer to representative of cluster */
+  NPDPoint            **points;            /* array of pointers to points */
+};
+
+typedef struct
+{
+  gint                  num_of_bones;
+  gint                  num_of_overlapping_points;
+  gboolean              ASAP;              /* don't change directly!
+                                            * use npd_set_deformation_type function */
+  gboolean              MLS_weights;       /* don't change directly!
+                                            * use npd_set_deformation_type function */
+  gfloat                MLS_weights_alpha;
+  NPDBone              *current_bones;     /* array of current bones */
+  NPDBone              *reference_bones;   /* array of reference bones */
+  NPDOverlappingPoints *list_of_overlapping_points; /* array of overlapping points */
+} NPDHiddenModel;
+
+typedef struct
+{
+  NPDPoint              point;
+  NPDOverlappingPoints *overlapping_points; /* pointer to overlapping points */
+} NPDControlPoint;
+
+typedef struct
+{
+  gfloat                control_point_radius;
+  gboolean              control_points_visible;
+  gboolean              mesh_visible;
+  gboolean              texture_visible;
+  gint                  mesh_square_size;
+  GArray               *control_points;     /* GArray of control points */
+  NPDHiddenModel       *hidden_model;
+  NPDImage             *reference_image;
+  NPDDisplay           *display;
+} NPDModel;
+
+#define npd_init(set_pixel, get_pixel,                                         \
+                 draw_line)                                                    \
+npd_set_pixel_color      = set_pixel;                                          \
+npd_get_pixel_color      = get_pixel;                                          \
+npd_draw_line            = draw_line;
+
+void             npd_init_model                 (NPDModel        *model);
+void             npd_destroy_hidden_model       (NPDHiddenModel  *model);
+void             npd_destroy_model              (NPDModel        *model);
+
+NPDControlPoint *npd_add_control_point          (NPDModel        *model,
+                                                 NPDPoint        *coord);
+void             npd_remove_control_point       (NPDModel        *model,
+                                                 NPDControlPoint *control_point);
+void             npd_remove_control_points      (NPDModel        *model,
+                                                 GList           *control_points);
+void             npd_remove_all_control_points  (NPDModel        *model);
+void             npd_set_control_point_weight   (NPDControlPoint *cp,
+                                                 gfloat           weight);
+gboolean         npd_equal_coordinates          (NPDPoint        *p1,
+                                                 NPDPoint        *p2);
+gboolean         npd_equal_coordinates_epsilon  (NPDPoint        *p1,
+                                                 NPDPoint        *p2,
+                                                 gfloat           epsilon);
+NPDControlPoint *npd_get_control_point_with_radius_at
+                                                (NPDModel        *model,
+                                                 NPDPoint        *coord,
+                                                 gfloat           control_point_radius);
+NPDControlPoint *npd_get_control_point_at       (NPDModel        *model,
+                                                 NPDPoint        *coord);
+void             npd_create_square              (NPDBone         *square,
+                                                 gint             x,
+                                                 gint             y,
+                                                 gint             width,
+                                                 gint             height);
+void             npd_create_list_of_overlapping_points
+                                                (NPDHiddenModel  *model);
+void             add_point_to_suitable_cluster  (GHashTable      *coords_to_cluster,
+                                                 NPDPoint        *point,
+                                                 GPtrArray       *list_of_overlapping_points);
+void             npd_set_overlapping_points_weight
+                                                (NPDOverlappingPoints *op,
+                                                 gfloat           weight);
+void             npd_set_point_coordinates      (NPDPoint        *target,
+                                                 NPDPoint        *source);
+void             npd_set_deformation_type       (NPDModel        *model,
+                                                 gboolean         ASAP,
+                                                 gboolean         MLS_weights);
+void             npd_compute_MLS_weights        (NPDModel        *model);
+void             npd_reset_weights              (NPDHiddenModel  *hidden_model);
+
+void             npd_print_model                (NPDModel        *model,
+                                                 gboolean         print_control_points);
+void             npd_print_hidden_model         (NPDHiddenModel  *hm,
+                                                 gboolean         print_bones,
+                                                 gboolean         print_overlapping_points);
+void             npd_print_bone                 (NPDBone         *bone);
+void             npd_print_point                (NPDPoint        *point,
+                                                 gboolean         print_details);
+void             npd_print_overlapping_points   (NPDOverlappingPoints *op);
+#endif /* __NPD_COMMON_H__ */
diff --git a/libs/npd/npd_gegl.c b/libs/npd/npd_gegl.c
new file mode 100644
index 0000000..1a8abaf
--- /dev/null
+++ b/libs/npd/npd_gegl.c
@@ -0,0 +1,68 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#include "npd_gegl.h"
+#include <glib.h>
+
+void
+npd_new_matrix (NPDMatrix **matrix)
+{
+  *matrix = g_new (NPDMatrix, 1);
+}
+
+void
+npd_destroy_matrix (NPDMatrix **matrix)
+{
+  g_free (*matrix);
+}
+
+void
+npd_compute_affinity (NPDPoint  *p11,
+                      NPDPoint  *p21,
+                      NPDPoint  *p31,
+                      NPDPoint  *p12,
+                      NPDPoint  *p22,
+                      NPDPoint  *p32,
+                      NPDMatrix *T)
+{
+  GeglMatrix3 Y, X;
+  
+  Y.coeff[0][0] = p12->x; Y.coeff[1][0] = p12->y; Y.coeff[2][0] = 1;
+  Y.coeff[0][1] = p22->x; Y.coeff[1][1] = p22->y; Y.coeff[2][1] = 1;
+  Y.coeff[0][2] = p32->x; Y.coeff[1][2] = p32->y; Y.coeff[2][2] = 1;
+  
+  X.coeff[0][0] = p11->x; X.coeff[1][0] = p11->y; X.coeff[2][0] = 1;
+  X.coeff[0][1] = p21->x; X.coeff[1][1] = p21->y; X.coeff[2][1] = 1;
+  X.coeff[0][2] = p31->x; X.coeff[1][2] = p31->y; X.coeff[2][2] = 1;
+  
+  gegl_matrix3_invert (&X);
+  gegl_matrix3_multiply (&Y, &X, &T->matrix);
+}
+
+void
+npd_apply_transformation (NPDMatrix *T,
+                          NPDPoint  *src,
+                          NPDPoint  *dest)
+{
+  gdouble x = src->x, y = src->y;
+  gegl_matrix3_transform_point (&T->matrix, &x, &y);
+  dest->x = x; dest->y = y;
+}
diff --git a/libs/npd/npd_gegl.h b/libs/npd/npd_gegl.h
new file mode 100644
index 0000000..7ce4790
--- /dev/null
+++ b/libs/npd/npd_gegl.h
@@ -0,0 +1,33 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#ifndef __NPD_GEGL_H__
+#define        __NPD_GEGL_H__
+
+#include "npd_math.h"
+#include <gegl-matrix.h>
+
+struct _NPDMatrix
+{
+  GeglMatrix3 matrix;
+};
+
+#endif /* __NPD_GEGL_H__ */
diff --git a/libs/npd/npd_math.c b/libs/npd/npd_math.c
new file mode 100644
index 0000000..b14f73f
--- /dev/null
+++ b/libs/npd/npd_math.c
@@ -0,0 +1,47 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#include "npd_math.h"
+#include <math.h>
+
+gboolean
+npd_equal_floats (gfloat a,
+                  gfloat b)
+{
+  return npd_equal_floats_epsilon (a, b, NPD_EPSILON);
+}
+
+gboolean
+npd_equal_floats_epsilon (gfloat a,
+                          gfloat b,
+                          gfloat epsilon)
+{
+  return fabs (a - b) < epsilon;
+}
+
+gfloat
+npd_SED (NPDPoint *p1,
+         NPDPoint *p2)
+{
+  gint dx = p1->x - p2->x;
+  gint dy = p1->y - p2->y;
+  return dx * dx + dy * dy;
+}
diff --git a/libs/npd/npd_math.h b/libs/npd/npd_math.h
new file mode 100644
index 0000000..411f4e9
--- /dev/null
+++ b/libs/npd/npd_math.h
@@ -0,0 +1,58 @@
+/*
+ * This file is part of N-point image deformation library.
+ *
+ * N-point image deformation library is free software: you can
+ * redistribute it and/or modify it under the terms of the
+ * GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 3 of the License,
+ * or (at your option) any later version.
+ *
+ * N-point image deformation library is distributed in the hope
+ * that it will be useful, but WITHOUT ANY WARRANTY; without even
+ * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ * PURPOSE.  See the GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with N-point image deformation library.
+ * If not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#ifndef __NPD_MATH_H__
+#define        __NPD_MATH_H__
+
+#include "npd_common.h"
+
+#define NPD_EPSILON 0.00001
+
+void        npd_compute_homography   (NPDPoint  *p11,
+                                      NPDPoint  *p21,
+                                      NPDPoint  *p31,
+                                      NPDPoint  *p41,
+                                      NPDPoint  *p12,
+                                      NPDPoint  *p22,
+                                      NPDPoint  *p32,
+                                      NPDPoint  *p42,
+                                      NPDMatrix *T);
+void        npd_compute_affinity     (NPDPoint  *p11,
+                                      NPDPoint  *p21,
+                                      NPDPoint  *p31,
+                                      NPDPoint  *p12,
+                                      NPDPoint  *p22,
+                                      NPDPoint  *p32,
+                                      NPDMatrix *T);
+void        npd_apply_transformation (NPDMatrix *T,
+                                      NPDPoint  *src,
+                                      NPDPoint  *dest);
+gboolean    npd_equal_floats_epsilon (gfloat a,
+                                      gfloat b,
+                                      gfloat epsilon);
+gboolean    npd_equal_floats         (gfloat a,
+                                      gfloat b);
+gfloat      npd_SED                  (NPDPoint *p1,
+                                      NPDPoint *p2);
+void        npd_new_matrix           (NPDMatrix **matrix);
+void        npd_destroy_matrix       (NPDMatrix **matrix);
+
+#endif /* __NPD_MATH_H__ */
diff --git a/operations/external/Makefile.am b/operations/external/Makefile.am
index b929744..6b2dc4f 100644
--- a/operations/external/Makefile.am
+++ b/operations/external/Makefile.am
@@ -158,5 +158,13 @@ rgbe_save_la_SOURCES = rgbe-save.c
 rgbe_save_la_CFLAGS = $(AM_CFLAGS) -I $(top_srcdir)/libs
 rgbe_save_la_LIBADD = $(op_libs) $(top_builddir)/libs/rgbe/librgbe.la
 
+# Dependencies are in our source tree
+if HAVE_CAIRO
+ops += npd.la
+npd_la_SOURCES = npd.c
+npd_la_CFLAGS = $(AM_CFLAGS) $(NPD_CFLAGS) $(CAIRO_CFLAGS)
+npd_la_LIBADD = $(op_libs) $(NPD_LIBS) $(CAIRO_LIBS)
+endif
+
 opdir = $(libdir)/gegl- GEGL_API_VERSION@
 op_LTLIBRARIES = $(ops)
diff --git a/operations/external/npd.c b/operations/external/npd.c
new file mode 100644
index 0000000..63bfc37
--- /dev/null
+++ b/operations/external/npd.c
@@ -0,0 +1,329 @@
+/* This file is an image processing operation for GEGL
+ *
+ * GEGL is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public
+ * License as published by the Free Software Foundation; either
+ * version 3 of the License, or (at your option) any later version.
+ *
+ * GEGL is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GEGL; if not, see <http://www.gnu.org/licenses/>.
+ *
+ * Copyright (C) 2013 Marek Dvoroznak <dvoromar gmail com>
+ */
+
+#include "config.h"
+#include <glib/gi18n-lib.h>
+
+#ifdef GEGL_CHANT_PROPERTIES
+gegl_chant_pointer (model,       _("model"),
+                    _("Model - basic element we operate on"))
+
+gegl_chant_int     (square_size, _("square size"),
+                    5,  1000,  20,
+                    _("Size of an edge of square the mesh consists of"))
+
+gegl_chant_int     (rigidity,    _("rigidity"),
+                    0, 10000, 100,
+                    _("The number of deformation iterations"))
+
+gegl_chant_boolean (ASAP_deformation, _("ASAP deformation"),
+                    FALSE,
+                    _("ASAP deformation is performed when TRUE, ARAP deformation otherwise"))
+
+gegl_chant_boolean (MLS_weights, _("MLS weights"),
+                    FALSE,
+                    _("Use MLS weights"))
+
+gegl_chant_double  (MLS_weights_alpha, _("MLS weights alpha"),
+                    0.1, 2.0, 1.0,
+                    _("Alpha parameter of MLS weights"))
+
+gegl_chant_boolean (mesh_visible, _("mesh visible"),
+                    TRUE,
+                    _("Should the mesh be visible?"))
+#else
+
+#define GEGL_CHANT_TYPE_FILTER
+#define GEGL_CHANT_C_FILE       "npd.c"
+
+#include "gegl-chant.h"
+#include <stdio.h>
+#include <math.h>
+#include <npd/npd.h>
+#include <npd/npd_gegl.h>
+#include <cairo.h>
+
+struct _NPDImage
+{
+  gint     width;
+  gint     height;
+  NPDPoint position;
+  guchar  *buffer;
+};
+
+struct _NPDDisplay
+{
+  NPDImage  image;
+  cairo_t  *cr;
+};
+
+typedef struct
+{
+  gboolean  first_run;
+  NPDModel  model;
+} NPDProperties;
+
+void npd_create_image         (NPDImage   *image,
+                               GeglBuffer *gegl_buffer,
+                               const Babl *format);
+
+void npd_set_pixel_color_impl (NPDImage *image,
+                               gint      x,
+                               gint      y,
+                               NPDColor *color);
+
+void npd_get_pixel_color_impl (NPDImage *image,
+                               gint      x,
+                               gint      y,
+                               NPDColor *color);
+
+void npd_draw_line_impl       (NPDDisplay *display,
+                               gfloat      x0,
+                               gfloat      y0,
+                               gfloat      x1,
+                               gfloat      y1);
+
+void npd_set_pixel_color_impl (NPDImage *image,
+                               gint      x,
+                               gint      y,
+                               NPDColor *color)
+{
+  if (x > -1 && x < image->width &&
+      y > -1 && y < image->height)
+    {
+      gint position = 4 * (y * image->width + x);
+
+      image->buffer[position + 0] = color->r;
+      image->buffer[position + 1] = color->g;
+      image->buffer[position + 2] = color->b;
+      image->buffer[position + 3] = color->a;
+    }
+}
+
+void
+npd_get_pixel_color_impl (NPDImage *image,
+                          gint      x,
+                          gint      y,
+                          NPDColor *color)
+{
+  if (x > -1 && x < image->width &&
+      y > -1 && y < image->height)
+    {
+      gint position = 4 * (y * image->width + x);
+
+      color->r = image->buffer[position + 0];
+      color->g = image->buffer[position + 1];
+      color->b = image->buffer[position + 2];
+      color->a = image->buffer[position + 3];
+    }
+  else
+    {
+      color->r = color->g = color->b = color->a = 0;
+    }
+}
+
+void npd_draw_line_impl (NPDDisplay *display,
+                         gfloat      x0,
+                         gfloat      y0,
+                         gfloat      x1,
+                         gfloat      y1)
+{
+  cairo_move_to (display->cr, x0, y0);
+  cairo_line_to (display->cr, x1, y1);
+}
+
+void
+npd_draw_model (NPDModel   *model,
+                NPDDisplay *display)
+{
+  NPDHiddenModel *hm = model->hidden_model;
+  NPDImage *image = model->reference_image;
+  gint i;
+
+  /* draw texture */
+  if (model->texture_visible)
+    {
+      for (i = 0; i < hm->num_of_bones; i++)
+        {
+          npd_texture_quadrilateral(&hm->reference_bones[i],
+                                    &hm->current_bones[i],
+                                     image,
+                                    &display->image,
+                                     NPD_BILINEAR_INTERPOLATION | NPD_ALPHA_BLENDING);
+        }
+    }
+  
+  /* draw mesh */
+  if (model->mesh_visible)
+    {
+      cairo_surface_t *surface;
+
+      surface = cairo_image_surface_create_for_data (display->image.buffer,
+                                                     CAIRO_FORMAT_ARGB32,
+                                                     display->image.width,
+                                                     display->image.height,
+                                                     display->image.width * 4);
+      display->cr = cairo_create (surface);
+      cairo_set_line_width (display->cr, 1);
+      cairo_set_source_rgba (display->cr, 0, 0, 0, 1);
+      npd_draw_mesh (model, display);
+      cairo_stroke (display->cr);
+    }
+}
+
+void
+npd_create_model_from_image (NPDModel *model,
+                             NPDImage *image,
+                             gint      square_size)
+{
+  npd_init_model(model);
+  model->reference_image = image;
+  model->mesh_square_size = square_size;
+    
+  npd_create_mesh_from_image (model, image->width, image->height, 0, 0);
+}
+
+void
+npd_create_image (NPDImage   *image,
+                  GeglBuffer *gegl_buffer,
+                  const Babl *format)
+{
+  guchar *buffer;
+  buffer = g_new0 (guchar, gegl_buffer_get_pixel_count (gegl_buffer) * 4);
+  gegl_buffer_get (gegl_buffer, NULL, 1.0, format,
+                   buffer, GEGL_AUTO_ROWSTRIDE, GEGL_ABYSS_NONE);
+
+  image->buffer = buffer;
+  image->width = gegl_buffer_get_width (gegl_buffer);
+  image->height = gegl_buffer_get_height (gegl_buffer);
+}
+
+void
+npd_destroy_image (NPDImage *image)
+{
+  g_free(image->buffer);
+}
+
+static void
+prepare (GeglOperation *operation)
+{
+  GeglChantO    *o      = GEGL_CHANT_PROPERTIES (operation);
+  NPDProperties *props;
+
+  if (o->chant_data == NULL)
+    {
+      props = g_new (NPDProperties, 1);
+      props->first_run = TRUE;
+      o->chant_data    = props;
+    }
+  
+  gegl_operation_set_format (operation, "input",
+                             babl_format ("RGBA float"));
+  gegl_operation_set_format (operation, "output",
+                             babl_format ("RGBA float"));
+}
+
+static gboolean
+process (GeglOperation       *operation,
+         GeglBuffer          *input,
+         GeglBuffer          *output,
+         const GeglRectangle *result,
+         gint                 level)
+{
+  GeglChantO *o      = GEGL_CHANT_PROPERTIES (operation);
+  const Babl *format = babl_format ("RGBA u8");
+  NPDProperties *props = o->chant_data;
+  NPDModel *model = &props->model;
+  NPDHiddenModel *hm;
+  guchar *output_buffer;
+  gint length = gegl_buffer_get_pixel_count (input) * 4;
+
+  if (props->first_run)
+    {
+      gint width, height;
+      NPDImage *input_image = g_new (NPDImage, 1);
+      NPDDisplay *display = g_new (NPDDisplay, 1);
+
+      npd_init (npd_set_pixel_color_impl,
+                npd_get_pixel_color_impl,
+                npd_draw_line_impl);
+
+      npd_create_image (input_image, input, format);
+      width = input_image->width;
+      height = input_image->height;
+
+      output_buffer = g_new0 (guchar, length);
+      display->image.width = width;
+      display->image.height = height;
+      display->image.buffer = output_buffer;
+      model->display = display;
+
+      npd_create_model_from_image (model, input_image, o->square_size);
+      hm = model->hidden_model;
+/*      npd_create_list_of_overlapping_points (hm);*/
+
+      o->model = model;
+
+      memcpy (output_buffer, input_image->buffer, length);
+
+      props->first_run = FALSE;
+    }
+  else
+    {
+      npd_set_deformation_type (model, o->ASAP_deformation, o->MLS_weights);
+
+      if (o->MLS_weights &&
+          model->hidden_model->MLS_weights_alpha != o->MLS_weights_alpha)
+        {
+          model->hidden_model->MLS_weights_alpha = o->MLS_weights_alpha;
+          npd_compute_MLS_weights (model);
+        }
+
+      model->mesh_visible = o->mesh_visible;
+
+      output_buffer = model->display->image.buffer;
+      memset (output_buffer, 0, length);
+      npd_deform_model (model, o->rigidity);
+      npd_draw_model (model, model->display);
+    }
+
+  gegl_buffer_set (output, NULL, 0, format, output_buffer, GEGL_AUTO_ROWSTRIDE);
+  
+  return TRUE;
+}
+
+static void
+gegl_chant_class_init (GeglChantClass *klass)
+{
+  GeglOperationClass       *operation_class;
+  GeglOperationFilterClass *filter_class;
+
+  operation_class = GEGL_OPERATION_CLASS (klass);
+  filter_class    = GEGL_OPERATION_FILTER_CLASS (klass);
+
+  filter_class->process    = process;
+  operation_class->prepare = prepare;
+
+  gegl_operation_class_set_keys (operation_class,
+    "categories"  , "transform",
+    "name"        , "gegl:npd",
+    "description" , _("Performs n-point image deformation"),
+    NULL);
+}
+
+#endif
\ No newline at end of file
diff --git a/po/POTFILES.in b/po/POTFILES.in
index 91377d0..8b672fa 100644
--- a/po/POTFILES.in
+++ b/po/POTFILES.in
@@ -117,6 +117,7 @@ operations/external/jpg-load.c
 operations/external/jpg-save.c
 operations/external/lcms-from-profile.c
 operations/external/matting-levin.c
+operations/external/npd.c
 operations/external/npy-save.c
 operations/external/path.c
 operations/external/pixbuf.c



[Date Prev][Date Next]   [Thread Prev][Thread Next]   [Thread Index] [Date Index] [Author Index]