Mercurial > hgrepos > Python2 > PyMuPDF

comparison mupdf-source/thirdparty/leptonica/src/maze.c @ 2:b50eed0cc0ef upstream

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
ADD: MuPDF v1.26.7: the MuPDF source as downloaded by a default build of PyMuPDF 1.26.4. The directory name has changed: no version number in the expanded directory now.
author Franz Glasner <fzglas.hg@dom66.de>
date Mon, 15 Sep 2025 11:43:07 +0200
parents
children
comparison
equal deleted inserted replaced
1:1d09e1dec1d9 2:b50eed0cc0ef
1 /*====================================================================*
2 - Copyright (C) 2001 Leptonica. All rights reserved.
3 -
4 - Redistribution and use in source and binary forms, with or without
5 - modification, are permitted provided that the following conditions
6 - are met:
7 - 1. Redistributions of source code must retain the above copyright
8 - notice, this list of conditions and the following disclaimer.
9 - 2. Redistributions in binary form must reproduce the above
10 - copyright notice, this list of conditions and the following
11 - disclaimer in the documentation and/or other materials
12 - provided with the distribution.
13 -
14 - THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
15 - ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
16 - LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
17 - A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL ANY
18 - CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
19 - EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
20 - PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
21 - PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY
22 - OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
23 - NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
24 - SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
25 *====================================================================*/
26
27
28 /*!
29 * \file maze.c
30 * <pre>
31 *
32 * This is a game with a pedagogical slant. A maze is represented
33 * by a binary image. The ON pixels (fg) are walls. The goal is
34 * to navigate on OFF pixels (bg), using Manhattan steps
35 * (N, S, E, W), between arbitrary start and end positions.
36 * The problem is thus to find the shortest route between two points
37 * in a binary image that are 4-connected in the bg. This is done
38 * with a breadth-first search, implemented with a queue.
39 * We also use a queue of pointers to generate the maze (image).
40 *
41 * PIX *generateBinaryMaze()
42 * static MAZEEL *mazeelCreate()
43 *
44 * PIX *pixSearchBinaryMaze()
45 * static l_int32 localSearchForBackground()
46 *
47 * Generalizing a maze to a grayscale image, the search is
48 * now for the "shortest" or least cost path, for some given
49 * cost function.
50 *
51 * PIX *pixSearchGrayMaze()
52 * </pre>
53 */
54
55 #ifdef HAVE_CONFIG_H
56 #include <config_auto.h>
57 #endif /* HAVE_CONFIG_H */
58
59 #include <string.h>
60 #ifdef _WIN32
61 #include <stdlib.h>
62 #include <time.h>
63 #endif /* _WIN32 */
64 #include "allheaders.h"
65
66 static const l_int32 MinMazeWidth = 50;
67 static const l_int32 MinMazeHeight = 50;
68
69 static const l_float32 DefaultWallProbability = 0.65f;
70 static const l_float32 DefaultAnisotropyRatio = 0.25f;
71
72 enum { /* direction from parent to newly created element */
73 START_LOC = 0,
74 DIR_NORTH = 1,
75 DIR_SOUTH = 2,
76 DIR_WEST = 3,
77 DIR_EAST = 4
78 };
79
80 struct MazeElement {
81 l_float32 distance;
82 l_int32 x;
83 l_int32 y;
84 l_uint32 val; /* value of maze pixel at this location */
85 l_int32 dir; /* direction from parent to child */
86 };
87 typedef struct MazeElement MAZEEL;
88
89
90 static MAZEEL *mazeelCreate(l_int32 x, l_int32 y, l_int32 dir);
91 static l_int32 localSearchForBackground(PIX *pix, l_int32 *px,
92 l_int32 *py, l_int32 maxrad);
93
94 #ifndef NO_CONSOLE_IO
95 #define DEBUG_PATH 0
96 #define DEBUG_MAZE 0
97 #endif /* ~NO_CONSOLE_IO */
98
99 /*---------------------------------------------------------------------*
100 * Binary maze generation as cellular automaton *
101 *---------------------------------------------------------------------*/
102 /*!
103 * \brief generateBinaryMaze()
104 *
105 * \param[in] w, h size of maze
106 * \param[in] xi, yi initial location
107 * \param[in] wallps probability that a pixel to the side is ON
108 * \param[in] ranis ratio of prob that pixel in forward direction
109 * is a wall to the probability that pixel in
110 * side directions is a wall
111 * \return pix, or NULL on error
112 *
113 * <pre>
114 * Notes:
115 * (1) We have two input probability factors that determine the
116 * density of walls and average length of straight passages.
117 * When ranis < 1.0, you are more likely to generate a wall
118 * to the side than going forward. Enter 0.0 for either if
119 * you want to use the default values.
120 * (2) This is a type of percolation problem, and exhibits
121 * different phases for different parameters wallps and ranis.
122 * For larger values of these parameters, regions in the maze
123 * are not explored because the maze generator walls them
124 * off and cannot get through. The boundary between the
125 * two phases in this two-dimensional parameter space goes
126 * near these values:
127 * wallps ranis
128 * 0.35 1.00
129 * 0.40 0.85
130 * 0.45 0.70
131 * 0.50 0.50
132 * 0.55 0.40
133 * 0.60 0.30
134 * 0.65 0.25
135 * 0.70 0.19
136 * 0.75 0.15
137 * 0.80 0.11
138 * (3) Because there is a considerable amount of overhead in calling
139 * pixGetPixel() and pixSetPixel(), this function can be sped
140 * up with little effort using raster line pointers and the
141 * GET_DATA* and SET_DATA* macros.
142 * </pre>
143 */
144 PIX *
145 generateBinaryMaze(l_int32 w,
146 l_int32 h,
147 l_int32 xi,
148 l_int32 yi,
149 l_float32 wallps,
150 l_float32 ranis)
151 {
152 l_int32 x, y, dir;
153 l_uint32 val;
154 l_float32 frand, wallpf, testp;
155 MAZEEL *el, *elp;
156 PIX *pixd; /* the destination maze */
157 PIX *pixm; /* for bookkeeping, to indicate pixels already visited */
158 L_QUEUE *lq;
159
160 /* On Windows, seeding is apparently necessary to get decent mazes.
161 * Windows rand() returns a value up to 2^15 - 1, whereas unix
162 * rand() returns a value up to 2^31 - 1. Therefore the generated
163 * mazes will differ on the two platforms. */
164 #ifdef _WIN32
165 srand(28*333);
166 #endif /* _WIN32 */
167
168 if (w < MinMazeWidth)
169 w = MinMazeWidth;
170 if (h < MinMazeHeight)
171 h = MinMazeHeight;
172 if (xi <= 0 || xi >= w)
173 xi = w / 6;
174 if (yi <= 0 || yi >= h)
175 yi = h / 5;
176 if (wallps < 0.05 || wallps > 0.95)
177 wallps = DefaultWallProbability;
178 if (ranis < 0.05 || ranis > 1.0)
179 ranis = DefaultAnisotropyRatio;
180 wallpf = wallps * ranis;
181
182 #if DEBUG_MAZE
183 lept_stderr("(w, h) = (%d, %d), (xi, yi) = (%d, %d)\n", w, h, xi, yi);
184 lept_stderr("Using: prob(wall) = %7.4f, anisotropy factor = %7.4f\n",
185 wallps, ranis);
186 #endif /* DEBUG_MAZE */
187
188 /* These are initialized to OFF */
189 pixd = pixCreate(w, h, 1);
190 pixm = pixCreate(w, h, 1);
191
192 lq = lqueueCreate(0);
193
194 /* Prime the queue with the first pixel; it is OFF */
195 el = mazeelCreate(xi, yi, START_LOC);
196 pixSetPixel(pixm, xi, yi, 1); /* mark visited */
197 lqueueAdd(lq, el);
198
199 /* While we're at it ... */
200 while (lqueueGetCount(lq) > 0) {
201 elp = (MAZEEL *)lqueueRemove(lq);
202 x = elp->x;
203 y = elp->y;
204 dir = elp->dir;
205 if (x > 0) { /* check west */
206 pixGetPixel(pixm, x - 1, y, &val);
207 if (val == 0) { /* not yet visited */
208 pixSetPixel(pixm, x - 1, y, 1); /* mark visited */
209 frand = (l_float32)rand() / (l_float32)RAND_MAX;
210 testp = wallps;
211 if (dir == DIR_WEST)
212 testp = wallpf;
213 if (frand <= testp) { /* make it a wall */
214 pixSetPixel(pixd, x - 1, y, 1);
215 } else { /* not a wall */
216 el = mazeelCreate(x - 1, y, DIR_WEST);
217 lqueueAdd(lq, el);
218 }
219 }
220 }
221 if (y > 0) { /* check north */
222 pixGetPixel(pixm, x, y - 1, &val);
223 if (val == 0) { /* not yet visited */
224 pixSetPixel(pixm, x, y - 1, 1); /* mark visited */
225 frand = (l_float32)rand() / (l_float32)RAND_MAX;
226 testp = wallps;
227 if (dir == DIR_NORTH)
228 testp = wallpf;
229 if (frand <= testp) { /* make it a wall */
230 pixSetPixel(pixd, x, y - 1, 1);
231 } else { /* not a wall */
232 el = mazeelCreate(x, y - 1, DIR_NORTH);
233 lqueueAdd(lq, el);
234 }
235 }
236 }
237 if (x < w - 1) { /* check east */
238 pixGetPixel(pixm, x + 1, y, &val);
239 if (val == 0) { /* not yet visited */
240 pixSetPixel(pixm, x + 1, y, 1); /* mark visited */
241 frand = (l_float32)rand() / (l_float32)RAND_MAX;
242 testp = wallps;
243 if (dir == DIR_EAST)
244 testp = wallpf;
245 if (frand <= testp) { /* make it a wall */
246 pixSetPixel(pixd, x + 1, y, 1);
247 } else { /* not a wall */
248 el = mazeelCreate(x + 1, y, DIR_EAST);
249 lqueueAdd(lq, el);
250 }
251 }
252 }
253 if (y < h - 1) { /* check south */
254 pixGetPixel(pixm, x, y + 1, &val);
255 if (val == 0) { /* not yet visited */
256 pixSetPixel(pixm, x, y + 1, 1); /* mark visited */
257 frand = (l_float32)rand() / (l_float32)RAND_MAX;
258 testp = wallps;
259 if (dir == DIR_SOUTH)
260 testp = wallpf;
261 if (frand <= testp) { /* make it a wall */
262 pixSetPixel(pixd, x, y + 1, 1);
263 } else { /* not a wall */
264 el = mazeelCreate(x, y + 1, DIR_SOUTH);
265 lqueueAdd(lq, el);
266 }
267 }
268 }
269 LEPT_FREE(elp);
270 }
271
272 lqueueDestroy(&lq, TRUE);
273 pixDestroy(&pixm);
274 return pixd;
275 }
276
277
278 static MAZEEL *
279 mazeelCreate(l_int32 x,
280 l_int32 y,
281 l_int32 dir)
282 {
283 MAZEEL *el;
284
285 el = (MAZEEL *)LEPT_CALLOC(1, sizeof(MAZEEL));
286 el->x = x;
287 el->y = y;
288 el->dir = dir;
289 return el;
290 }
291
292
293 /*---------------------------------------------------------------------*
294 * Binary maze search *
295 *---------------------------------------------------------------------*/
296 /*!
297 * \brief pixSearchBinaryMaze()
298 *
299 * \param[in] pixs 1 bpp, maze
300 * \param[in] xi, yi beginning point; use same initial point
301 * that was used to generate the maze
302 * \param[in] xf, yf end point, or close to it
303 * \param[out] ppixd [optional] maze with path illustrated, or
304 * if no path possible, the part of the maze
305 * that was searched
306 * \return pta shortest path, or NULL if either no path
307 * exists or on error
308 *
309 * <pre>
310 * Notes:
311 * (1) Because of the overhead in calling pixGetPixel() and
312 * pixSetPixel(), we have used raster line pointers and the
313 * GET_DATA* and SET_DATA* macros for many of the pix accesses.
314 * (2) Commentary:
315 * The goal is to find the shortest path between beginning and
316 * end points, without going through walls, and there are many
317 * ways to solve this problem.
318 * We use a queue to implement a breadth-first search. Two auxiliary
319 * "image" data structures can be used: one to mark the visited
320 * pixels and one to give the direction to the parent for each
321 * visited pixel. The first structure is used to avoid putting
322 * pixels on the queue more than once, and the second is used
323 * for retracing back to the origin, like the breadcrumbs in
324 * Hansel and Gretel. Each pixel taken off the queue is destroyed
325 * after it is used to locate the allowed neighbors. In fact,
326 * only one distance image is required, if you initialize it
327 * to some value that signifies "not yet visited." (We use
328 * a binary image for marking visited pixels because it is clearer.)
329 * This method for a simple search of a binary maze is implemented in
330 * pixSearchBinaryMaze().
331 * An alternative method would store the (manhattan) distance
332 * from the start point with each pixel on the queue. The children
333 * of each pixel get a distance one larger than the parent. These
334 * values can be stored in an auxiliary distance map image
335 * that is constructed simultaneously with the search. Once the
336 * end point is reached, the distance map is used to backtrack
337 * along a minimum path. There may be several equal length
338 * minimum paths, any one of which can be chosen this way.
339 * </pre>
340 */
341 PTA *
342 pixSearchBinaryMaze(PIX *pixs,
343 l_int32 xi,
344 l_int32 yi,
345 l_int32 xf,
346 l_int32 yf,
347 PIX **ppixd)
348 {
349 l_int32 i, j, x, y, w, h, d, found;
350 l_uint32 val, rpixel, gpixel, bpixel;
351 void **lines1, **linem1, **linep8, **lined32;
352 MAZEEL *el, *elp;
353 PIX *pixd; /* the shortest path written on the maze image */
354 PIX *pixm; /* for bookkeeping, to indicate pixels already visited */
355 PIX *pixp; /* for bookkeeping, to indicate direction to parent */
356 L_QUEUE *lq;
357 PTA *pta;
358
359 if (ppixd) *ppixd = NULL;
360 if (!pixs)
361 return (PTA *)ERROR_PTR("pixs not defined", __func__, NULL);
362 pixGetDimensions(pixs, &w, &h, &d);
363 if (d != 1)
364 return (PTA *)ERROR_PTR("pixs not 1 bpp", __func__, NULL);
365 if (xi <= 0 || xi >= w)
366 return (PTA *)ERROR_PTR("xi not valid", __func__, NULL);
367 if (yi <= 0 || yi >= h)
368 return (PTA *)ERROR_PTR("yi not valid", __func__, NULL);
369 pixGetPixel(pixs, xi, yi, &val);
370 if (val != 0)
371 return (PTA *)ERROR_PTR("(xi,yi) not bg pixel", __func__, NULL);
372 pixd = NULL;
373 pta = NULL;
374
375 /* Find a bg pixel near input point (xf, yf) */
376 localSearchForBackground(pixs, &xf, &yf, 5);
377
378 #if DEBUG_MAZE
379 lept_stderr("(xi, yi) = (%d, %d), (xf, yf) = (%d, %d)\n",
380 xi, yi, xf, yf);
381 #endif /* DEBUG_MAZE */
382
383 pixm = pixCreate(w, h, 1); /* initialized to OFF */
384 pixp = pixCreate(w, h, 8); /* direction to parent stored as enum val */
385 lines1 = pixGetLinePtrs(pixs, NULL);
386 linem1 = pixGetLinePtrs(pixm, NULL);
387 linep8 = pixGetLinePtrs(pixp, NULL);
388
389 lq = lqueueCreate(0);
390
391 /* Prime the queue with the first pixel; it is OFF */
392 el = mazeelCreate(xi, yi, 0); /* don't need direction here */
393 pixSetPixel(pixm, xi, yi, 1); /* mark visited */
394 lqueueAdd(lq, el);
395
396 /* Fill up the pix storing directions to parents,
397 * stopping when we hit the point (xf, yf) */
398 found = FALSE;
399 while (lqueueGetCount(lq) > 0) {
400 elp = (MAZEEL *)lqueueRemove(lq);
401 x = elp->x;
402 y = elp->y;
403 if (x == xf && y == yf) {
404 found = TRUE;
405 LEPT_FREE(elp);
406 break;
407 }
408
409 if (x > 0) { /* check to west */
410 val = GET_DATA_BIT(linem1[y], x - 1);
411 if (val == 0) { /* not yet visited */
412 SET_DATA_BIT(linem1[y], x - 1); /* mark visited */
413 val = GET_DATA_BIT(lines1[y], x - 1);
414 if (val == 0) { /* bg, not a wall */
415 SET_DATA_BYTE(linep8[y], x - 1, DIR_EAST); /* parent E */
416 el = mazeelCreate(x - 1, y, 0);
417 lqueueAdd(lq, el);
418 }
419 }
420 }
421 if (y > 0) { /* check north */
422 val = GET_DATA_BIT(linem1[y - 1], x);
423 if (val == 0) { /* not yet visited */
424 SET_DATA_BIT(linem1[y - 1], x); /* mark visited */
425 val = GET_DATA_BIT(lines1[y - 1], x);
426 if (val == 0) { /* bg, not a wall */
427 SET_DATA_BYTE(linep8[y - 1], x, DIR_SOUTH); /* parent S */
428 el = mazeelCreate(x, y - 1, 0);
429 lqueueAdd(lq, el);
430 }
431 }
432 }
433 if (x < w - 1) { /* check east */
434 val = GET_DATA_BIT(linem1[y], x + 1);
435 if (val == 0) { /* not yet visited */
436 SET_DATA_BIT(linem1[y], x + 1); /* mark visited */
437 val = GET_DATA_BIT(lines1[y], x + 1);
438 if (val == 0) { /* bg, not a wall */
439 SET_DATA_BYTE(linep8[y], x + 1, DIR_WEST); /* parent W */
440 el = mazeelCreate(x + 1, y, 0);
441 lqueueAdd(lq, el);
442 }
443 }
444 }
445 if (y < h - 1) { /* check south */
446 val = GET_DATA_BIT(linem1[y + 1], x);
447 if (val == 0) { /* not yet visited */
448 SET_DATA_BIT(linem1[y + 1], x); /* mark visited */
449 val = GET_DATA_BIT(lines1[y + 1], x);
450 if (val == 0) { /* bg, not a wall */
451 SET_DATA_BYTE(linep8[y + 1], x, DIR_NORTH); /* parent N */
452 el = mazeelCreate(x, y + 1, 0);
453 lqueueAdd(lq, el);
454 }
455 }
456 }
457 LEPT_FREE(elp);
458 }
459
460 lqueueDestroy(&lq, TRUE);
461 pixDestroy(&pixm);
462 LEPT_FREE(linem1);
463
464 if (ppixd) {
465 pixd = pixUnpackBinary(pixs, 32, 1);
466 *ppixd = pixd;
467 }
468 composeRGBPixel(255, 0, 0, &rpixel); /* start point */
469 composeRGBPixel(0, 255, 0, &gpixel);
470 composeRGBPixel(0, 0, 255, &bpixel); /* end point */
471
472 if (found) {
473 L_INFO(" Path found\n", __func__);
474 pta = ptaCreate(0);
475 x = xf;
476 y = yf;
477 while (1) {
478 ptaAddPt(pta, x, y);
479 if (x == xi && y == yi)
480 break;
481 if (pixd) /* write 'gpixel' onto the path */
482 pixSetPixel(pixd, x, y, gpixel);
483 pixGetPixel(pixp, x, y, &val);
484 if (val == DIR_NORTH)
485 y--;
486 else if (val == DIR_SOUTH)
487 y++;
488 else if (val == DIR_EAST)
489 x++;
490 else if (val == DIR_WEST)
491 x--;
492 }
493 } else {
494 L_INFO(" No path found\n", __func__);
495 if (pixd) { /* paint all visited locations */
496 lined32 = pixGetLinePtrs(pixd, NULL);
497 for (i = 0; i < h; i++) {
498 for (j = 0; j < w; j++) {
499 if (GET_DATA_BYTE(linep8[i], j) != 0)
500 SET_DATA_FOUR_BYTES(lined32[i], j, gpixel);
501 }
502 }
503 LEPT_FREE(lined32);
504 }
505 }
506 if (pixd) {
507 pixSetPixel(pixd, xi, yi, rpixel);
508 pixSetPixel(pixd, xf, yf, bpixel);
509 }
510
511 pixDestroy(&pixp);
512 LEPT_FREE(lines1);
513 LEPT_FREE(linep8);
514 return pta;
515 }
516
517
518 /*!
519 * \brief localSearchForBackground()
520 *
521 * \param[in] pix
522 * \param[out] px, py starting position for search; return found position
523 * \param[in] maxrad max distance to search from starting location
524 * \return 0 if bg pixel found; 1 if not found
525 */
526 static l_int32
527 localSearchForBackground(PIX *pix,
528 l_int32 *px,
529 l_int32 *py,
530 l_int32 maxrad)
531 {
532 l_int32 x, y, w, h, r, i, j;
533 l_uint32 val;
534
535 x = *px;
536 y = *py;
537 pixGetPixel(pix, x, y, &val);
538 if (val == 0) return 0;
539
540 /* For each value of r, restrict the search to the boundary
541 * pixels in a square centered on (x,y), clipping to the
542 * image boundaries if necessary. */
543 pixGetDimensions(pix, &w, &h, NULL);
544 for (r = 1; r < maxrad; r++) {
545 for (i = -r; i <= r; i++) {
546 if (y + i < 0 || y + i >= h)
547 continue;
548 for (j = -r; j <= r; j++) {
549 if (x + j < 0 || x + j >= w)
550 continue;
551 if (L_ABS(i) != r && L_ABS(j) != r) /* not on "r ring" */
552 continue;
553 pixGetPixel(pix, x + j, y + i, &val);
554 if (val == 0) {
555 *px = x + j;
556 *py = y + i;
557 return 0;
558 }
559 }
560 }
561 }
562 return 1;
563 }
564
565
566
567 /*---------------------------------------------------------------------*
568 * Gray maze search *
569 *---------------------------------------------------------------------*/
570 /*!
571 * \brief pixSearchGrayMaze()
572 *
573 * \param[in] pixs 1 bpp maze; w and h must be >= 50
574 * \param[in] xi, yi beginning point; use same initial point
575 * that was used to generate the maze
576 * \param[in] xf, yf end point, or close to it
577 * \param[out] ppixd [optional] maze with path illustrated, or
578 * if no path possible, the part of the maze
579 * that was searched
580 * \return pta shortest path, or NULL if either no path exists or on error
581 *
582 * <pre>
583 * Commentary:
584 * Consider first a slight generalization of the binary maze
585 * search problem. Suppose that you can go through walls,
586 * but the cost is higher say, an increment of 3 to go into
587 * a wall pixel rather than 1? You're still trying to find
588 * the shortest path. One way to do this is with an ordered
589 * queue, and a simple way to visualize an ordered queue is as
590 * a set of stacks, each stack being marked with the distance
591 * of each pixel in the stack from the start. We place the
592 * start pixel in stack 0, pop it, and process its 4 children.
593 * Each pixel is given a distance that is incremented from that
594 * of its parent 0 in this case, depending on if it is a wall
595 * pixel or not. That value may be recorded on a distance map,
596 * according to the algorithm below. For children of the first
597 * pixel, those not on a wall go in stack 1, and wall
598 * children go in stack 3. Stack 0 being emptied, the process
599 * then continues with pixels being popped from stack 1.
600 * Here is the algorithm for each child pixel. The pixel's
601 * distance value, were it to be placed on a stack, is compared
602 * with the value for it that is on the distance map. There
603 * are three possible cases:
604 * 1 If the pixel has not yet been registered, it is pushed
605 * on its stack and the distance is written to the map.
606 * 2 If it has previously been registered with a higher distance,
607 * the distance on the map is relaxed to that of the
608 * current pixel, which is then placed on its stack.
609 * 3 If it has previously been registered with an equal
610 * or lower value, the pixel is discarded.
611 * The pixels are popped and processed successively from
612 * stack 1, and when stack 1 is empty, popping starts on stack 2.
613 * This continues until the destination pixel is popped off
614 * a stack. The minimum path is then derived from the distance map,
615 * going back from the end point as before. This is just Dijkstra's
616 * algorithm for a directed graph; here, the underlying graph
617 * consisting of the pixels and four edges connecting each pixel
618 * to its 4-neighbor is a special case of a directed graph, where
619 * each edge is bi-directional. The implementation of this generalized
620 * maze search is left as an exercise to the reader.
621 *
622 * Let's generalize a bit further. Suppose the "maze" is just
623 * a grayscale image -- think of it as an elevation map. The cost
624 * of moving on this surface depends on the height, or the gradient,
625 * or whatever you want. All that is required is that the cost
626 * is specified and non-negative on each link between adjacent
627 * pixels. Now the problem becomes: find the least cost path
628 * moving on this surface between two specified end points.
629 * For example, if the cost across an edge between two pixels
630 * depends on the "gradient", you can use:
631 * cost = 1 + L_ABSdeltaV
632 * where deltaV is the difference in value between two adjacent
633 * pixels. If the costs are all integers, we can still use an array
634 * of stacks to avoid ordering the queue e.g., by using a heap sort.
635 * This is a neat problem, because you don't even have to build a
636 * maze -- you can can use it on any grayscale image!
637 *
638 * Rather than using an array of stacks, a more practical
639 * approach is to implement with a priority queue, which is
640 * a queue that is sorted so that the elements with the largest
641 * or smallest key values always come off first. The
642 * priority queue is efficiently implemented as a heap, and
643 * this is how we do it. Suppose you run the algorithm
644 * using a priority queue, doing the bookkeeping with an
645 * auxiliary image data structure that saves the distance of
646 * each pixel put on the queue as before, according to the method
647 * described above. We implement it as a 2-way choice by
648 * initializing the distance array to a large value and putting
649 * a pixel on the queue if its distance is less than the value
650 * found on the array. When you finally pop the end pixel from
651 * the queue, you're done, and you can trace the path backward,
652 * either always going downhill or using an auxiliary image to
653 * give you the direction to go at each step. This is implemented
654 * here in searchGrayMaze.
655 *
656 * Do we really have to use a sorted queue? Can we solve this
657 * generalized maze with an unsorted queue of pixels? Or even
658 * an unsorted stack, doing a depth-first search (DFS)?
659 * Consider a different algorithm for this generalized maze, where
660 * we travel again breadth first, but this time use a single,
661 * unsorted queue. An auxiliary image is used as before to
662 * store the distances and to determine if pixels get pushed
663 * on the stack or dropped. As before, we must allow pixels
664 * to be revisited, with relaxation of the distance if a shorter
665 * path arrives later. As a result, we will in general have
666 * multiple instances of the same pixel on the stack with different
667 * distances. However, because the queue is not ordered, some of
668 * these pixels will be popped when another instance with a lower
669 * distance is still on the stack. Here, we're just popping them
670 * in the order they go on, rather than setting up a priority
671 * based on minimum distance. Thus, unlike the priority queue,
672 * when a pixel is popped we have to check the distance map to
673 * see if a pixel with a lower distance has been put on the queue,
674 * and, if so, we discard the pixel we just popped. So the
675 * "while" loop looks like this:
676 * ~ pop a pixel from the queue
677 * ~ check its distance against the distance stored in the
678 * distance map; if larger, discard
679 * ~ otherwise, for each of its neighbors:
680 * ~ compute its distance from the start pixel
681 * ~ compare this distance with that on the distance map:
682 * ~ if the distance map value higher, relax the distance
683 * and push the pixel on the queue
684 * ~ if the distance map value is lower, discard the pixel
685 *
686 * How does this loop terminate? Before, with an ordered queue,
687 * it terminates when you pop the end pixel. But with an unordered
688 * queue or stack, the first time you hit the end pixel, the
689 * distance is not guaranteed to be correct, because the pixels
690 * along the shortest path may not have yet been visited and relaxed.
691 * Because the shortest path can theoretically go anywhere,
692 * we must keep going. How do we know when to stop? Dijkstra
693 * uses an ordered queue to systematically remove nodes from
694 * further consideration. Each time a pixel is popped, we're
695 * done with it; it's "finalized" in the Dijkstra sense because
696 * we know the shortest path to it. However, with an unordered
697 * queue, the brute force answer is: stop when the queue
698 * or stack is empty, because then every pixel in the image
699 * has been assigned its minimum "distance" from the start pixel.
700 *
701 * This is similar to the situation when you use a stack for the
702 * simpler uniform-step problem: with breadth-first search BFS
703 * the pixels on the queue are automatically ordered, so you are
704 * done when you locate the end pixel as a neighbor of a popped pixel;
705 * whereas depth-first search DFS, using a stack, requires,
706 * in general, a search of every accessible pixel. Further, if
707 * a pixel is revisited with a smaller distance, that distance is
708 * recorded and the pixel is put on the stack again.
709 *
710 * But surely, you ask, can't we stop sooner? What if the
711 * start and end pixels are very close to each other?
712 * OK, suppose they are, and you have very high walls and a
713 * long snaking level path that is actually the minimum cost.
714 * That long path can wind back and forth across the entire
715 * maze many times before ending up at the end point, which
716 * could be just over a wall from the start. With the unordered
717 * queue, you very quickly get a high distance for the end
718 * pixel, which will be relaxed to the minimum distance only
719 * after all the pixels of the path have been visited and placed
720 * on the queue, multiple times for many of them. So that's the
721 * price for not ordering the queue!
722 * </pre>
723 */
724 PTA *
725 pixSearchGrayMaze(PIX *pixs,
726 l_int32 xi,
727 l_int32 yi,
728 l_int32 xf,
729 l_int32 yf,
730 PIX **ppixd)
731 {
732 l_int32 x, y, w, h, d;
733 l_uint32 val, valr, vals, rpixel, gpixel, bpixel;
734 void **lines8, **liner32, **linep8;
735 l_int32 cost, dist, distparent, sival, sivals;
736 MAZEEL *el, *elp;
737 PIX *pixd; /* optionally plot the path on this RGB version of pixs */
738 PIX *pixr; /* for bookkeeping, to indicate the minimum distance */
739 /* to pixels already visited */
740 PIX *pixp; /* for bookkeeping, to indicate direction to parent */
741 L_HEAP *lh;
742 PTA *pta;
743
744 if (ppixd) *ppixd = NULL;
745 if (!pixs)
746 return (PTA *)ERROR_PTR("pixs not defined", __func__, NULL);
747 pixGetDimensions(pixs, &w, &h, &d);
748 if (w < 50 || h < 50)
749 return (PTA *)ERROR_PTR("too small: w and h not >= 50", __func__, NULL);
750 if (d != 8)
751 return (PTA *)ERROR_PTR("pixs not 8 bpp", __func__, NULL);
752 if (xi <= 0 || xi >= w)
753 return (PTA *)ERROR_PTR("xi not valid", __func__, NULL);
754 if (yi <= 0 || yi >= h)
755 return (PTA *)ERROR_PTR("yi not valid", __func__, NULL);
756 pixd = NULL;
757 pta = NULL;
758
759 /* Allocate stuff */
760 pixr = pixCreate(w, h, 32);
761 pixSetAll(pixr); /* initialize to max value */
762 pixp = pixCreate(w, h, 8); /* direction to parent stored as enum val */
763 lines8 = pixGetLinePtrs(pixs, NULL);
764 linep8 = pixGetLinePtrs(pixp, NULL);
765 liner32 = pixGetLinePtrs(pixr, NULL);
766 lh = lheapCreate(0, L_SORT_INCREASING); /* always remove closest pixels */
767
768 /* Prime the heap with the first pixel */
769 pixGetPixel(pixs, xi, yi, &val);
770 el = mazeelCreate(xi, yi, 0); /* don't need direction here */
771 el->distance = 0;
772 pixGetPixel(pixs, xi, yi, &val);
773 el->val = val;
774 pixSetPixel(pixr, xi, yi, 0); /* distance is 0 */
775 lheapAdd(lh, el);
776
777 /* Breadth-first search with priority queue (implemented by
778 a heap), labeling direction to parents in pixp and minimum
779 distance to visited pixels in pixr. Stop when we pull
780 the destination point (xf, yf) off the queue. */
781 while (lheapGetCount(lh) > 0) {
782 elp = (MAZEEL *)lheapRemove(lh);
783 if (!elp) {
784 L_ERROR("heap broken!!\n", __func__);
785 goto cleanup_stuff;
786 }
787 x = elp->x;
788 y = elp->y;
789 if (x == xf && y == yf) { /* exit condition */
790 LEPT_FREE(elp);
791 break;
792 }
793 distparent = (l_int32)elp->distance;
794 val = elp->val;
795 sival = val;
796
797 if (x > 0) { /* check to west */
798 vals = GET_DATA_BYTE(lines8[y], x - 1);
799 valr = GET_DATA_FOUR_BYTES(liner32[y], x - 1);
800 sivals = (l_int32)vals;
801 cost = 1 + L_ABS(sivals - sival); /* cost to move to this pixel */
802 dist = distparent + cost;
803 if (dist < valr) { /* shortest path so far to this pixel */
804 SET_DATA_FOUR_BYTES(liner32[y], x - 1, dist); /* new dist */
805 SET_DATA_BYTE(linep8[y], x - 1, DIR_EAST); /* parent to E */
806 el = mazeelCreate(x - 1, y, 0);
807 el->val = vals;
808 el->distance = dist;
809 lheapAdd(lh, el);
810 }
811 }
812 if (y > 0) { /* check north */
813 vals = GET_DATA_BYTE(lines8[y - 1], x);
814 valr = GET_DATA_FOUR_BYTES(liner32[y - 1], x);
815 sivals = (l_int32)vals;
816 cost = 1 + L_ABS(sivals - sival); /* cost to move to this pixel */
817 dist = distparent + cost;
818 if (dist < valr) { /* shortest path so far to this pixel */
819 SET_DATA_FOUR_BYTES(liner32[y - 1], x, dist); /* new dist */
820 SET_DATA_BYTE(linep8[y - 1], x, DIR_SOUTH); /* parent to S */
821 el = mazeelCreate(x, y - 1, 0);
822 el->val = vals;
823 el->distance = dist;
824 lheapAdd(lh, el);
825 }
826 }
827 if (x < w - 1) { /* check east */
828 vals = GET_DATA_BYTE(lines8[y], x + 1);
829 valr = GET_DATA_FOUR_BYTES(liner32[y], x + 1);
830 sivals = (l_int32)vals;
831 cost = 1 + L_ABS(sivals - sival); /* cost to move to this pixel */
832 dist = distparent + cost;
833 if (dist < valr) { /* shortest path so far to this pixel */
834 SET_DATA_FOUR_BYTES(liner32[y], x + 1, dist); /* new dist */
835 SET_DATA_BYTE(linep8[y], x + 1, DIR_WEST); /* parent to W */
836 el = mazeelCreate(x + 1, y, 0);
837 el->val = vals;
838 el->distance = dist;
839 lheapAdd(lh, el);
840 }
841 }
842 if (y < h - 1) { /* check south */
843 vals = GET_DATA_BYTE(lines8[y + 1], x);
844 valr = GET_DATA_FOUR_BYTES(liner32[y + 1], x);
845 sivals = (l_int32)vals;
846 cost = 1 + L_ABS(sivals - sival); /* cost to move to this pixel */
847 dist = distparent + cost;
848 if (dist < valr) { /* shortest path so far to this pixel */
849 SET_DATA_FOUR_BYTES(liner32[y + 1], x, dist); /* new dist */
850 SET_DATA_BYTE(linep8[y + 1], x, DIR_NORTH); /* parent to N */
851 el = mazeelCreate(x, y + 1, 0);
852 el->val = vals;
853 el->distance = dist;
854 lheapAdd(lh, el);
855 }
856 }
857 LEPT_FREE(elp);
858 }
859
860 lheapDestroy(&lh, TRUE);
861
862 if (ppixd) {
863 pixd = pixConvert8To32(pixs);
864 *ppixd = pixd;
865 }
866 composeRGBPixel(255, 0, 0, &rpixel); /* start point */
867 composeRGBPixel(0, 255, 0, &gpixel);
868 composeRGBPixel(0, 0, 255, &bpixel); /* end point */
869
870 x = xf;
871 y = yf;
872 pta = ptaCreate(0);
873 while (1) { /* write path onto pixd */
874 ptaAddPt(pta, x, y);
875 if (x == xi && y == yi)
876 break;
877 if (pixd)
878 pixSetPixel(pixd, x, y, gpixel);
879 pixGetPixel(pixp, x, y, &val);
880 if (val == DIR_NORTH)
881 y--;
882 else if (val == DIR_SOUTH)
883 y++;
884 else if (val == DIR_EAST)
885 x++;
886 else if (val == DIR_WEST)
887 x--;
888 pixGetPixel(pixr, x, y, &val);
889
890 #if DEBUG_PATH
891 lept_stderr("(x,y) = (%d, %d); dist = %d\n", x, y, val);
892 #endif /* DEBUG_PATH */
893
894 }
895 if (pixd) {
896 pixSetPixel(pixd, xi, yi, rpixel);
897 pixSetPixel(pixd, xf, yf, bpixel);
898 }
899
900 cleanup_stuff:
901 lheapDestroy(&lh, TRUE);
902 pixDestroy(&pixp);
903 pixDestroy(&pixr);
904 LEPT_FREE(lines8);
905 LEPT_FREE(linep8);
906 LEPT_FREE(liner32);
907 return pta;
908 }