From 7a49a338e41da142b1b9c5253efd79d3663d750e Mon Sep 17 00:00:00 2001
From: sam <sam@92316355-f0b4-4df1-b90c-862c8a59935f>
Date: Sat, 30 May 2009 16:01:09 +0000
Subject: [PATCH] img2twit: prevent two points from having the same
 coordinates, since this is forbidden by the image format; simplify operation
 computations by not packing coordinates at runtime; add a version information
 in the file format for future extension.

git-svn-id: file:///srv/caca.zoy.org/var/lib/svn/libpipi/trunk@3537 92316355-f0b4-4df1-b90c-862c8a59935f
---
 examples/img2twit.cpp | 303 +++++++++++++++++++++---------------------
 1 file changed, 155 insertions(+), 148 deletions(-)

diff --git a/examples/img2twit.cpp b/examples/img2twit.cpp
index add22e3..e97328a 100644
--- a/examples/img2twit.cpp
+++ b/examples/img2twit.cpp
@@ -66,9 +66,10 @@ static const uint32_t unichars_symbols[] = { RANGE_SYMBOLS, RANGE_END };
 /* The Unicode characters at disposal */
 static const uint32_t *unichars;
 
-/* The maximum image size we want to support */
-#define RANGE_W 4000
-#define RANGE_H 4000
+/* The maximum image size we want to support, and the version range */
+#define RANGE_W 2000
+#define RANGE_H 2000
+#define RANGE_V 10
 
 /* How does the algorithm work: one point per cell, or two? XXX: changing
  * this value breaks compatibility with other images. */
@@ -117,13 +118,9 @@ static int NUM_CHARACTERS;
 static int MAX_ITERATIONS;
 static unsigned int TOTAL_CELLS;
 
-#define RANGE_XY (RANGE_Y*RANGE_X)
 #define RANGE_XY2 (RANGE_Y*RANGE_X*(RANGE_Y*RANGE_X+1)/2)
-#define RANGE_RXY (RANGE_Y*RANGE_X*RANGE_R)
-#define RANGE_GRXY (RANGE_Y*RANGE_X*RANGE_R*RANGE_G)
-#define RANGE_BGRXY (RANGE_Y*RANGE_X*RANGE_R*RANGE_G*RANGE_B)
-#define RANGE_SBGRXY (RANGE_Y*RANGE_X*RANGE_R*RANGE_G*RANGE_B*RANGE_S)
 #define RANGE_SBGR (RANGE_R*RANGE_G*RANGE_B*RANGE_S)
+#define RANGE_SBGRXY (RANGE_Y*RANGE_X*RANGE_R*RANGE_G*RANGE_B*RANGE_S)
 
 struct K : CGAL::Exact_predicates_inexact_constructions_kernel {};
 typedef CGAL::Delaunay_triangulation_2<K> Delaunay_triangulation;
@@ -132,8 +129,16 @@ typedef std::vector<std::pair<K::Point_2, K::FT> > Point_coordinate_vector;
 /* Global aspect ratio */
 static unsigned int dw, dh;
 
+/* Algorithm version */
+static unsigned int version;
+
 /* Global point encoding */
-static uint32_t points[4096]; /* FIXME: allocate this dynamically */
+typedef struct point
+{
+    uint8_t x, y, r, g, b, s;
+}
+point_t;
+static point_t points[4096]; /* FIXME: allocate this dynamically */
 static int npoints = 0;
 
 /* Global triangulation */
@@ -146,7 +151,7 @@ static Delaunay_triangulation dt;
 void compute_ranges(int width, int height)
 {
     TOTAL_BITS = MAX_MSG_LEN * logf(NUM_CHARACTERS) / logf(2);
-    HEADER_BITS = logf(RANGE_W * RANGE_H) / logf(2);
+    HEADER_BITS = logf(RANGE_W * RANGE_H * RANGE_V) / logf(2);
     DATA_BITS = TOTAL_BITS - HEADER_BITS;
 #if POINTS_PER_CELL == 2
     CELL_BITS = (2 * logf(RANGE_SBGR) + logf(RANGE_XY2)) / logf(2);
@@ -507,47 +512,43 @@ static inline void set_point(int index, float x, float y, float r,
     float fx = (x - dx * RANGE_X) / RANGE_X;
     float fy = (y - dy * RANGE_Y) / RANGE_Y;
 
-    int iy = range2int(fy, RANGE_Y);
-    int ix = range2int(fx, RANGE_X);
+    points[index].x = range2int(fx, RANGE_X);
+    points[index].y = range2int(fy, RANGE_Y);
 
-    int ir = range2int(r, RANGE_R);
-    int ig = range2int(g, RANGE_G);
-    int ib = range2int(b, RANGE_B);
+    points[index].r = range2int(r, RANGE_R);
+    points[index].g = range2int(g, RANGE_G);
+    points[index].b = range2int(b, RANGE_B);
 
-    int is = range2int(s, RANGE_S);
-
-    points[index] = iy + RANGE_Y * (ix + RANGE_X * (ir + RANGE_R *
-                      (ig + (RANGE_G * ib + (RANGE_B * is)))));
+    points[index].s = range2int(s, RANGE_S);
 }
 
 static inline void get_point(int index, float *x, float *y, float *r,
                              float *g, float *b, float *s, bool final = false)
 {
-    uint32_t pt = points[index];
     int dx, dy;
 
     index2cell(index, &dx, &dy);
 
-    float fy = int2midrange(pt % RANGE_Y, RANGE_Y); pt /= RANGE_Y;
-    float fx = int2midrange(pt % RANGE_X, RANGE_X); pt /= RANGE_X;
+    float fx = int2midrange(points[index].x, RANGE_X);
+    float fy = int2midrange(points[index].y, RANGE_Y);
 
     *y = (fy + dy) * RANGE_Y /*+ 0.5 * (index & 1)*/;
     *x = (fx + dx) * RANGE_X /*+ 0.5 * (index & 1)*/;
 
     if(final)
     {
-        *r = int2fullrange(pt % RANGE_R, RANGE_R); pt /= RANGE_R;
-        *g = int2fullrange(pt % RANGE_G, RANGE_G); pt /= RANGE_G;
-        *b = int2fullrange(pt % RANGE_B, RANGE_B); pt /= RANGE_B;
+        *r = int2fullrange(points[index].r, RANGE_R);
+        *g = int2fullrange(points[index].g, RANGE_G);
+        *b = int2fullrange(points[index].b, RANGE_B);
     }
     else
     {
-        *r = int2midrange(pt % RANGE_R, RANGE_R); pt /= RANGE_R;
-        *g = int2midrange(pt % RANGE_G, RANGE_G); pt /= RANGE_G;
-        *b = int2midrange(pt % RANGE_B, RANGE_B); pt /= RANGE_B;
+        *r = int2midrange(points[index].r, RANGE_R);
+        *g = int2midrange(points[index].g, RANGE_G);
+        *b = int2midrange(points[index].b, RANGE_B);
     }
 
-    *s = int2fullrange(pt % RANGE_S, RANGE_S); pt /= RANGE_S;
+    *s = int2fullrange(points[index].s, RANGE_S);
 }
 
 static void add_point(float x, float y, float r, float g, float b, float s)
@@ -589,7 +590,12 @@ static void unpack_coords(uint32_t pack, int *x1, int *y1, int *x2, int *y2)
 #if RANDOM_START == 1
 static void add_random_point()
 {
-    points[npoints] = det_rand(RANGE_SBGRXY);
+    points[npoints].x = det_rand(RANGE_X);
+    points[npoints].y = det_rand(RANGE_Y);
+    points[npoints].r = det_rand(RANGE_R);
+    points[npoints].g = det_rand(RANGE_G);
+    points[npoints].b = det_rand(RANGE_B);
+    points[npoints].s = det_rand(RANGE_S);
     npoints++;
 }
 #endif
@@ -613,104 +619,75 @@ static uint8_t rand_op(void)
     return x;
 }
 
-static uint32_t apply_op(uint8_t op, uint32_t val)
+static void apply_op(uint8_t op, point_t *val)
 {
-    uint32_t rem, ext;
-
     switch(op)
     {
     case 0: /* Flip strength value */
     case 1:
         /* Statistics show that this helps often, but does not reduce
          * the error significantly. */
-        rem = val % RANGE_BGRXY;
-        ext = val / RANGE_BGRXY;
-        ext ^= 1;
-        return ext * RANGE_BGRXY + rem;
+        val->s ^= 1; break;
     case 2: /* Move up; if impossible, down */
-        ext = val % RANGE_Y;
-        ext = ext > 0 ? ext - 1 : ext + 1;
-        return val / RANGE_Y * RANGE_Y + ext;
+        val->y = val->y > 0 ? val->y - 1 : val->y + 1; break;
     case 3: /* Move down; if impossible, up */
-        ext = val % RANGE_Y;
-        ext = ext < RANGE_Y - 1 ? ext + 1 : ext - 1;
-        return val / RANGE_Y * RANGE_Y + ext;
+        val->y = val->y + 1U < RANGE_Y ? val->y + 1 : val->y - 1; break;
     case 4: /* Move left; if impossible, right */
-        rem = val % RANGE_Y;
-        ext = (val / RANGE_Y) % RANGE_X;
-        ext = ext > 0 ? ext - 1 : ext + 1;
-        return (val / RANGE_XY * RANGE_X + ext) * RANGE_Y + rem;
-    case 5: /* Move left; if impossible, right */
-        rem = val % RANGE_Y;
-        ext = (val / RANGE_Y) % RANGE_X;
-        ext = ext < RANGE_X - 1 ? ext + 1 : ext - 1;
-        return (val / RANGE_XY * RANGE_X + ext) * RANGE_Y + rem;
+        val->x = val->x > 0 ? val->x - 1 : val->x + 1; break;
+    case 5: /* Move right; if impossible, left */
+        val->x = val->x + 1U < RANGE_X ? val->x + 1 : val->x - 1; break;
     case 6: /* Corner 1 */
-        return apply_op(2, apply_op(4, val));
+        val->y = val->y > 0 ? val->y - 1 : val->y + 1;
+        val->x = val->x > 0 ? val->x - 1 : val->x + 1; break;
     case 7: /* Corner 2 */
-        return apply_op(2, apply_op(5, val));
+        val->y = val->y > 0 ? val->y - 1 : val->y + 1;
+        val->x = val->x + 1U < RANGE_X ? val->x + 1 : val->x - 1; break;
     case 8: /* Corner 3 */
-        return apply_op(3, apply_op(5, val));
+        val->y = val->y + 1U < RANGE_Y ? val->y + 1 : val->y - 1;
+        val->x = val->x + 1U < RANGE_X ? val->x + 1 : val->x - 1; break;
     case 9: /* Corner 4 */
-        return apply_op(3, apply_op(4, val));
+        val->y = val->y + 1U < RANGE_Y ? val->y + 1 : val->y - 1;
+        val->x = val->x > 0 ? val->x - 1 : val->x + 1; break;
     case 16: /* Double up */
-        return apply_op(2, apply_op(2, val));
+        val->y = val->y > 1 ? val->y - 2 : val->y + 2; break;
     case 17: /* Double down */
-        return apply_op(3, apply_op(3, val));
+        val->y = val->y + 2U < RANGE_Y ? val->y + 2 : val->y - 2; break;
     case 18: /* Double left */
-        return apply_op(4, apply_op(4, val));
+        val->x = val->x > 1 ? val->x - 2 : val->x + 2; break;
     case 19: /* Double right */
-        return apply_op(5, apply_op(5, val));
+        val->x = val->x + 2U < RANGE_X ? val->x + 2 : val->x - 2; break;
     case 10: /* R-- (or R++) */
-        rem = val % RANGE_XY;
-        ext = (val / RANGE_XY) % RANGE_R;
-        ext = ext > 0 ? ext - 1 : ext + 1;
-        return (val / RANGE_RXY * RANGE_R + ext) * RANGE_XY + rem;
+        val->r = val->r > 0 ? val->r - 1 : val->r + 1; break;
     case 11: /* R++ (or R--) */
-        rem = val % RANGE_XY;
-        ext = (val / RANGE_XY) % RANGE_R;
-        ext = ext < RANGE_R - 1 ? ext + 1 : ext - 1;
-        return (val / RANGE_RXY * RANGE_R + ext) * RANGE_XY + rem;
+        val->r = val->r + 1U < RANGE_R ? val->r + 1 : val->r - 1; break;
     case 12: /* G-- (or G++) */
-        rem = val % RANGE_RXY;
-        ext = (val / RANGE_RXY) % RANGE_G;
-        ext = ext > 0 ? ext - 1 : ext + 1;
-        return (val / RANGE_GRXY * RANGE_G + ext) * RANGE_RXY + rem;
+        val->g = val->g > 0 ? val->g - 1 : val->g + 1; break;
     case 13: /* G++ (or G--) */
-        rem = val % RANGE_RXY;
-        ext = (val / RANGE_RXY) % RANGE_G;
-        ext = ext < RANGE_G - 1 ? ext + 1 : ext - 1;
-        return (val / RANGE_GRXY * RANGE_G + ext) * RANGE_RXY + rem;
+        val->g = val->g + 1U < RANGE_G ? val->g + 1 : val->g - 1; break;
     case 14: /* B-- (or B++) */
-        rem = val % RANGE_GRXY;
-        ext = (val / RANGE_GRXY) % RANGE_B;
-        ext = ext > 0 ? ext - 1 : ext + 1;
-        return (val / RANGE_BGRXY * RANGE_B + ext) * RANGE_GRXY + rem;
+        val->b = val->b > 0 ? val->g - 1 : val->b + 1; break;
     case 15: /* B++ (or B--) */
-        rem = val % RANGE_GRXY;
-        ext = (val / RANGE_GRXY) % RANGE_B;
-        ext = ext < RANGE_B - 1 ? ext + 1 : ext - 1;
-        return (val / RANGE_BGRXY * RANGE_B + ext) * RANGE_GRXY + rem;
+        val->b = val->b + 1U < RANGE_B ? val->b + 1 : val->b - 1; break;
 #if 0
     case 15: /* Brightness-- */
-        return apply_op(9, apply_op(11, apply_op(13, val)));
+        apply_op(9, val); apply_op(11, val); apply_op(13, val); break;
     case 16: /* Brightness++ */
-        return apply_op(10, apply_op(12, apply_op(14, val)));
+        apply_op(10, val); apply_op(12, val); apply_op(14, val); break;
     case 17: /* RG-- */
-        return apply_op(9, apply_op(11, val));
+        apply_op(9, val); apply_op(11, val); break;
     case 18: /* RG++ */
-        return apply_op(10, apply_op(12, val));
+        apply_op(10, val); apply_op(12, val); break;
     case 19: /* GB-- */
-        return apply_op(11, apply_op(13, val));
+        apply_op(11, val); apply_op(13, val); break;
     case 20: /* GB++ */
-        return apply_op(12, apply_op(14, val));
+        apply_op(12, val); apply_op(14, val); break;
     case 21: /* RB-- */
-        return apply_op(9, apply_op(13, val));
+        apply_op(9, val); apply_op(13, val); break;
     case 22: /* RB++ */
-        return apply_op(10, apply_op(14, val));
+        apply_op(10, val); apply_op(14, val); break;
 #endif
     default:
-        return val;
+        break;
     }
 }
 
@@ -1072,10 +1049,13 @@ int main(int argc, char *argv[])
         for(int i = MAX_MSG_LEN; i--; )
             b.push(uni2index(unicode_data[i]), NUM_CHARACTERS);
 
+        /* The first thing we pop from the stream is the version information */
+        version = b.pop(RANGE_V);
+
         /* Read width and height from bitstream */
-        src = NULL;
         width = b.pop(RANGE_W) + 1;
         height = b.pop(RANGE_H) + 1;
+        src = NULL;
     }
     else
     {
@@ -1088,6 +1068,7 @@ int main(int argc, char *argv[])
             return EXIT_FAILURE;
         }
 
+        version = 1;
         width = pipi_get_image_width(src);
         height = pipi_get_image_height(src);
     }
@@ -1099,6 +1080,12 @@ int main(int argc, char *argv[])
         return EXIT_FAILURE;
     }
 
+    if(version > 1)
+    {
+        fprintf(stderr, "Error: unsupported algorithm version %i\n", version);
+        return EXIT_FAILURE;
+    }
+
     compute_ranges(width, height);
 
     /* Try to cram some more information into our points as long as it
@@ -1141,6 +1128,7 @@ int main(int argc, char *argv[])
         fprintf(stderr, "Available characters: %i\n", NUM_CHARACTERS);
         fprintf(stderr, "Available bits: %f\n", TOTAL_BITS);
         fprintf(stderr, "Width/Height ranges: %ix%i\n", RANGE_W, RANGE_H);
+        fprintf(stderr, "Algorithm version: %i\n", RANGE_V);
         fprintf(stderr, "Image resolution: %ix%i\n", width, height);
         fprintf(stderr, "Header bits: %f\n", HEADER_BITS);
         fprintf(stderr, "Bits available for data: %f\n", DATA_BITS);
@@ -1194,7 +1182,7 @@ int main(int argc, char *argv[])
 
             /* Choose a point at random */
             int pt = det_rand(npoints);
-            uint32_t oldval = points[pt];
+            point_t oldpt = points[pt];
 
             /* Compute the affected image zone */
             float fx, fy, fr, fg, fb, fs;
@@ -1212,67 +1200,61 @@ int main(int argc, char *argv[])
             uint8_t op1 = rand_op();
             //uint8_t op2 = rand_op();
 
-            uint32_t candidates[3];
-            double besterr = error + 1.0;
-            int bestop = -1;
-            candidates[0] = apply_op(op1, oldval);
-            //candidates[1] = apply_op(op2, oldval);
-            //candidates[2] = apply_op(op1, apply_op(op2, oldval));
+            apply_op(op1, &points[pt]);
 
-            for(int i = 0; i < 1; i++)
-            //for(int i = 0; i < 3; i++)
+#if POINTS_PER_CELL == 2
+            /* Check that two points don't fall at the same place */
+            while(points[pt].x == points[pt ^ 1].x
+                   && points[pt].y == points[pt ^ 1].y)
             {
-                if(oldval == candidates[i])
-                    continue;
+                points[pt] = oldpt;
+                op1 = rand_op();
+                apply_op(op1, &points[pt]);
+            }
+#endif
 
-                points[pt] = candidates[i];
+            render(scrap, zonex * RANGE_X, zoney * RANGE_Y,
+                   zonew * RANGE_X, zoneh * RANGE_Y, false);
 
-                render(scrap, zonex * RANGE_X, zoney * RANGE_Y,
-                       zonew * RANGE_X, zoneh * RANGE_Y, false);
-
-                double newerr = pipi_measure_rmsd(src, scrap);
-                if(newerr < besterr)
-                {
-                    besterr = newerr;
-                    bestop = i;
-                }
-            }
+            double newerr = pipi_measure_rmsd(src, scrap);
 
             opstats[op1 * 2]++;
             //opstats[op2 * 2]++;
 
-            if(besterr < error)
+            if(newerr < error)
             {
-                points[pt] = candidates[bestop];
-                /* Redraw image if the last check wasn't the best one */
-                if(bestop != 0)
-                    render(scrap, zonex * RANGE_X, zoney * RANGE_Y,
-                           zonew * RANGE_X, zoneh * RANGE_Y, false);
-
                 pipi_free(tmp);
+
+#if 0
+                /* Save image! */
+                if((success % 10) == 0)
+                {
+                    char buf[128];
+                    sprintf(buf, "twit%08i.bmp", success);
+                    tmp = pipi_new(width, height);
+                    render(tmp, 0, 0, width, height, true);
+                    pipi_save(tmp, buf);
+                    pipi_free(tmp);
+                }
+#endif
+
                 tmp = scrap;
 
                 if(DEBUG_MODE)
-                    fprintf(stderr, "%08i -.%08i %2.010g after op%i(%i)\n",
-                            iter, (int)((error - besterr) * 100000000), error,
-                            op1, pt);
+                    fprintf(stderr, "%08i -0.%010i %2.010g after op%i(%i)\n",
+                            iter, (int)((error - newerr) * 10000000000L),
+                            error, op1, pt);
 
-                error = besterr;
+                error = newerr;
                 opstats[op1 * 2 + 1]++;
                 //opstats[op2 * 2 + 1]++;
                 failures = 0;
                 success++;
-
-                /* Save image! */
-                char buf[128];
-                sprintf(buf, "twit%08i.bmp", success);
-                if((success % 10) == 0)
-                    pipi_save(tmp, buf);
             }
             else
             {
                 pipi_free(scrap);
-                points[pt] = oldval;
+                points[pt] = oldpt;
                 failures++;
             }
         }
@@ -1298,7 +1280,7 @@ int main(int argc, char *argv[])
         else
             fprintf(stderr, "\r                    \r");
 
-#if 0
+#if 1
         dst = pipi_resize(tmp, width, height);
         pipi_free(tmp);
 
@@ -1312,23 +1294,35 @@ int main(int argc, char *argv[])
         {
 #if POINTS_PER_CELL == 2
             int x1, y1, x2, y2;
-            x1 = (points[i] / RANGE_Y) % RANGE_X;
-            y1 = points[i] % RANGE_Y;
-            x2 = (points[i + 1] / RANGE_Y) % RANGE_X;
-            y2 = points[i + 1] % RANGE_Y;
+            x1 = points[i].x;
+            y1 = points[i].y;
+            x2 = points[i + 1].x;
+            y2 = points[i + 1].y;
 
             bool swap;
             uint32_t pack = pack_coords(x1, y1, x2, y2, &swap);
 
-            b.push(points[i + (swap ? 1 : 0)] / RANGE_XY, RANGE_SBGR);
-            b.push(points[i + (swap ? 0 : 1)] / RANGE_XY, RANGE_SBGR);
+            b.push(points[i + (swap ? 1 : 0)].s, RANGE_S);
+            b.push(points[i + (swap ? 1 : 0)].b, RANGE_B);
+            b.push(points[i + (swap ? 1 : 0)].g, RANGE_G);
+            b.push(points[i + (swap ? 1 : 0)].r, RANGE_R);
+            b.push(points[i + (swap ? 0 : 1)].s, RANGE_S);
+            b.push(points[i + (swap ? 0 : 1)].b, RANGE_B);
+            b.push(points[i + (swap ? 0 : 1)].g, RANGE_G);
+            b.push(points[i + (swap ? 0 : 1)].r, RANGE_R);
             b.push(pack, RANGE_XY2);
 #else
-            b.push(points[i], RANGE_SBGRXY);
+            b.push(points[i].s, RANGE_S);
+            b.push(points[i].b, RANGE_B);
+            b.push(points[i].g, RANGE_G);
+            b.push(points[i].r, RANGE_R);
+            b.push(points[i].x, RANGE_X);
+            b.push(points[i].y, RANGE_Y);
 #endif
         }
         b.push(height - 1, RANGE_H);
         b.push(width - 1, RANGE_W);
+        b.push(version, RANGE_V);
 
         /* Pop Unicode characters from the bitstream and print them */
         for(int i = 0; i < MAX_MSG_LEN; i++)
@@ -1342,15 +1336,28 @@ int main(int argc, char *argv[])
         {
 #if POINTS_PER_CELL == 2
             uint32_t pack = b.pop(RANGE_XY2);
-            uint32_t p2 = b.pop(RANGE_SBGR);
-            uint32_t p1 = b.pop(RANGE_SBGR);
             int x1, y1, x2, y2;
-
             unpack_coords(pack, &x1, &y1, &x2, &y2);
-            points[i * 2] = p1 * RANGE_XY + x1 * RANGE_Y + y1;
-            points[i * 2 + 1] = p2 * RANGE_XY + x2 * RANGE_Y + y2;
+
+            points[i * 2 + 1].y = y2;
+            points[i * 2 + 1].x = x2;
+            points[i * 2 + 1].r = b.pop(RANGE_R);
+            points[i * 2 + 1].g = b.pop(RANGE_G);
+            points[i * 2 + 1].b = b.pop(RANGE_B);
+            points[i * 2 + 1].s = b.pop(RANGE_S);
+            points[i * 2].y = y1;
+            points[i * 2].x = x1;
+            points[i * 2].r = b.pop(RANGE_R);
+            points[i * 2].g = b.pop(RANGE_G);
+            points[i * 2].b = b.pop(RANGE_B);
+            points[i * 2].s = b.pop(RANGE_S);
 #else
-            points[i] = b.pop(RANGE_SBGRXY);
+            points[i].y = b.pop(RANGE_Y);
+            points[i].x = b.pop(RANGE_X);
+            points[i].r = b.pop(RANGE_R);
+            points[i].g = b.pop(RANGE_G);
+            points[i].b = b.pop(RANGE_B);
+            points[i].s = b.pop(RANGE_S);
 #endif
         }
         npoints = dw * dh * POINTS_PER_CELL;