source: libpipi/trunk/examples/img2twit.cpp @ 3523

Last change on this file since 3523 was 3523, checked in by Sam Hocevar, 11 years ago

Make img2twit message length configurable at runtime, improve the initial
guess by smoothing hard colour gradients, and fixed the bitstack copy
constructor to avoid corruption on exit.

  • Property svn:keywords set to Id
File size: 34.0 KB
Line 
1/*
2 *  img2twit      Image to short text message encoder/decoder
3 *  Copyright (c) 2009 Sam Hocevar <sam@hocevar.net>
4 *                All Rights Reserved
5 *
6 *  This program is free software. It comes without any warranty, to
7 *  the extent permitted by applicable law. You can redistribute it
8 *  and/or modify it under the terms of the Do What The Fuck You Want
9 *  To Public License, Version 2, as published by Sam Hocevar. See
10 *  http://sam.zoy.org/wtfpl/COPYING for more details.
11 */
12
13#include "config.h"
14
15#include <stdio.h>
16#include <stdlib.h>
17#include <string.h>
18#include <math.h>
19
20#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
21#include <CGAL/Delaunay_triangulation_2.h>
22#include <CGAL/natural_neighbor_coordinates_2.h>
23
24#include <pipi.h>
25
26#include "../genethumb/mygetopt.h"
27
28/*
29 * User-definable settings.
30 */
31
32/* The Unicode characters at disposal - XXX: must be _ordered_ */
33static const uint32_t unichars[] =
34{
35    /* Printable ASCII (except space) */
36    //0x0021, 0x007f,
37
38    /* Stupid symbols and Dingbats shit */
39    //0x25a0, 0x2600, /* Geometric Shapes */
40    //0x2600, 0x269e, 0x26a0, 0x26bd, 0x26c0, 0x26c4, /* Misc. Symbols */
41    //0x2701, 0x2705, 0x2706, 0x270a, 0x270c, 0x2728, 0x2729, 0x274c,
42    //  0x274d, 0x274e, 0x274f, 0x2753, 0x2756, 0x2757, 0x2758, 0x275f,
43    //  0x2761, 0x2795, 0x2798, 0x27b0, 0x27b1, 0x27bf, /* Dingbats */
44
45    /* Chinese-looking stuff */
46    //0x2e80, 0x2e9a, 0x2e9b, 0x2ef4, /* CJK Radicals Supplement */
47    //0x2f00, 0x2fd6, /* Kangxi Radicals */
48    //0x3400, 0x4db6, /* CJK Unified Ideographs Extension A */
49    0x4e00, 0x9fa6, /* CJK Unified Ideographs */
50
51    /* Korean - most people don't know the difference anyway */
52    //0xac00, 0xd7a4, /* Hangul Syllables */
53
54    /* More Chinese */
55    //0xf900, 0xfa2e, 0xfa30, 0xfa6b, 0xfa70, 0xfada, /* CJK Compat. Idgphs. */
56
57    /* TODO: there's also the U+20000 and U+2f800 planes, but they're
58     * not supported by the Twitter Javascript filter (yet?). */
59
60    /* End of list marker - XXX: don't remove! */
61    0x0000, 0x0000
62};
63
64/* The maximum image size we want to support */
65#define MAX_W 4000
66#define MAX_H 4000
67
68/* How does the algorithm work: one point per cell, or two */
69#define POINTS_PER_CELL 2
70
71/*
72 * These values can be overwritten at runtime
73 */
74
75/* Debug mode */
76static bool DEBUG_MODE = false;
77
78/* The maximum message length */
79static int MAX_MSG_LEN = 140;
80
81/* Iterations per point -- larger means slower but nicer */
82static int ITERATIONS_PER_POINT = 50;
83
84/* The range value for point parameters: X Y, red/green/blue, "strength"
85 * Tested values (on Mona Lisa) are:
86 *  16 16 5 5 5 2 -> 0.06511725914
87 *  16 16 6 7 6 1 -> 0.05731491348 *
88 *  16 16 7 6 6 1 -> 0.06450513783
89 *  14 14 7 7 6 1 -> 0.0637207893
90 *  19 19 6 6 5 1 -> 0.06801999094 */
91static unsigned int RANGE_X = 16;
92static unsigned int RANGE_Y = 16;
93static unsigned int RANGE_R = 6;
94static unsigned int RANGE_G = 6;
95static unsigned int RANGE_B = 6;
96static unsigned int RANGE_S = 1;
97
98/*
99 * These values are computed at runtime
100 */
101
102static float TOTAL_BITS;
103static float HEADER_BITS;
104static float DATA_BITS;
105static float CELL_BITS;
106
107static int NUM_CHARACTERS;
108static int MAX_ITERATIONS;
109static unsigned int TOTAL_CELLS;
110
111#define RANGE_SY (RANGE_S*RANGE_Y)
112#define RANGE_SYX (RANGE_S*RANGE_Y*RANGE_X)
113#define RANGE_SYXR (RANGE_S*RANGE_Y*RANGE_X*RANGE_R)
114#define RANGE_SYXRG (RANGE_S*RANGE_Y*RANGE_X*RANGE_R*RANGE_G)
115#define RANGE_SYXRGB (RANGE_S*RANGE_Y*RANGE_X*RANGE_R*RANGE_G*RANGE_B)
116
117struct K : CGAL::Exact_predicates_inexact_constructions_kernel {};
118typedef CGAL::Delaunay_triangulation_2<K> Delaunay_triangulation;
119typedef std::vector<std::pair<K::Point_2, K::FT> > Point_coordinate_vector;
120
121/* Global aspect ratio */
122static unsigned int dw, dh;
123
124/* Global point encoding */
125static uint32_t points[4096]; /* FIXME: allocate this dynamically */
126static int npoints = 0;
127
128/* Global triangulation */
129static Delaunay_triangulation dt;
130
131/*
132 * Unicode stuff handling
133 */
134
135/* Return the number of chars in the unichars table */
136static int count_unichars(void)
137{
138    int ret = 0;
139
140    for(int u = 0; unichars[u] != unichars[u + 1]; u += 2)
141        ret += unichars[u + 1] - unichars[u];
142
143    return ret;
144}
145
146/* Get the ith Unicode character in our list */
147static uint32_t index2uni(uint32_t i)
148{
149    for(int u = 0; unichars[u] != unichars[u + 1]; u += 2)
150        if(i < unichars[u + 1] - unichars[u])
151            return unichars[u] + i;
152        else
153            i -= unichars[u + 1] - unichars[u];
154
155    return 0; /* Should not happen! */
156}
157
158/* Convert a Unicode character to its position in the compact list */
159static uint32_t uni2index(uint32_t x)
160{
161    uint32_t ret = 0;
162
163    for(int u = 0; unichars[u] != unichars[u + 1]; u += 2)
164        if(x < unichars[u + 1])
165            return ret + x - unichars[u];
166        else
167            ret += unichars[u + 1] - unichars[u];
168
169    return ret; /* Should not happen! */
170}
171
172static uint8_t const utf8_trailing[256] =
173{
174    0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
175    0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
176    0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
177    0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
178    0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
179    0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
180    1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1, 1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,
181    2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2, 3,3,3,3,3,3,3,3,4,4,4,4,5,5,5,5
182};
183
184static uint32_t const utf8_offsets[6] =
185{
186    0x00000000UL, 0x00003080UL, 0x000E2080UL,
187    0x03C82080UL, 0xFA082080UL, 0x82082080UL
188};
189
190static uint32_t fread_utf8(FILE *f)
191{
192    int ch, i = 0, todo = -1;
193    uint32_t ret = 0;
194
195    for(;;)
196    {
197        ch = fgetc(f);
198        if(!ch)
199            return 0;
200        if(todo == -1)
201            todo = utf8_trailing[ch];
202        ret += ((uint32_t)ch) << (6 * (todo - i));
203        if(todo == i++)
204            return ret - utf8_offsets[todo];
205    }
206}
207
208static void fwrite_utf8(FILE *f, uint32_t x)
209{
210    static const uint8_t mark[7] =
211    {
212        0x00, 0x00, 0xC0, 0xE0, 0xF0, 0xF8, 0xFC
213    };
214
215    char buf[8];
216    char *parser = buf;
217    size_t bytes;
218
219    if(x < 0x80)
220    {
221        fprintf(f, "%c", x);
222        return;
223    }
224
225    bytes = (x < 0x800) ? 2 : (x < 0x10000) ? 3 : 4;
226    parser += bytes;
227    *parser = '\0';
228
229    switch(bytes)
230    {
231        case 4: *--parser = (x | 0x80) & 0xbf; x >>= 6;
232        case 3: *--parser = (x | 0x80) & 0xbf; x >>= 6;
233        case 2: *--parser = (x | 0x80) & 0xbf; x >>= 6;
234    }
235    *--parser = x | mark[bytes];
236
237    fprintf(f, "%s", buf);
238}
239
240/*
241 * Our nifty non-power-of-two bitstack handling
242 */
243
244class bitstack
245{
246public:
247    bitstack() { alloc(); init(0); }
248
249    ~bitstack() { delete[] digits; delete[] str; }
250
251    char const *tostring()
252    {
253        int pos = sprintf(str, "0x%x", digits[msb]);
254
255        for(int i = msb - 1; i >= 0; i--)
256            pos += sprintf(str + pos, "%08x", digits[i]);
257
258        return str;
259    }
260
261    void push(uint32_t val, uint32_t range)
262    {
263        if(!range)
264            return;
265
266        mul(range);
267        add(val % range);
268    }
269
270    uint32_t pop(uint32_t range)
271    {
272        if(!range)
273            return 0;
274
275        return div(range);
276    }
277
278    bool isempty()
279    {
280        for(int i = msb; i >= 0; i--)
281            if(digits[i])
282                return false;
283
284        return true;
285    }
286
287private:
288    bitstack(uint32_t i) { alloc(); init(i); }
289
290    bitstack(bitstack &b)
291    {
292        alloc();
293        msb = b.msb;
294        memcpy(digits, b.digits, (MAX_MSG_LEN + 1) * sizeof(uint32_t));
295    }
296
297    bitstack(bitstack const &b)
298    {
299        alloc();
300        msb = b.msb;
301        memcpy(digits, b.digits, (MAX_MSG_LEN + 1) * sizeof(uint32_t));
302    }
303
304    void alloc()
305    {
306        digits = new uint32_t[MAX_MSG_LEN + 1];
307        str = new char[(MAX_MSG_LEN + 1) * 8 + 1];
308    }
309
310    void init(uint32_t i)
311    {
312        msb = 0;
313        memset(digits, 0, (MAX_MSG_LEN + 1) * sizeof(uint32_t));
314        digits[0] = i;
315    }
316
317    /* Could be done much faster, but we don't care! */
318    void add(uint32_t x) { add(bitstack(x)); }
319    void sub(uint32_t x) { sub(bitstack(x)); }
320
321    void add(bitstack const &_b)
322    {
323        /* Copy the operand in case we get added to ourselves */
324        bitstack b(_b);
325        uint64_t x = 0;
326
327        if(msb < b.msb)
328            msb = b.msb;
329
330        for(int i = 0; i <= msb; i++)
331        {
332            uint64_t tmp = (uint64_t)digits[i] + (uint64_t)b.digits[i] + x;
333            digits[i] = tmp;
334            if((uint64_t)digits[i] == tmp)
335                x = 0;
336            else
337            {
338                x = 1;
339                if(i == msb)
340                    msb++;
341            }
342        }
343    }
344
345    void sub(bitstack const &_b)
346    {
347        /* Copy the operand in case we get substracted from ourselves */
348        bitstack b(_b);
349        uint64_t x = 0;
350
351        /* We cannot substract a larger number! */
352        if(msb < b.msb)
353        {
354            init(0);
355            return;
356        }
357
358        for(int i = 0; i <= msb; i++)
359        {
360            uint64_t tmp = (uint64_t)digits[i] - (uint64_t)b.digits[i] - x;
361            digits[i] = tmp;
362            if((uint64_t)digits[i] == tmp)
363                x = 0;
364            else
365            {
366                x = 1;
367                if(i == msb)
368                {
369                    /* Error: carry into MSB! */
370                    init(0);
371                    return;
372                }
373            }
374        }
375
376        while(msb > 0 && digits[msb] == 0) msb--;
377    }
378
379    void mul(uint32_t x)
380    {
381        bitstack b(*this);
382        init(0);
383
384        while(x)
385        {
386            if(x & 1)
387                add(b);
388            x /= 2;
389            b.add(b);
390        }
391    }
392
393    uint32_t div(uint32_t x)
394    {
395        bitstack b(*this);
396
397        for(int i = msb; i >= 0; i--)
398        {
399            uint64_t tmp = b.digits[i] + (((uint64_t)b.digits[i + 1]) << 32);
400            uint32_t res = tmp / x;
401            uint32_t rem = tmp % x;
402            digits[i]= res;
403            b.digits[i + 1] = 0;
404            b.digits[i] = rem;
405        }
406
407        while(msb > 0 && digits[msb] == 0) msb--;
408
409        return b.digits[0];
410    }
411
412    int msb;
413    uint32_t *digits;
414    char *str;
415};
416
417/*
418 * Point handling
419 */
420
421static unsigned int det_rand(unsigned int mod)
422{
423    static unsigned long next = 1;
424    next = next * 1103515245 + 12345;
425    return ((unsigned)(next / 65536) % 32768) % mod;
426}
427
428static inline int range2int(float val, int range)
429{
430    int ret = (int)(val * ((float)range - 0.0001));
431    return ret < 0 ? 0 : ret > range - 1 ? range - 1 : ret;
432}
433
434static inline float int2midrange(int val, int range)
435{
436    return (float)(1 + 2 * val) / (float)(2 * range);
437}
438
439static inline float int2fullrange(int val, int range)
440{
441    return range > 1 ? (float)val / (float)(range - 1) : 0.0;
442}
443
444static inline void set_point(int index, float x, float y, float r,
445                             float g, float b, float s)
446{
447    int dx = (index / POINTS_PER_CELL) % dw;
448    int dy = (index / POINTS_PER_CELL) / dw;
449
450    float fx = (x - dx * RANGE_X) / RANGE_X;
451    float fy = (y - dy * RANGE_Y) / RANGE_Y;
452
453    int is = range2int(s, RANGE_S);
454
455    int ix = range2int(fx, RANGE_X);
456    int iy = range2int(fy, RANGE_Y);
457
458    int ir = range2int(r, RANGE_R);
459    int ig = range2int(g, RANGE_G);
460    int ib = range2int(b, RANGE_B);
461
462    points[index] = is + RANGE_S * (iy + RANGE_Y * (ix + RANGE_X *
463                               (ib + RANGE_B * (ig + (RANGE_R * ir)))));
464}
465
466static inline void get_point(int index, float *x, float *y, float *r,
467                             float *g, float *b, float *s)
468{
469    uint32_t pt = points[index];
470
471    unsigned int dx = (index / POINTS_PER_CELL) % dw;
472    unsigned int dy = (index / POINTS_PER_CELL) / dw;
473
474    *s = int2fullrange(pt % RANGE_S, RANGE_S); pt /= RANGE_S;
475
476    float fy = int2midrange(pt % RANGE_Y, RANGE_Y); pt /= RANGE_Y;
477    float fx = int2midrange(pt % RANGE_X, RANGE_X); pt /= RANGE_X;
478
479    *x = (fx + dx) * RANGE_X /*+ 0.5 * (index & 1)*/;
480    *y = (fy + dy) * RANGE_Y /*+ 0.5 * (index & 1)*/;
481
482    *b = int2midrange(pt % RANGE_R, RANGE_R); pt /= RANGE_R;
483    *g = int2midrange(pt % RANGE_G, RANGE_G); pt /= RANGE_G;
484    *r = int2midrange(pt % RANGE_B, RANGE_B); pt /= RANGE_B;
485}
486
487static inline float clip(float x, int modulo)
488{
489    float mul = (float)modulo + 0.9999;
490    int round = (int)(x * mul);
491    return (float)round / (float)modulo;
492}
493
494static void add_point(float x, float y, float r, float g, float b, float s)
495{
496    set_point(npoints, x, y, r, g, b, s);
497    npoints++;
498}
499
500#if 0
501static void add_random_point()
502{
503    points[npoints] = det_rand(RANGE_SYXRGB);
504    npoints++;
505}
506#endif
507
508#define NB_OPS 20
509
510static uint8_t rand_op(void)
511{
512    uint8_t x = det_rand(NB_OPS);
513
514    /* Randomly ignore statistically less efficient ops */
515    if(x == 0)
516        return rand_op();
517    if(x == 1 && (RANGE_S == 1 || det_rand(2)))
518        return rand_op();
519    if(x <= 5 && det_rand(2))
520        return rand_op();
521    //if((x < 10 || x > 15) && !det_rand(4)) /* Favour colour changes */
522    //    return rand_op();
523
524    return x;
525}
526
527static uint32_t apply_op(uint8_t op, uint32_t val)
528{
529    uint32_t rem, ext;
530
531    switch(op)
532    {
533    case 0: /* Flip strength value */
534    case 1:
535        /* Statistics show that this helps often, but does not reduce
536         * the error significantly. */
537        return val ^ 1;
538    case 2: /* Move up; if impossible, down */
539        rem = val % RANGE_S;
540        ext = (val / RANGE_S) % RANGE_Y;
541        ext = ext > 0 ? ext - 1 : ext + 1;
542        return (val / RANGE_SY * RANGE_Y + ext) * RANGE_S + rem;
543    case 3: /* Move down; if impossible, up */
544        rem = val % RANGE_S;
545        ext = (val / RANGE_S) % RANGE_Y;
546        ext = ext < RANGE_Y - 1 ? ext + 1 : ext - 1;
547        return (val / RANGE_SY * RANGE_Y + ext) * RANGE_S + rem;
548    case 4: /* Move left; if impossible, right */
549        rem = val % RANGE_SY;
550        ext = (val / RANGE_SY) % RANGE_X;
551        ext = ext > 0 ? ext - 1 : ext + 1;
552        return (val / RANGE_SYX * RANGE_X + ext) * RANGE_SY + rem;
553    case 5: /* Move left; if impossible, right */
554        rem = val % RANGE_SY;
555        ext = (val / RANGE_SY) % RANGE_X;
556        ext = ext < RANGE_X - 1 ? ext + 1 : ext - 1;
557        return (val / RANGE_SYX * RANGE_X + ext) * RANGE_SY + rem;
558    case 6: /* Corner 1 */
559        return apply_op(2, apply_op(4, val));
560    case 7: /* Corner 2 */
561        return apply_op(2, apply_op(5, val));
562    case 8: /* Corner 3 */
563        return apply_op(3, apply_op(5, val));
564    case 9: /* Corner 4 */
565        return apply_op(3, apply_op(4, val));
566    case 16: /* Double up */
567        return apply_op(2, apply_op(2, val));
568    case 17: /* Double down */
569        return apply_op(3, apply_op(3, val));
570    case 18: /* Double left */
571        return apply_op(4, apply_op(4, val));
572    case 19: /* Double right */
573        return apply_op(5, apply_op(5, val));
574    case 10: /* R-- (or R++) */
575        rem = val % RANGE_SYX;
576        ext = (val / RANGE_SYX) % RANGE_R;
577        ext = ext > 0 ? ext - 1 : ext + 1;
578        return (val / RANGE_SYXR * RANGE_R + ext) * RANGE_SYX + rem;
579    case 11: /* R++ (or R--) */
580        rem = val % RANGE_SYX;
581        ext = (val / RANGE_SYX) % RANGE_R;
582        ext = ext < RANGE_R - 1 ? ext + 1 : ext - 1;
583        return (val / RANGE_SYXR * RANGE_R + ext) * RANGE_SYX + rem;
584    case 12: /* G-- (or G++) */
585        rem = val % RANGE_SYXR;
586        ext = (val / RANGE_SYXR) % RANGE_G;
587        ext = ext > 0 ? ext - 1 : ext + 1;
588        return (val / RANGE_SYXRG * RANGE_G + ext) * RANGE_SYXR + rem;
589    case 13: /* G++ (or G--) */
590        rem = val % RANGE_SYXR;
591        ext = (val / RANGE_SYXR) % RANGE_G;
592        ext = ext < RANGE_G - 1 ? ext + 1 : ext - 1;
593        return (val / RANGE_SYXRG * RANGE_G + ext) * RANGE_SYXR + rem;
594    case 14: /* B-- (or B++) */
595        rem = val % RANGE_SYXRG;
596        ext = (val / RANGE_SYXRG) % RANGE_B;
597        ext = ext > 0 ? ext - 1 : ext + 1;
598        return ext * RANGE_SYXRG + rem;
599    case 15: /* B++ (or B--) */
600        rem = val % RANGE_SYXRG;
601        ext = (val / RANGE_SYXRG) % RANGE_B;
602        ext = ext < RANGE_B - 1 ? ext + 1 : ext - 1;
603        return ext * RANGE_SYXRG + rem;
604#if 0
605    case 15: /* Brightness-- */
606        return apply_op(9, apply_op(11, apply_op(13, val)));
607    case 16: /* Brightness++ */
608        return apply_op(10, apply_op(12, apply_op(14, val)));
609    case 17: /* RG-- */
610        return apply_op(9, apply_op(11, val));
611    case 18: /* RG++ */
612        return apply_op(10, apply_op(12, val));
613    case 19: /* GB-- */
614        return apply_op(11, apply_op(13, val));
615    case 20: /* GB++ */
616        return apply_op(12, apply_op(14, val));
617    case 21: /* RB-- */
618        return apply_op(9, apply_op(13, val));
619    case 22: /* RB++ */
620        return apply_op(10, apply_op(14, val));
621#endif
622    default:
623        return val;
624    }
625}
626
627static void render(pipi_image_t *dst, int rx, int ry, int rw, int rh)
628{
629    int lookup[dw * RANGE_X * 2 * dh * RANGE_Y * 2];
630    pipi_pixels_t *p = pipi_get_pixels(dst, PIPI_PIXELS_RGBA_F32);
631    float *data = (float *)p->pixels;
632    int x, y;
633
634    memset(lookup, 0, sizeof(lookup));
635    dt.clear();
636    for(int i = 0; i < npoints; i++)
637    {
638        float fx, fy, fr, fg, fb, fs;
639        get_point(i, &fx, &fy, &fr, &fg, &fb, &fs);
640        dt.insert(K::Point_2(fx, fy));
641        /* Keep link to point */
642        lookup[(int)(fx * 2) + dw * RANGE_X * 2 * (int)(fy * 2)] = i;
643    }
644
645    /* Add fake points to close the triangulation */
646    dt.insert(K::Point_2(-p->w, -p->h));
647    dt.insert(K::Point_2(2 * p->w, -p->h));
648    dt.insert(K::Point_2(-p->w, 2 * p->h));
649    dt.insert(K::Point_2(2 * p->w, 2 * p->h));
650
651    for(y = ry; y < ry + rh; y++)
652    {
653        for(x = rx; x < rx + rw; x++)
654        {
655            K::Point_2 m(x, y);
656            Point_coordinate_vector coords;
657            CGAL::Triple<
658              std::back_insert_iterator<Point_coordinate_vector>,
659              K::FT, bool> result =
660              CGAL::natural_neighbor_coordinates_2(dt, m,
661                                                   std::back_inserter(coords));
662
663            float r = 0.0f, g = 0.0f, b = 0.0f, norm = 0.000000000000001f;
664
665            Point_coordinate_vector::iterator it;
666            for(it = coords.begin(); it != coords.end(); ++it)
667            {
668                float fx, fy, fr, fg, fb, fs;
669
670                fx = (*it).first.x();
671                fy = (*it).first.y();
672
673                if(fx < 0 || fy < 0 || fx > p->w - 1 || fy > p->h - 1)
674                    continue;
675
676                int index = lookup[(int)(fx * 2)
677                                    + dw * RANGE_X * 2 * (int)(fy * 2)];
678
679                get_point(index, &fx, &fy, &fr, &fg, &fb, &fs);
680
681                //float k = pow((*it).second * (1.0 + fs), 1.2);
682                float k = (*it).second * (1.00f + fs);
683                //float k = (*it).second * (0.60f + fs);
684                //float k = pow((*it).second, (1.0f + fs));
685
686                r += k * fr;
687                g += k * fg;
688                b += k * fb;
689                norm += k;
690            }
691
692            data[4 * (x + y * p->w) + 0] = r / norm;
693            data[4 * (x + y * p->w) + 1] = g / norm;
694            data[4 * (x + y * p->w) + 2] = b / norm;
695            data[4 * (x + y * p->w) + 3] = 0.0;
696        }
697    }
698
699    pipi_release_pixels(dst, p);
700}
701
702static void analyse(pipi_image_t *src)
703{
704    pipi_pixels_t *p = pipi_get_pixels(src, PIPI_PIXELS_RGBA_F32);
705    float *data = (float *)p->pixels;
706
707    for(unsigned int dy = 0; dy < dh; dy++)
708        for(unsigned int dx = 0; dx < dw; dx++)
709        {
710            float min = 1.1f, max = -0.1f, mr = 0.0f, mg = 0.0f, mb = 0.0f;
711            float total = 0.0;
712            int xmin = 0, xmax = 0, ymin = 0, ymax = 0;
713            int npixels = 0;
714
715            for(unsigned int iy = RANGE_Y * dy; iy < RANGE_Y * (dy + 1); iy++)
716                for(unsigned int ix = RANGE_X * dx; ix < RANGE_X * (dx + 1); ix++)
717                {
718                    float lum = 0.0f;
719
720                    lum += data[4 * (ix + iy * p->w) + 0];
721                    lum += data[4 * (ix + iy * p->w) + 1];
722                    lum += data[4 * (ix + iy * p->w) + 2];
723                    lum /= 3;
724
725                    mr += data[4 * (ix + iy * p->w) + 0];
726                    mg += data[4 * (ix + iy * p->w) + 1];
727                    mb += data[4 * (ix + iy * p->w) + 2];
728
729                    if(lum < min)
730                    {
731                        min = lum;
732                        xmin = ix;
733                        ymin = iy;
734                    }
735
736                    if(lum > max)
737                    {
738                        max = lum;
739                        xmax = ix;
740                        ymax = iy;
741                    }
742
743                    total += lum;
744                    npixels++;
745                }
746
747            total /= npixels;
748            mr /= npixels;
749            mg /= npixels;
750            mb /= npixels;
751
752            float wmin, wmax;
753
754            if(total < min + (max - min) / 4)
755                wmin = 1.0, wmax = 0.0;
756            else if(total < min + (max - min) / 4 * 3)
757                wmin = 0.0, wmax = 0.0;
758            else
759                wmin = 0.0, wmax = 1.0;
760
761#if 0
762            add_random_point();
763            add_random_point();
764#else
765            /* 0.80 and 0.20 were chosen empirically, it gives a 10% better
766             * initial distance. Definitely worth it. */
767#if POINTS_PER_CELL == 1
768            if(total < min + (max - min) / 2)
769            {
770#endif
771            add_point(xmin, ymin,
772                      data[4 * (xmin + ymin * p->w) + 0] * 0.80 + mr * 0.20,
773                      data[4 * (xmin + ymin * p->w) + 1] * 0.80 + mg * 0.20,
774                      data[4 * (xmin + ymin * p->w) + 2] * 0.80 + mb * 0.20,
775                      wmin);
776#if POINTS_PER_CELL == 1
777            }
778            else
779            {
780#endif
781            add_point(xmax, ymax,
782                      data[4 * (xmax + ymax * p->w) + 0] * 0.80 + mr * 0.20,
783                      data[4 * (xmax + ymax * p->w) + 1] * 0.80 + mg * 0.20,
784                      data[4 * (xmax + ymax * p->w) + 2] * 0.80 + mb * 0.20,
785                      wmax);
786#if POINTS_PER_CELL == 1
787            }
788#endif
789#endif
790        }
791}
792
793#define MOREINFO "Try `%s --help' for more information.\n"
794
795int main(int argc, char *argv[])
796{
797    int opstats[2 * NB_OPS];
798    char const *srcname = NULL, *dstname = NULL;
799    pipi_image_t *src, *tmp, *dst;
800    double error = 1.0;
801    int width, height, ret = 0;
802
803    /* Parse command-line options */
804    for(;;)
805    {
806        int option_index = 0;
807        static struct myoption long_options[] =
808        {
809            { "output",      1, NULL, 'o' },
810            { "length",      1, NULL, 'l' },
811            { "quality",     1, NULL, 'q' },
812            { "debug",       0, NULL, 'd' },
813            { "help",        0, NULL, 'h' },
814            { NULL,          0, NULL, 0   },
815        };
816        int c = mygetopt(argc, argv, "o:l:q:dh", long_options, &option_index);
817
818        if(c == -1)
819            break;
820
821        switch(c)
822        {
823        case 'o':
824            dstname = myoptarg;
825            break;
826        case 'l':
827            MAX_MSG_LEN = atoi(myoptarg);
828            if(MAX_MSG_LEN < 16)
829            {
830                fprintf(stderr, "Warning: rounding minimum message length to 16\n");
831                MAX_MSG_LEN = 16;
832            }
833            break;
834        case 'q':
835            ITERATIONS_PER_POINT = 10 * atoi(myoptarg);
836            if(ITERATIONS_PER_POINT < 0)
837                ITERATIONS_PER_POINT = 0;
838            else if(ITERATIONS_PER_POINT > 100)
839                ITERATIONS_PER_POINT = 100;
840            break;
841        case 'd':
842            DEBUG_MODE = true;
843            break;
844        case 'h':
845            printf("Usage: img2twit [OPTIONS] SOURCE\n");
846            printf("       img2twit [OPTIONS] -o DESTINATION\n");
847            printf("Encode SOURCE image to stdout or decode stdin to DESTINATION.\n");
848            printf("\n");
849            printf("Mandatory arguments to long options are mandatory for short options too.\n");
850            printf("  -o, --output <filename>   output resulting image to filename\n");
851            printf("  -l, --length <size>       message length in characters (default 140)\n");
852            printf("  -q, --quality <rate>      set image quality (0 - 10) (default 5)\n");
853            printf("  -d, --debug               print debug information\n");
854            printf("  -h, --help                display this help and exit\n");
855            printf("\n");
856            printf("Written by Sam Hocevar. Report bugs to <sam@hocevar.net>.\n");
857            return EXIT_SUCCESS;
858        default:
859            fprintf(stderr, "%s: invalid option -- %c\n", argv[0], c);
860            printf(MOREINFO, argv[0]);
861            return EXIT_FAILURE;
862        }
863    }
864
865    if(myoptind == argc && !dstname)
866    {
867        fprintf(stderr, "%s: too few arguments\n", argv[0]);
868        printf(MOREINFO, argv[0]);
869        return EXIT_FAILURE;
870    }
871
872    if((myoptind == argc - 1 && dstname) || myoptind < argc - 1)
873    {
874        fprintf(stderr, "%s: too many arguments\n", argv[0]);
875        printf(MOREINFO, argv[0]);
876        return EXIT_FAILURE;
877    }
878
879    if(myoptind == argc - 1)
880        srcname = argv[myoptind];
881
882    pipi_set_gamma(1.0);
883
884    /* Precompute bit allocation */
885    NUM_CHARACTERS = count_unichars();
886    TOTAL_BITS = MAX_MSG_LEN * logf(NUM_CHARACTERS) / logf(2);
887    HEADER_BITS = logf(MAX_W * MAX_H) / logf(2);
888    DATA_BITS = TOTAL_BITS - HEADER_BITS;
889#if POINTS_PER_CELL == 1
890    CELL_BITS = logf(RANGE_SYXRGB) / logf(2);
891#else
892    // TODO: implement the following shit
893    //float coord_bits = logf((RANGE_Y * RANGE_X) * (RANGE_Y * RANGE_X + 1) / 2);
894    //float other_bits = logf(RANGE_R * RANGE_G * RANGE_B * RANGE_S);
895    //CELL_BITS = (coord_bits + 2 * other_bits) / logf(2);
896    CELL_BITS = 2 * logf(RANGE_SYXRGB) / logf(2);
897#endif
898    TOTAL_CELLS = (int)(DATA_BITS / CELL_BITS);
899    MAX_ITERATIONS = ITERATIONS_PER_POINT * POINTS_PER_CELL * TOTAL_CELLS;
900
901    bitstack b; /* We cannot declare this before, because MAX_MSG_LEN
902                 * wouldn't be defined. */
903
904    if(dstname)
905    {
906        /* Decoding mode: read UTF-8 text from stdin, find each
907         * character's index in our character list, and push it to our
908         * wonderful custom bitstream. */
909        uint32_t data[MAX_MSG_LEN];
910        for(int i = 0; i < MAX_MSG_LEN; i++)
911            data[i] = uni2index(fread_utf8(stdin));
912        for(int i = MAX_MSG_LEN; i--; )
913            b.push(data[i], NUM_CHARACTERS);
914
915        /* Read width and height from bitstream */
916        src = NULL;
917        width = b.pop(MAX_W);
918        height = b.pop(MAX_H);
919    }
920    else
921    {
922        /* Argument given: open image for encoding */
923        src = pipi_load(srcname);
924
925        if(!src)
926        {
927            fprintf(stderr, "Error loading %s\n", srcname);
928            return EXIT_FAILURE;
929        }
930
931        width = pipi_get_image_width(src);
932        height = pipi_get_image_height(src);
933    }
934
935    /* Compute best w/h ratio */
936    dw = 1; dh = TOTAL_CELLS;
937    for(unsigned int i = 1; i <= TOTAL_CELLS; i++)
938    {
939        int j = TOTAL_CELLS / i;
940
941        float r = (float)width / (float)height;
942        float ir = (float)i / (float)j;
943        float dwr = (float)dw / (float)dh;
944
945        if(fabs(logf(r / ir)) < fabs(logf(r / dwr)))
946        {
947            dw = i;
948            dh = TOTAL_CELLS / dw;
949        }
950    }
951    while((dh + 1) * dw <= TOTAL_CELLS) dh++;
952    while(dw * (dh + 1) <= TOTAL_CELLS) dw++;
953
954    /* Print debug information */
955    if(DEBUG_MODE)
956    {
957        fprintf(stderr, "Maximum message size: %i\n", MAX_MSG_LEN);
958        fprintf(stderr, "Available characters: %i\n", NUM_CHARACTERS);
959        fprintf(stderr, "Available bits: %f\n", TOTAL_BITS);
960        fprintf(stderr, "Maximum image resolution: %ix%i\n", MAX_W, MAX_H);
961        fprintf(stderr, "Image resolution: %ix%i\n", width, height);
962        fprintf(stderr, "Header bits: %f\n", HEADER_BITS);
963        fprintf(stderr, "Bits available for data: %f\n", DATA_BITS);
964        fprintf(stderr, "Cell bits: %f\n", CELL_BITS);
965        fprintf(stderr, "Available cells: %i\n", TOTAL_CELLS);
966        fprintf(stderr, "Wasted bits: %f\n",
967                DATA_BITS - CELL_BITS * TOTAL_CELLS);
968        fprintf(stderr, "Chosen image ratio: %i:%i (wasting %i point cells)\n",
969                dw, dh, TOTAL_CELLS - dw * dh);
970        fprintf(stderr, "Total wasted bits: %f\n",
971                DATA_BITS - CELL_BITS * dw * dh);
972    }
973
974    if(srcname)
975    {
976        /* Resize and filter image to better state */
977        tmp = pipi_resize(src, dw * RANGE_X, dh * RANGE_Y);
978        pipi_free(src);
979        src = pipi_median_ext(tmp, 1, 1);
980        pipi_free(tmp);
981
982        /* Analyse image */
983        analyse(src);
984
985        /* Render what we just computed */
986        tmp = pipi_new(dw * RANGE_X, dh * RANGE_Y);
987        render(tmp, 0, 0, dw * RANGE_X, dh * RANGE_Y);
988        error = pipi_measure_rmsd(src, tmp);
989
990        if(DEBUG_MODE)
991            fprintf(stderr, "Initial distance: %2.10g\n", error);
992
993        memset(opstats, 0, sizeof(opstats));
994        for(int iter = 0, stuck = 0, failures = 0, success = 0;
995            iter < MAX_ITERATIONS /* && stuck < 5 && */;
996            iter++)
997        {
998            if(failures > 500)
999            {
1000                stuck++;
1001                failures = 0;
1002            }
1003
1004            if(!DEBUG_MODE && !(iter % 16))
1005                fprintf(stderr, "\rEncoding... %i%%",
1006                        iter * 100 / MAX_ITERATIONS);
1007
1008            pipi_image_t *scrap = pipi_copy(tmp);
1009
1010            /* Choose a point at random */
1011            int pt = det_rand(npoints);
1012            uint32_t oldval = points[pt];
1013
1014            /* Compute the affected image zone */
1015            float fx, fy, fr, fg, fb, fs;
1016            get_point(pt, &fx, &fy, &fr, &fg, &fb, &fs);
1017            int zonex = (int)fx / RANGE_X - 1;
1018            int zoney = (int)fy / RANGE_Y - 1;
1019            int zonew = 3;
1020            int zoneh = 3;
1021            if(zonex < 0) { zonex = 0; zonew--; }
1022            if(zoney < 0) { zoney = 0; zoneh--; }
1023            if(zonex + zonew >= (int)dw) { zonew--; }
1024            if(zoney + zoneh >= (int)dh) { zoneh--; }
1025
1026            /* Choose random operations and measure their effect */
1027            uint8_t op1 = rand_op();
1028            //uint8_t op2 = rand_op();
1029
1030            uint32_t candidates[3];
1031            double besterr = error + 1.0;
1032            int bestop = -1;
1033            candidates[0] = apply_op(op1, oldval);
1034            //candidates[1] = apply_op(op2, oldval);
1035            //candidates[2] = apply_op(op1, apply_op(op2, oldval));
1036
1037            for(int i = 0; i < 1; i++)
1038            //for(int i = 0; i < 3; i++)
1039            {
1040                if(oldval == candidates[i])
1041                    continue;
1042
1043                points[pt] = candidates[i];
1044
1045                render(scrap, zonex * RANGE_X, zoney * RANGE_Y,
1046                       zonew * RANGE_X, zoneh * RANGE_Y);
1047
1048                double newerr = pipi_measure_rmsd(src, scrap);
1049                if(newerr < besterr)
1050                {
1051                    besterr = newerr;
1052                    bestop = i;
1053                }
1054            }
1055
1056            opstats[op1 * 2]++;
1057            //opstats[op2 * 2]++;
1058
1059            if(besterr < error)
1060            {
1061                points[pt] = candidates[bestop];
1062                /* Redraw image if the last check wasn't the best one */
1063                if(bestop != 2)
1064                    render(scrap, zonex * RANGE_X, zoney * RANGE_Y,
1065                           zonew * RANGE_X, zoneh * RANGE_Y);
1066
1067                pipi_free(tmp);
1068                tmp = scrap;
1069
1070                if(DEBUG_MODE)
1071                    fprintf(stderr, "%08i -.%08i %2.010g after op%i(%i)\n",
1072                            iter, (int)((error - besterr) * 100000000), error,
1073                            op1, pt);
1074
1075                error = besterr;
1076                opstats[op1 * 2 + 1]++;
1077                //opstats[op2 * 2 + 1]++;
1078                failures = 0;
1079                success++;
1080
1081                /* Save image! */
1082                //char buf[128];
1083                //sprintf(buf, "twit%08i.bmp", success);
1084                //if((success % 10) == 0)
1085                //    pipi_save(tmp, buf);
1086            }
1087            else
1088            {
1089                pipi_free(scrap);
1090                points[pt] = oldval;
1091                failures++;
1092            }
1093        }
1094
1095        if(DEBUG_MODE)
1096        {
1097            for(int j = 0; j < 2; j++)
1098            {
1099                fprintf(stderr,   "operation: ");
1100                for(int i = NB_OPS / 2 * j; i < NB_OPS / 2 * (j + 1); i++)
1101                    fprintf(stderr, "%4i ", i);
1102                fprintf(stderr, "\nattempts:  ");
1103                for(int i = NB_OPS / 2 * j; i < NB_OPS / 2 * (j + 1); i++)
1104                    fprintf(stderr, "%4i ", opstats[i * 2]);
1105                fprintf(stderr, "\nsuccesses: ");
1106                for(int i = NB_OPS / 2 * j; i < NB_OPS / 2 * (j + 1); i++)
1107                    fprintf(stderr, "%4i ", opstats[i * 2 + 1]);
1108                fprintf(stderr, "\n");
1109            }
1110
1111            fprintf(stderr, "Distance: %2.10g\n", error);
1112        }
1113        else
1114            fprintf(stderr, "\r                    \r");
1115
1116#if 0
1117        dst = pipi_resize(tmp, width, height);
1118        pipi_free(tmp);
1119
1120        /* Save image and bail out */
1121        pipi_save(dst, "lol.bmp");
1122        pipi_free(dst);
1123#endif
1124
1125        /* Push our points to the bitstream */
1126        for(int i = 0; i < npoints; i++)
1127            b.push(points[i], RANGE_SYXRGB);
1128        b.push(height, MAX_H);
1129        b.push(width, MAX_W);
1130
1131        /* Pop Unicode characters from the bitstream and print them */
1132        for(int i = 0; i < MAX_MSG_LEN; i++)
1133            fwrite_utf8(stdout, index2uni(b.pop(NUM_CHARACTERS)));
1134        fprintf(stdout, "\n");
1135    }
1136    else
1137    {
1138        /* Pop points from the bitstream */
1139        for(int i = dw * dh; i--; )
1140        {
1141#if POINTS_PER_CELL == 2
1142            points[i * 2 + 1] = b.pop(RANGE_SYXRGB);
1143            points[i * 2] = b.pop(RANGE_SYXRGB);
1144#else
1145            points[i] = b.pop(RANGE_SYXRGB);
1146#endif
1147        }
1148        npoints = dw * dh * POINTS_PER_CELL;
1149
1150        /* Render these points to a new image */
1151        tmp = pipi_new(dw * RANGE_X, dh * RANGE_Y);
1152        render(tmp, 0, 0, dw * RANGE_X, dh * RANGE_Y);
1153
1154        /* TODO: render directly to the final image; scaling sucks */
1155        dst = pipi_resize(tmp, width, height);
1156        pipi_free(tmp);
1157
1158        /* Save image and bail out */
1159        pipi_save(dst, dstname);
1160        pipi_free(dst);
1161    }
1162
1163    return ret;
1164}
1165
Note: See TracBrowser for help on using the repository browser.