• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 /*====================================================================*
2  -  Copyright (C) 2001 Leptonica.  All rights reserved.
3  -  This software is distributed in the hope that it will be
4  -  useful, but with NO WARRANTY OF ANY KIND.
5  -  No author or distributor accepts responsibility to anyone for the
6  -  consequences of using this software, or for whether it serves any
7  -  particular purpose or works at all, unless he or she says so in
8  -  writing.  Everyone is granted permission to copy, modify and
9  -  redistribute this source code, for commercial or non-commercial
10  -  purposes, with the following restrictions: (1) the origin of this
11  -  source code must not be misrepresented; (2) modified versions must
12  -  be plainly marked as such; and (3) this notice may not be removed
13  -  or altered from any source or modified source distribution.
14  *====================================================================*/
15 
16 /*
17  *  pix4.c
18  *
19  *    This file has these operations:
20  *
21  *      (1) Pixel histograms
22  *      (2) Foreground/background estimation
23  *      (3) Rectangle extraction
24  *
25  *    Pixel histogram, rank val, averaging and min/max
26  *           NUMA       *pixGetGrayHistogram()
27  *           NUMA       *pixGetGrayHistogramMasked()
28  *           l_int32     pixGetColorHistogram()
29  *           l_int32     pixGetColorHistogramMasked()
30  *           l_int32     pixGetRankValueMaskedRGB()
31  *           l_int32     pixGetRankValueMasked()
32  *           l_int32     pixGetAverageMaskedRGB()
33  *           l_int32     pixGetAverageMasked()
34  *           l_int32     pixGetAverageTiledRGB()
35  *           PIX        *pixGetAverageTiled()
36  *           l_int32     pixGetExtremeValue()
37  *           l_int32     pixGetMaxValueInRect()
38  *
39  *    Pixelwise aligned statistics
40  *           PIX        *pixaGetAlignedStats()
41  *           l_int32     pixaExtractColumnFromEachPix()
42  *           l_int32     pixGetRowStats()
43  *           l_int32     pixGetColumnStats()
44  *           l_int32     pixSetPixelColumn()
45  *
46  *    Foreground/background estimation
47  *           l_int32     pixThresholdForFgBg()
48  *           l_int32     pixSplitDistributionFgBg()
49  *
50  *    Measurement of properties
51  *           l_int32     pixaFindDimensions()
52  *           NUMA        pixaFindAreaPerimRatio()
53  *           l_int32     pixFindAreaPerimRatio()
54  *           NUMA        pixaFindPerimSizeRatio()
55  *           l_int32     pixFindPerimSizeRatio()
56  *           NUMA        pixaFindAreaFraction()
57  *           l_int32     pixFindAreaFraction()
58  *           NUMA        pixaFindWidthHeightRatio()
59  *           NUMA        pixaFindWidthHeightProduct()
60  *
61  *    Extract rectangle
62  *           PIX        *pixClipRectangle()
63  *           PIX        *pixClipMasked()
64  *
65  *    Clip to foreground
66  *           PIX        *pixClipToForeground()
67  *           l_int32     pixClipBoxToForeground()
68  *           l_int32     pixScanForForeground()
69  *           l_int32     pixClipBoxToEdges()
70  *           l_int32     pixScanForEdge()
71  */
72 
73 #include <stdio.h>
74 #include <stdlib.h>
75 #include <string.h>
76 #include <math.h>
77 #include "allheaders.h"
78 
79 static const l_uint32 rmask32[] = {0x0,
80     0x00000001, 0x00000003, 0x00000007, 0x0000000f,
81     0x0000001f, 0x0000003f, 0x0000007f, 0x000000ff,
82     0x000001ff, 0x000003ff, 0x000007ff, 0x00000fff,
83     0x00001fff, 0x00003fff, 0x00007fff, 0x0000ffff,
84     0x0001ffff, 0x0003ffff, 0x0007ffff, 0x000fffff,
85     0x001fffff, 0x003fffff, 0x007fffff, 0x00ffffff,
86     0x01ffffff, 0x03ffffff, 0x07ffffff, 0x0fffffff,
87     0x1fffffff, 0x3fffffff, 0x7fffffff, 0xffffffff};
88 
89 #ifndef  NO_CONSOLE_IO
90 #define  DEBUG_EDGES         0
91 #endif  /* ~NO_CONSOLE_IO */
92 
93 
94 /*------------------------------------------------------------------*
95  *                  Pixel histogram and averaging                   *
96  *------------------------------------------------------------------*/
97 /*!
98  *  pixGetGrayHistogram()
99  *
100  *      Input:  pixs (1, 2, 4, 8, 16 bpp; can be colormapped)
101  *              factor (subsampling factor; integer >= 1)
102  *      Return: na (histogram), or null on error
103  *
104  *  Notes:
105  *      (1) This generates a histogram of gray or cmapped pixels.
106  *          The image must not be rgb.
107  *      (2) If pixs does not have a colormap, the output histogram is
108  *          of size 2^d, where d is the depth of pixs.
109  *      (3) If pixs has has a colormap with color entries, the histogram
110  *          generated is of the colormap indices, and is of size 2^d.
111  *      (4) If pixs has a gray (r=g=b) colormap, a temporary 8 bpp image
112  *          without the colormap is used to construct a histogram of
113  *          size 256.
114  *      (5) Set the subsampling factor > 1 to reduce the amount of computation.
115  */
116 NUMA *
pixGetGrayHistogram(PIX * pixs,l_int32 factor)117 pixGetGrayHistogram(PIX     *pixs,
118                     l_int32  factor)
119 {
120 l_int32     i, j, w, h, d, wpl, val, size, count, colorfound;
121 l_uint32   *data, *line;
122 l_float32  *array;
123 NUMA       *na;
124 PIX        *pixg;
125 PIXCMAP    *cmap;
126 
127     PROCNAME("pixGetGrayHistogram");
128 
129     if (!pixs)
130         return (NUMA *)ERROR_PTR("pixs not defined", procName, NULL);
131     d = pixGetDepth(pixs);
132     if (d > 16)
133         return (NUMA *)ERROR_PTR("depth not in {1,2,4,8,16}", procName, NULL);
134     if (factor < 1)
135         return (NUMA *)ERROR_PTR("sampling factor < 1", procName, NULL);
136 
137     if ((cmap = pixGetColormap(pixs)) != NULL)
138         pixcmapHasColor(cmap, &colorfound);
139     if (cmap && !colorfound)
140         pixg = pixRemoveColormap(pixs, REMOVE_CMAP_TO_GRAYSCALE);
141     else
142         pixg = pixClone(pixs);
143 
144     pixGetDimensions(pixg, &w, &h, &d);
145     size = 1 << d;
146     if ((na = numaCreate(size)) == NULL)
147         return (NUMA *)ERROR_PTR("na not made", procName, NULL);
148     numaSetCount(na, size);  /* all initialized to 0.0 */
149     array = numaGetFArray(na, L_NOCOPY);
150 
151     if (d == 1) {  /* special case */
152         pixCountPixels(pixg, &count, NULL);
153         array[0] = w * h - count;
154         array[1] = count;
155         pixDestroy(&pixg);
156         return na;
157     }
158 
159     wpl = pixGetWpl(pixg);
160     data = pixGetData(pixg);
161     for (i = 0; i < h; i += factor) {
162         line = data + i * wpl;
163         switch (d)
164         {
165         case 2:
166             for (j = 0; j < w; j += factor) {
167                 val = GET_DATA_DIBIT(line, j);
168                 array[val] += 1.0;
169             }
170             break;
171         case 4:
172             for (j = 0; j < w; j += factor) {
173                 val = GET_DATA_QBIT(line, j);
174                 array[val] += 1.0;
175             }
176             break;
177         case 8:
178             for (j = 0; j < w; j += factor) {
179                 val = GET_DATA_BYTE(line, j);
180                 array[val] += 1.0;
181             }
182             break;
183         case 16:
184             for (j = 0; j < w; j += factor) {
185                 val = GET_DATA_TWO_BYTES(line, j);
186                 array[val] += 1.0;
187             }
188             break;
189         default:
190             numaDestroy(&na);
191             return (NUMA *)ERROR_PTR("illegal depth", procName, NULL);
192         }
193     }
194 
195     pixDestroy(&pixg);
196     return na;
197 }
198 
199 
200 /*!
201  *  pixGetGrayHistogramMasked()
202  *
203  *      Input:  pixs (8 bpp, or colormapped)
204  *              pixm (<optional> 1 bpp mask over which histogram is
205  *                    to be computed; use use all pixels if null)
206  *              x, y (UL corner of pixm relative to the UL corner of pixs;
207  *                    can be < 0; these values are ignored if pixm is null)
208  *              factor (subsampling factor; integer >= 1)
209  *      Return: na (histogram), or null on error
210  *
211  *  Notes:
212  *      (1) If pixs is cmapped, it is converted to 8 bpp gray.
213  *      (2) This always returns a 256-value histogram of pixel values.
214  *      (3) Set the subsampling factor > 1 to reduce the amount of computation.
215  *      (4) Clipping of pixm (if it exists) to pixs is done in the inner loop.
216  *      (5) Input x,y are ignored unless pixm exists.
217  */
218 NUMA *
pixGetGrayHistogramMasked(PIX * pixs,PIX * pixm,l_int32 x,l_int32 y,l_int32 factor)219 pixGetGrayHistogramMasked(PIX        *pixs,
220                           PIX        *pixm,
221                           l_int32     x,
222                           l_int32     y,
223                           l_int32     factor)
224 {
225 l_int32     i, j, w, h, wm, hm, dm, wplg, wplm, val;
226 l_uint32   *datag, *datam, *lineg, *linem;
227 l_float32  *array;
228 NUMA       *na;
229 PIX        *pixg;
230 
231     PROCNAME("pixGetGrayHistogramMasked");
232 
233     if (!pixm)
234         return pixGetGrayHistogram(pixs, factor);
235 
236     if (!pixs)
237         return (NUMA *)ERROR_PTR("pixs not defined", procName, NULL);
238     if (pixGetDepth(pixs) != 8 && !pixGetColormap(pixs))
239         return (NUMA *)ERROR_PTR("pixs neither 8 bpp nor colormapped",
240                                  procName, NULL);
241     pixGetDimensions(pixm, &wm, &hm, &dm);
242     if (dm != 1)
243         return (NUMA *)ERROR_PTR("pixm not 1 bpp", procName, NULL);
244     if (factor < 1)
245         return (NUMA *)ERROR_PTR("sampling factor < 1", procName, NULL);
246 
247     if ((na = numaCreate(256)) == NULL)
248         return (NUMA *)ERROR_PTR("na not made", procName, NULL);
249     numaSetCount(na, 256);  /* all initialized to 0.0 */
250     array = numaGetFArray(na, L_NOCOPY);
251 
252     if (pixGetColormap(pixs))
253         pixg = pixRemoveColormap(pixs, REMOVE_CMAP_TO_GRAYSCALE);
254     else
255         pixg = pixClone(pixs);
256     pixGetDimensions(pixg, &w, &h, NULL);
257     datag = pixGetData(pixg);
258     wplg = pixGetWpl(pixg);
259     datam = pixGetData(pixm);
260     wplm = pixGetWpl(pixm);
261 
262         /* Generate the histogram */
263     for (i = 0; i < hm; i += factor) {
264         if (y + i < 0 || y + i >= h) continue;
265         lineg = datag + (y + i) * wplg;
266         linem = datam + i * wplm;
267         for (j = 0; j < wm; j += factor) {
268             if (x + j < 0 || x + j >= w) continue;
269             if (GET_DATA_BIT(linem, j)) {
270                 val = GET_DATA_BYTE(lineg, x + j);
271                 array[val] += 1.0;
272             }
273         }
274     }
275 
276     pixDestroy(&pixg);
277     return na;
278 }
279 
280 
281 /*!
282  *  pixGetColorHistogram()
283  *
284  *      Input:  pixs (rgb or colormapped)
285  *              factor (subsampling factor; integer >= 1)
286  *              &nar (<return> red histogram)
287  *              &nag (<return> green histogram)
288  *              &nab (<return> blue histogram)
289  *      Return: 0 if OK, 1 on error
290  *
291  *  Notes:
292  *      (1) This generates a set of three 256 entry histograms,
293  *          one for each color component (r,g,b).
294  *      (2) Set the subsampling factor > 1 to reduce the amount of computation.
295  */
296 l_int32
pixGetColorHistogram(PIX * pixs,l_int32 factor,NUMA ** pnar,NUMA ** pnag,NUMA ** pnab)297 pixGetColorHistogram(PIX     *pixs,
298                      l_int32  factor,
299                      NUMA   **pnar,
300                      NUMA   **pnag,
301                      NUMA   **pnab)
302 {
303 l_int32     i, j, w, h, d, wpl, index, rval, gval, bval;
304 l_uint32   *data, *line;
305 l_float32  *rarray, *garray, *barray;
306 NUMA       *nar, *nag, *nab;
307 PIXCMAP    *cmap;
308 
309     PROCNAME("pixGetColorHistogram");
310 
311     if (!pnar || !pnag || !pnab)
312         return ERROR_INT("&nar, &nag, &nab not all defined", procName, 1);
313     *pnar = *pnag = *pnab = NULL;
314     if (!pixs)
315         return ERROR_INT("pixs not defined", procName, 1);
316     pixGetDimensions(pixs, &w, &h, &d);
317     cmap = pixGetColormap(pixs);
318     if (cmap && (d != 2 && d != 4 && d != 8))
319         return ERROR_INT("colormap and not 2, 4, or 8 bpp", procName, 1);
320     if (!cmap && d != 32)
321         return ERROR_INT("no colormap and not rgb", procName, 1);
322     if (factor < 1)
323         return ERROR_INT("sampling factor < 1", procName, 1);
324 
325         /* Set up the histogram arrays */
326     nar = numaCreate(256);
327     nag = numaCreate(256);
328     nab = numaCreate(256);
329     numaSetCount(nar, 256);
330     numaSetCount(nag, 256);
331     numaSetCount(nab, 256);
332     rarray = numaGetFArray(nar, L_NOCOPY);
333     garray = numaGetFArray(nag, L_NOCOPY);
334     barray = numaGetFArray(nab, L_NOCOPY);
335     *pnar = nar;
336     *pnag = nag;
337     *pnab = nab;
338 
339         /* Generate the color histograms */
340     data = pixGetData(pixs);
341     wpl = pixGetWpl(pixs);
342     if (cmap) {
343         for (i = 0; i < h; i += factor) {
344             line = data + i * wpl;
345             for (j = 0; j < w; j += factor) {
346                 if (d == 8)
347                     index = GET_DATA_BYTE(line, j);
348                 else if (d == 4)
349                     index = GET_DATA_QBIT(line, j);
350                 else   /* 2 bpp */
351                     index = GET_DATA_DIBIT(line, j);
352                 pixcmapGetColor(cmap, index, &rval, &gval, &bval);
353                 rarray[rval] += 1.0;
354                 garray[gval] += 1.0;
355                 barray[bval] += 1.0;
356             }
357         }
358     }
359     else {  /* 32 bpp rgb */
360         for (i = 0; i < h; i += factor) {
361             line = data + i * wpl;
362             for (j = 0; j < w; j += factor) {
363                 extractRGBValues(line[j], &rval, &gval, &bval);
364                 rarray[rval] += 1.0;
365                 garray[gval] += 1.0;
366                 barray[bval] += 1.0;
367             }
368         }
369     }
370 
371     return 0;
372 }
373 
374 
375 /*!
376  *  pixGetColorHistogramMasked()
377  *
378  *      Input:  pixs (32 bpp rgb, or colormapped)
379  *              pixm (<optional> 1 bpp mask over which histogram is
380  *                    to be computed; use use all pixels if null)
381  *              x, y (UL corner of pixm relative to the UL corner of pixs;
382  *                    can be < 0; these values are ignored if pixm is null)
383  *              factor (subsampling factor; integer >= 1)
384  *              &nar (<return> red histogram)
385  *              &nag (<return> green histogram)
386  *              &nab (<return> blue histogram)
387  *      Return: 0 if OK, 1 on error
388  *
389  *  Notes:
390  *      (1) This generates a set of three 256 entry histograms,
391  *      (2) Set the subsampling factor > 1 to reduce the amount of computation.
392  *      (3) Clipping of pixm (if it exists) to pixs is done in the inner loop.
393  *      (4) Input x,y are ignored unless pixm exists.
394  */
395 l_int32
pixGetColorHistogramMasked(PIX * pixs,PIX * pixm,l_int32 x,l_int32 y,l_int32 factor,NUMA ** pnar,NUMA ** pnag,NUMA ** pnab)396 pixGetColorHistogramMasked(PIX        *pixs,
397                            PIX        *pixm,
398                            l_int32     x,
399                            l_int32     y,
400                            l_int32     factor,
401                            NUMA      **pnar,
402                            NUMA      **pnag,
403                            NUMA      **pnab)
404 {
405 l_int32     i, j, w, h, d, wm, hm, dm, wpls, wplm, index, rval, gval, bval;
406 l_uint32   *datas, *datam, *lines, *linem;
407 l_float32  *rarray, *garray, *barray;
408 NUMA       *nar, *nag, *nab;
409 PIXCMAP    *cmap;
410 
411     PROCNAME("pixGetColorHistogramMasked");
412 
413     if (!pixm)
414         return pixGetColorHistogram(pixs, factor, pnar, pnag, pnab);
415 
416     if (!pnar || !pnag || !pnab)
417         return ERROR_INT("&nar, &nag, &nab not all defined", procName, 1);
418     *pnar = *pnag = *pnab = NULL;
419     if (!pixs)
420         return ERROR_INT("pixs not defined", procName, 1);
421     pixGetDimensions(pixs, &w, &h, &d);
422     cmap = pixGetColormap(pixs);
423     if (cmap && (d != 2 && d != 4 && d != 8))
424         return ERROR_INT("colormap and not 2, 4, or 8 bpp", procName, 1);
425     if (!cmap && d != 32)
426         return ERROR_INT("no colormap and not rgb", procName, 1);
427     pixGetDimensions(pixm, &wm, &hm, &dm);
428     if (dm != 1)
429         return ERROR_INT("pixm not 1 bpp", procName, 1);
430     if (factor < 1)
431         return ERROR_INT("sampling factor < 1", procName, 1);
432 
433         /* Set up the histogram arrays */
434     nar = numaCreate(256);
435     nag = numaCreate(256);
436     nab = numaCreate(256);
437     numaSetCount(nar, 256);
438     numaSetCount(nag, 256);
439     numaSetCount(nab, 256);
440     rarray = numaGetFArray(nar, L_NOCOPY);
441     garray = numaGetFArray(nag, L_NOCOPY);
442     barray = numaGetFArray(nab, L_NOCOPY);
443     *pnar = nar;
444     *pnag = nag;
445     *pnab = nab;
446 
447         /* Generate the color histograms */
448     datas = pixGetData(pixs);
449     wpls = pixGetWpl(pixs);
450     datam = pixGetData(pixm);
451     wplm = pixGetWpl(pixm);
452     if (cmap) {
453         for (i = 0; i < hm; i += factor) {
454             if (y + i < 0 || y + i >= h) continue;
455             lines = datas + (y + i) * wpls;
456             linem = datam + i * wplm;
457             for (j = 0; j < wm; j += factor) {
458                 if (x + j < 0 || x + j >= w) continue;
459                 if (GET_DATA_BIT(linem, j)) {
460                     if (d == 8)
461                         index = GET_DATA_BYTE(lines, x + j);
462                     else if (d == 4)
463                         index = GET_DATA_QBIT(lines, x + j);
464                     else   /* 2 bpp */
465                         index = GET_DATA_DIBIT(lines, x + j);
466                     pixcmapGetColor(cmap, index, &rval, &gval, &bval);
467                     rarray[rval] += 1.0;
468                     garray[gval] += 1.0;
469                     barray[bval] += 1.0;
470                 }
471             }
472         }
473     }
474     else {  /* 32 bpp rgb */
475         for (i = 0; i < hm; i += factor) {
476             if (y + i < 0 || y + i >= h) continue;
477             lines = datas + (y + i) * wpls;
478             linem = datam + i * wplm;
479             for (j = 0; j < wm; j += factor) {
480                 if (x + j < 0 || x + j >= w) continue;
481                 if (GET_DATA_BIT(linem, j)) {
482                     extractRGBValues(lines[x + j], &rval, &gval, &bval);
483                     rarray[rval] += 1.0;
484                     garray[gval] += 1.0;
485                     barray[bval] += 1.0;
486                 }
487             }
488         }
489     }
490 
491     return 0;
492 }
493 
494 
495 /*!
496  *  pixGetRankValueMaskedRGB()
497  *
498  *      Input:  pixs (32 bpp)
499  *              pixm (<optional> 1 bpp mask over which rank val is to be taken;
500  *                    use all pixels if null)
501  *              x, y (UL corner of pixm relative to the UL corner of pixs;
502  *                    can be < 0; these values are ignored if pixm is null)
503  *              factor (subsampling factor; integer >= 1)
504  *              rank (between 0.0 and 1.0; 1.0 is brightest, 0.0 is darkest)
505  *              &rval (<optional return> red component val for to input rank)
506  *              &gval (<optional return> green component val for to input rank)
507  *              &bval (<optional return> blue component val for to input rank)
508  *      Return: 0 if OK, 1 on error
509  *
510  *  Notes:
511  *      (1) Computes the rank component values of pixels in pixs that
512  *          are under the fg of the optional mask.  If the mask is null, it
513  *          computes the average of the pixels in pixs.
514  *      (2) Set the subsampling factor > 1 to reduce the amount of
515  *          computation.
516  *      (4) Input x,y are ignored unless pixm exists.
517  *      (5) The rank must be in [0.0 ... 1.0], where the brightest pixel
518  *          has rank 1.0.  For the median pixel value, use 0.5.
519  */
520 l_int32
pixGetRankValueMaskedRGB(PIX * pixs,PIX * pixm,l_int32 x,l_int32 y,l_int32 factor,l_float32 rank,l_float32 * prval,l_float32 * pgval,l_float32 * pbval)521 pixGetRankValueMaskedRGB(PIX        *pixs,
522                          PIX        *pixm,
523                          l_int32     x,
524                          l_int32     y,
525                          l_int32     factor,
526                          l_float32   rank,
527                          l_float32  *prval,
528                          l_float32  *pgval,
529                          l_float32  *pbval)
530 {
531 l_float32  scale;
532 PIX       *pixmt, *pixt;
533 
534     PROCNAME("pixGetRankValueMaskedRGB");
535 
536     if (!pixs)
537         return ERROR_INT("pixs not defined", procName, 1);
538     if (pixGetDepth(pixs) != 32)
539         return ERROR_INT("pixs not 32 bpp", procName, 1);
540     if (pixm && pixGetDepth(pixm) != 1)
541         return ERROR_INT("pixm not 1 bpp", procName, 1);
542     if (factor < 1)
543         return ERROR_INT("sampling factor < 1", procName, 1);
544     if (rank < 0.0 || rank > 1.0)
545         return ERROR_INT("rank not in [0.0 ... 1.0]", procName, 1);
546     if (!prval && !pgval && !pbval)
547         return ERROR_INT("no results requested", procName, 1);
548 
549     pixmt = NULL;
550     if (pixm) {
551         scale = 1.0 / (l_float32)factor;
552         pixmt = pixScale(pixm, scale, scale);
553     }
554     if (prval) {
555         pixt = pixScaleRGBToGrayFast(pixs, factor, COLOR_RED);
556         pixGetRankValueMasked(pixt, pixmt, x / factor, y / factor,
557                               factor, rank, prval, NULL);
558         pixDestroy(&pixt);
559     }
560     if (pgval) {
561         pixt = pixScaleRGBToGrayFast(pixs, factor, COLOR_GREEN);
562         pixGetRankValueMasked(pixt, pixmt, x / factor, y / factor,
563                               factor, rank, pgval, NULL);
564         pixDestroy(&pixt);
565     }
566     if (pbval) {
567         pixt = pixScaleRGBToGrayFast(pixs, factor, COLOR_BLUE);
568         pixGetRankValueMasked(pixt, pixmt, x / factor, y / factor,
569                               factor, rank, pbval, NULL);
570         pixDestroy(&pixt);
571     }
572     pixDestroy(&pixmt);
573     return 0;
574 }
575 
576 
577 /*!
578  *  pixGetRankValueMasked()
579  *
580  *      Input:  pixs (8 bpp, or colormapped)
581  *              pixm (<optional> 1 bpp mask over which rank val is to be taken;
582  *                    use all pixels if null)
583  *              x, y (UL corner of pixm relative to the UL corner of pixs;
584  *                    can be < 0; these values are ignored if pixm is null)
585  *              factor (subsampling factor; integer >= 1)
586  *              rank (between 0.0 and 1.0; 1.0 is brightest, 0.0 is darkest)
587  *              &val (<return> pixel value corresponding to input rank)
588  *              &na (<optional return> of histogram)
589  *      Return: 0 if OK, 1 on error
590  *
591  *  Notes:
592  *      (1) Computes the rank value of pixels in pixs that are under
593  *          the fg of the optional mask.  If the mask is null, it
594  *          computes the average of the pixels in pixs.
595  *      (2) Set the subsampling factor > 1 to reduce the amount of
596  *          computation.
597  *      (3) Clipping of pixm (if it exists) to pixs is done in the inner loop.
598  *      (4) Input x,y are ignored unless pixm exists.
599  *      (5) The rank must be in [0.0 ... 1.0], where the brightest pixel
600  *          has rank 1.0.  For the median pixel value, use 0.5.
601  *      (6) The histogram can optionally be returned, so that other rank
602  *          values can be extracted without recomputing the histogram.
603  *          In that case, just use
604  *              numaHistogramGetValFromRank(na, rank, &val);
605  *          on the returned Numa for additional rank values.
606  */
607 l_int32
pixGetRankValueMasked(PIX * pixs,PIX * pixm,l_int32 x,l_int32 y,l_int32 factor,l_float32 rank,l_float32 * pval,NUMA ** pna)608 pixGetRankValueMasked(PIX        *pixs,
609                       PIX        *pixm,
610                       l_int32     x,
611                       l_int32     y,
612                       l_int32     factor,
613                       l_float32   rank,
614                       l_float32  *pval,
615                       NUMA      **pna)
616 {
617 NUMA  *na;
618 
619     PROCNAME("pixGetRankValueMasked");
620 
621     if (pna)
622         *pna = NULL;
623     if (!pixs)
624         return ERROR_INT("pixs not defined", procName, 1);
625     if (pixGetDepth(pixs) != 8 && !pixGetColormap(pixs))
626         return ERROR_INT("pixs neither 8 bpp nor colormapped", procName, 1);
627     if (pixm && pixGetDepth(pixm) != 1)
628         return ERROR_INT("pixm not 1 bpp", procName, 1);
629     if (factor < 1)
630         return ERROR_INT("sampling factor < 1", procName, 1);
631     if (rank < 0.0 || rank > 1.0)
632         return ERROR_INT("rank not in [0.0 ... 1.0]", procName, 1);
633     if (!pval)
634         return ERROR_INT("&val not defined", procName, 1);
635     *pval = 0.0;  /* init */
636 
637     if ((na = pixGetGrayHistogramMasked(pixs, pixm, x, y, factor)) == NULL)
638         return ERROR_INT("na not made", procName, 1);
639     numaHistogramGetValFromRank(na, rank, pval);
640     if (pna)
641         *pna = na;
642     else
643         numaDestroy(&na);
644 
645     return 0;
646 }
647 
648 
649 /*!
650  *  pixGetAverageMaskedRGB()
651  *
652  *      Input:  pixs (32 bpp, or colormapped)
653  *              pixm (<optional> 1 bpp mask over which average is to be taken;
654  *                    use all pixels if null)
655  *              x, y (UL corner of pixm relative to the UL corner of pixs;
656  *                    can be < 0)
657  *              factor (subsampling factor; >= 1)
658  *              type (L_MEAN_ABSVAL, L_ROOT_MEAN_SQUARE,
659  *                    L_STANDARD_DEVIATION, L_VARIANCE)
660  *              &rval (<return optional> measured red value of given 'type')
661  *              &gval (<return optional> measured green value of given 'type')
662  *              &bval (<return optional> measured blue value of given 'type')
663  *      Return: 0 if OK, 1 on error
664  *
665  *  Notes:
666  *      (1) For usage, see pixGetAverageMasked().
667  *      (2) If there is a colormap, it is removed before the 8 bpp
668  *          component images are extracted.
669  */
670 l_int32
pixGetAverageMaskedRGB(PIX * pixs,PIX * pixm,l_int32 x,l_int32 y,l_int32 factor,l_int32 type,l_float32 * prval,l_float32 * pgval,l_float32 * pbval)671 pixGetAverageMaskedRGB(PIX        *pixs,
672                        PIX        *pixm,
673                        l_int32     x,
674                        l_int32     y,
675                        l_int32     factor,
676                        l_int32     type,
677                        l_float32  *prval,
678                        l_float32  *pgval,
679                        l_float32  *pbval)
680 {
681 l_int32   color;
682 PIX      *pixt;
683 PIXCMAP  *cmap;
684 
685     PROCNAME("pixGetAverageMaskedRGB");
686 
687     if (!pixs)
688         return ERROR_INT("pixs not defined", procName, 1);
689     color = 0;
690     cmap = pixGetColormap(pixs);
691     if (pixGetDepth(pixs) != 32 && !cmap)
692         return ERROR_INT("pixs neither 32 bpp nor colormapped", procName, 1);
693     if (pixm && pixGetDepth(pixm) != 1)
694         return ERROR_INT("pixm not 1 bpp", procName, 1);
695     if (factor < 1)
696         return ERROR_INT("subsampling factor < 1", procName, 1);
697     if (type != L_MEAN_ABSVAL && type != L_ROOT_MEAN_SQUARE &&
698         type != L_STANDARD_DEVIATION && type != L_VARIANCE)
699         return ERROR_INT("invalid measure type", procName, 1);
700     if (!prval && !pgval && !pbval)
701         return ERROR_INT("no values requested", procName, 1);
702 
703     if (prval) {
704         if (cmap)
705             pixt = pixGetRGBComponentCmap(pixs, COLOR_RED);
706         else
707             pixt = pixGetRGBComponent(pixs, COLOR_RED);
708         pixGetAverageMasked(pixt, pixm, x, y, factor, type, prval);
709         pixDestroy(&pixt);
710     }
711     if (pgval) {
712         if (cmap)
713             pixt = pixGetRGBComponentCmap(pixs, COLOR_GREEN);
714         else
715             pixt = pixGetRGBComponent(pixs, COLOR_GREEN);
716         pixGetAverageMasked(pixt, pixm, x, y, factor, type, pgval);
717         pixDestroy(&pixt);
718     }
719     if (pbval) {
720         if (cmap)
721             pixt = pixGetRGBComponentCmap(pixs, COLOR_BLUE);
722         else
723             pixt = pixGetRGBComponent(pixs, COLOR_BLUE);
724         pixGetAverageMasked(pixt, pixm, x, y, factor, type, pbval);
725         pixDestroy(&pixt);
726     }
727 
728     return 0;
729 }
730 
731 
732 /*!
733  *  pixGetAverageMasked()
734  *
735  *      Input:  pixs (8 or 16 bpp, or colormapped)
736  *              pixm (<optional> 1 bpp mask over which average is to be taken;
737  *                    use all pixels if null)
738  *              x, y (UL corner of pixm relative to the UL corner of pixs;
739  *                    can be < 0)
740  *              factor (subsampling factor; >= 1)
741  *              type (L_MEAN_ABSVAL, L_ROOT_MEAN_SQUARE,
742  *                    L_STANDARD_DEVIATION, L_VARIANCE)
743  *              &val (<return> measured value of given 'type')
744  *      Return: 0 if OK, 1 on error
745  *
746  *  Notes:
747  *      (1) Use L_MEAN_ABSVAL to get the average value of pixels in pixs
748  *          that are under the fg of the optional mask.  If the mask
749  *          is null, it finds the average of the pixels in pixs.
750  *      (2) Likewise, use L_ROOT_MEAN_SQUARE to get the rms value of
751  *          pixels in pixs, either masked or not; L_STANDARD_DEVIATION
752  *          to get the standard deviation from the mean of the pixels;
753  *          L_VARIANCE to get the average squared difference from the
754  *          expected value.  The variance is the square of the stdev.
755  *          For the standard deviation, we use
756  *              sqrt(<(<x> - x)>^2) = sqrt(<x^2> - <x>^2)
757  *      (3) Set the subsampling factor > 1 to reduce the amount of
758  *          computation.
759  *      (4) Clipping of pixm (if it exists) to pixs is done in the inner loop.
760  *      (5) Input x,y are ignored unless pixm exists.
761  */
762 l_int32
pixGetAverageMasked(PIX * pixs,PIX * pixm,l_int32 x,l_int32 y,l_int32 factor,l_int32 type,l_float32 * pval)763 pixGetAverageMasked(PIX        *pixs,
764                     PIX        *pixm,
765                     l_int32     x,
766                     l_int32     y,
767                     l_int32     factor,
768                     l_int32     type,
769                     l_float32  *pval)
770 {
771 l_int32    i, j, w, h, d, wm, hm, wplg, wplm, val, count;
772 l_uint32  *datag, *datam, *lineg, *linem;
773 l_float64  sumave, summs, ave, meansq, var;
774 PIX       *pixg;
775 
776     PROCNAME("pixGetAverageMasked");
777 
778     if (!pixs)
779         return ERROR_INT("pixs not defined", procName, 1);
780     d = pixGetDepth(pixs);
781     if (d != 8 && d != 16 && !pixGetColormap(pixs))
782         return ERROR_INT("pixs not 8 or 16 bpp or colormapped", procName, 1);
783     if (pixm && pixGetDepth(pixm) != 1)
784         return ERROR_INT("pixm not 1 bpp", procName, 1);
785     if (factor < 1)
786         return ERROR_INT("subsampling factor < 1", procName, 1);
787     if (type != L_MEAN_ABSVAL && type != L_ROOT_MEAN_SQUARE &&
788         type != L_STANDARD_DEVIATION && type != L_VARIANCE)
789         return ERROR_INT("invalid measure type", procName, 1);
790     if (!pval)
791         return ERROR_INT("&val not defined", procName, 1);
792     *pval = 0.0;  /* init */
793 
794     if (pixGetColormap(pixs))
795         pixg = pixRemoveColormap(pixs, REMOVE_CMAP_TO_GRAYSCALE);
796     else
797         pixg = pixClone(pixs);
798     pixGetDimensions(pixg, &w, &h, &d);
799     datag = pixGetData(pixg);
800     wplg = pixGetWpl(pixg);
801 
802     sumave = summs = 0.0;
803     count = 0;
804     if (!pixm) {
805         for (i = 0; i < h; i += factor) {
806             lineg = datag + i * wplg;
807             for (j = 0; j < w; j += factor) {
808                 if (d == 8)
809                     val = GET_DATA_BYTE(lineg, j);
810                 else  /* d == 16 */
811                     val = GET_DATA_TWO_BYTES(lineg, j);
812                 if (type != L_ROOT_MEAN_SQUARE)
813                     sumave += val;
814                 if (type != L_MEAN_ABSVAL)
815                     summs += val * val;
816                 count++;
817             }
818         }
819     }
820     else {
821         pixGetDimensions(pixm, &wm, &hm, NULL);
822         datam = pixGetData(pixm);
823         wplm = pixGetWpl(pixm);
824         for (i = 0; i < hm; i += factor) {
825             if (y + i < 0 || y + i >= h) continue;
826             lineg = datag + (y + i) * wplg;
827             linem = datam + i * wplm;
828             for (j = 0; j < wm; j += factor) {
829                 if (x + j < 0 || x + j >= w) continue;
830                 if (GET_DATA_BIT(linem, j)) {
831                     if (d == 8)
832                         val = GET_DATA_BYTE(lineg, x + j);
833                     else  /* d == 16 */
834                         val = GET_DATA_TWO_BYTES(lineg, x + j);
835                     if (type != L_ROOT_MEAN_SQUARE)
836                         sumave += val;
837                     if (type != L_MEAN_ABSVAL)
838                         summs += val * val;
839                     count++;
840                 }
841             }
842         }
843     }
844 
845     pixDestroy(&pixg);
846     if (count == 0)
847         return ERROR_INT("no pixels sampled", procName, 1);
848     ave = sumave / (l_float64)count;
849     meansq = summs / (l_float64)count;
850     var = meansq - ave * ave;
851     if (type == L_MEAN_ABSVAL)
852         *pval = (l_float32)ave;
853     else if (type == L_ROOT_MEAN_SQUARE)
854         *pval = (l_float32)sqrt(meansq);
855     else if (type == L_STANDARD_DEVIATION)
856         *pval = (l_float32)sqrt(var);
857     else  /* type == L_VARIANCE */
858         *pval = (l_float32)var;
859 
860     return 0;
861 }
862 
863 
864 /*!
865  *  pixGetAverageTiledRGB()
866  *
867  *      Input:  pixs (32 bpp, or colormapped)
868  *              sx, sy (tile size; must be at least 2 x 2)
869  *              type (L_MEAN_ABSVAL, L_ROOT_MEAN_SQUARE, L_STANDARD_DEVIATION)
870  *              &pixr (<optional return> tiled 'average' of red component)
871  *              &pixg (<optional return> tiled 'average' of green component)
872  *              &pixb (<optional return> tiled 'average' of blue component)
873  *      Return: 0 if OK, 1 on error
874  *
875  *  Notes:
876  *      (1) For usage, see pixGetAverageTiled().
877  *      (2) If there is a colormap, it is removed before the 8 bpp
878  *          component images are extracted.
879  */
880 l_int32
pixGetAverageTiledRGB(PIX * pixs,l_int32 sx,l_int32 sy,l_int32 type,PIX ** ppixr,PIX ** ppixg,PIX ** ppixb)881 pixGetAverageTiledRGB(PIX     *pixs,
882                       l_int32  sx,
883                       l_int32  sy,
884                       l_int32  type,
885                       PIX    **ppixr,
886                       PIX    **ppixg,
887                       PIX    **ppixb)
888 {
889 PIX      *pixt;
890 PIXCMAP  *cmap;
891 
892     PROCNAME("pixGetAverageTiledRGB");
893 
894     if (!pixs)
895         return ERROR_INT("pixs not defined", procName, 1);
896     cmap = pixGetColormap(pixs);
897     if (pixGetDepth(pixs) != 32 && !cmap)
898         return ERROR_INT("pixs neither 32 bpp nor colormapped", procName, 1);
899     if (sx < 2 || sy < 2)
900         return ERROR_INT("sx and sy not both > 1", procName, 1);
901     if (type != L_MEAN_ABSVAL && type != L_ROOT_MEAN_SQUARE &&
902         type != L_STANDARD_DEVIATION)
903         return ERROR_INT("invalid measure type", procName, 1);
904     if (!ppixr && !ppixg && !ppixb)
905         return ERROR_INT("no returned data requested", procName, 1);
906 
907     if (ppixr) {
908         if (cmap)
909             pixt = pixGetRGBComponentCmap(pixs, COLOR_RED);
910         else
911             pixt = pixGetRGBComponent(pixs, COLOR_RED);
912         *ppixr = pixGetAverageTiled(pixt, sx, sy, type);
913         pixDestroy(&pixt);
914     }
915     if (ppixg) {
916         if (cmap)
917             pixt = pixGetRGBComponentCmap(pixs, COLOR_GREEN);
918         else
919             pixt = pixGetRGBComponent(pixs, COLOR_GREEN);
920         *ppixg = pixGetAverageTiled(pixt, sx, sy, type);
921         pixDestroy(&pixt);
922     }
923     if (ppixb) {
924         if (cmap)
925             pixt = pixGetRGBComponentCmap(pixs, COLOR_BLUE);
926         else
927             pixt = pixGetRGBComponent(pixs, COLOR_BLUE);
928         *ppixb = pixGetAverageTiled(pixt, sx, sy, type);
929         pixDestroy(&pixt);
930     }
931 
932     return 0;
933 }
934 
935 
936 /*!
937  *  pixGetAverageTiled()
938  *
939  *      Input:  pixs (8 bpp, or colormapped)
940  *              sx, sy (tile size; must be at least 2 x 2)
941  *              type (L_MEAN_ABSVAL, L_ROOT_MEAN_SQUARE, L_STANDARD_DEVIATION)
942  *      Return: pixd (average values in each tile), or null on error
943  *
944  *  Notes:
945  *      (1) Only computes for tiles that are entirely contained in pixs.
946  *      (2) Use L_MEAN_ABSVAL to get the average abs value within the tile;
947  *          L_ROOT_MEAN_SQUARE to get the rms value within each tile;
948  *          L_STANDARD_DEVIATION to get the standard dev. from the average
949  *          within each tile.
950  *      (3) If colormapped, converts to 8 bpp gray.
951  */
952 PIX *
pixGetAverageTiled(PIX * pixs,l_int32 sx,l_int32 sy,l_int32 type)953 pixGetAverageTiled(PIX     *pixs,
954                    l_int32  sx,
955                    l_int32  sy,
956                    l_int32  type)
957 {
958 l_int32    i, j, k, m, w, h, wd, hd, d, pos, wplt, wpld, valt;
959 l_uint32  *datat, *datad, *linet, *lined, *startt;
960 l_float64  sumave, summs, ave, meansq, normfact;
961 PIX       *pixt, *pixd;
962 
963     PROCNAME("pixGetAverageTiled");
964 
965     if (!pixs)
966         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
967     pixGetDimensions(pixs, &w, &h, &d);
968     if (d != 8 && !pixGetColormap(pixs))
969         return (PIX *)ERROR_PTR("pixs not 8 bpp or cmapped", procName, NULL);
970     if (sx < 2 || sy < 2)
971         return (PIX *)ERROR_PTR("sx and sy not both > 1", procName, NULL);
972     wd = w / sx;
973     hd = h / sy;
974     if (wd < 1 || hd < 1)
975         return (PIX *)ERROR_PTR("wd or hd == 0", procName, NULL);
976     if (type != L_MEAN_ABSVAL && type != L_ROOT_MEAN_SQUARE &&
977         type != L_STANDARD_DEVIATION)
978         return (PIX *)ERROR_PTR("invalid measure type", procName, NULL);
979 
980     pixt = pixRemoveColormap(pixs, REMOVE_CMAP_TO_GRAYSCALE);
981     pixd = pixCreate(wd, hd, 8);
982     datat = pixGetData(pixt);
983     wplt = pixGetWpl(pixt);
984     datad = pixGetData(pixd);
985     wpld = pixGetWpl(pixd);
986     normfact = 1. / (l_float64)(sx * sy);
987     for (i = 0; i < hd; i++) {
988         lined = datad + i * wpld;
989         linet = datat + i * sy * wplt;
990         for (j = 0; j < wd; j++) {
991             if (type == L_MEAN_ABSVAL || type == L_STANDARD_DEVIATION) {
992                 sumave = 0.0;
993                 for (k = 0; k < sy; k++) {
994                     startt = linet + k * wplt;
995                     for (m = 0; m < sx; m++) {
996                         pos = j * sx + m;
997                         valt = GET_DATA_BYTE(startt, pos);
998                         sumave += valt;
999                     }
1000                 }
1001                 ave = normfact * sumave;
1002             }
1003             if (type == L_ROOT_MEAN_SQUARE || type == L_STANDARD_DEVIATION) {
1004                 summs = 0.0;
1005                 for (k = 0; k < sy; k++) {
1006                     startt = linet + k * wplt;
1007                     for (m = 0; m < sx; m++) {
1008                         pos = j * sx + m;
1009                         valt = GET_DATA_BYTE(startt, pos);
1010                         summs += valt * valt;
1011                     }
1012                 }
1013                 meansq = normfact * summs;
1014             }
1015             if (type == L_MEAN_ABSVAL)
1016                 valt = (l_int32)(ave + 0.5);
1017             else if (type == L_ROOT_MEAN_SQUARE)
1018                 valt = (l_int32)(sqrt(meansq) + 0.5);
1019             else  /* type == L_STANDARD_DEVIATION */
1020                 valt = (l_int32)(sqrt(meansq - ave * ave) + 0.5);
1021             SET_DATA_BYTE(lined, j, valt);
1022         }
1023     }
1024 
1025     pixDestroy(&pixt);
1026     return pixd;
1027 }
1028 
1029 
1030 /*!
1031  *  pixGetExtremeValue()
1032  *
1033  *      Input:  pixs (8 bpp grayscale, 32 bpp rgb, or colormapped)
1034  *              factor (subsampling factor; >= 1; ignored if colormapped)
1035  *              type (L_CHOOSE_MIN or L_CHOOSE_MAX)
1036  *              &rval (<optional return> red component)
1037  *              &gval (<optional return> green component)
1038  *              &bval (<optional return> blue component)
1039  *              &grayval (<optional return> min gray value)
1040  *      Return: 0 if OK, 1 on error
1041  *
1042  *  Notes:
1043  *      (1) If pixs is grayscale, the result is returned in &grayval.
1044  *          Otherwise, if there is a colormap or d == 32,
1045  *          each requested color component is returned.  At least
1046  *          one color component (address) must be input.
1047  */
1048 l_int32
pixGetExtremeValue(PIX * pixs,l_int32 factor,l_int32 type,l_int32 * prval,l_int32 * pgval,l_int32 * pbval,l_int32 * pgrayval)1049 pixGetExtremeValue(PIX      *pixs,
1050                    l_int32   factor,
1051                    l_int32   type,
1052                    l_int32  *prval,
1053                    l_int32  *pgval,
1054                    l_int32  *pbval,
1055                    l_int32  *pgrayval)
1056 {
1057 l_int32    i, j, w, h, d, wpl;
1058 l_int32    val, extval, rval, gval, bval, extrval, extgval, extbval;
1059 l_uint32   pixel;
1060 l_uint32  *data, *line;
1061 PIXCMAP   *cmap;
1062 
1063     PROCNAME("pixGetExtremeValue");
1064 
1065     if (!pixs)
1066         return ERROR_INT("pixs not defined", procName, 1);
1067 
1068     cmap = pixGetColormap(pixs);
1069     if (cmap)
1070         return pixcmapGetExtremeValue(cmap, type, prval, pgval, pbval);
1071 
1072     pixGetDimensions(pixs, &w, &h, &d);
1073     if (type != L_CHOOSE_MIN && type != L_CHOOSE_MAX)
1074         return ERROR_INT("invalid type", procName, 1);
1075     if (factor < 1)
1076         return ERROR_INT("subsampling factor < 1", procName, 1);
1077     if (d != 8 && d != 32)
1078         return ERROR_INT("pixs not 8 or 32 bpp", procName, 1);
1079     if (d == 8 && !pgrayval)
1080         return ERROR_INT("can't return result in grayval", procName, 1);
1081     if (d == 32 && !prval && !pgval && !pbval)
1082         return ERROR_INT("can't return result in r/g/b-val", procName, 1);
1083 
1084     data = pixGetData(pixs);
1085     wpl = pixGetWpl(pixs);
1086     if (d == 8) {
1087         if (type == L_CHOOSE_MIN)
1088             extval = 100000;
1089         else  /* get max */
1090             extval = 0;
1091 
1092         for (i = 0; i < h; i += factor) {
1093             line = data + i * wpl;
1094             for (j = 0; j < w; j += factor) {
1095                 val = GET_DATA_BYTE(line, j);
1096                 if ((type == L_CHOOSE_MIN && val < extval) ||
1097                     (type == L_CHOOSE_MAX && val > extval))
1098                     extval = val;
1099             }
1100         }
1101         *pgrayval = extval;
1102         return 0;
1103     }
1104 
1105         /* 32 bpp rgb */
1106     if (type == L_CHOOSE_MIN) {
1107         extrval = 100000;
1108         extgval = 100000;
1109         extbval = 100000;
1110     }
1111     else {
1112         extrval = 0;
1113         extgval = 0;
1114         extbval = 0;
1115     }
1116     for (i = 0; i < h; i += factor) {
1117         line = data + i * wpl;
1118         for (j = 0; j < w; j += factor) {
1119             pixel = line[j];
1120             if (prval) {
1121                 rval = (pixel >> L_RED_SHIFT) & 0xff;
1122                 if ((type == L_CHOOSE_MIN && rval < extrval) ||
1123                     (type == L_CHOOSE_MAX && rval > extrval))
1124                     extrval = rval;
1125             }
1126             if (pgval) {
1127                 gval = (pixel >> L_GREEN_SHIFT) & 0xff;
1128                 if ((type == L_CHOOSE_MIN && gval < extgval) ||
1129                     (type == L_CHOOSE_MAX && gval > extgval))
1130                     extgval = gval;
1131             }
1132             if (pbval) {
1133                 bval = (pixel >> L_BLUE_SHIFT) & 0xff;
1134                 if ((type == L_CHOOSE_MIN && bval < extbval) ||
1135                     (type == L_CHOOSE_MAX && bval > extbval))
1136                     extbval = bval;
1137             }
1138         }
1139     }
1140     if (prval) *prval = extrval;
1141     if (pgval) *pgval = extgval;
1142     if (pbval) *pbval = extbval;
1143     return 0;
1144 }
1145 
1146 
1147 /*!
1148  *  pixGetMaxValueInRect()
1149  *
1150  *      Input:  pixs (8 bpp or 32 bpp grayscale; no color space components)
1151  *              box (<optional> region; set box = NULL to use entire pixs)
1152  *              &maxval (<optional return> max value in region)
1153  *              &xmax (<optional return> x location of max value)
1154  *              &ymax (<optional return> y location of max value)
1155  *      Return: 0 if OK, 1 on error
1156  *
1157  *  Notes:
1158  *      (1) This can be used to find the maximum and its location
1159  *          in a 2-dimensional histogram, where the x and y directions
1160  *          represent two color components (e.g., saturation and hue).
1161  *      (2) Note that here a 32 bpp pixs has pixel values that are simply
1162  *          numbers.  They are not 8 bpp components in a colorspace.
1163  */
1164 l_int32
pixGetMaxValueInRect(PIX * pixs,BOX * box,l_uint32 * pmaxval,l_int32 * pxmax,l_int32 * pymax)1165 pixGetMaxValueInRect(PIX       *pixs,
1166                      BOX       *box,
1167                      l_uint32  *pmaxval,
1168                      l_int32   *pxmax,
1169                      l_int32   *pymax)
1170 {
1171 l_int32    i, j, w, h, d, wpl, bw, bh;
1172 l_int32    xstart, ystart, xend, yend, xmax, ymax;
1173 l_uint32   val, maxval;
1174 l_uint32  *data, *line;
1175 
1176     PROCNAME("pixGetMaxValueInRect");
1177 
1178     if (!pmaxval && !pxmax && !pymax)
1179         return ERROR_INT("nothing to do", procName, 1);
1180     if (pmaxval) *pmaxval = 0;
1181     if (pxmax) *pxmax = 0;
1182     if (pymax) *pymax = 0;
1183     if (!pixs)
1184         return ERROR_INT("pixs not defined", procName, 1);
1185     if (pixGetColormap(pixs) != NULL)
1186         return ERROR_INT("pixs has colormap", procName, 1);
1187     pixGetDimensions(pixs, &w, &h, &d);
1188     if (d != 8 && d != 32)
1189         return ERROR_INT("pixs not 8 or 32 bpp", procName, 1);
1190 
1191     xstart = ystart = 0;
1192     xend = w - 1;
1193     yend = h - 1;
1194     if (box) {
1195         boxGetGeometry(box, &xstart, &ystart, &bw, &bh);
1196         xend = xstart + bw - 1;
1197         yend = ystart + bh - 1;
1198     }
1199 
1200     data = pixGetData(pixs);
1201     wpl = pixGetWpl(pixs);
1202     maxval = 0;
1203     xmax = ymax = 0;
1204     for (i = ystart; i <= yend; i++) {
1205         line = data + i * wpl;
1206         for (j = xstart; j <= xend; j++) {
1207             if (d == 8)
1208                 val = GET_DATA_BYTE(line, j);
1209             else  /* d == 32 */
1210                 val = line[j];
1211             if (val > maxval) {
1212                 maxval = val;
1213                 xmax = j;
1214                 ymax = i;
1215             }
1216         }
1217     }
1218     if (maxval == 0) {  /* no counts; pick the center of the rectangle */
1219         xmax = (xstart + xend) / 2;
1220         ymax = (ystart + yend) / 2;
1221     }
1222 
1223     if (pmaxval) *pmaxval = maxval;
1224     if (pxmax) *pxmax = xmax;
1225     if (pymax) *pymax = ymax;
1226     return 0;
1227 }
1228 
1229 
1230 /*-------------------------------------------------------------*
1231  *                 Pixelwise aligned statistics                *
1232  *-------------------------------------------------------------*/
1233 /*!
1234  *  pixaGetAlignedStats()
1235  *
1236  *      Input:  pixa (of identically sized, 8 bpp pix; not cmapped)
1237  *              type (L_MEAN_ABSVAL, L_MEDIAN_VAL, L_MODE_VAL, L_MODE_COUNT)
1238  *              nbins (of histogram for median and mode; ignored for mean)
1239  *              thresh (on histogram for mode val; ignored for all other types)
1240  *      Return: pix (with pixelwise aligned stats), or null on error.
1241  *
1242  *  Notes:
1243  *      (1) Each pixel in the returned pix represents an average
1244  *          (or median, or mode) over the corresponding pixels in each
1245  *          pix in the pixa.
1246  *      (2) The @thresh parameter works with L_MODE_VAL only, and
1247  *          sets a minimum occupancy of the mode bin.
1248  *          If the occupancy of the mode bin is less than @thresh, the
1249  *          mode value is returned as 0.  To always return the actual
1250  *          mode value, set @thresh = 0.  See pixGetRowStats().
1251  */
1252 PIX *
pixaGetAlignedStats(PIXA * pixa,l_int32 type,l_int32 nbins,l_int32 thresh)1253 pixaGetAlignedStats(PIXA     *pixa,
1254                     l_int32   type,
1255                     l_int32   nbins,
1256                     l_int32   thresh)
1257 {
1258 l_int32     j, n, w, h, d;
1259 l_float32  *colvect;
1260 PIX        *pixt, *pixd;
1261 
1262     PROCNAME("pixaGetAlignedStats");
1263 
1264     if (!pixa)
1265         return (PIX *)ERROR_PTR("pixa not defined", procName, NULL);
1266     if (type != L_MEAN_ABSVAL && type != L_MEDIAN_VAL &&
1267         type != L_MODE_VAL && type != L_MODE_COUNT)
1268         return (PIX *)ERROR_PTR("invalid type", procName, NULL);
1269     n = pixaGetCount(pixa);
1270     if (n == 0)
1271         return (PIX *)ERROR_PTR("no pix in pixa", procName, NULL);
1272     pixaGetPixDimensions(pixa, 0, &w, &h, &d);
1273     if (d != 8)
1274         return (PIX *)ERROR_PTR("pix not 8 bpp", procName, NULL);
1275 
1276     pixd = pixCreate(w, h, 8);
1277     pixt = pixCreate(n, h, 8);
1278     colvect = (l_float32 *)CALLOC(h, sizeof(l_float32));
1279     for (j = 0; j < w; j++) {
1280         pixaExtractColumnFromEachPix(pixa, j, pixt);
1281         pixGetRowStats(pixt, type, nbins, thresh, colvect);
1282         pixSetPixelColumn(pixd, j, colvect);
1283     }
1284 
1285     FREE(colvect);
1286     pixDestroy(&pixt);
1287     return pixd;
1288 }
1289 
1290 
1291 /*!
1292  *  pixaExtractColumnFromEachPix()
1293  *
1294  *      Input:  pixa (of identically sized, 8 bpp; not cmapped)
1295  *              col (column index)
1296  *              pixd (pix into which each column is inserted)
1297  *      Return: 0 if OK, 1 on error
1298  */
1299 l_int32
pixaExtractColumnFromEachPix(PIXA * pixa,l_int32 col,PIX * pixd)1300 pixaExtractColumnFromEachPix(PIXA    *pixa,
1301                              l_int32  col,
1302                              PIX     *pixd)
1303 {
1304 l_int32    i, k, n, w, h, ht, val, wplt, wpld;
1305 l_uint32  *datad, *datat;
1306 PIX       *pixt;
1307 
1308     PROCNAME("pixaExtractColumnFromEachPix");
1309 
1310     if (!pixa)
1311         return ERROR_INT("pixa not defined", procName, 1);
1312     if (!pixd || pixGetDepth(pixd) != 8)
1313         return ERROR_INT("pixa not defined or not 8 bpp", procName, 1);
1314     n = pixaGetCount(pixa);
1315     pixGetDimensions(pixd, &w, &h, NULL);
1316     if (n != w)
1317         return ERROR_INT("pix width != n", procName, 1);
1318     pixt = pixaGetPix(pixa, 0, L_CLONE);
1319     wplt = pixGetWpl(pixt);
1320     pixGetDimensions(pixt, NULL, &ht, NULL);
1321     pixDestroy(&pixt);
1322     if (h != ht)
1323         return ERROR_INT("pixd height != column height", procName, 1);
1324 
1325     datad = pixGetData(pixd);
1326     wpld = pixGetWpl(pixd);
1327     for (k = 0; k < n; k++) {
1328         pixt = pixaGetPix(pixa, k, L_CLONE);
1329         datat = pixGetData(pixt);
1330         for (i = 0; i < h; i++) {
1331             val = GET_DATA_BYTE(datat, col);
1332             SET_DATA_BYTE(datad + i * wpld, k, val);
1333             datat += wplt;
1334         }
1335         pixDestroy(&pixt);
1336     }
1337 
1338     return 0;
1339 }
1340 
1341 
1342 /*!
1343  *  pixGetRowStats()
1344  *
1345  *      Input:  pixs (8 bpp; not cmapped)
1346  *              type (L_MEAN_ABSVAL, L_MEDIAN_VAL, L_MODE_VAL, L_MODE_COUNT)
1347  *              nbins (of histogram for median and mode; ignored for mean)
1348  *              thresh (on histogram for mode; ignored for mean and median)
1349  *              colvect (vector of results gathered across the rows of pixs)
1350  *      Return: 0 if OK, 1 on error
1351  *
1352  *  Notes:
1353  *      (1) This computes a column vector of statistics using each
1354  *          row of a Pix.  The result is put in @colvect.
1355  *      (2) The @thresh parameter works with L_MODE_VAL only, and
1356  *          sets a minimum occupancy of the mode bin.
1357  *          If the occupancy of the mode bin is less than @thresh, the
1358  *          mode value is returned as 0.  To always return the actual
1359  *          mode value, set @thresh = 0.
1360  *      (3) What is the meaning of this @thresh parameter?
1361  *          For each row, the total count in the histogram is w, the
1362  *          image width.  So @thresh, relative to w, gives a measure
1363  *          of the ratio of the bin width to the width of the distribution.
1364  *          The larger @thresh, the narrower the distribution must be
1365  *          for the mode value to be returned (instead of returning 0).
1366  *      (4) If the Pix consists of a set of corresponding columns,
1367  *          one for each Pix in a Pixa, the width of the Pix is the
1368  *          number of Pix in the Pixa and the column vector can
1369  *          be stored as a column in a Pix of the same size as
1370  *          each Pix in the Pixa.
1371  */
1372 l_int32
pixGetRowStats(PIX * pixs,l_int32 type,l_int32 nbins,l_int32 thresh,l_float32 * colvect)1373 pixGetRowStats(PIX        *pixs,
1374                l_int32     type,
1375                l_int32     nbins,
1376                l_int32     thresh,
1377                l_float32  *colvect)
1378 {
1379 l_int32    i, j, k, w, h, val, wpls, sum, target, max, modeval;
1380 l_int32   *histo, *gray2bin, *bin2gray;
1381 l_uint32  *lines, *datas;
1382 
1383     PROCNAME("pixGetRowStats");
1384 
1385     if (!pixs || pixGetDepth(pixs) != 8)
1386         return ERROR_INT("pixs not defined or not 8 bpp", procName, 1);
1387     if (!colvect)
1388         return ERROR_INT("colvect not defined", procName, 1);
1389     if (type != L_MEAN_ABSVAL && type != L_MEDIAN_VAL &&
1390         type != L_MODE_VAL && type != L_MODE_COUNT)
1391         return ERROR_INT("invalid type", procName, 1);
1392     if (type != L_MEAN_ABSVAL && (nbins < 1 || nbins > 256))
1393         return ERROR_INT("invalid nbins", procName, 1);
1394     pixGetDimensions(pixs, &w, &h, NULL);
1395 
1396     datas = pixGetData(pixs);
1397     wpls = pixGetWpl(pixs);
1398     if (type == L_MEAN_ABSVAL) {
1399         for (i = 0; i < h; i++) {
1400             sum = 0;
1401             lines = datas + i * wpls;
1402             for (j = 0; j < w; j++)
1403                 sum += GET_DATA_BYTE(lines, j);
1404             colvect[i] = (l_float32)sum / (l_float32)w;
1405         }
1406         return 0;
1407     }
1408 
1409         /* We need a histogram; binwidth ~ 256 / nbins */
1410     histo = (l_int32 *)CALLOC(nbins, sizeof(l_int32));
1411     gray2bin = (l_int32 *)CALLOC(256, sizeof(l_int32));
1412     bin2gray = (l_int32 *)CALLOC(nbins, sizeof(l_int32));
1413     for (i = 0; i < 256; i++)  /* gray value --> histo bin */
1414         gray2bin[i] = (i * nbins) / 256;
1415     for (i = 0; i < nbins; i++)  /* histo bin --> gray value */
1416         bin2gray[i] = (i * 255 + 128) / nbins;
1417 
1418     for (i = 0; i < h; i++) {
1419         lines = datas + i * wpls;
1420         for (k = 0; k < nbins; k++)
1421             histo[k] = 0;
1422         for (j = 0; j < w; j++) {
1423             val = GET_DATA_BYTE(lines, j);
1424             histo[gray2bin[val]]++;
1425         }
1426 
1427         if (type == L_MEDIAN_VAL) {
1428             sum = 0;
1429             target = w / 2;
1430             for (k = 0; k < nbins; k++) {
1431                 sum += histo[k];
1432                 if (sum >= target) {
1433                     colvect[i] = bin2gray[k];
1434                     break;
1435                 }
1436             }
1437         }
1438         else if (type == L_MODE_VAL) {
1439             max = 0;
1440             modeval = 0;
1441             for (k = 0; k < nbins; k++) {
1442                 if (histo[k] > max) {
1443                     max = histo[k];
1444                     modeval = k;
1445                 }
1446             }
1447             if (max < thresh)
1448                 colvect[i] = 0;
1449             else
1450                 colvect[i] = bin2gray[modeval];
1451         }
1452         else {  /* type == L_MODE_COUNT */
1453             max = 0;
1454             modeval = 0;
1455             for (k = 0; k < nbins; k++) {
1456                 if (histo[k] > max) {
1457                     max = histo[k];
1458                     modeval = k;
1459                 }
1460             }
1461             colvect[i] = max;
1462         }
1463     }
1464 
1465     FREE(histo);
1466     FREE(gray2bin);
1467     FREE(bin2gray);
1468     return 0;
1469 }
1470 
1471 
1472 /*!
1473  *  pixGetColumnStats()
1474  *
1475  *      Input:  pixs (8 bpp; not cmapped)
1476  *              type (L_MEAN_ABSVAL, L_MEDIAN_VAL, L_MODE_VAL, L_MODE_COUNT)
1477  *              nbins (of histogram for median and mode; ignored for mean)
1478  *              thresh (on histogram for mode val; ignored for all other types)
1479  *              rowvect (vector of results gathered down the columns of pixs)
1480  *      Return: 0 if OK, 1 on error
1481  *
1482  *  Notes:
1483  *      (1) This computes a row vector of statistics using each
1484  *          column of a Pix.  The result is put in @rowvect.
1485  *      (2) The @thresh parameter works with L_MODE_VAL only, and
1486  *          sets a minimum occupancy of the mode bin.
1487  *          If the occupancy of the mode bin is less than @thresh, the
1488  *          mode value is returned as 0.  To always return the actual
1489  *          mode value, set @thresh = 0.
1490  *      (3) What is the meaning of this @thresh parameter?
1491  *          For each column, the total count in the histogram is h, the
1492  *          image height.  So @thresh, relative to h, gives a measure
1493  *          of the ratio of the bin width to the width of the distribution.
1494  *          The larger @thresh, the narrower the distribution must be
1495  *          for the mode value to be returned (instead of returning 0).
1496  */
1497 l_int32
pixGetColumnStats(PIX * pixs,l_int32 type,l_int32 nbins,l_int32 thresh,l_float32 * rowvect)1498 pixGetColumnStats(PIX        *pixs,
1499                   l_int32     type,
1500                   l_int32     nbins,
1501                   l_int32     thresh,
1502                   l_float32  *rowvect)
1503 {
1504 l_int32    i, j, k, w, h, val, wpls, sum, target, max, modeval;
1505 l_int32   *histo, *gray2bin, *bin2gray;
1506 l_uint32  *datas;
1507 
1508     PROCNAME("pixGetColumnStats");
1509 
1510     if (!pixs || pixGetDepth(pixs) != 8)
1511         return ERROR_INT("pixs not defined or not 8 bpp", procName, 1);
1512     if (!rowvect)
1513         return ERROR_INT("rowvect not defined", procName, 1);
1514     if (type != L_MEAN_ABSVAL && type != L_MEDIAN_VAL &&
1515         type != L_MODE_VAL && type != L_MODE_COUNT)
1516         return ERROR_INT("invalid type", procName, 1);
1517     if (type != L_MEAN_ABSVAL && (nbins < 1 || nbins > 256))
1518         return ERROR_INT("invalid nbins", procName, 1);
1519     pixGetDimensions(pixs, &w, &h, NULL);
1520 
1521     datas = pixGetData(pixs);
1522     wpls = pixGetWpl(pixs);
1523     if (type == L_MEAN_ABSVAL) {
1524         for (j = 0; j < w; j++) {
1525             sum = 0;
1526             for (i = 0; i < h; i++)
1527                 sum += GET_DATA_BYTE(datas + i * wpls, j);
1528             rowvect[j] = (l_float32)sum / (l_float32)h;
1529         }
1530         return 0;
1531     }
1532 
1533         /* We need a histogram; binwidth ~ 256 / nbins */
1534     histo = (l_int32 *)CALLOC(nbins, sizeof(l_int32));
1535     gray2bin = (l_int32 *)CALLOC(256, sizeof(l_int32));
1536     bin2gray = (l_int32 *)CALLOC(nbins, sizeof(l_int32));
1537     for (i = 0; i < 256; i++)  /* gray value --> histo bin */
1538         gray2bin[i] = (i * nbins) / 256;
1539     for (i = 0; i < nbins; i++)  /* histo bin --> gray value */
1540         bin2gray[i] = (i * 255 + 128) / nbins;
1541 
1542     for (j = 0; j < w; j++) {
1543         for (i = 0; i < h; i++) {
1544             val = GET_DATA_BYTE(datas + i * wpls, j);
1545             histo[gray2bin[val]]++;
1546         }
1547 
1548         if (type == L_MEDIAN_VAL) {
1549             sum = 0;
1550             target = h / 2;
1551             for (k = 0; k < nbins; k++) {
1552                 sum += histo[k];
1553                 if (sum >= target) {
1554                     rowvect[j] = bin2gray[k];
1555                     break;
1556                 }
1557             }
1558         }
1559         else if (type == L_MODE_VAL) {
1560             max = 0;
1561             modeval = 0;
1562             for (k = 0; k < nbins; k++) {
1563                 if (histo[k] > max) {
1564                     max = histo[k];
1565                     modeval = k;
1566                 }
1567             }
1568             if (max < thresh)
1569                 rowvect[j] = 0;
1570             else
1571                 rowvect[j] = bin2gray[modeval];
1572         }
1573         else {  /* type == L_MODE_COUNT */
1574             max = 0;
1575             modeval = 0;
1576             for (k = 0; k < nbins; k++) {
1577                 if (histo[k] > max) {
1578                     max = histo[k];
1579                     modeval = k;
1580                 }
1581             }
1582             rowvect[j] = max;
1583         }
1584         for (k = 0; k < nbins; k++)
1585             histo[k] = 0;
1586     }
1587 
1588     FREE(histo);
1589     FREE(gray2bin);
1590     FREE(bin2gray);
1591     return 0;
1592 }
1593 
1594 
1595 /*!
1596  *  pixSetPixelColumn()
1597  *
1598  *      Input:  pix (8 bpp; not cmapped)
1599  *              col (column index)
1600  *              colvect (vector of floats)
1601  *      Return: 0 if OK, 1 on error
1602  */
1603 l_int32
pixSetPixelColumn(PIX * pix,l_int32 col,l_float32 * colvect)1604 pixSetPixelColumn(PIX        *pix,
1605                   l_int32     col,
1606                   l_float32  *colvect)
1607 {
1608 l_int32    i, w, h, wpl;
1609 l_uint32  *data;
1610 
1611     PROCNAME("pixSetCPixelColumn");
1612 
1613     if (!pix || pixGetDepth(pix) != 8)
1614         return ERROR_INT("pix not defined or not 8 bpp", procName, 1);
1615     if (!colvect)
1616         return ERROR_INT("colvect not defined", procName, 1);
1617     pixGetDimensions(pix, &w, &h, NULL);
1618     if (col < 0 || col > w)
1619         return ERROR_INT("invalid col", procName, 1);
1620 
1621     data = pixGetData(pix);
1622     wpl = pixGetWpl(pix);
1623     for (i = 0; i < h; i++)
1624         SET_DATA_BYTE(data + i * wpl, col, (l_int32)colvect[i]);
1625 
1626     return 0;
1627 }
1628 
1629 
1630 /*-------------------------------------------------------------*
1631  *              Foreground/background estimation               *
1632  *-------------------------------------------------------------*/
1633 /*!
1634  *  pixThresholdForFgBg()
1635  *
1636  *      Input:  pixs (any depth; cmapped ok)
1637  *              factor (subsampling factor; integer >= 1)
1638  *              thresh (threshold for generating foreground mask)
1639  *              &fgval (<optional return> average foreground value)
1640  *              &bgval (<optional return> average background value)
1641  *      Return: 0 if OK, 1 on error
1642  */
1643 l_int32
pixThresholdForFgBg(PIX * pixs,l_int32 factor,l_int32 thresh,l_int32 * pfgval,l_int32 * pbgval)1644 pixThresholdForFgBg(PIX      *pixs,
1645                     l_int32   factor,
1646                     l_int32   thresh,
1647                     l_int32  *pfgval,
1648                     l_int32  *pbgval)
1649 {
1650 l_float32  fval;
1651 PIX       *pixg, *pixm;
1652 
1653     PROCNAME("pixThresholdForFgBg");
1654 
1655     if (!pixs)
1656         return ERROR_INT("pixs not defined", procName, 1);
1657 
1658         /* Generate a subsampled 8 bpp version and a mask over the fg */
1659     pixg = pixConvertTo8BySampling(pixs, factor, 0);
1660     pixm = pixThresholdToBinary(pixg, thresh);
1661 
1662     if (pfgval) {
1663         pixGetAverageMasked(pixg, pixm, 0, 0, 1, L_MEAN_ABSVAL, &fval);
1664         *pfgval = (l_int32)(fval + 0.5);
1665     }
1666 
1667     if (pbgval) {
1668         pixInvert(pixm, pixm);
1669         pixGetAverageMasked(pixg, pixm, 0, 0, 1, L_MEAN_ABSVAL, &fval);
1670         *pbgval = (l_int32)(fval + 0.5);
1671     }
1672 
1673     pixDestroy(&pixg);
1674     pixDestroy(&pixm);
1675     return 0;
1676 }
1677 
1678 
1679 /*!
1680  *  pixSplitDistributionFgBg()
1681  *
1682  *      Input:  pixs (any depth; cmapped ok)
1683  *              scorefract (fraction of the max score, used to determine
1684  *                          the range over which the histogram min is searched)
1685  *              factor (subsampling factor; integer >= 1)
1686  *              &thresh (<optional return> best threshold for separating)
1687  *              &fgval (<optional return> average foreground value)
1688  *              &bgval (<optional return> average background value)
1689  *              debugflag (1 for plotting of distribution and split point)
1690  *      Return: 0 if OK, 1 on error
1691  *
1692  *  Notes:
1693  *      (1) See numaSplitDistribution() for details on the underlying
1694  *          method of choosing a threshold.
1695  */
1696 l_int32
pixSplitDistributionFgBg(PIX * pixs,l_float32 scorefract,l_int32 factor,l_int32 * pthresh,l_int32 * pfgval,l_int32 * pbgval,l_int32 debugflag)1697 pixSplitDistributionFgBg(PIX       *pixs,
1698                          l_float32  scorefract,
1699                          l_int32    factor,
1700                          l_int32   *pthresh,
1701                          l_int32   *pfgval,
1702                          l_int32   *pbgval,
1703                          l_int32    debugflag)
1704 {
1705 l_int32    thresh;
1706 l_float32  avefg, avebg, maxnum;
1707 GPLOT     *gplot;
1708 NUMA      *na, *nascore, *nax, *nay;
1709 PIX       *pixg;
1710 
1711     PROCNAME("pixSplitDistributionFgBg");
1712 
1713     if (pthresh) *pthresh = 0;
1714     if (pfgval) *pfgval = 0;
1715     if (pbgval) *pbgval = 0;
1716     if (!pixs)
1717         return ERROR_INT("pixs not defined", procName, 1);
1718 
1719         /* Generate a subsampled 8 bpp version */
1720     pixg = pixConvertTo8BySampling(pixs, factor, 0);
1721 
1722         /* Make the fg/bg estimates */
1723     na = pixGetGrayHistogram(pixg, 1);
1724     if (debugflag) {
1725         numaSplitDistribution(na, scorefract, &thresh, &avefg, &avebg,
1726                               NULL, NULL, &nascore);
1727         numaDestroy(&nascore);
1728     }
1729     else
1730         numaSplitDistribution(na, scorefract, &thresh, &avefg, &avebg,
1731                               NULL, NULL, NULL);
1732 
1733     if (pthresh) *pthresh = thresh;
1734     if (pfgval) *pfgval = (l_int32)(avefg + 0.5);
1735     if (pbgval) *pbgval = (l_int32)(avebg + 0.5);
1736 
1737     if (debugflag) {
1738         gplot = gplotCreate("junk_histplot", GPLOT_X11, "Histogram",
1739                             "Grayscale value", "Number of pixels");
1740         gplotAddPlot(gplot, NULL, na, GPLOT_LINES, NULL);
1741         nax = numaMakeConstant(thresh, 2);
1742         numaGetMax(na, &maxnum, NULL);
1743         nay = numaMakeConstant(0, 2);
1744         numaReplaceNumber(nay, 1, (l_int32)(0.5 * maxnum));
1745         gplotAddPlot(gplot, nax, nay, GPLOT_LINES, NULL);
1746         gplotMakeOutput(gplot);
1747         gplotDestroy(&gplot);
1748         numaDestroy(&nax);
1749         numaDestroy(&nay);
1750     }
1751 
1752     pixDestroy(&pixg);
1753     numaDestroy(&na);
1754     return 0;
1755 }
1756 
1757 
1758 /*-------------------------------------------------------------*
1759  *                 Measurement of properties                   *
1760  *-------------------------------------------------------------*/
1761 /*!
1762  *  pixaFindDimensions()
1763  *
1764  *      Input:  pixa
1765  *              &naw (<optional return> numa of pix widths)
1766  *              &nah (<optional return> numa of pix heights)
1767  *      Return: 0 if OK, 1 on error
1768  */
1769 l_int32
pixaFindDimensions(PIXA * pixa,NUMA ** pnaw,NUMA ** pnah)1770 pixaFindDimensions(PIXA   *pixa,
1771                    NUMA  **pnaw,
1772                    NUMA  **pnah)
1773 {
1774 l_int32  i, n, w, h;
1775 PIX     *pixt;
1776 
1777     PROCNAME("pixaFindDimensions");
1778 
1779     if (!pixa)
1780         return ERROR_INT("pixa not defined", procName, 1);
1781     if (!pnaw && !pnah)
1782         return 0;
1783 
1784     n = pixaGetCount(pixa);
1785     if (pnaw) *pnaw = numaCreate(n);
1786     if (pnah) *pnah = numaCreate(n);
1787     for (i = 0; i < n; i++) {
1788         pixt = pixaGetPix(pixa, i, L_CLONE);
1789         pixGetDimensions(pixt, &w, &h, NULL);
1790         if (pnaw)
1791             numaAddNumber(*pnaw, w);
1792         if (pnah)
1793             numaAddNumber(*pnah, h);
1794         pixDestroy(&pixt);
1795     }
1796     return 0;
1797 }
1798 
1799 
1800 /*!
1801  *  pixaFindAreaPerimRatio()
1802  *
1803  *      Input:  pixa (of 1 bpp pix)
1804  *      Return: na (of area/perimeter ratio for each pix), or null on error
1805  *
1806  *  Notes:
1807  *      (1) This is typically used for a pixa consisting of
1808  *          1 bpp connected components.
1809  */
1810 NUMA *
pixaFindAreaPerimRatio(PIXA * pixa)1811 pixaFindAreaPerimRatio(PIXA  *pixa)
1812 {
1813 l_int32    i, n;
1814 l_int32   *tab;
1815 l_float32  fract;
1816 NUMA      *na;
1817 PIX       *pixt;
1818 
1819     PROCNAME("pixaFindAreaPerimRatio");
1820 
1821     if (!pixa)
1822         return (NUMA *)ERROR_PTR("pixa not defined", procName, NULL);
1823 
1824     n = pixaGetCount(pixa);
1825     na = numaCreate(n);
1826     tab = makePixelSumTab8();
1827     for (i = 0; i < n; i++) {
1828         pixt = pixaGetPix(pixa, i, L_CLONE);
1829         pixFindAreaPerimRatio(pixt, tab, &fract);
1830         numaAddNumber(na, fract);
1831         pixDestroy(&pixt);
1832     }
1833     FREE(tab);
1834     return na;
1835 }
1836 
1837 
1838 /*!
1839  *  pixFindAreaPerimRatio()
1840  *
1841  *      Input:  pixs (1 bpp)
1842  *              tab (<optional> pixel sum table, can be NULL)
1843  *              &fract (<return> area/perimeter ratio)
1844  *      Return: 0 if OK, 1 on error
1845  *
1846  *  Notes:
1847  *      (1) The area is the number of fg pixels that are not on the
1848  *          boundary (i.e., not 8-connected to a bg pixel), and the
1849  *          perimeter is the number of boundary fg pixels.
1850  *      (2) This is typically used for a pixa consisting of
1851  *          1 bpp connected components.
1852  */
1853 l_int32
pixFindAreaPerimRatio(PIX * pixs,l_int32 * tab,l_float32 * pfract)1854 pixFindAreaPerimRatio(PIX        *pixs,
1855                       l_int32    *tab,
1856                       l_float32  *pfract)
1857 {
1858 l_int32  *tab8;
1859 l_int32   nin, nbound;
1860 PIX      *pixt;
1861 
1862     PROCNAME("pixFindAreaPerimRatio");
1863 
1864     if (!pfract)
1865         return ERROR_INT("&fract not defined", procName, 1);
1866     *pfract = 0.0;
1867     if (!pixs || pixGetDepth(pixs) != 1)
1868         return ERROR_INT("pixs not defined or not 1 bpp", procName, 1);
1869 
1870     if (!tab)
1871         tab8 = makePixelSumTab8();
1872     else
1873         tab8 = tab;
1874 
1875     pixt = pixErodeBrick(NULL, pixs, 3, 3);
1876     pixCountPixels(pixt, &nin, tab8);
1877     pixXor(pixt, pixt, pixs);
1878     pixCountPixels(pixt, &nbound, tab8);
1879     *pfract = (l_float32)nin / (l_float32)nbound;
1880 
1881     if (!tab)
1882         FREE(tab8);
1883     pixDestroy(&pixt);
1884     return 0;
1885 }
1886 
1887 
1888 /*!
1889  *  pixaFindPerimSizeRatio()
1890  *
1891  *      Input:  pixa (of 1 bpp pix)
1892  *      Return: na (of fg perimeter/(w*h) ratio for each pix), or null on error
1893  *
1894  *  Notes:
1895  *      (1) This is typically used for a pixa consisting of
1896  *          1 bpp connected components.
1897  */
1898 NUMA *
pixaFindPerimSizeRatio(PIXA * pixa)1899 pixaFindPerimSizeRatio(PIXA  *pixa)
1900 {
1901 l_int32    i, n;
1902 l_int32   *tab;
1903 l_float32  ratio;
1904 NUMA      *na;
1905 PIX       *pixt;
1906 
1907     PROCNAME("pixaFindPerimSizeRatio");
1908 
1909     if (!pixa)
1910         return (NUMA *)ERROR_PTR("pixa not defined", procName, NULL);
1911 
1912     n = pixaGetCount(pixa);
1913     na = numaCreate(n);
1914     tab = makePixelSumTab8();
1915     for (i = 0; i < n; i++) {
1916         pixt = pixaGetPix(pixa, i, L_CLONE);
1917         pixFindPerimSizeRatio(pixt, tab, &ratio);
1918         numaAddNumber(na, ratio);
1919         pixDestroy(&pixt);
1920     }
1921     FREE(tab);
1922     return na;
1923 }
1924 
1925 
1926 /*!
1927  *  pixFindPerimSizeRatio()
1928  *
1929  *      Input:  pixs (1 bpp)
1930  *              tab (<optional> pixel sum table, can be NULL)
1931  *              &ratio (<return> perimeter/size ratio)
1932  *      Return: 0 if OK, 1 on error
1933  *
1934  *  Notes:
1935  *      (1) The size is the sum of the width and height of the pix,
1936  *          and the perimeter is the number of boundary fg pixels.
1937  *      (2) This has a large value for dendritic, fractal-like components
1938  *          with highly irregular boundaries.
1939  *      (3) This is typically used for a single connected component.
1940  */
1941 l_int32
pixFindPerimSizeRatio(PIX * pixs,l_int32 * tab,l_float32 * pratio)1942 pixFindPerimSizeRatio(PIX        *pixs,
1943                       l_int32    *tab,
1944                       l_float32  *pratio)
1945 {
1946 l_int32  *tab8;
1947 l_int32   w, h, nbound;
1948 PIX      *pixt;
1949 
1950     PROCNAME("pixFindPerimSizeRatio");
1951 
1952     if (!pratio)
1953         return ERROR_INT("&ratio not defined", procName, 1);
1954     *pratio = 0.0;
1955     if (!pixs || pixGetDepth(pixs) != 1)
1956         return ERROR_INT("pixs not defined or not 1 bpp", procName, 1);
1957 
1958     if (!tab)
1959         tab8 = makePixelSumTab8();
1960     else
1961         tab8 = tab;
1962 
1963     pixt = pixErodeBrick(NULL, pixs, 3, 3);
1964     pixXor(pixt, pixt, pixs);
1965     pixCountPixels(pixt, &nbound, tab8);
1966     pixGetDimensions(pixs, &w, &h, NULL);
1967     *pratio = (l_float32)nbound / (l_float32)(w + h);
1968 
1969     if (!tab)
1970         FREE(tab8);
1971     pixDestroy(&pixt);
1972     return 0;
1973 }
1974 
1975 
1976 /*!
1977  *  pixaFindAreaFraction()
1978  *
1979  *      Input:  pixa (of 1 bpp pix)
1980  *      Return: na (of area fractions for each pix), or null on error
1981  *
1982  *  Notes:
1983  *      (1) This is typically used for a pixa consisting of
1984  *          1 bpp connected components.
1985  */
1986 NUMA *
pixaFindAreaFraction(PIXA * pixa)1987 pixaFindAreaFraction(PIXA  *pixa)
1988 {
1989 l_int32    i, n;
1990 l_int32   *tab;
1991 l_float32  fract;
1992 NUMA      *na;
1993 PIX       *pixt;
1994 
1995     PROCNAME("pixaFindAreaFraction");
1996 
1997     if (!pixa)
1998         return (NUMA *)ERROR_PTR("pixa not defined", procName, NULL);
1999 
2000     n = pixaGetCount(pixa);
2001     na = numaCreate(n);
2002     tab = makePixelSumTab8();
2003     for (i = 0; i < n; i++) {
2004         pixt = pixaGetPix(pixa, i, L_CLONE);
2005         pixFindAreaFraction(pixt, tab, &fract);
2006         numaAddNumber(na, fract);
2007         pixDestroy(&pixt);
2008     }
2009     FREE(tab);
2010     return na;
2011 }
2012 
2013 
2014 /*!
2015  *  pixFindAreaFraction()
2016  *
2017  *      Input:  pixs (1 bpp)
2018  *              tab (<optional> pixel sum table, can be NULL)
2019  *              &fract (<return> fg area/size ratio)
2020  *      Return: 0 if OK, 1 on error
2021  *
2022  *  Notes:
2023  *      (1) This finds the ratio of the number of fg pixels to the
2024  *          size of the pix (w * h).  It is typically used for a
2025  *          single connected component.
2026  */
2027 l_int32
pixFindAreaFraction(PIX * pixs,l_int32 * tab,l_float32 * pfract)2028 pixFindAreaFraction(PIX        *pixs,
2029                     l_int32    *tab,
2030                     l_float32  *pfract)
2031 {
2032 l_int32   w, h, d, sum;
2033 l_int32  *tab8;
2034 
2035     PROCNAME("pixFindAreaFraction");
2036 
2037     if (!pfract)
2038         return ERROR_INT("&fract not defined", procName, 1);
2039     *pfract = 0.0;
2040     pixGetDimensions(pixs, &w, &h, &d);
2041     if (!pixs || d != 1)
2042         return ERROR_INT("pixs not defined or not 1 bpp", procName, 1);
2043 
2044     if (!tab)
2045         tab8 = makePixelSumTab8();
2046     else
2047         tab8 = tab;
2048 
2049     pixCountPixels(pixs, &sum, tab8);
2050     *pfract = (l_float32)sum / (l_float32)(w * h);
2051 
2052     if (!tab)
2053         FREE(tab8);
2054     return 0;
2055 }
2056 
2057 
2058 /*!
2059  *  pixaFindWidthHeightRatio()
2060  *
2061  *      Input:  pixa (of 1 bpp pix)
2062  *      Return: na (of width/height ratios for each pix), or null on error
2063  *
2064  *  Notes:
2065  *      (1) This is typically used for a pixa consisting of
2066  *          1 bpp connected components.
2067  */
2068 NUMA *
pixaFindWidthHeightRatio(PIXA * pixa)2069 pixaFindWidthHeightRatio(PIXA  *pixa)
2070 {
2071 l_int32  i, n, w, h;
2072 NUMA    *na;
2073 PIX     *pixt;
2074 
2075     PROCNAME("pixaFindWidthHeightRatio");
2076 
2077     if (!pixa)
2078         return (NUMA *)ERROR_PTR("pixa not defined", procName, NULL);
2079 
2080     n = pixaGetCount(pixa);
2081     na = numaCreate(n);
2082     for (i = 0; i < n; i++) {
2083         pixt = pixaGetPix(pixa, i, L_CLONE);
2084         pixGetDimensions(pixt, &w, &h, NULL);
2085         numaAddNumber(na, (l_float32)w / (l_float32)h);
2086         pixDestroy(&pixt);
2087     }
2088     return na;
2089 }
2090 
2091 
2092 /*!
2093  *  pixaFindWidthHeightProduct()
2094  *
2095  *      Input:  pixa (of 1 bpp pix)
2096  *      Return: na (of width*height products for each pix), or null on error
2097  *
2098  *  Notes:
2099  *      (1) This is typically used for a pixa consisting of
2100  *          1 bpp connected components.
2101  */
2102 NUMA *
pixaFindWidthHeightProduct(PIXA * pixa)2103 pixaFindWidthHeightProduct(PIXA  *pixa)
2104 {
2105 l_int32  i, n, w, h;
2106 NUMA    *na;
2107 PIX     *pixt;
2108 
2109     PROCNAME("pixaFindWidthHeightProduct");
2110 
2111     if (!pixa)
2112         return (NUMA *)ERROR_PTR("pixa not defined", procName, NULL);
2113 
2114     n = pixaGetCount(pixa);
2115     na = numaCreate(n);
2116     for (i = 0; i < n; i++) {
2117         pixt = pixaGetPix(pixa, i, L_CLONE);
2118         pixGetDimensions(pixt, &w, &h, NULL);
2119         numaAddNumber(na, w * h);
2120         pixDestroy(&pixt);
2121     }
2122     return na;
2123 }
2124 
2125 
2126 /*-------------------------------------------------------------*
2127  *                Extract rectangular region                   *
2128  *-------------------------------------------------------------*/
2129 /*!
2130  *  pixClipRectangle()
2131  *
2132  *      Input:  pixs
2133  *              box  (requested clipping region; const)
2134  *              &boxc (<optional return> actual box of clipped region)
2135  *      Return: clipped pix, or null on error or if rectangle
2136  *              doesn't intersect pixs
2137  *
2138  *  Notes:
2139  *
2140  *  This should be simple, but there are choices to be made.
2141  *  The box is defined relative to the pix coordinates.  However,
2142  *  if the box is not contained within the pix, we have two choices:
2143  *
2144  *      (1) clip the box to the pix
2145  *      (2) make a new pix equal to the full box dimensions,
2146  *          but let rasterop do the clipping and positioning
2147  *          of the src with respect to the dest
2148  *
2149  *  Choice (2) immediately brings up the problem of what pixel values
2150  *  to use that were not taken from the src.  For example, on a grayscale
2151  *  image, do you want the pixels not taken from the src to be black
2152  *  or white or something else?  To implement choice 2, one needs to
2153  *  specify the color of these extra pixels.
2154  *
2155  *  So we adopt (1), and clip the box first, if necessary,
2156  *  before making the dest pix and doing the rasterop.  But there
2157  *  is another issue to consider.  If you want to paste the
2158  *  clipped pix back into pixs, it must be properly aligned, and
2159  *  it is necessary to use the clipped box for alignment.
2160  *  Accordingly, this function has a third (optional) argument, which is
2161  *  the input box clipped to the src pix.
2162  */
2163 PIX *
pixClipRectangle(PIX * pixs,BOX * box,BOX ** pboxc)2164 pixClipRectangle(PIX   *pixs,
2165                  BOX   *box,
2166                  BOX  **pboxc)
2167 {
2168 l_int32  w, h, d, bx, by, bw, bh;
2169 BOX     *boxc;
2170 PIX     *pixd;
2171 
2172     PROCNAME("pixClipRectangle");
2173 
2174     if (pboxc)
2175         *pboxc = NULL;
2176     if (!pixs)
2177         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
2178     if (!box)
2179         return (PIX *)ERROR_PTR("box not defined", procName, NULL);
2180 
2181         /* Clip the input box to the pix */
2182     pixGetDimensions(pixs, &w, &h, &d);
2183     if ((boxc = boxClipToRectangle(box, w, h)) == NULL) {
2184         L_WARNING("box doesn't overlap pix", procName);
2185         return NULL;
2186     }
2187     boxGetGeometry(boxc, &bx, &by, &bw, &bh);
2188 
2189         /* Extract the block */
2190     if ((pixd = pixCreate(bw, bh, d)) == NULL)
2191         return (PIX *)ERROR_PTR("pixd not made", procName, NULL);
2192     pixCopyResolution(pixd, pixs);
2193     pixCopyColormap(pixd, pixs);
2194     pixRasterop(pixd, 0, 0, bw, bh, PIX_SRC, pixs, bx, by);
2195 
2196     if (pboxc)
2197         *pboxc = boxc;
2198     else
2199         boxDestroy(&boxc);
2200 
2201     return pixd;
2202 }
2203 
2204 
2205 /*!
2206  *  pixClipMasked()
2207  *
2208  *      Input:  pixs (1, 2, 4, 8, 16, 32 bpp; colormap ok)
2209  *              pixm  (clipping mask, 1 bpp)
2210  *              x, y (origin of clipping mask relative to pixs)
2211  *              outval (val to use for pixels that are outside the mask)
2212  *      Return: pixd, (clipped pix) or null on error or if pixm doesn't
2213  *              intersect pixs
2214  *
2215  *  Notes:
2216  *      (1) If pixs has a colormap, it is preserved in pixd.
2217  *      (2) The depth of pixd is the same as that of pixs.
2218  *      (3) If the depth of pixs is 1, use @outval = 0 for white background
2219  *          and 1 for black; otherwise, use the max value for white
2220  *          and 0 for black.  If pixs has a colormap, the max value for
2221  *          @outval is 0xffffffff; otherwise, it is 2^d - 1.
2222  *      (4) When using 1 bpp pixs, this is a simple clip and
2223  *          blend operation.  For example, if both pix1 and pix2 are
2224  *          black text on white background, and you want to OR the
2225  *          fg on the two images, let pixm be the inverse of pix2.
2226  *          Then the operation takes all of pix1 that's in the bg of
2227  *          pix2, and for the remainder (which are the pixels
2228  *          corresponding to the fg of the pix2), paint them black
2229  *          (1) in pix1.  The function call looks like
2230  *             pixClipMasked(pix2, pixInvert(pix1, pix1), x, y, 1);
2231  */
2232 PIX *
pixClipMasked(PIX * pixs,PIX * pixm,l_int32 x,l_int32 y,l_uint32 outval)2233 pixClipMasked(PIX      *pixs,
2234               PIX      *pixm,
2235               l_int32   x,
2236               l_int32   y,
2237               l_uint32  outval)
2238 {
2239 l_int32   wm, hm, d, index, rval, gval, bval;
2240 l_uint32  pixel;
2241 BOX      *box;
2242 PIX      *pixmi, *pixd;
2243 PIXCMAP  *cmap;
2244 
2245     PROCNAME("pixClipMasked");
2246 
2247     if (!pixs)
2248         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
2249     if (!pixm || pixGetDepth(pixm) != 1)
2250         return (PIX *)ERROR_PTR("pixm undefined or not 1 bpp", procName, NULL);
2251 
2252         /* Clip out the region specified by pixm and (x,y) */
2253     pixGetDimensions(pixm, &wm, &hm, NULL);
2254     box = boxCreate(x, y, wm, hm);
2255     pixd = pixClipRectangle(pixs, box, NULL);
2256 
2257         /* Paint 'outval' (or something close to it if cmapped) through
2258          * the pixels not masked by pixm */
2259     cmap = pixGetColormap(pixd);
2260     pixmi = pixInvert(NULL, pixm);
2261     d = pixGetDepth(pixd);
2262     if (cmap) {
2263         extractRGBValues(outval, &rval, &gval, &bval);
2264         pixcmapGetNearestIndex(cmap, rval, gval, bval, &index);
2265         pixcmapGetColor(cmap, index, &rval, &gval, &bval);
2266         composeRGBPixel(rval, gval, bval, &pixel);
2267         pixPaintThroughMask(pixd, pixmi, 0, 0, pixel);
2268     }
2269     else
2270         pixPaintThroughMask(pixd, pixmi, 0, 0, outval);
2271 
2272     boxDestroy(&box);
2273     pixDestroy(&pixmi);
2274     return pixd;
2275 }
2276 
2277 
2278 /*---------------------------------------------------------------------*
2279  *                          Clipping to Foreground                     *
2280  *---------------------------------------------------------------------*/
2281 /*!
2282  *  pixClipToForeground()
2283  *
2284  *      Input:  pixs (1 bpp)
2285  *              &pixd  (<optional return> clipped pix returned)
2286  *              &box   (<optional return> bounding box)
2287  *      Return: 0 if OK; 1 on error or if there are no fg pixels
2288  *
2289  *  Notes:
2290  *      (1) At least one of {&pixd, &box} must be specified.
2291  *      (2) If there are no fg pixels, the returned ptrs are null.
2292  */
2293 l_int32
pixClipToForeground(PIX * pixs,PIX ** ppixd,BOX ** pbox)2294 pixClipToForeground(PIX   *pixs,
2295                     PIX  **ppixd,
2296                     BOX  **pbox)
2297 {
2298 l_int32    w, h, wpl, nfullwords, extra, i, j;
2299 l_int32    minx, miny, maxx, maxy;
2300 l_uint32   result, mask;
2301 l_uint32  *data, *line;
2302 BOX       *box;
2303 
2304     PROCNAME("pixClipToForeground");
2305 
2306     if (!ppixd && !pbox)
2307         return ERROR_INT("neither &pixd nor &box defined", procName, 1);
2308     if (ppixd)
2309         *ppixd = NULL;
2310     if (pbox)
2311         *pbox = NULL;
2312     if (!pixs || (pixGetDepth(pixs) != 1))
2313         return ERROR_INT("pixs not defined or not 1 bpp", procName, 1);
2314 
2315     pixGetDimensions(pixs, &w, &h, NULL);
2316     nfullwords = w / 32;
2317     extra = w & 31;
2318     mask = ~rmask32[32 - extra];
2319     wpl = pixGetWpl(pixs);
2320     data = pixGetData(pixs);
2321 
2322     result = 0;
2323     for (i = 0, miny = 0; i < h; i++, miny++) {
2324         line = data + i * wpl;
2325         for (j = 0; j < nfullwords; j++)
2326             result |= line[j];
2327         if (extra)
2328             result |= (line[j] & mask);
2329         if (result)
2330             break;
2331     }
2332     if (miny == h)  /* no ON pixels */
2333         return 1;
2334 
2335     result = 0;
2336     for (i = h - 1, maxy = h - 1; i >= 0; i--, maxy--) {
2337         line = data + i * wpl;
2338         for (j = 0; j < nfullwords; j++)
2339             result |= line[j];
2340         if (extra)
2341             result |= (line[j] & mask);
2342         if (result)
2343             break;
2344     }
2345 
2346     minx = 0;
2347     for (j = 0, minx = 0; j < w; j++, minx++) {
2348         for (i = 0; i < h; i++) {
2349             line = data + i * wpl;
2350             if (GET_DATA_BIT(line, j))
2351                 goto minx_found;
2352         }
2353     }
2354 
2355 minx_found:
2356     for (j = w - 1, maxx = w - 1; j >= 0; j--, maxx--) {
2357         for (i = 0; i < h; i++) {
2358             line = data + i * wpl;
2359             if (GET_DATA_BIT(line, j))
2360                 goto maxx_found;
2361         }
2362     }
2363 
2364 maxx_found:
2365     box = boxCreate(minx, miny, maxx - minx + 1, maxy - miny + 1);
2366 
2367     if (ppixd)
2368         *ppixd = pixClipRectangle(pixs, box, NULL);
2369     if (pbox)
2370         *pbox = box;
2371     else
2372         boxDestroy(&box);
2373 
2374     return 0;
2375 }
2376 
2377 
2378 /*!
2379  *  pixClipBoxToForeground()
2380  *
2381  *      Input:  pixs (1 bpp)
2382  *              boxs  (<optional> ; use full image if null)
2383  *              &pixd  (<optional return> clipped pix returned)
2384  *              &boxd  (<optional return> bounding box)
2385  *      Return: 0 if OK; 1 on error or if there are no fg pixels
2386  *
2387  *  Notes:
2388  *      (1) At least one of {&pixd, &boxd} must be specified.
2389  *      (2) If there are no fg pixels, the returned ptrs are null.
2390  *      (3) Do not use &pixs for the 3rd arg or &boxs for the 4th arg;
2391  *          this will leak memory.
2392  */
2393 l_int32
pixClipBoxToForeground(PIX * pixs,BOX * boxs,PIX ** ppixd,BOX ** pboxd)2394 pixClipBoxToForeground(PIX   *pixs,
2395                        BOX   *boxs,
2396                        PIX  **ppixd,
2397                        BOX  **pboxd)
2398 {
2399 l_int32  w, h, bx, by, bw, bh, cbw, cbh, left, right, top, bottom;
2400 BOX     *boxt, *boxd;
2401 
2402     PROCNAME("pixClipBoxToForeground");
2403 
2404     if (!ppixd && !pboxd)
2405         return ERROR_INT("neither &pixd nor &boxd defined", procName, 1);
2406     if (ppixd) *ppixd = NULL;
2407     if (pboxd) *pboxd = NULL;
2408     if (!pixs || (pixGetDepth(pixs) != 1))
2409         return ERROR_INT("pixs not defined or not 1 bpp", procName, 1);
2410 
2411     if (!boxs)
2412         return pixClipToForeground(pixs, ppixd, pboxd);
2413 
2414     pixGetDimensions(pixs, &w, &h, NULL);
2415     if (boxs) {
2416         boxGetGeometry(boxs, &bx, &by, &bw, &bh);
2417         cbw = L_MIN(bw, w - bx);
2418         cbh = L_MIN(bh, h - by);
2419         if (cbw < 0 || cbh < 0)
2420             return ERROR_INT("box not within image", procName, 1);
2421         boxt = boxCreate(bx, by, cbw, cbh);
2422     }
2423     else
2424         boxt = boxCreate(0, 0, w, h);
2425 
2426 
2427     if (pixScanForForeground(pixs, boxt, L_FROM_LEFT, &left)) {
2428         boxDestroy(&boxt);
2429         return 1;
2430     }
2431     pixScanForForeground(pixs, boxt, L_FROM_RIGHT, &right);
2432     pixScanForForeground(pixs, boxt, L_FROM_TOP, &top);
2433     pixScanForForeground(pixs, boxt, L_FROM_BOTTOM, &bottom);
2434 
2435     boxd = boxCreate(left, top, right - left + 1, bottom - top + 1);
2436     if (ppixd)
2437         *ppixd = pixClipRectangle(pixs, boxd, NULL);
2438     if (pboxd)
2439         *pboxd = boxd;
2440     else
2441         boxDestroy(&boxd);
2442 
2443     boxDestroy(&boxt);
2444     return 0;
2445 }
2446 
2447 
2448 /*!
2449  *  pixScanForForeground()
2450  *
2451  *      Input:  pixs (1 bpp)
2452  *              box  (<optional> within which the search is conducted)
2453  *              scanflag (direction of scan; e.g., L_FROM_LEFT)
2454  *              &loc (location in scan direction of first black pixel)
2455  *      Return: 0 if OK; 1 on error or if no fg pixels are found
2456  *
2457  *  Notes:
2458  *      (1) If there are no fg pixels, the position is set to 0.
2459  *          Caller must check the return value!
2460  *      (2) Use @box == NULL to scan from edge of pixs
2461  */
2462 l_int32
pixScanForForeground(PIX * pixs,BOX * box,l_int32 scanflag,l_int32 * ploc)2463 pixScanForForeground(PIX      *pixs,
2464                      BOX      *box,
2465                      l_int32   scanflag,
2466                      l_int32  *ploc)
2467 {
2468 l_int32    bx, by, bw, bh, x, xstart, xend, y, ystart, yend, wpl;
2469 l_uint32  *data, *line;
2470 BOX       *boxt;
2471 
2472     PROCNAME("pixScanForForeground");
2473 
2474     if (!ploc)
2475         return ERROR_INT("&ploc not defined", procName, 1);
2476     *ploc = 0;
2477     if (!pixs || (pixGetDepth(pixs) != 1))
2478         return ERROR_INT("pixs not defined or not 1 bpp", procName, 1);
2479 
2480         /* Clip box to pixs if it exists */
2481     pixGetDimensions(pixs, &bw, &bh, NULL);
2482     if (box) {
2483         if ((boxt = boxClipToRectangle(box, bw, bh)) == NULL)
2484             return ERROR_INT("invalid box", procName, 1);
2485         boxGetGeometry(boxt, &bx, &by, &bw, &bh);
2486         boxDestroy(&boxt);
2487     }
2488     else
2489         bx = by = 0;
2490     xstart = bx;
2491     ystart = by;
2492     xend = bx + bw - 1;
2493     yend = by + bh - 1;
2494 
2495     data = pixGetData(pixs);
2496     wpl = pixGetWpl(pixs);
2497     if (scanflag == L_FROM_LEFT) {
2498         for (x = xstart; x <= xend; x++) {
2499             for (y = ystart; y <= yend; y++) {
2500                 line = data + y * wpl;
2501                 if (GET_DATA_BIT(line, x)) {
2502                     *ploc = x;
2503                     return 0;
2504                 }
2505             }
2506         }
2507     }
2508     else if (scanflag == L_FROM_RIGHT) {
2509         for (x = xend; x >= xstart; x--) {
2510             for (y = ystart; y <= yend; y++) {
2511                 line = data + y * wpl;
2512                 if (GET_DATA_BIT(line, x)) {
2513                     *ploc = x;
2514                     return 0;
2515                 }
2516             }
2517         }
2518     }
2519     else if (scanflag == L_FROM_TOP) {
2520         for (y = ystart; y <= yend; y++) {
2521             line = data + y * wpl;
2522             for (x = xstart; x <= xend; x++) {
2523                 if (GET_DATA_BIT(line, x)) {
2524                     *ploc = y;
2525                     return 0;
2526                 }
2527             }
2528         }
2529     }
2530     else if (scanflag == L_FROM_BOTTOM) {
2531         for (y = yend; y >= ystart; y--) {
2532             line = data + y * wpl;
2533             for (x = xstart; x <= xend; x++) {
2534                 if (GET_DATA_BIT(line, x)) {
2535                     *ploc = y;
2536                     return 0;
2537                 }
2538             }
2539         }
2540     }
2541     else
2542         return ERROR_INT("invalid scanflag", procName, 1);
2543 
2544     return 1;  /* no fg found */
2545 }
2546 
2547 
2548 /*!
2549  *  pixClipBoxToEdges()
2550  *
2551  *      Input:  pixs (1 bpp)
2552  *              boxs  (<optional> ; use full image if null)
2553  *              lowthresh (threshold to choose clipping location)
2554  *              highthresh (threshold required to find an edge)
2555  *              maxwidth (max allowed width between low and high thresh locs)
2556  *              factor (sampling factor along pixel counting direction)
2557  *              &pixd  (<optional return> clipped pix returned)
2558  *              &boxd  (<optional return> bounding box)
2559  *      Return: 0 if OK; 1 on error or if a fg edge is not found from
2560  *              all four sides.
2561  *
2562  *  Notes:
2563  *      (1) At least one of {&pixd, &boxd} must be specified.
2564  *      (2) If there are no fg pixels, the returned ptrs are null.
2565  *      (3) This function attempts to locate rectangular "image" regions
2566  *          of high-density fg pixels, that have well-defined edges
2567  *          on the four sides.
2568  *      (4) Edges are searched for on each side, iterating in order
2569  *          from left, right, top and bottom.  As each new edge is
2570  *          found, the search box is resized to use that location.
2571  *          Once an edge is found, it is held.  If no more edges
2572  *          are found in one iteration, the search fails.
2573  *      (5) See pixScanForEdge() for usage of the thresholds and @maxwidth.
2574  *      (6) The thresholds must be at least 1, and the low threshold
2575  *          cannot be larger than the high threshold.
2576  *      (7) If the low and high thresholds are both 1, this is equivalent
2577  *          to pixClipBoxToForeground().
2578  */
2579 l_int32
pixClipBoxToEdges(PIX * pixs,BOX * boxs,l_int32 lowthresh,l_int32 highthresh,l_int32 maxwidth,l_int32 factor,PIX ** ppixd,BOX ** pboxd)2580 pixClipBoxToEdges(PIX     *pixs,
2581                   BOX     *boxs,
2582                   l_int32  lowthresh,
2583                   l_int32  highthresh,
2584                   l_int32  maxwidth,
2585                   l_int32  factor,
2586                   PIX    **ppixd,
2587                   BOX    **pboxd)
2588 {
2589 l_int32  w, h, bx, by, bw, bh, cbw, cbh, left, right, top, bottom;
2590 l_int32  lfound, rfound, tfound, bfound, change;
2591 BOX     *boxt, *boxd;
2592 
2593     PROCNAME("pixClipBoxToEdges");
2594 
2595     if (!ppixd && !pboxd)
2596         return ERROR_INT("neither &pixd nor &boxd defined", procName, 1);
2597     if (ppixd) *ppixd = NULL;
2598     if (pboxd) *pboxd = NULL;
2599     if (!pixs || (pixGetDepth(pixs) != 1))
2600         return ERROR_INT("pixs not defined or not 1 bpp", procName, 1);
2601     if (lowthresh < 1 || highthresh < 1 ||
2602         lowthresh > highthresh || maxwidth < 1)
2603         return ERROR_INT("invalid thresholds", procName, 1);
2604     factor = L_MIN(1, factor);
2605 
2606     if (lowthresh == 1 && highthresh == 1)
2607         return pixClipBoxToForeground(pixs, boxs, ppixd, pboxd);
2608 
2609     pixGetDimensions(pixs, &w, &h, NULL);
2610     if (boxs) {
2611         boxGetGeometry(boxs, &bx, &by, &bw, &bh);
2612         cbw = L_MIN(bw, w - bx);
2613         cbh = L_MIN(bh, h - by);
2614         if (cbw < 0 || cbh < 0)
2615             return ERROR_INT("box not within image", procName, 1);
2616         boxt = boxCreate(bx, by, cbw, cbh);
2617     }
2618     else
2619         boxt = boxCreate(0, 0, w, h);
2620 
2621     lfound = rfound = tfound = bfound = 0;
2622     while (!lfound || !rfound || !tfound || !bfound) {
2623         change = 0;
2624         if (!lfound) {
2625             if (!pixScanForEdge(pixs, boxt, lowthresh, highthresh, maxwidth,
2626                                 factor, L_FROM_LEFT, &left)) {
2627                 lfound = 1;
2628                 change = 1;
2629                 boxRelocateOneSide(boxt, boxt, left, L_FROM_LEFT);
2630             }
2631         }
2632         if (!rfound) {
2633             if (!pixScanForEdge(pixs, boxt, lowthresh, highthresh, maxwidth,
2634                                 factor, L_FROM_RIGHT, &right)) {
2635                 rfound = 1;
2636                 change = 1;
2637                 boxRelocateOneSide(boxt, boxt, right, L_FROM_RIGHT);
2638             }
2639         }
2640         if (!tfound) {
2641             if (!pixScanForEdge(pixs, boxt, lowthresh, highthresh, maxwidth,
2642                                 factor, L_FROM_TOP, &top)) {
2643                 tfound = 1;
2644                 change = 1;
2645                 boxRelocateOneSide(boxt, boxt, top, L_FROM_TOP);
2646             }
2647         }
2648         if (!bfound) {
2649             if (!pixScanForEdge(pixs, boxt, lowthresh, highthresh, maxwidth,
2650                                 factor, L_FROM_BOTTOM, &bottom)) {
2651                 bfound = 1;
2652                 change = 1;
2653                 boxRelocateOneSide(boxt, boxt, bottom, L_FROM_BOTTOM);
2654             }
2655         }
2656 
2657 #if DEBUG_EDGES
2658         fprintf(stderr, "iter: %d %d %d %d\n", lfound, rfound, tfound, bfound);
2659 #endif  /* DEBUG_EDGES */
2660 
2661         if (change == 0) break;
2662     }
2663     boxDestroy(&boxt);
2664 
2665     if (change == 0)
2666         return ERROR_INT("not all edges found", procName, 1);
2667 
2668     boxd = boxCreate(left, top, right - left + 1, bottom - top + 1);
2669     if (ppixd)
2670         *ppixd = pixClipRectangle(pixs, boxd, NULL);
2671     if (pboxd)
2672         *pboxd = boxd;
2673     else
2674         boxDestroy(&boxd);
2675 
2676     return 0;
2677 }
2678 
2679 
2680 /*!
2681  *  pixScanForEdge()
2682  *
2683  *      Input:  pixs (1 bpp)
2684  *              box  (<optional> within which the search is conducted)
2685  *              lowthresh (threshold to choose clipping location)
2686  *              highthresh (threshold required to find an edge)
2687  *              maxwidth (max allowed width between low and high thresh locs)
2688  *              factor (sampling factor along pixel counting direction)
2689  *              scanflag (direction of scan; e.g., L_FROM_LEFT)
2690  *              &loc (location in scan direction of first black pixel)
2691  *      Return: 0 if OK; 1 on error or if the edge is not found
2692  *
2693  *  Notes:
2694  *      (1) If there are no fg pixels, the position is set to 0.
2695  *          Caller must check the return value!
2696  *      (2) Use @box == NULL to scan from edge of pixs
2697  *      (3) As the scan progresses, the location where the sum of
2698  *          pixels equals or excees @lowthresh is noted (loc).  The
2699  *          scan is stopped when the sum of pixels equals or exceeds
2700  *          @highthresh.  If the scan distance between loc and that
2701  *          point does not exceed @maxwidth, an edge is found and
2702  *          its position is taken to be loc.  @maxwidth implicitly
2703  *          sets a minimum on the required gradient of the edge.
2704  *      (4) The thresholds must be at least 1, and the low threshold
2705  *          cannot be larger than the high threshold.
2706  */
2707 l_int32
pixScanForEdge(PIX * pixs,BOX * box,l_int32 lowthresh,l_int32 highthresh,l_int32 maxwidth,l_int32 factor,l_int32 scanflag,l_int32 * ploc)2708 pixScanForEdge(PIX      *pixs,
2709                BOX      *box,
2710                l_int32   lowthresh,
2711                l_int32   highthresh,
2712                l_int32   maxwidth,
2713                l_int32   factor,
2714                l_int32   scanflag,
2715                l_int32  *ploc)
2716 {
2717 l_int32    bx, by, bw, bh, foundmin, loc, sum, wpl;
2718 l_int32    x, xstart, xend, y, ystart, yend;
2719 l_uint32  *data, *line;
2720 BOX       *boxt;
2721 
2722     PROCNAME("pixScanForEdge");
2723 
2724     if (!ploc)
2725         return ERROR_INT("&ploc not defined", procName, 1);
2726     *ploc = 0;
2727     if (!pixs || (pixGetDepth(pixs) != 1))
2728         return ERROR_INT("pixs not defined or not 1 bpp", procName, 1);
2729     if (lowthresh < 1 || highthresh < 1 ||
2730         lowthresh > highthresh || maxwidth < 1)
2731         return ERROR_INT("invalid thresholds", procName, 1);
2732     factor = L_MIN(1, factor);
2733 
2734         /* Clip box to pixs if it exists */
2735     pixGetDimensions(pixs, &bw, &bh, NULL);
2736     if (box) {
2737         if ((boxt = boxClipToRectangle(box, bw, bh)) == NULL)
2738             return ERROR_INT("invalid box", procName, 1);
2739         boxGetGeometry(boxt, &bx, &by, &bw, &bh);
2740         boxDestroy(&boxt);
2741     }
2742     else
2743         bx = by = 0;
2744     xstart = bx;
2745     ystart = by;
2746     xend = bx + bw - 1;
2747     yend = by + bh - 1;
2748 
2749     data = pixGetData(pixs);
2750     wpl = pixGetWpl(pixs);
2751     foundmin = 0;
2752     if (scanflag == L_FROM_LEFT) {
2753         for (x = xstart; x <= xend; x++) {
2754             sum = 0;
2755             for (y = ystart; y <= yend; y += factor) {
2756                 line = data + y * wpl;
2757                 if (GET_DATA_BIT(line, x))
2758                     sum++;
2759             }
2760             if (!foundmin && sum < lowthresh)
2761                 continue;
2762             if (!foundmin) {  /* save the loc of the beginning of the edge */
2763                 foundmin = 1;
2764                 loc = x;
2765             }
2766             if (sum >= highthresh) {
2767 #if DEBUG_EDGES
2768                 fprintf(stderr, "Left: x = %d, loc = %d\n", x, loc);
2769 #endif  /* DEBUG_EDGES */
2770                 if (x - loc < maxwidth) {
2771                     *ploc = loc;
2772                     return 0;
2773                 }
2774                 else return 1;
2775             }
2776         }
2777     }
2778     else if (scanflag == L_FROM_RIGHT) {
2779         for (x = xend; x >= xstart; x--) {
2780             sum = 0;
2781             for (y = ystart; y <= yend; y += factor) {
2782                 line = data + y * wpl;
2783                 if (GET_DATA_BIT(line, x))
2784                     sum++;
2785             }
2786             if (!foundmin && sum < lowthresh)
2787                 continue;
2788             if (!foundmin) {
2789                 foundmin = 1;
2790                 loc = x;
2791             }
2792             if (sum >= highthresh) {
2793 #if DEBUG_EDGES
2794                 fprintf(stderr, "Right: x = %d, loc = %d\n", x, loc);
2795 #endif  /* DEBUG_EDGES */
2796                 if (loc - x < maxwidth) {
2797                     *ploc = loc;
2798                     return 0;
2799                 }
2800                 else return 1;
2801             }
2802         }
2803     }
2804     else if (scanflag == L_FROM_TOP) {
2805         for (y = ystart; y <= yend; y++) {
2806             sum = 0;
2807             line = data + y * wpl;
2808             for (x = xstart; x <= xend; x += factor) {
2809                 if (GET_DATA_BIT(line, x))
2810                     sum++;
2811             }
2812             if (!foundmin && sum < lowthresh)
2813                 continue;
2814             if (!foundmin) {
2815                 foundmin = 1;
2816                 loc = y;
2817             }
2818             if (sum >= highthresh) {
2819 #if DEBUG_EDGES
2820                 fprintf(stderr, "Top: y = %d, loc = %d\n", y, loc);
2821 #endif  /* DEBUG_EDGES */
2822                 if (y - loc < maxwidth) {
2823                     *ploc = loc;
2824                     return 0;
2825                 }
2826                 else return 1;
2827             }
2828         }
2829     }
2830     else if (scanflag == L_FROM_BOTTOM) {
2831         for (y = yend; y >= ystart; y--) {
2832             sum = 0;
2833             line = data + y * wpl;
2834             for (x = xstart; x <= xend; x += factor) {
2835                 if (GET_DATA_BIT(line, x))
2836                     sum++;
2837             }
2838             if (!foundmin && sum < lowthresh)
2839                 continue;
2840             if (!foundmin) {
2841                 foundmin = 1;
2842                 loc = y;
2843             }
2844             if (sum >= highthresh) {
2845 #if DEBUG_EDGES
2846                 fprintf(stderr, "Bottom: y = %d, loc = %d\n", y, loc);
2847 #endif  /* DEBUG_EDGES */
2848                 if (loc - y < maxwidth) {
2849                     *ploc = loc;
2850                     return 0;
2851                 }
2852                 else return 1;
2853             }
2854         }
2855     }
2856     else
2857         return ERROR_INT("invalid scanflag", procName, 1);
2858 
2859     return 1;  /* edge not found */
2860 }
2861 
2862