1 /* image-proc.c: image processing routines */
2
3 #ifdef HAVE_CONFIG_H
4 #include "config.h"
5 #endif /* Def: HAVE_CONFIG_H */
6
7 #include <assert.h>
8 #include <math.h>
9 #include "xstd.h"
10 #include "image-proc.h"
11
12 #define BLACK 0
13 #define WHITE 0xff
14 #ifndef M_SQRT2
15 #define M_SQRT2 1.41421356237
16 #endif
17
18 #if 0
19 struct etyp { int t00, t11, t01, t01s; };
20
21
22 static at_bool get_edge(bitmap_type, int y, int x, struct etyp *t);
23 static void check(int v1, int v2, int v3, struct etyp *t);
24 #endif
25
26
27 /* Allocate storage for a new distance map with the same dimensions
28 as BITMAP and initialize it so that pixels in BITMAP with value
29 TARGET_VALUE are at distance zero and all other pixels are at
30 distance infinity. Then compute the gray-weighted distance from
31 every non-target point to the nearest target point. */
32
33 distance_map_type
new_distance_map(bitmap_type bitmap,unsigned char target_value,at_bool padded,at_exception_type * exp)34 new_distance_map(bitmap_type bitmap, unsigned char target_value, at_bool padded, at_exception_type * exp)
35 {
36 signed x, y;
37 float d, min;
38 distance_map_type dist;
39 unsigned char *b = BITMAP_BITS(bitmap);
40 unsigned w = BITMAP_WIDTH(bitmap);
41 unsigned h = BITMAP_HEIGHT(bitmap);
42 unsigned spp = BITMAP_PLANES(bitmap);
43
44 dist.height = h; dist.width = w;
45 XMALLOC(dist.d, h * sizeof(float*));
46 XMALLOC(dist.weight, h * sizeof(float*));
47 for (y = 0; y < (signed) h; y++)
48 {
49 XCALLOC(dist.d[y], w * sizeof(float));
50 XMALLOC(dist.weight[y], w * sizeof(float));
51 }
52
53 if (spp == 3)
54 {
55 for (y = 0; y < (signed) h; y++)
56 {
57 for (x = 0; x < (signed) w; x++, b += spp)
58 {
59 int gray; float fgray;
60 gray = (int)LUMINANCE(b[0], b[1], b[2]);
61 dist.d[y][x] = (gray == target_value ? 0.0F : 1.0e10F);
62 fgray = gray * 0.0039215686F; /* = gray / 255.0F */
63 dist.weight[y][x] = 1.0F - fgray;
64 /* dist.weight[y][x] = 1.0F - (fgray * fgray);*/
65 /* dist.weight[y][x] = (fgray < 0.5F ? 1.0F - fgray : -2.0F * fgray * (fgray - 1.0F));*/
66 }
67 }
68 }
69 else
70 {
71 for (y = 0; y < (signed) h; y++)
72 {
73 for (x = 0; x < (signed) w; x++, b += spp)
74 {
75 int gray; float fgray;
76 gray = b[0];
77 dist.d[y][x] = (gray == target_value ? 0.0F : 1.0e10F);
78 fgray = gray * 0.0039215686F; /* = gray / 255.0F */
79 dist.weight[y][x] = 1.0F - fgray;
80 /* dist.weight[y][x] = 1.0F - (fgray * fgray);*/
81 /* dist.weight[y][x] = (fgray < 0.5F ? 1.0F - fgray : -2.0F * fgray * (fgray - 1.0F)); */
82 }
83 }
84 }
85
86 /* If the image is padded then border points are all at most
87 one unit away from the nearest target point. */
88 if (padded)
89 {
90 for (y = 0; y < (signed) h; y++)
91 {
92 if (dist.d[y][0] > dist.weight[y][0])
93 dist.d[y][0] = dist.weight[y][0];
94 if (dist.d[y][w - 1] > dist.weight[y][w - 1])
95 dist.d[y][w - 1] = dist.weight[y][w - 1];
96 }
97 for (x = 0; x < (signed) w; x++)
98 {
99 if (dist.d[0][x] > dist.weight[0][x])
100 dist.d[0][x] = dist.weight[0][x];
101 if (dist.d[h - 1][x] > dist.weight[h - 1][x])
102 dist.d[h - 1][x] = dist.weight[h - 1][x];
103 }
104 }
105
106 /* Scan the image from left to right, top to bottom.
107 Examine the already-visited neighbors of each point (those
108 situated above or to the left of it). Each neighbor knows
109 the distance to its nearest target point; add to this distance
110 the distance from the central point to the neighbor (either
111 sqrt(2) or one) multiplied by the central point's weight
112 (derived from its gray level). Replace the distance already
113 stored at the central point if the new distance is smaller. */
114 for (y = 1; y < (signed) h; y++)
115 {
116 for (x = 1; x < (signed) w; x++)
117 {
118 if (dist.d[y][x] == 0.0F) continue;
119
120 min = dist.d[y][x];
121
122 /* upper-left neighbor */
123 d = dist.d[y - 1][x - 1] + (float) M_SQRT2 * dist.weight[y][x];
124 if (d < min) min = dist.d[y][x] = d;
125
126 /* upper neighbor */
127 d = dist.d[y - 1][x] + dist.weight[y][x];
128 if (d < min) min = dist.d[y][x] = d;
129
130 /* left neighbor */
131 d = dist.d[y][x - 1] + dist.weight[y][x];
132 if (d < min) min = dist.d[y][x] = d;
133
134 /* upper-right neighbor (except at the last column) */
135 if (x + 1 < (signed) w)
136 {
137 d = dist.d[y - 1][x + 1] + (float) M_SQRT2 * dist.weight[y][x];
138 if (d < min) min = dist.d[y][x] = d;
139 }
140 }
141 }
142
143 /* Same as above, but now scanning right to left, bottom to top. */
144 for (y = h - 2; y >= 0; y--)
145 {
146 for (x = w - 2; x >= 0; x--)
147 {
148 min = dist.d[y][x];
149
150 /* lower-right neighbor */
151 d = dist.d[y + 1][x + 1] + (float) M_SQRT2 * dist.weight[y][x];
152 if (d < min) min = dist.d[y][x] = d;
153
154 /* lower neighbor */
155 d = dist.d[y + 1][x] + dist.weight[y][x];
156 if (d < min) min = dist.d[y][x] = d;
157
158 /* right neighbor */
159 d = dist.d[y][x + 1] + dist.weight[y][x];
160 if (d < min) min = dist.d[y][x] = d;
161
162 /* lower-left neighbor (except at the first column) */
163 if (x - 1 >= 0)
164 {
165 d = dist.d[y + 1][x - 1] + (float) M_SQRT2 * dist.weight[y][x];
166 if (d < min) min = dist.d[y][x] = d;
167 }
168 }
169 }
170 return dist;
171 }
172
173
174 /* Free the dynamically-allocated storage associated with a distance map. */
175
176 void
free_distance_map(distance_map_type * dist)177 free_distance_map(distance_map_type *dist)
178 {
179 unsigned y, h;
180
181 if (!dist) return;
182
183 h = BITMAP_HEIGHT(*dist);
184
185 if (dist->d != NULL)
186 {
187 for (y = 0; y < h; y++)
188 free((at_address*)dist->d[y]);
189 free((at_address*)dist->d);
190 }
191 if (dist->weight != NULL)
192 {
193 for (y = 0; y < h; y++)
194 free((at_address*)dist->weight[y]);
195 free((at_address*)dist->weight);
196 }
197 }
198
199
200 #if 0
201 void
202 medial_axis(bitmap_type *bitmap, distance_map_type *dist,
203 const color_type *bg_color)
204 {
205 unsigned x, y, test;
206 unsigned w, h;
207 unsigned char *b;
208 float **d, f;
209 color_type bg;
210
211 assert(bitmap != NULL);
212
213 assert(BITMAP_PLANES(*bitmap) == 1);
214
215 b = BITMAP_BITS(*bitmap);
216 assert(b != NULL);
217 assert(dist != NULL);
218 d = dist->d;
219 assert(d != NULL);
220
221 h = BITMAP_HEIGHT(*dist);
222 w = BITMAP_WIDTH(*dist);
223 assert(BITMAP_WIDTH(*bitmap) == w && BITMAP_HEIGHT(*bitmap) == h);
224
225 if (bg_color) bg = *bg_color;
226 else bg.r = bg.g = bg.b = 255;
227
228 f = d[0][0] + 0.5;
229 test = (f < d[1][0]) + (f < d[1][1]) + (f < d[0][1]);
230 if (test > 1) b[0] = bg.r;
231
232 f = d[0][w-1] + 0.5;
233 test = (f < d[1][w-1]) + (f < d[1][w-2]) + (f < d[0][w-2]);
234 if (test > 1) b[w-1] = bg.r;
235
236 for (x = 1; x < w - 1; x++)
237 {
238 f = d[0][x] + 0.5;
239 test = (f < d[0][x-1]) + (f < d[0][x+1]) + (f < d[1][x-1])
240 + (f < d[1][x]) + (f < d[1][x+1]);
241 if (test > 1) b[x] = bg.r;
242 }
243 b += w;
244
245 for (y = 1; y < h - 1; y++)
246 {
247 f = d[y][0] + 0.5;
248 test = (f < d[y-1][0]) + (f < d[y-1][1]) + (f < d[y][1])
249 + (f < d[y+1][0]) + (f < d[y+1][1]);
250 if (test > 1) b[0] = bg.r;
251
252 for (x = 1; x < w - 1; x++)
253 {
254 f = d[y][x] + 0.5;
255 test = (f < d[y-1][x-1]) + (f < d[y-1][x]) + (f < d[y-1][x+1])
256 + (f < d[y][x-1]) + (f < d[y][x+1])
257 + (f < d[y+1][x-1]) + (f < d[y+1][x]) + (f < d[y+1][x+1]);
258 if (test > 1) b[x] = bg.r;
259 }
260
261 f = d[y][w-1] + 0.5;
262 test = (f < d[y-1][w-1]) + (f < d[y-1][w-2]) + (f < d[y][w-2])
263 + (f < d[y+1][w-1]) + (f < d[y+1][w-2]);
264 if (test > 1)
265 b[w-1] = bg.r;
266
267 b += w;
268 }
269
270 for (x = 1; x < w - 1; x++)
271 {
272 f = d[h-1][x] + 0.5;
273 test = (f < d[h-1][x-1]) + (f < d[h-1][x+1])
274 + (f < d[h-2][x-1]) + (f < d[h-2][x]) + (f < d[h-2][x+1]);
275 if (test > 1) b[x] = bg.r;
276 }
277
278 f = d[h-1][0] + 0.5;
279 test = (f < d[h-2][0]) + (f < d[h-2][1]) + (f < d[h-1][1]);
280 if (test > 1) b[0] = bg.r;
281
282 f = d[h-1][w-1] + 0.5;
283 test = (f < d[h-2][w-1]) + (f < d[h-2][w-2]) + (f < d[h-1][w-2]);
284 if (test > 1) b[w-1] = bg.r;
285 }
286 #endif
287
288
289 /* Binarize a grayscale or color image. */
290
291 void
binarize(bitmap_type * bitmap)292 binarize(bitmap_type *bitmap)
293 {
294 unsigned i, npixels, spp;
295 unsigned char *b;
296
297 assert(bitmap != NULL);
298 assert(BITMAP_BITS(*bitmap) != NULL);
299
300 b = BITMAP_BITS(*bitmap);
301 spp = BITMAP_PLANES(*bitmap);
302 npixels = BITMAP_WIDTH(*bitmap) * BITMAP_HEIGHT(*bitmap);
303
304 if (spp == 1)
305 {
306 for (i = 0; i < npixels; i++)
307 b[i] = (b[i] > GRAY_THRESHOLD ? WHITE : BLACK);
308 }
309 else if (spp == 3)
310 {
311 unsigned char *rgb = b;
312 for (i = 0; i < npixels; i++, rgb += 3)
313 {
314 b[i] = (LUMINANCE(rgb[0], rgb[1], rgb[2]) > GRAY_THRESHOLD
315 ? WHITE : BLACK);
316 }
317 XREALLOC(BITMAP_BITS(*bitmap), npixels);
318 BITMAP_PLANES(*bitmap) = 1;
319 }
320 else
321 {
322 WARNING1("binarize: %u-plane images are not supported", spp);
323 }
324 }
325
326
327 #if 0
328 /* Thin a binary image, replacing the original image with the thinned one. */
329
330 bitmap_type
331 ip_thin(bitmap_type input_b)
332 {
333 unsigned y, x, i;
334 at_bool k, again;
335 struct etyp t;
336 unsigned w = BITMAP_WIDTH(input_b);
337 unsigned h = BITMAP_HEIGHT(input_b);
338 size_t num_bytes = w * h;
339 bitmap_type b = input_b;
340
341 if (BITMAP_PLANES(input_b) != 1)
342 {
343 FATAL1("thin: single-plane image required; "
344 "%u-plane images cannot be thinned", BITMAP_PLANES(input_b));
345 return b;
346 }
347
348 /* Process and return a copy of the input image. */
349 XMALLOC(b.bitmap, num_bytes);
350 memcpy(b.bitmap, input_b.bitmap, num_bytes);
351
352 /* Set background pixels to zero, foreground pixels to one. */
353 for (i = 0; i < num_bytes; i++)
354 b.bitmap[i] = (b.bitmap[i] == BLACK ? 1 : 0);
355
356 again = true;
357 while (again)
358 {
359 again = false;
360
361 for (y = 1; y < h - 1; y++)
362 {
363 for (x = 1; x < w - 1; x++)
364 {
365 /* During processing, pixels are used to store edge
366 type codes, so we can't just test for WHITE or BLACK. */
367 if (*BITMAP_PIXEL(b, y, x) == 0) continue;
368
369 k = (!get_edge(b, y, x, &t)
370 || (get_edge(b, y, x+1, &t) && *BITMAP_PIXEL(b, y-1, x)
371 && *BITMAP_PIXEL(b, y+1, x))
372 || (get_edge(b, y+1, x, &t) && *BITMAP_PIXEL(b, y, x-1)
373 && *BITMAP_PIXEL(b, y, x+1))
374 || (get_edge(b, y, x+1, &t) && get_edge(b, y+1, x+1, &t)
375 && get_edge(b, y+1, x, &t)));
376 if (k) continue;
377
378 get_edge(b, y, x, &t);
379 if (t.t01) *BITMAP_PIXEL(b, y, x) |= 4;
380 *BITMAP_PIXEL(b, y, x) |= 2;
381 again = true;
382 }
383 }
384
385 for (y = 0; y < h; y++)
386 for (x = 0; x < w; x++)
387 if (*BITMAP_PIXEL(b, y, x) & 02) *BITMAP_PIXEL(b, y, x) = 0;
388
389 for (y = 1; y < h - 1; y++)
390 {
391 for (x = 1; x < w - 1; x++)
392 {
393 if (*BITMAP_PIXEL(b, y, x) == 0) continue;
394
395 k = (!get_edge(b, y, x, &t)
396 || ((*BITMAP_PIXEL(b, y, x) & 04) == 0)
397 || (get_edge(b, y+1, x, &t) && (*BITMAP_PIXEL(b, y, x-1))
398 && *BITMAP_PIXEL(b, y, x+1))
399 || (get_edge(b, y, x+1, &t) && *BITMAP_PIXEL(b, y-1, x)
400 && *BITMAP_PIXEL(b, y+1, x))
401 || (get_edge(b, y+1, x, &t) && get_edge(b, y, x+1, &t)
402 && get_edge(b, y+1, x+1, &t)));
403 if (k) continue;
404
405 *BITMAP_PIXEL(b, y, x) |= 02;
406 again = true;
407 }
408 }
409
410 for (y = 0; y < h; y++)
411 {
412 for (x = 0; x < w; x++)
413 {
414 if (*BITMAP_PIXEL(b, y, x) & 02) *BITMAP_PIXEL(b, y, x) = 0;
415 else if (*BITMAP_PIXEL(b, y, x) > 0) *BITMAP_PIXEL(b, y, x) = 1;
416 }
417 }
418 }
419
420 /* Staircase removal; northward bias. */
421 for (y = 1; y < h - 1; y++)
422 {
423 for (x = 1; x < w - 1; x++)
424 {
425 if (*BITMAP_PIXEL(b, y, x) == 0) continue;
426
427 k = !(*BITMAP_PIXEL(b, y-1, x)
428 && ((*BITMAP_PIXEL(b, y, x+1) && !*BITMAP_PIXEL(b, y-1, x+1)
429 && !*BITMAP_PIXEL(b, y+1, x-1)
430 && (!*BITMAP_PIXEL(b, y, x-1) || !*BITMAP_PIXEL(b, y+1, x)))
431 || (*BITMAP_PIXEL(b, y, x-1) && !*BITMAP_PIXEL(b, y-1, x-1)
432 && !*BITMAP_PIXEL(b, y+1, x+1) &&
433 (!*BITMAP_PIXEL(b, y, x+1) || !*BITMAP_PIXEL(b, y+1, x)))));
434 if (k) continue;
435
436 *BITMAP_PIXEL(b, y, x) |= 02;
437 }
438 }
439 for (y = 0; y < h; y++)
440 {
441 for (x = 0; x < w; x++)
442 {
443 if (*BITMAP_PIXEL(b, y, x) & 02) *BITMAP_PIXEL(b, y, x) = 0;
444 else if (*BITMAP_PIXEL(b, y, x) > 0) *BITMAP_PIXEL(b, y, x) = 1;
445 }
446 }
447
448 /* Southward bias */
449 for (y = 1; y < h - 1; y++)
450 {
451 for (x = 1; x < w - 1; x++)
452 {
453 if (*BITMAP_PIXEL(b, y, x) == 0) continue;
454
455 k = !(*BITMAP_PIXEL(b, y+1, x)
456 && ((*BITMAP_PIXEL(b, y, x+1) && !*BITMAP_PIXEL(b, y+1, x+1)
457 && !*BITMAP_PIXEL(b, y-1, x-1) && (!*BITMAP_PIXEL(b, y, x-1)
458 || !*BITMAP_PIXEL(b, y-1, x))) || (*BITMAP_PIXEL(b, y, x-1)
459 && !*BITMAP_PIXEL(b, y+1, x-1) && !*BITMAP_PIXEL(b, y-1, x+1)
460 && (!*BITMAP_PIXEL(b, y, x+1) || !*BITMAP_PIXEL(b, y-1, x)) )));
461 if (k) continue;
462
463 *BITMAP_PIXEL(b, y, x) |= 02;
464 }
465 }
466 for (y = 0; y < h; y++)
467 {
468 for (x = 0; x < w; x++)
469 {
470 if (*BITMAP_PIXEL(b, y, x) & 02) *BITMAP_PIXEL(b, y, x) = 0;
471 else if (*BITMAP_PIXEL(b, y, x) > 0) *BITMAP_PIXEL(b, y, x) = 1;
472 }
473 }
474
475 /* Set background pixels to WHITE, foreground pixels to BLACK. */
476 for (i = 0; i < num_bytes; i++)
477 b.bitmap[i] = (b.bitmap[i] == 0 ? WHITE : BLACK);
478 return b;
479 }
480
481
482 at_bool get_edge(bitmap_type b, int y, int x, struct etyp *t)
483 {
484 t->t00 = 0; t->t01 = 0; t->t01s = 0; t->t11 = 0;
485 check(*BITMAP_PIXEL(b, y - 1, x - 1), *BITMAP_PIXEL(b, y - 1, x),
486 *BITMAP_PIXEL(b, y - 1, x + 1), t);
487 check(*BITMAP_PIXEL(b, y - 1, x + 1), *BITMAP_PIXEL(b, y, x + 1),
488 *BITMAP_PIXEL(b, y + 1, x + 1), t);
489 check(*BITMAP_PIXEL(b, y + 1, x + 1), *BITMAP_PIXEL(b, y + 1, x),
490 *BITMAP_PIXEL(b, y + 1, x - 1), t);
491 check(*BITMAP_PIXEL(b, y + 1, x - 1), *BITMAP_PIXEL(b, y, x - 1),
492 *BITMAP_PIXEL(b, y - 1, x - 1), t);
493 return *BITMAP_PIXEL(b, y, x) && t->t00 && t->t11 && !t->t01s;
494 }
495
496
497 void check(int v1, int v2, int v3, struct etyp *t)
498 {
499 if (!v2 && (!v1 || !v3)) t->t00 = 1;
500 if (v2 && (v1 || v3)) t->t11 = 1;
501 if ((!v1 && v2) || (!v2 && v3)) { t->t01s = t->t01; t->t01 = 1; }
502 }
503 #endif
504