]> mj.ucw.cz Git - libucw.git/blob - images/sig-cmp-gen.h
Merge with git+ssh://git.ucw.cz/projects/sherlock/GIT/sherlock.git
[libucw.git] / images / sig-cmp-gen.h
1 #ifdef EXPLAIN
2 #  define MSG(x...) do{ line += sprintf(line, x); }while(0)
3 #  define LINE do{ line = buf; msg(line, param); }while(0)
4
5 static void
6 explain_signature(struct image_signature *sig, void (*msg)(byte *text, void *param), void *param)
7 {
8   byte buf[1024], *line = buf;
9   MSG("signature: flags=0x%x df=%u dh=%u f=(%u", sig->flags, sig->df, sig->dh, sig->vec.f[0]);
10   for (uns i = 1; i < IMAGE_VEC_F; i++)
11     MSG(" %u", sig->vec.f[i]);
12   MSG(")");
13   LINE;
14   for (uns j = 0; j < sig->len; j++)
15     {
16       struct image_region *reg = sig->reg + j;
17       MSG("region %u: wa=%u wb=%u f=(%u", j, reg->wa, reg->wb, reg->f[0]);
18       for (uns i = 1; i < IMAGE_VEC_F; i++)
19         MSG(" %u", reg->f[i]);
20       MSG(") h=(%u", reg->h[0]);
21       for (uns i = 1; i < IMAGE_REG_H; i++)
22         MSG(" %u", reg->h[i]);
23       MSG(")");
24       LINE;
25     }
26 }
27
28 #else
29 #  define MSG(x...) do{}while(0)
30 #  define LINE do{}while(0)
31 #endif
32
33 #define MSGL(x...) do{ MSG(x); LINE; }while(0)
34
35 #ifndef EXPLAIN
36 static uns
37 image_signatures_dist_integrated(struct image_signature *sig1, struct image_signature *sig2)
38 #else
39 static uns
40 image_signatures_dist_integrated_explain(struct image_signature *sig1, struct image_signature *sig2, void (*msg)(byte *text, void *param), void *param)
41 #endif
42 {
43   uns dist[IMAGE_REG_MAX * IMAGE_REG_MAX], p[IMAGE_REG_MAX], q[IMAGE_REG_MAX];
44   uns n, i, j, k, l, s, d;
45   struct image_region *reg1, *reg2;
46 #ifdef EXPLAIN
47   byte buf[1024], *line = buf;
48   MSGL("Integrated matching");
49   explain_signature(sig1, msg, param);
50   explain_signature(sig2, msg, param);
51 #endif
52
53   /* FIXME: do not mux textured and non-textured images (should be split in clusters tree) */
54   if ((sig1->flags ^ sig2->flags) & IMAGE_SIG_TEXTURED)
55     {
56       MSGL("Textured vs non-textured");
57       return ~0U;
58     }
59
60   /* Compute distance matrix */
61   n = 0;
62   MSGL("Distance matrix:");
63   /* ... for non-textured images */
64   if (!((sig1->flags | sig2->flags) & IMAGE_SIG_TEXTURED))
65     for (j = 0, reg2 = sig2->reg; j < sig2->len; j++, reg2++)
66       for (i = 0, reg1 = sig1->reg; i < sig1->len; i++, reg1++)
67         {
68           uns dt = 0, ds = 0, dp = 0, d;
69           for (uns i = 0; i < IMAGE_VEC_F; i++)
70             dt += image_sig_cmp_features_weights[i] * isqr((int)reg1->f[i] - (int)reg2->f[i]);
71           for (uns i = 0; i < 3; i++)
72             ds += image_sig_cmp_features_weights[IMAGE_VEC_F + i] * isqr((int)reg1->h[i] - (int)reg2->h[i]);
73           for (uns i = 3; i < 5; i++)
74             dp += image_sig_cmp_features_weights[IMAGE_VEC_F + i] * isqr((int)reg1->h[i] - (int)reg2->h[i]);
75 #if 0
76           int x1, y1, x2, y2;
77           if (sig1->cols > sig1->rows)
78             {
79               x1 = reg1->h[3];
80               y1 = ((int)reg1->h[4] - 64) * (int)sig1->rows / (int)sig1->cols + 64;
81             }
82           else
83             {
84               y1 = reg1->h[4];
85               x1 = ((int)reg1->h[3] - 64) * (int)sig1->cols / (int)sig1->rows + 64;
86             }
87           if (sig2->cols > sig2->rows)
88             {
89               x2 = reg2->h[3];
90               y2 = ((int)reg2->h[4] - 64) * (int)sig2->rows / (int)sig2->cols + 64;
91             }
92           else
93             {
94               y2 = reg2->h[4];
95               x2 = ((int)reg2->h[3] - 64) * (int)sig2->cols / (int)sig2->rows + 64;
96             }
97           MSGL("%d %d %d %d", x1, y1, x2, y2);
98           dp = image_sig_cmp_features_weights[IMAGE_VEC_F + 3] * isqr(x1 - x2) +
99                image_sig_cmp_features_weights[IMAGE_VEC_F + 4] * isqr(y1 - y2);
100 #endif
101 #if 0
102           d = dt * (4 + MIN(8, (ds >> 12))) * (4 + MIN(8, (dp >> 10))) + (ds >> 11) + (dp >> 10);
103           MSG("[%u, %u] d=%u=(%u * %u * %u + %u + %u) dt=%u ds=%u dp=%u df=(%d", i, j, d,
104               dt, 4 + MIN(8, (ds >> 12)), 4 + MIN(8, dp >> 10), ds >> 11, dp >> 10, dt, ds, dp, (int)reg1->f[0] - (int)reg2->f[0]);
105 #endif
106 #if 1
107           d = dt;
108           if (ds < 1000)
109             d = d * 4;
110           else if (ds < 4000)
111             d = d * 6 + 8;
112           else if (ds < 10000)
113             d = d * 8 + 20;
114           else if (ds < 50000)
115             d = d * 10 + 50;
116           else
117             d = d * 12 + 100;
118           if (dp < 1000)
119             d = d * 2;
120           else if (dp < 4000)
121             d = d * 3 + 100;
122           else if (dp < 10000)
123             d = d * 4 + 800;
124           else
125             d = d * 5 + 3000;
126 #endif
127           dist[n++] = (d << 8) + i + (j << 4);
128           MSG("[%u, %u] d=%u dt=%u ds=%u dp=%u df=(%d", i, j, d, dt, ds, dp, (int)reg1->f[0] - (int)reg2->f[0]);
129 #ifdef EXPLAIN
130           for (uns i = 1; i < IMAGE_VEC_F; i++)
131             MSG(" %d", (int)reg1->f[i] - (int)reg2->f[i]);
132           MSG(") dh=(%d", (int)reg1->h[0] - (int)reg2->h[0]);
133           for (uns i = 1; i < IMAGE_REG_H; i++)
134             MSG(" %d", (int)reg1->h[i] - (int)reg2->h[i]);
135           MSGL(")");
136 #endif
137         }
138   /* ... for textured images (ignore shape properties) */
139   else
140     for (j = 0, reg2 = sig2->reg; j < sig2->len; j++, reg2++)
141       for (i = 0, reg1 = sig1->reg; i < sig1->len; i++, reg1++)
142         {
143           uns dt = 0;
144           for (uns i = 0; i < IMAGE_VEC_F; i++)
145             dt += image_sig_cmp_features_weights[i] * isqr((int)reg1->f[i] - (int)reg2->f[i]);
146           dist[n++] = (dt << 12) + i + (j << 4);
147 #ifdef EXPLAIN
148           MSG("[%u, %u] dt=%u df=(%d", i, j, dt, (int)reg1->f[0] - (int)reg2->f[0]);
149           for (uns i = 1; i < IMAGE_VEC_F; i++)
150             MSG(" %d", (int)reg1->f[i] - (int)reg2->f[i]);
151           MSGL(")");
152 #endif
153         }
154
155   /* One or both signatures have no regions */
156   if (!n)
157     return ~0U;
158
159   /* Get percentages */
160   for (i = 0, reg1 = sig1->reg; i < sig1->len; i++, reg1++)
161     p[i] = reg1->wb;
162   for (i = 0, reg2 = sig2->reg; i < sig2->len; i++, reg2++)
163     q[i] = reg2->wb;
164
165   /* Sort entries in distance matrix */
166   image_signatures_dist_integrated_sort(n, dist);
167
168   /* Compute significance matrix and resulting distance */
169   uns sum = 0;
170   MSGL("Significance matrix:");
171   for (k = 0, l = 128; l; k++)
172     {
173       i = dist[k] & 15;
174       j = (dist[k] >> 4) & 15;
175       d = dist[k] >> 8;
176       if (p[i] <= q[j])
177         {
178           s = p[i];
179           q[j] -= p[i];
180           p[i] = 0;
181         }
182       else
183         {
184           s = q[j];
185           p[i] -= q[j];
186           q[j] = 0;
187         }
188       l -= s;
189       sum += s * d;
190 #ifdef EXPLAIN
191       reg1 = sig1->reg + i;
192       reg2 = sig2->reg + j;
193       MSG("[%u, %u] s=%u d=%u df=(%d", i, j, s, d, (int)reg1->f[0] - (int)reg2->f[0]);
194       for (uns i = 1; i < IMAGE_VEC_F; i++)
195         MSG(" %d", (int)reg1->f[i] - (int)reg2->f[i]);
196       if (!((sig1->flags | sig2->flags) & IMAGE_SIG_TEXTURED))
197         {
198           MSG(") dh=(%d", (int)reg1->h[0] - (int)reg2->h[0]);
199           for (uns i = 1; i < IMAGE_REG_H; i++)
200             MSG(" %d", (int)reg1->h[i] - (int)reg2->h[i]);
201         }
202       MSGL(")");
203 #endif
204     }
205
206   d = sum / 32;
207
208   uns a = sig1->cols * sig2->rows;
209   uns b = sig1->rows * sig2->cols;
210   if (a < 2 * b && b < 2 * a)
211     d = d * 2;
212   else if (a < 4 * b && b < 4 * a)
213     d = d * 3;
214   else
215     d = d * 5;
216
217   a = sig1->cols * sig1->rows;
218   b = sig2->cols * sig2->rows;
219
220   if ((a < 1000 && b > 5000) || (b < 1000 && a > 5000))
221     d = d * 2;
222   else if ((a < 5000 && b > 20000) || (b < 5000 && a > 20000))
223     d = d * 3 / 2;
224
225   return d;
226 }
227
228 #ifndef EXPLAIN
229 static uns
230 image_signatures_dist_fuzzy(struct image_signature *sig1, struct image_signature *sig2)
231 #else
232 static uns
233 image_signatures_dist_fuzzy_explain(struct image_signature *sig1, struct image_signature *sig2, void (*msg)(byte *text, void *param), void *param)
234 #endif
235 {
236 #ifdef EXPLAIN
237   byte buf[1024], *line = buf;
238   MSGL("Fuzzy matching");
239   explain_signature(sig1, msg, param);
240   explain_signature(sig2, msg, param);
241 #endif
242
243   /* FIXME: do not mux textured and non-textured images (should be split in clusters tree) */
244   if ((sig1->flags ^ sig2->flags) & IMAGE_SIG_TEXTURED)
245     {
246       MSGL("Textured vs non-textured");
247       return ~0U;
248     }
249
250   uns cnt1 = sig1->len;
251   uns cnt2 = sig2->len;
252   struct image_region *reg1 = sig1->reg;
253   struct image_region *reg2 = sig2->reg;
254   uns mf[IMAGE_REG_MAX][IMAGE_REG_MAX], mh[IMAGE_REG_MAX][IMAGE_REG_MAX];
255   uns lf[IMAGE_REG_MAX * 2], lh[IMAGE_REG_MAX * 2];
256   uns df = sig1->df + sig2->df, dh = sig1->dh + sig2->dh;
257
258   /* Compute distance matrix */
259   for (uns i = 0; i < cnt1; i++)
260     for (uns j = 0; j < cnt2; j++)
261       {
262         uns d = 0;
263         for (uns k = 0; k < IMAGE_VEC_F; k++)
264           {
265             int dif = reg1[i].f[k] - reg2[j].f[k];
266             d += image_sig_cmp_features_weights[k] * dif * dif;
267           }
268         mf[i][j] = d;
269         d = 0;
270         for (uns k = 0; k < IMAGE_REG_H; k++)
271           {
272             int dif = reg1[i].h[k] - reg2[j].h[k];
273             d += image_sig_cmp_features_weights[k + IMAGE_VEC_F] * dif * dif;
274           }
275         mh[i][j] = d;
276       }
277
278   uns lfs = 0, lhs = 0;
279   for (uns i = 0; i < cnt1; i++)
280     {
281       uns f = mf[i][0], h = mh[i][0];
282       for (uns j = 1; j < cnt2; j++)
283         {
284           f = MIN(f, mf[i][j]);
285           h = MIN(h, mh[i][j]);
286         }
287       lf[i] = (df * 0x10000) / (df + fast_sqrt_u32(f));
288       lh[i] = (dh * 0x10000) / (dh + fast_sqrt_u32(h));
289       lfs += lf[i] * (6 * reg1[i].wa + 2 * reg1[i].wb);
290       lhs += lh[i] * reg1[i].wa;
291     }
292   for (uns i = 0; i < cnt2; i++)
293     {
294       uns f = mf[0][i], h = mh[0][i];
295       for (uns j = 1; j < cnt1; j++)
296         {
297           f = MIN(f, mf[j][i]);
298           h = MIN(h, mh[j][i]);
299         }
300       lf[i + cnt1] = (df * 0x10000) / (df + fast_sqrt_u32(f));
301       lh[i + cnt1] = (dh * 0x10000) / (dh + fast_sqrt_u32(h));
302       lfs += lf[i] * (6 * reg2[i].wa + 2 * reg2[i].wb);
303       lhs += lh[i] * reg2[i].wa;
304     }
305
306   uns measure = lfs * 6 + lhs * 2 * 8;
307
308 #ifdef EXPLAIN
309   /* Display similarity vectors */
310   MSG("Lf=(");
311   for (uns i = 0; i < cnt1 + cnt2; i++)
312     {
313       if (i)
314         MSG(" ");
315       if (i == cnt1)
316         MSG("~ ");
317       MSG("%.4f", (double)lf[i] / 0x10000);
318     }
319   MSGL(")");
320   MSG("Lh=(");
321   for (uns i = 0; i < cnt1 + cnt2; i++)
322     {
323       if (i)
324         MSG(" ");
325       if (i == cnt1)
326         MSG("~ ");
327       MSG("%.4f", (double)lh[i] / 0x10000);
328     }
329   MSGL(")");
330   MSGL("Lfm=%.4f", lfs / (double)(1 << (3 + 8 + 16)));
331   MSGL("Lhm=%.4f", lhs / (double)(1 << (8 + 16)));
332   MSGL("measure=%.4f", measure / (double)(1 << (3 + 3 + 8 + 16)));
333 #endif
334
335   return (1 << (3 + 3 + 8 + 16)) - measure;
336 }
337
338 #ifndef EXPLAIN
339 static uns
340 image_signatures_dist_average(struct image_signature *sig1, struct image_signature *sig2)
341 #else
342 static uns
343 image_signatures_dist_average_explain(struct image_signature *sig1, struct image_signature *sig2, void (*msg)(byte *text, void *param), void *param)
344 #endif
345 {
346 #ifdef EXPLAIN
347   byte buf[1024], *line = buf;
348   MSGL("Average matching");
349 #endif
350
351   uns dist = 0;
352   for (uns i = 0; i < IMAGE_VEC_F; i++)
353     {
354       uns d = image_sig_cmp_features_weights[0] * isqr((int)sig1->vec.f[i] - (int)sig2->vec.f[i]); 
355       MSGL("feature %u: d=%u (%u %u)", i, d, sig1->vec.f[i], sig2->vec.f[i]);
356       dist += d;
357     }
358
359   MSGL("dist=%u", dist);
360   return dist;
361 }
362
363 #ifndef EXPLAIN
364 #define CALL(x) image_signatures_dist_##x(sig1, sig2)
365 uns
366 image_signatures_dist(struct image_signature *sig1, struct image_signature *sig2)
367 #else
368 #define CALL(x) image_signatures_dist_##x##_explain(sig1, sig2, msg, param)
369 uns
370 image_signatures_dist_explain(struct image_signature *sig1, struct image_signature *sig2, void (*msg)(byte *text, void *param), void *param)
371 #endif
372 {
373   if (!sig1->len)
374     return CALL(average);
375   else
376     switch (image_sig_compare_method)
377       {
378         case 0:
379           return CALL(integrated);
380         case 1:
381           return CALL(fuzzy);
382         case 2:
383           return CALL(average);
384         default:
385           ASSERT(0);
386       }
387 }
388 #undef CALL
389
390 #undef EXPLAIN
391 #undef MSG
392 #undef LINE
393 #undef MSGL