• 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  *  adaptmap.c
18  *
19  *  ===================================================================
20  *  Image binarization algorithms are found in:
21  *     grayquant.c:   standard, simple, general grayscale quantization
22  *     adaptmap.c:    local adaptive; mostly gray-to-gray in preparation
23  *                    for binarization
24  *     binarize.c:    special binarization methods, locally adaptive.
25  *  ===================================================================
26  *
27  *      Adaptive background normalization (top-level functions)
28  *          PIX       *pixBackgroundNormSimple()     8 and 32 bpp
29  *          PIX       *pixBackgroundNorm()           8 and 32 bpp
30  *          PIX       *pixBackgroundNormMorph()      8 and 32 bpp
31  *
32  *      Arrays of inverted background values for normalization (16 bpp)
33  *          l_int32    pixBackgroundNormGrayArray()   8 bpp input
34  *          l_int32    pixBackgroundNormRGBArrays()   32 bpp input
35  *          l_int32    pixBackgroundNormGrayArrayMorph()   8 bpp input
36  *          l_int32    pixBackgroundNormRGBArraysMorph()   32 bpp input
37  *
38  *      Measurement of local background
39  *          l_int32    pixGetBackgroundGrayMap()        8 bpp
40  *          l_int32    pixGetBackgroundRGBMap()         32 bpp
41  *          l_int32    pixGetBackgroundGrayMapMorph()   8 bpp
42  *          l_int32    pixGetBackgroundRGBMapMorph()    32 bpp
43  *          l_int32    pixFillMapHoles()
44  *          PIX       *pixExtendByReplication()         8 bpp
45  *          l_int32    pixSmoothConnectedRegions()      8 bpp
46  *
47  *      Measurement of local foreground
48  *          l_int32    pixGetForegroundGrayMap()        8 bpp
49  *
50  *      Generate inverted background map for each component
51  *          PIX       *pixGetInvBackgroundMap()   16 bpp
52  *
53  *      Apply inverse background map to image
54  *          PIX       *pixApplyInvBackgroundGrayMap()   8 bpp
55  *          PIX       *pixApplyInvBackgroundRGBMap()    32 bpp
56  *
57  *      Apply variable map
58  *          PIX       *pixApplyVariableGrayMap()        8 bpp
59  *
60  *      Non-adaptive (global) mapping
61  *          PIX       *pixGlobalNormRGB()               32 bpp or cmapped
62  *          PIX       *pixGlobalNormNoSatRGB()          32 bpp
63  *
64  *      Adaptive threshold spread normalization
65  *          l_int32    pixThresholdSpreadNorm()         8 bpp
66  *
67  *      Adaptive background normalization (flexible adaptaption)
68  *          PIX       *pixBackgroundNormFlex()          8 bpp
69  *
70  *      Adaptive contrast normalization
71  *          PIX             *pixContrastNorm()          8 bpp
72  *          l_int32          pixMinMaxTiles()
73  *          l_int32          pixSetLowContrast()
74  *          PIX             *pixLinearTRCTiled()
75  *          static l_int32  *iaaGetLinearTRC()
76  *
77  *  Background normalization is done by generating a reduced map (or set
78  *  of maps) representing the estimated background value of the
79  *  input image, and using this to shift the pixel values so that
80  *  this background value is set to some constant value.
81  *
82  *  Specifically, normalization has 3 steps:
83  *    (1) Generate a background map at a reduced scale.
84  *    (2) Make the array of inverted background values by inverting
85  *        the map.  The result is an array of local multiplicative factors.
86  *    (3) Apply this inverse background map to the image
87  *
88  *  The inverse background arrays can be generated in two different ways here:
89  *    (1) Remove the 'foreground' pixels and average over the remaining
90  *        pixels in each tile.  Propagate values into tiles where
91  *        values have not been assigned, either because there was not
92  *        enough background in the tile or because the tile is covered
93  *        by a foreground region described by an image mask.
94  *        After the background map is made, the inverse map is generated by
95  *        smoothing over some number of adjacent tiles
96  *        (block convolution) and then inverting.
97  *    (2) Remove the foreground pixels using a morphological closing
98  *        on a subsampled version of the image.  Propagate values
99  *        into pixels covered by an optional image mask.  Invert the
100  *        background map without preconditioning by convolutional smoothing.
101  *
102  *  Note: Several of these functions make an implicit assumption about RGB
103  *        component ordering.
104  *
105  *  Other methods for adaptively normalizing the image are also given here.
106  *
107  *  (1) pixThresholdSpreadNorm() computes a local threshold over the image
108  *      and normalizes the input pixel values so that this computed threshold
109  *      is a constant across the entire image.
110  *
111  *  (2) pixContrastNorm() computes and applies a local TRC so that the
112  *      local dynamic range is expanded to the full 8 bits, where the
113  *      darkest pixels are mapped to 0 and the lightest to 255.  This is
114  *      useful for improving the appearance of pages with very light
115  *      foreground or very dark background, and where the local TRC
116  *      function doesn't change rapidly with position.
117  */
118 
119 #include <stdio.h>
120 #include <stdlib.h>
121 #include "allheaders.h"
122 
123 
124     /* Default input parameters for pixBackgroundNormSimple()
125      * Note:
126      *    (1) mincount must never exceed the tile area (width * height)
127      *    (2) bgval must be sufficiently below 255 to avoid accidental
128      *        saturation; otherwise it should be large to avoid
129      *        shrinking the dynamic range
130      *    (3) results should otherwise not be sensitive to these values
131      */
132 static const l_int32  DEFAULT_TILE_WIDTH = 10;
133 static const l_int32  DEFAULT_TILE_HEIGHT = 15;
134 static const l_int32  DEFAULT_FG_THRESHOLD = 60;
135 static const l_int32  DEFAULT_MIN_COUNT = 40;
136 static const l_int32  DEFAULT_BG_VAL = 200;
137 static const l_int32  DEFAULT_X_SMOOTH_SIZE = 2;
138 static const l_int32  DEFAULT_Y_SMOOTH_SIZE = 1;
139 
140 static l_int32 *iaaGetLinearTRC(l_int32 **iaa, l_int32 diff);
141 
142 #ifndef  NO_CONSOLE_IO
143 #define  DEBUG_GLOBAL    0
144 #endif  /* ~NO_CONSOLE_IO */
145 
146 
147 
148 /*------------------------------------------------------------------*
149  *                Adaptive background normalization                 *
150  *------------------------------------------------------------------*/
151 /*!
152  *  pixBackgroundNormSimple()
153  *
154  *      Input:  pixs (8 bpp grayscale or 32 bpp rgb)
155  *              pixim (<optional> 1 bpp 'image' mask; can be null)
156  *              pixg (<optional> 8 bpp grayscale version; can be null)
157  *      Return: pixd (8 bpp or 32 bpp rgb), or null on error
158  *
159  *  Notes:
160  *    (1) This is a simplified interface to pixBackgroundNorm(),
161  *        where seven parameters are defaulted.
162  *    (2) The input image is either grayscale or rgb.
163  *    (3) See pixBackgroundNorm() for usage and function.
164  */
165 PIX *
pixBackgroundNormSimple(PIX * pixs,PIX * pixim,PIX * pixg)166 pixBackgroundNormSimple(PIX  *pixs,
167                         PIX  *pixim,
168                         PIX  *pixg)
169 {
170     return pixBackgroundNorm(pixs, pixim, pixg,
171                              DEFAULT_TILE_WIDTH, DEFAULT_TILE_HEIGHT,
172                              DEFAULT_FG_THRESHOLD, DEFAULT_MIN_COUNT,
173                              DEFAULT_BG_VAL, DEFAULT_X_SMOOTH_SIZE,
174                              DEFAULT_Y_SMOOTH_SIZE);
175 }
176 
177 
178 /*!
179  *  pixBackgroundNorm()
180  *
181  *      Input:  pixs (8 bpp grayscale or 32 bpp rgb)
182  *              pixim (<optional> 1 bpp 'image' mask; can be null)
183  *              pixg (<optional> 8 bpp grayscale version; can be null)
184  *              sx, sy (tile size in pixels)
185  *              thresh (threshold for determining foreground)
186  *              mincount (min threshold on counts in a tile)
187  *              bgval (target bg val; typ. > 128)
188  *              smoothx (half-width of block convolution kernel width)
189  *              smoothy (half-width of block convolution kernel height)
190  *      Return: pixd (8 bpp or 32 bpp rgb), or null on error
191  *
192  *  Notes:
193  *    (1) This is a top-level interface for normalizing the image intensity
194  *        by mapping the image so that the background is near the input
195  *        value 'bgval'.
196  *    (2) The input image is either grayscale or rgb.
197  *    (3) For each component in the input image, the background value
198  *        in each tile is estimated using the values in the tile that
199  *        are not part of the foreground, where the foreground is
200  *        determined by the input 'thresh' argument.
201  *    (4) An optional binary mask can be specified, with the foreground
202  *        pixels typically over image regions.  The resulting background
203  *        map values will be determined by surrounding pixels that are
204  *        not under the mask foreground.  The origin (0,0) of this mask
205  *        is assumed to be aligned with the origin of the input image.
206  *        This binary mask must not fully cover pixs, because then there
207  *        will be no pixels in the input image available to compute
208  *        the background.
209  *    (5) An optional grayscale version of the input pixs can be supplied.
210  *        The only reason to do this is if the input is RGB and this
211  *        grayscale version can be used elsewhere.  If the input is RGB
212  *        and this is not supplied, it is made internally using only
213  *        the green component, and destroyed after use.
214  *    (6) The dimensions of the pixel tile (sx, sy) give the amount by
215  *        by which the map is reduced in size from the input image.
216  *    (7) The threshold is used to binarize the input image, in order to
217  *        locate the foreground components.  If this is set too low,
218  *        some actual foreground may be used to determine the maps;
219  *        if set too high, there may not be enough background
220  *        to determine the map values accurately.  Typically, it's
221  *        better to err by setting the threshold too high.
222  *    (8) A 'mincount' threshold is a minimum count of pixels in a
223  *        tile for which a background reading is made, in order for that
224  *        pixel in the map to be valid.  This number should perhaps be
225  *        at least 1/3 the size of the tile.
226  *    (9) A 'bgval' target background value for the normalized image.  This
227  *        should be at least 128.  If set too close to 255, some
228  *        clipping will occur in the result.
229  *    (10) Two factors, 'smoothx' and 'smoothy', are input for smoothing
230  *        the map.  Each low-pass filter kernel dimension is
231  *        is 2 * (smoothing factor) + 1, so a
232  *        value of 0 means no smoothing. A value of 1 or 2 is recommended.
233  */
234 PIX *
pixBackgroundNorm(PIX * pixs,PIX * pixim,PIX * pixg,l_int32 sx,l_int32 sy,l_int32 thresh,l_int32 mincount,l_int32 bgval,l_int32 smoothx,l_int32 smoothy)235 pixBackgroundNorm(PIX     *pixs,
236                   PIX     *pixim,
237                   PIX     *pixg,
238                   l_int32  sx,
239                   l_int32  sy,
240                   l_int32  thresh,
241                   l_int32  mincount,
242                   l_int32  bgval,
243                   l_int32  smoothx,
244                   l_int32  smoothy)
245 {
246 l_int32  d, allfg;
247 PIX     *pixm, *pixmi, *pixd;
248 PIX     *pixmr, *pixmg, *pixmb, *pixmri, *pixmgi, *pixmbi;
249 
250     PROCNAME("pixBackgroundNorm");
251 
252     if (!pixs)
253         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
254     d = pixGetDepth(pixs);
255     if (d != 8 && d != 32)
256         return (PIX *)ERROR_PTR("pixs not 8 or 32 bpp", procName, NULL);
257     if (sx < 4 || sy < 4)
258         return (PIX *)ERROR_PTR("sx and sy must be >= 4", procName, NULL);
259     if (mincount > sx * sy) {
260         L_WARNING("mincount too large for tile size", procName);
261         mincount = (sx * sy) / 3;
262     }
263 
264         /* If pixim exists, verify that it is not all foreground. */
265     if (pixim) {
266         pixInvert(pixim, pixim);
267         pixZero(pixim, &allfg);
268         pixInvert(pixim, pixim);
269         if (allfg)
270             return (PIX *)ERROR_PTR("pixim all foreground", procName, NULL);
271     }
272 
273     pixd = NULL;
274     if (d == 8) {
275         pixm = NULL;
276         pixGetBackgroundGrayMap(pixs, pixim, sx, sy, thresh, mincount, &pixm);
277         if (!pixm) {
278             L_WARNING("map not made; returning a copy of the source", procName);
279             return pixCopy(NULL, pixs);
280         }
281 
282         pixmi = pixGetInvBackgroundMap(pixm, bgval, smoothx, smoothy);
283         if (!pixmi)
284             ERROR_PTR("pixmi not made", procName, NULL);
285         else
286             pixd = pixApplyInvBackgroundGrayMap(pixs, pixmi, sx, sy);
287 
288         pixDestroy(&pixm);
289         pixDestroy(&pixmi);
290     }
291     else {
292         pixmr = pixmg = pixmb = NULL;
293         pixGetBackgroundRGBMap(pixs, pixim, pixg, sx, sy, thresh,
294                                mincount, &pixmr, &pixmg, &pixmb);
295         if (!pixmr || !pixmg || !pixmb) {
296             pixDestroy(&pixmr);
297             pixDestroy(&pixmg);
298             pixDestroy(&pixmb);
299             L_WARNING("map not made; returning a copy of the source", procName);
300             return pixCopy(NULL, pixs);
301         }
302 
303         pixmri = pixGetInvBackgroundMap(pixmr, bgval, smoothx, smoothy);
304         pixmgi = pixGetInvBackgroundMap(pixmg, bgval, smoothx, smoothy);
305         pixmbi = pixGetInvBackgroundMap(pixmb, bgval, smoothx, smoothy);
306         if (!pixmri || !pixmgi || !pixmbi)
307             ERROR_PTR("not all pixm*i are made", procName, NULL);
308         else
309             pixd = pixApplyInvBackgroundRGBMap(pixs, pixmri, pixmgi, pixmbi,
310                                                sx, sy);
311 
312         pixDestroy(&pixmr);
313         pixDestroy(&pixmg);
314         pixDestroy(&pixmb);
315         pixDestroy(&pixmri);
316         pixDestroy(&pixmgi);
317         pixDestroy(&pixmbi);
318     }
319 
320     if (!pixd)
321         ERROR_PTR("pixd not made", procName, NULL);
322     return pixd;
323 }
324 
325 
326 /*!
327  *  pixBackgroundNormMorph()
328  *
329  *      Input:  pixs (8 bpp grayscale or 32 bpp rgb)
330  *              pixim (<optional> 1 bpp 'image' mask; can be null)
331  *              reduction (at which morph closings are done; between 2 and 16)
332  *              size (of square Sel for the closing; use an odd number)
333  *              bgval (target bg val; typ. > 128)
334  *      Return: pixd (8 bpp), or null on error
335  *
336  *  Notes:
337  *    (1) This is a top-level interface for normalizing the image intensity
338  *        by mapping the image so that the background is near the input
339  *        value 'bgval'.
340  *    (2) The input image is either grayscale or rgb.
341  *    (3) For each component in the input image, the background value
342  *        is estimated using a grayscale closing; hence the 'Morph'
343  *        in the function name.
344  *    (4) An optional binary mask can be specified, with the foreground
345  *        pixels typically over image regions.  The resulting background
346  *        map values will be determined by surrounding pixels that are
347  *        not under the mask foreground.  The origin (0,0) of this mask
348  *        is assumed to be aligned with the origin of the input image.
349  *        This binary mask must not fully cover pixs, because then there
350  *        will be no pixels in the input image available to compute
351  *        the background.
352  *    (5) The map is computed at reduced size (given by 'reduction')
353  *        from the input pixs and optional pixim.  At this scale,
354  *        pixs is closed to remove the background, using a square Sel
355  *        of odd dimension.  The product of reduction * size should be
356  *        large enough to remove most of the text foreground.
357  *    (6) No convolutional smoothing needs to be done on the map before
358  *        inverting it.
359  *    (7) A 'bgval' target background value for the normalized image.  This
360  *        should be at least 128.  If set too close to 255, some
361  *        clipping will occur in the result.
362  */
363 PIX *
pixBackgroundNormMorph(PIX * pixs,PIX * pixim,l_int32 reduction,l_int32 size,l_int32 bgval)364 pixBackgroundNormMorph(PIX     *pixs,
365                        PIX     *pixim,
366                        l_int32  reduction,
367                        l_int32  size,
368                        l_int32  bgval)
369 {
370 l_int32    d, allfg;
371 PIX       *pixm, *pixmi, *pixd;
372 PIX       *pixmr, *pixmg, *pixmb, *pixmri, *pixmgi, *pixmbi;
373 
374     PROCNAME("pixBackgroundNormMorph");
375 
376     if (!pixs)
377         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
378     d = pixGetDepth(pixs);
379     if (d != 8 && d != 32)
380         return (PIX *)ERROR_PTR("pixs not 8 or 32 bpp", procName, NULL);
381     if (reduction < 2 || reduction > 16)
382         return (PIX *)ERROR_PTR("reduction must be between 2 and 16",
383                                 procName, NULL);
384 
385         /* If pixim exists, verify that it is not all foreground. */
386     if (pixim) {
387         pixInvert(pixim, pixim);
388         pixZero(pixim, &allfg);
389         pixInvert(pixim, pixim);
390         if (allfg)
391             return (PIX *)ERROR_PTR("pixim all foreground", procName, NULL);
392     }
393 
394     pixd = NULL;
395     if (d == 8) {
396         pixGetBackgroundGrayMapMorph(pixs, pixim, reduction, size, &pixm);
397         if (!pixm)
398             return (PIX *)ERROR_PTR("pixm not made", procName, NULL);
399         pixmi = pixGetInvBackgroundMap(pixm, bgval, 0, 0);
400         if (!pixmi)
401             ERROR_PTR("pixmi not made", procName, NULL);
402         else
403             pixd = pixApplyInvBackgroundGrayMap(pixs, pixmi,
404                                                 reduction, reduction);
405         pixDestroy(&pixm);
406         pixDestroy(&pixmi);
407     }
408     else {  /* d == 32 */
409         pixmr = pixmg = pixmb = NULL;
410         pixGetBackgroundRGBMapMorph(pixs, pixim, reduction, size,
411                                     &pixmr, &pixmg, &pixmb);
412         if (!pixmr || !pixmg || !pixmb) {
413             pixDestroy(&pixmr);
414             pixDestroy(&pixmg);
415             pixDestroy(&pixmb);
416             return (PIX *)ERROR_PTR("not all pixm*", procName, NULL);
417         }
418 
419         pixmri = pixGetInvBackgroundMap(pixmr, bgval, 0, 0);
420         pixmgi = pixGetInvBackgroundMap(pixmg, bgval, 0, 0);
421         pixmbi = pixGetInvBackgroundMap(pixmb, bgval, 0, 0);
422         if (!pixmri || !pixmgi || !pixmbi)
423             ERROR_PTR("not all pixm*i are made", procName, NULL);
424         else
425             pixd = pixApplyInvBackgroundRGBMap(pixs, pixmri, pixmgi, pixmbi,
426                                                reduction, reduction);
427 
428         pixDestroy(&pixmr);
429         pixDestroy(&pixmg);
430         pixDestroy(&pixmb);
431         pixDestroy(&pixmri);
432         pixDestroy(&pixmgi);
433         pixDestroy(&pixmbi);
434     }
435 
436     if (!pixd)
437         ERROR_PTR("pixd not made", procName, NULL);
438     return pixd;
439 }
440 
441 
442 /*-------------------------------------------------------------------------*
443  *      Arrays of inverted background values for normalization             *
444  *-------------------------------------------------------------------------*
445  *  Notes for these four functions:                                        *
446  *      (1) They are useful if you need to save the actual mapping array.  *
447  *      (2) They could be used in the top-level functions but are          *
448  *          not because their use makes those functions less clear.        *
449  *      (3) Each component in the input pixs generates a 16 bpp pix array. *
450  *-------------------------------------------------------------------------*/
451 /*!
452  *  pixBackgroundNormGrayArray()
453  *
454  *      Input:  pixs (8 bpp grayscale)
455  *              pixim (<optional> 1 bpp 'image' mask; can be null)
456  *              sx, sy (tile size in pixels)
457  *              thresh (threshold for determining foreground)
458  *              mincount (min threshold on counts in a tile)
459  *              bgval (target bg val; typ. > 128)
460  *              smoothx (half-width of block convolution kernel width)
461  *              smoothy (half-width of block convolution kernel height)
462  *              &pixd (<return> 16 bpp array of inverted background value)
463  *      Return: 0 if OK, 1 on error
464  *
465  *  Notes:
466  *    (1) See notes in pixBackgroundNorm().
467  *    (2) This returns a 16 bpp pix that can be used by
468  *        pixApplyInvBackgroundGrayMap() to generate a normalized version
469  *        of the input pixs.
470  */
471 l_int32
pixBackgroundNormGrayArray(PIX * pixs,PIX * pixim,l_int32 sx,l_int32 sy,l_int32 thresh,l_int32 mincount,l_int32 bgval,l_int32 smoothx,l_int32 smoothy,PIX ** ppixd)472 pixBackgroundNormGrayArray(PIX     *pixs,
473                            PIX     *pixim,
474                            l_int32  sx,
475                            l_int32  sy,
476                            l_int32  thresh,
477                            l_int32  mincount,
478                            l_int32  bgval,
479                            l_int32  smoothx,
480                            l_int32  smoothy,
481                            PIX    **ppixd)
482 {
483 l_int32  allfg;
484 PIX     *pixm;
485 
486     PROCNAME("pixBackgroundNormGrayArray");
487 
488     if (!ppixd)
489         return ERROR_INT("&pixd not defined", procName, 1);
490     *ppixd = NULL;
491     if (!pixs)
492         return ERROR_INT("pixs not defined", procName, 1);
493     if (pixGetDepth(pixs) != 8)
494         return ERROR_INT("pixs not 8 bpp", procName, 1);
495     if (pixim && pixGetDepth(pixim) != 1)
496         return ERROR_INT("pixim not 1 bpp", procName, 1);
497     if (sx < 4 || sy < 4)
498         return ERROR_INT("sx and sy must be >= 4", procName, 1);
499     if (mincount > sx * sy) {
500         L_WARNING("mincount too large for tile size", procName);
501         mincount = (sx * sy) / 3;
502     }
503 
504         /* If pixim exists, verify that it is not all foreground. */
505     if (pixim) {
506         pixInvert(pixim, pixim);
507         pixZero(pixim, &allfg);
508         pixInvert(pixim, pixim);
509         if (allfg)
510             return ERROR_INT("pixim all foreground", procName, 1);
511     }
512 
513     pixGetBackgroundGrayMap(pixs, pixim, sx, sy, thresh, mincount, &pixm);
514     if (!pixm)
515         return ERROR_INT("pixm not made", procName, 1);
516     *ppixd = pixGetInvBackgroundMap(pixm, bgval, smoothx, smoothy);
517     pixDestroy(&pixm);
518     return 0;
519 }
520 
521 
522 /*!
523  *  pixBackgroundNormRGBArrays()
524  *
525  *      Input:  pixs (32 bpp rgb)
526  *              pixim (<optional> 1 bpp 'image' mask; can be null)
527  *              pixg (<optional> 8 bpp grayscale version; can be null)
528  *              sx, sy (tile size in pixels)
529  *              thresh (threshold for determining foreground)
530  *              mincount (min threshold on counts in a tile)
531  *              bgval (target bg val; typ. > 128)
532  *              smoothx (half-width of block convolution kernel width)
533  *              smoothy (half-width of block convolution kernel height)
534  *              &pixr (<return> 16 bpp array of inverted R background value)
535  *              &pixg (<return> 16 bpp array of inverted G background value)
536  *              &pixb (<return> 16 bpp array of inverted B background value)
537  *      Return: 0 if OK, 1 on error
538  *
539  *  Notes:
540  *    (1) See notes in pixBackgroundNorm().
541  *    (2) This returns a set of three 16 bpp pix that can be used by
542  *        pixApplyInvBackgroundGrayMap() to generate a normalized version
543  *        of each component of the input pixs.
544  */
545 l_int32
pixBackgroundNormRGBArrays(PIX * pixs,PIX * pixim,PIX * pixg,l_int32 sx,l_int32 sy,l_int32 thresh,l_int32 mincount,l_int32 bgval,l_int32 smoothx,l_int32 smoothy,PIX ** ppixr,PIX ** ppixg,PIX ** ppixb)546 pixBackgroundNormRGBArrays(PIX     *pixs,
547                            PIX     *pixim,
548                            PIX     *pixg,
549                            l_int32  sx,
550                            l_int32  sy,
551                            l_int32  thresh,
552                            l_int32  mincount,
553                            l_int32  bgval,
554                            l_int32  smoothx,
555                            l_int32  smoothy,
556                            PIX    **ppixr,
557                            PIX    **ppixg,
558                            PIX    **ppixb)
559 {
560 l_int32  allfg;
561 PIX     *pixmr, *pixmg, *pixmb;
562 
563     PROCNAME("pixBackgroundNormRGBArrays");
564 
565     if (!ppixr || !ppixg || !ppixb)
566         return ERROR_INT("&pixr, &pixg, &pixb not all defined", procName, 1);
567     *ppixr = *ppixg = *ppixb = NULL;
568     if (!pixs)
569         return ERROR_INT("pixs not defined", procName, 1);
570     if (pixGetDepth(pixs) != 32)
571         return ERROR_INT("pixs not 32 bpp", procName, 1);
572     if (pixim && pixGetDepth(pixim) != 1)
573         return ERROR_INT("pixim not 1 bpp", procName, 1);
574     if (sx < 4 || sy < 4)
575         return ERROR_INT("sx and sy must be >= 4", procName, 1);
576     if (mincount > sx * sy) {
577         L_WARNING("mincount too large for tile size", procName);
578         mincount = (sx * sy) / 3;
579     }
580 
581         /* If pixim exists, verify that it is not all foreground. */
582     if (pixim) {
583         pixInvert(pixim, pixim);
584         pixZero(pixim, &allfg);
585         pixInvert(pixim, pixim);
586         if (allfg)
587             return ERROR_INT("pixim all foreground", procName, 1);
588     }
589 
590     pixGetBackgroundRGBMap(pixs, pixim, pixg, sx, sy, thresh, mincount,
591                            &pixmr, &pixmg, &pixmb);
592     if (!pixmr || !pixmg || !pixmb) {
593         pixDestroy(&pixmr);
594         pixDestroy(&pixmg);
595         pixDestroy(&pixmb);
596         return ERROR_INT("not all pixm* made", procName, 1);
597     }
598 
599     *ppixr = pixGetInvBackgroundMap(pixmr, bgval, smoothx, smoothy);
600     *ppixg = pixGetInvBackgroundMap(pixmg, bgval, smoothx, smoothy);
601     *ppixb = pixGetInvBackgroundMap(pixmb, bgval, smoothx, smoothy);
602     pixDestroy(&pixmr);
603     pixDestroy(&pixmg);
604     pixDestroy(&pixmb);
605     return 0;
606 }
607 
608 
609 /*!
610  *  pixBackgroundNormGrayArrayMorph()
611  *
612  *      Input:  pixs (8 bpp grayscale)
613  *              pixim (<optional> 1 bpp 'image' mask; can be null)
614  *              reduction (at which morph closings are done; between 2 and 16)
615  *              size (of square Sel for the closing; use an odd number)
616  *              bgval (target bg val; typ. > 128)
617  *              &pixd (<return> 16 bpp array of inverted background value)
618  *      Return: 0 if OK, 1 on error
619  *
620  *  Notes:
621  *    (1) See notes in pixBackgroundNormMorph().
622  *    (2) This returns a 16 bpp pix that can be used by
623  *        pixApplyInvBackgroundGrayMap() to generate a normalized version
624  *        of the input pixs.
625  */
626 l_int32
pixBackgroundNormGrayArrayMorph(PIX * pixs,PIX * pixim,l_int32 reduction,l_int32 size,l_int32 bgval,PIX ** ppixd)627 pixBackgroundNormGrayArrayMorph(PIX     *pixs,
628                                 PIX     *pixim,
629                                 l_int32  reduction,
630                                 l_int32  size,
631                                 l_int32  bgval,
632                                 PIX    **ppixd)
633 {
634 l_int32  allfg;
635 PIX     *pixm;
636 
637     PROCNAME("pixBackgroundNormGrayArrayMorph");
638 
639     if (!ppixd)
640         return ERROR_INT("&pixd not defined", procName, 1);
641     *ppixd = NULL;
642     if (!pixs)
643         return ERROR_INT("pixs not defined", procName, 1);
644     if (pixGetDepth(pixs) != 8)
645         return ERROR_INT("pixs not 8 bpp", procName, 1);
646     if (pixim && pixGetDepth(pixim) != 1)
647         return ERROR_INT("pixim not 1 bpp", procName, 1);
648     if (reduction < 2 || reduction > 16)
649         return ERROR_INT("reduction must be between 2 and 16", procName, 1);
650 
651         /* If pixim exists, verify that it is not all foreground. */
652     if (pixim) {
653         pixInvert(pixim, pixim);
654         pixZero(pixim, &allfg);
655         pixInvert(pixim, pixim);
656         if (allfg)
657             return ERROR_INT("pixim all foreground", procName, 1);
658     }
659 
660     pixGetBackgroundGrayMapMorph(pixs, pixim, reduction, size, &pixm);
661     if (!pixm)
662         return ERROR_INT("pixm not made", procName, 1);
663     *ppixd = pixGetInvBackgroundMap(pixm, bgval, 0, 0);
664     pixDestroy(&pixm);
665     return 0;
666 }
667 
668 
669 /*!
670  *  pixBackgroundNormRGBArraysMorph()
671  *
672  *      Input:  pixs (32 bpp rgb)
673  *              pixim (<optional> 1 bpp 'image' mask; can be null)
674  *              reduction (at which morph closings are done; between 2 and 16)
675  *              size (of square Sel for the closing; use an odd number)
676  *              bgval (target bg val; typ. > 128)
677  *              &pixr (<return> 16 bpp array of inverted R background value)
678  *              &pixg (<return> 16 bpp array of inverted G background value)
679  *              &pixb (<return> 16 bpp array of inverted B background value)
680  *      Return: 0 if OK, 1 on error
681  *
682  *  Notes:
683  *    (1) See notes in pixBackgroundNormMorph().
684  *    (2) This returns a set of three 16 bpp pix that can be used by
685  *        pixApplyInvBackgroundGrayMap() to generate a normalized version
686  *        of each component of the input pixs.
687  */
688 l_int32
pixBackgroundNormRGBArraysMorph(PIX * pixs,PIX * pixim,l_int32 reduction,l_int32 size,l_int32 bgval,PIX ** ppixr,PIX ** ppixg,PIX ** ppixb)689 pixBackgroundNormRGBArraysMorph(PIX     *pixs,
690                                 PIX     *pixim,
691                                 l_int32  reduction,
692                                 l_int32  size,
693                                 l_int32  bgval,
694                                 PIX    **ppixr,
695                                 PIX    **ppixg,
696                                 PIX    **ppixb)
697 {
698 l_int32  allfg;
699 PIX     *pixmr, *pixmg, *pixmb;
700 
701     PROCNAME("pixBackgroundNormRGBArraysMorph");
702 
703     if (!ppixr || !ppixg || !ppixb)
704         return ERROR_INT("&pixr, &pixg, &pixb not all defined", procName, 1);
705     *ppixr = *ppixg = *ppixb = NULL;
706     if (!pixs)
707         return ERROR_INT("pixs not defined", procName, 1);
708     if (pixGetDepth(pixs) != 32)
709         return ERROR_INT("pixs not 32 bpp", procName, 1);
710     if (pixim && pixGetDepth(pixim) != 1)
711         return ERROR_INT("pixim not 1 bpp", procName, 1);
712     if (reduction < 2 || reduction > 16)
713         return ERROR_INT("reduction must be between 2 and 16", procName, 1);
714 
715         /* If pixim exists, verify that it is not all foreground. */
716     if (pixim) {
717         pixInvert(pixim, pixim);
718         pixZero(pixim, &allfg);
719         pixInvert(pixim, pixim);
720         if (allfg)
721             return ERROR_INT("pixim all foreground", procName, 1);
722     }
723 
724     pixGetBackgroundRGBMapMorph(pixs, pixim, reduction, size,
725                                 &pixmr, &pixmg, &pixmb);
726     if (!pixmr || !pixmg || !pixmb) {
727         pixDestroy(&pixmr);
728         pixDestroy(&pixmg);
729         pixDestroy(&pixmb);
730         return ERROR_INT("not all pixm* made", procName, 1);
731     }
732 
733     *ppixr = pixGetInvBackgroundMap(pixmr, bgval, 0, 0);
734     *ppixg = pixGetInvBackgroundMap(pixmg, bgval, 0, 0);
735     *ppixb = pixGetInvBackgroundMap(pixmb, bgval, 0, 0);
736     pixDestroy(&pixmr);
737     pixDestroy(&pixmg);
738     pixDestroy(&pixmb);
739     return 0;
740 }
741 
742 
743 /*------------------------------------------------------------------*
744  *                 Measurement of local background                  *
745  *------------------------------------------------------------------*/
746 /*!
747  *  pixGetBackgroundGrayMap()
748  *
749  *      Input:  pixs (8 bpp)
750  *              pixim (<optional> 1 bpp 'image' mask; can be null; it
751  *                     should not have all foreground pixels)
752  *              sx, sy (tile size in pixels)
753  *              thresh (threshold for determining foreground)
754  *              mincount (min threshold on counts in a tile)
755  *              &pixd (<return> 8 bpp grayscale map)
756  *      Return: 0 if OK, 1 on error
757  *
758  *  Notes:
759  *      (1) The background is measured in regions that don't have
760  *          images.  It is then propagated into the image regions,
761  *          and finally smoothed in each image region.
762  */
763 l_int32
pixGetBackgroundGrayMap(PIX * pixs,PIX * pixim,l_int32 sx,l_int32 sy,l_int32 thresh,l_int32 mincount,PIX ** ppixd)764 pixGetBackgroundGrayMap(PIX     *pixs,
765                         PIX     *pixim,
766                         l_int32  sx,
767                         l_int32  sy,
768                         l_int32  thresh,
769                         l_int32  mincount,
770                         PIX    **ppixd)
771 {
772 l_int32    w, h, wd, hd, wim, him, wpls, wplim, wpld, wplf;
773 l_int32    xim, yim, delx, nx, ny, i, j, k, m;
774 l_int32    count, sum, val8;
775 l_int32    empty, fgpixels;
776 l_uint32  *datas, *dataim, *datad, *dataf, *lines, *lineim, *lined, *linef;
777 l_float32  scalex, scaley;
778 PIX       *pixd, *piximi, *pixb, *pixf, *pixims;
779 
780     PROCNAME("pixGetBackgroundGrayMap");
781 
782     if (!ppixd)
783         return ERROR_INT("&pixd not defined", procName, 1);
784     *ppixd = NULL;
785     if (!pixs)
786         return ERROR_INT("pixs not defined", procName, 1);
787     if (pixGetDepth(pixs) != 8)
788         return ERROR_INT("pixs not 8 bpp", procName, 1);
789     if (pixim && pixGetDepth(pixim) != 1)
790         return ERROR_INT("pixim not 1 bpp", procName, 1);
791     if (sx < 4 || sy < 4)
792         return ERROR_INT("sx and sy must be >= 4", procName, 1);
793     if (mincount > sx * sy) {
794         L_WARNING("mincount too large for tile size", procName);
795         mincount = (sx * sy) / 3;
796     }
797 
798         /* Evaluate the 'image' mask, pixim, and make sure
799          * it is not all fg. */
800     fgpixels = 0;  /* boolean for existence of fg pixels in the image mask. */
801     if (pixim) {
802         piximi = pixInvert(NULL, pixim);  /* set non-'image' pixels to 1 */
803         pixZero(piximi, &empty);
804         pixDestroy(&piximi);
805         if (empty)
806             return ERROR_INT("pixim all fg; no background", procName, 1);
807         pixZero(pixim, &empty);
808         if (!empty)  /* there are fg pixels in pixim */
809             fgpixels = 1;
810     }
811 
812         /* Generate the foreground mask, pixf, which is at
813          * full resolution.  These pixels will be ignored when
814          * computing the background values. */
815     pixb = pixThresholdToBinary(pixs, thresh);
816     pixf = pixMorphSequence(pixb, "d7.1 + d1.7", 0);
817     pixDestroy(&pixb);
818 
819 
820     /* ------------- Set up the output map pixd --------------- */
821         /* Generate pixd, which is reduced by the factors (sx, sy). */
822     w = pixGetWidth(pixs);
823     h = pixGetHeight(pixs);
824     wd = (w + sx - 1) / sx;
825     hd = (h + sy - 1) / sy;
826     pixd = pixCreate(wd, hd, 8);
827 
828         /* Note: we only compute map values in tiles that are complete.
829          * In general, tiles at right and bottom edges will not be
830          * complete, and we must fill them in later. */
831     nx = w / sx;
832     ny = h / sy;
833     wpls = pixGetWpl(pixs);
834     datas = pixGetData(pixs);
835     wpld = pixGetWpl(pixd);
836     datad = pixGetData(pixd);
837     wplf = pixGetWpl(pixf);
838     dataf = pixGetData(pixf);
839     for (i = 0; i < ny; i++) {
840         lines = datas + sy * i * wpls;
841         linef = dataf + sy * i * wplf;
842         lined = datad + i * wpld;
843         for (j = 0; j < nx; j++) {
844             delx = j * sx;
845             sum = 0;
846             count = 0;
847             for (k = 0; k < sy; k++) {
848                 for (m = 0; m < sx; m++) {
849                     if (GET_DATA_BIT(linef + k * wplf, delx + m) == 0) {
850                         sum += GET_DATA_BYTE(lines + k * wpls, delx + m);
851                         count++;
852                     }
853                 }
854             }
855             if (count >= mincount) {
856                 val8 = sum / count;
857                 SET_DATA_BYTE(lined, j, val8);
858             }
859         }
860     }
861     pixDestroy(&pixf);
862 
863         /* If there is an optional mask with fg pixels, erase the previous
864          * calculation for the corresponding map pixels, setting the
865          * map values to 0.   Then, when all the map holes are filled,
866          * these erased pixels will be set by the surrounding map values.
867          *
868          * The calculation here is relatively efficient: for each pixel
869          * in pixd (which corresponds to a tile of mask pixels in pixim)
870          * we look only at the pixel in pixim that is at the center
871          * of the tile.  If the mask pixel is ON, we reset the map
872          * pixel in pixd to 0, so that it can later be filled in. */
873     pixims = NULL;
874     if (pixim && fgpixels) {
875         wim = pixGetWidth(pixim);
876         him = pixGetHeight(pixim);
877         dataim = pixGetData(pixim);
878         wplim = pixGetWpl(pixim);
879         for (i = 0; i < ny; i++) {
880             yim = i * sy + sy / 2;
881             if (yim >= him)
882                 break;
883             lineim = dataim + yim * wplim;
884             for (j = 0; j < nx; j++) {
885                 xim = j * sx + sx / 2;
886                 if (xim >= wim)
887                     break;
888                 if (GET_DATA_BIT(lineim, xim))
889                     pixSetPixel(pixd, j, i, 0);
890             }
891         }
892     }
893 
894         /* Fill all the holes in the map. */
895     if (pixFillMapHoles(pixd, nx, ny, L_FILL_BLACK)) {
896         pixDestroy(&pixd);
897         L_WARNING("can't make the map", procName);
898         return 1;
899     }
900 
901         /* Finally, for each connected region corresponding to the
902          * 'image' mask, reset all pixels to their average value.
903          * Each of these components represents an image (or part of one)
904          * in the input, and this smooths the background values
905          * in each of these regions. */
906     if (pixim && fgpixels) {
907         scalex = 1. / (l_float32)sx;
908         scaley = 1. / (l_float32)sy;
909         pixims = pixScaleBySampling(pixim, scalex, scaley);
910         pixSmoothConnectedRegions(pixd, pixims, 2);
911         pixDestroy(&pixims);
912     }
913 
914     *ppixd = pixd;
915     return 0;
916 }
917 
918 
919 /*!
920  *  pixGetBackgroundRGBMap()
921  *
922  *      Input:  pixs (32 bpp rgb)
923  *              pixim (<optional> 1 bpp 'image' mask; can be null; it
924  *                     should not have all foreground pixels)
925  *              pixg (<optional> 8 bpp grayscale version; can be null)
926  *              sx, sy (tile size in pixels)
927  *              thresh (threshold for determining foreground)
928  *              mincount (min threshold on counts in a tile)
929  *              &pixmr, &pixmg, &pixmb (<return> rgb maps)
930  *      Return: 0 if OK, 1 on error
931  *
932  *  Notes:
933  *      (1) If pixg, which is a grayscale version of pixs, is provided,
934  *          use this internally to generate the foreground mask.
935  *          Otherwise, a grayscale version of pixs will be generated
936  *          from the green component only, used, and destroyed.
937  */
938 l_int32
pixGetBackgroundRGBMap(PIX * pixs,PIX * pixim,PIX * pixg,l_int32 sx,l_int32 sy,l_int32 thresh,l_int32 mincount,PIX ** ppixmr,PIX ** ppixmg,PIX ** ppixmb)939 pixGetBackgroundRGBMap(PIX     *pixs,
940                        PIX     *pixim,
941                        PIX     *pixg,
942                        l_int32  sx,
943                        l_int32  sy,
944                        l_int32  thresh,
945                        l_int32  mincount,
946                        PIX    **ppixmr,
947                        PIX    **ppixmg,
948                        PIX    **ppixmb)
949 {
950 l_int32    w, h, wm, hm, wim, him, wpls, wplim, wplf;
951 l_int32    xim, yim, delx, nx, ny, i, j, k, m;
952 l_int32    count, rsum, gsum, bsum, rval, gval, bval;
953 l_int32    empty, fgpixels;
954 l_uint32   pixel;
955 l_uint32  *datas, *dataim, *dataf, *lines, *lineim, *linef;
956 l_float32  scalex, scaley;
957 PIX       *piximi, *pixgc, *pixb, *pixf, *pixims;
958 PIX       *pixmr, *pixmg, *pixmb;
959 
960     PROCNAME("pixGetBackgroundRGBMap");
961 
962     if (!ppixmr || !ppixmg || !ppixmb)
963         return ERROR_INT("&pixm* not all defined", procName, 1);
964     *ppixmr = *ppixmg = *ppixmb = NULL;
965     if (!pixs)
966         return ERROR_INT("pixs not defined", procName, 1);
967     if (pixGetDepth(pixs) != 32)
968         return ERROR_INT("pixs not 32 bpp", procName, 1);
969     if (pixim && pixGetDepth(pixim) != 1)
970         return ERROR_INT("pixim not 1 bpp", procName, 1);
971     if (sx < 4 || sy < 4)
972         return ERROR_INT("sx and sy must be >= 4", procName, 1);
973     if (mincount > sx * sy) {
974         L_WARNING("mincount too large for tile size", procName);
975         mincount = (sx * sy) / 3;
976     }
977 
978         /* Evaluate the mask pixim and make sure it is not all foreground */
979     fgpixels = 0;  /* boolean for existence of fg mask pixels */
980     if (pixim) {
981         piximi = pixInvert(NULL, pixim);  /* set non-'image' pixels to 1 */
982         pixZero(piximi, &empty);
983         pixDestroy(&piximi);
984         if (empty)
985             return ERROR_INT("pixim all fg; no background", procName, 1);
986         pixZero(pixim, &empty);
987         if (!empty)  /* there are fg pixels in pixim */
988             fgpixels = 1;
989     }
990 
991         /* Generate the foreground mask.  These pixels will be
992          * ignored when computing the background values. */
993     if (pixg)  /* use the input grayscale version if it is provided */
994         pixgc = pixClone(pixg);
995     else
996         pixgc = pixConvertRGBToGrayFast(pixs);
997     pixb = pixThresholdToBinary(pixgc, thresh);
998     pixf = pixMorphSequence(pixb, "d7.1 + d1.7", 0);
999     pixDestroy(&pixgc);
1000     pixDestroy(&pixb);
1001 
1002         /* Generate the output mask images */
1003     w = pixGetWidth(pixs);
1004     h = pixGetHeight(pixs);
1005     wm = (w + sx - 1) / sx;
1006     hm = (h + sy - 1) / sy;
1007     pixmr = pixCreate(wm, hm, 8);
1008     pixmg = pixCreate(wm, hm, 8);
1009     pixmb = pixCreate(wm, hm, 8);
1010 
1011     /* ------------- Set up the mapping images --------------- */
1012         /* Note: we only compute map values in tiles that are complete.
1013          * In general, tiles at right and bottom edges will not be
1014          * complete, and we must fill them in later. */
1015     nx = w / sx;
1016     ny = h / sy;
1017     wpls = pixGetWpl(pixs);
1018     datas = pixGetData(pixs);
1019     wplf = pixGetWpl(pixf);
1020     dataf = pixGetData(pixf);
1021     for (i = 0; i < ny; i++) {
1022         lines = datas + sy * i * wpls;
1023         linef = dataf + sy * i * wplf;
1024         for (j = 0; j < nx; j++) {
1025             delx = j * sx;
1026             rsum = gsum = bsum = 0;
1027             count = 0;
1028             for (k = 0; k < sy; k++) {
1029                 for (m = 0; m < sx; m++) {
1030                     if (GET_DATA_BIT(linef + k * wplf, delx + m) == 0) {
1031                         pixel = *(lines + k * wpls + delx + m);
1032                         rsum += (pixel >> 24);
1033                         gsum += ((pixel >> 16) & 0xff);
1034                         bsum += ((pixel >> 8) & 0xff);
1035                         count++;
1036                     }
1037                 }
1038             }
1039             if (count >= mincount) {
1040                 rval = rsum / count;
1041                 gval = gsum / count;
1042                 bval = bsum / count;
1043                 pixSetPixel(pixmr, j, i, rval);
1044                 pixSetPixel(pixmg, j, i, gval);
1045                 pixSetPixel(pixmb, j, i, bval);
1046             }
1047         }
1048     }
1049     pixDestroy(&pixf);
1050 
1051         /* If there is an optional mask with fg pixels, erase the previous
1052          * calculation for the corresponding map pixels, setting the
1053          * map values in each of the 3 color maps to 0.   Then, when
1054          * all the map holes are filled, these erased pixels will
1055          * be set by the surrounding map values. */
1056     if (pixim) {
1057         wim = pixGetWidth(pixim);
1058         him = pixGetHeight(pixim);
1059         dataim = pixGetData(pixim);
1060         wplim = pixGetWpl(pixim);
1061         for (i = 0; i < ny; i++) {
1062             yim = i * sy + sy / 2;
1063             if (yim >= him)
1064                 break;
1065             lineim = dataim + yim * wplim;
1066             for (j = 0; j < nx; j++) {
1067                 xim = j * sx + sx / 2;
1068                 if (xim >= wim)
1069                     break;
1070                 if (GET_DATA_BIT(lineim, xim)) {
1071                     pixSetPixel(pixmr, j, i, 0);
1072                     pixSetPixel(pixmg, j, i, 0);
1073                     pixSetPixel(pixmb, j, i, 0);
1074                 }
1075             }
1076         }
1077     }
1078 
1079     /* ----------------- Now fill in the holes ----------------------- */
1080     if (pixFillMapHoles(pixmr, nx, ny, L_FILL_BLACK) ||
1081         pixFillMapHoles(pixmg, nx, ny, L_FILL_BLACK) ||
1082         pixFillMapHoles(pixmb, nx, ny, L_FILL_BLACK)) {
1083         pixDestroy(&pixmr);
1084         pixDestroy(&pixmg);
1085         pixDestroy(&pixmb);
1086         L_WARNING("can't make the maps", procName);
1087         return 1;
1088     }
1089 
1090         /* Finally, for each connected region corresponding to the
1091          * fg mask, reset all pixels to their average value. */
1092     if (pixim && fgpixels) {
1093         scalex = 1. / (l_float32)sx;
1094         scaley = 1. / (l_float32)sy;
1095         pixims = pixScaleBySampling(pixim, scalex, scaley);
1096         pixSmoothConnectedRegions(pixmr, pixims, 2);
1097         pixSmoothConnectedRegions(pixmg, pixims, 2);
1098         pixSmoothConnectedRegions(pixmb, pixims, 2);
1099         pixDestroy(&pixims);
1100     }
1101 
1102     *ppixmr = pixmr;
1103     *ppixmg = pixmg;
1104     *ppixmb = pixmb;
1105     return 0;
1106 }
1107 
1108 
1109 /*!
1110  *  pixGetBackgroundGrayMapMorph()
1111  *
1112  *      Input:  pixs (8 bpp)
1113  *              pixim (<optional> 1 bpp 'image' mask; can be null; it
1114  *                     should not have all foreground pixels)
1115  *              reduction (factor at which closing is performed)
1116  *              size (of square Sel for the closing; use an odd number)
1117  *              &pixm (<return> grayscale map)
1118  *      Return: 0 if OK, 1 on error
1119  */
1120 l_int32
pixGetBackgroundGrayMapMorph(PIX * pixs,PIX * pixim,l_int32 reduction,l_int32 size,PIX ** ppixm)1121 pixGetBackgroundGrayMapMorph(PIX     *pixs,
1122                              PIX     *pixim,
1123                              l_int32  reduction,
1124                              l_int32  size,
1125                              PIX    **ppixm)
1126 {
1127 l_int32    nx, ny, empty, fgpixels;
1128 l_float32  scale;
1129 PIX       *pixm, *pixt1, *pixt2, *pixt3, *pixims;
1130 
1131     PROCNAME("pixGetBackgroundGrayMapMorph");
1132 
1133     if (!ppixm)
1134         return ERROR_INT("&pixm not defined", procName, 1);
1135     *ppixm = NULL;
1136     if (!pixs)
1137         return ERROR_INT("pixs not defined", procName, 1);
1138     if (pixGetDepth(pixs) != 8)
1139         return ERROR_INT("pixs not 8 bpp", procName, 1);
1140     if (pixim && pixGetDepth(pixim) != 1)
1141         return ERROR_INT("pixim not 1 bpp", procName, 1);
1142 
1143         /* Evaluate the mask pixim and make sure it is not all foreground. */
1144     fgpixels = 0;  /* boolean for existence of fg mask pixels */
1145     if (pixim) {
1146         pixInvert(pixim, pixim);  /* set background pixels to 1 */
1147         pixZero(pixim, &empty);
1148         if (empty)
1149             return ERROR_INT("pixim all fg; no background", procName, 1);
1150         pixInvert(pixim, pixim);  /* revert to original mask */
1151         pixZero(pixim, &empty);
1152         if (!empty)  /* there are fg pixels in pixim */
1153             fgpixels = 1;
1154     }
1155 
1156         /* Downscale as requested and do the closing to get the background. */
1157     scale = 1. / (l_float32)reduction;
1158     pixt1 = pixScaleBySampling(pixs, scale, scale);
1159     pixt2 = pixCloseGray(pixt1, size, size);
1160     pixt3 = pixExtendByReplication(pixt2, 1, 1);
1161 
1162         /* Downscale the image mask, if any, and remove it from the
1163          * background.  These pixels will be filled in (twice). */
1164     pixims = NULL;
1165     if (pixim) {
1166         pixims = pixScale(pixim, scale, scale);
1167         pixm = pixConvertTo8(pixims, FALSE);
1168         pixAnd(pixm, pixm, pixt3);
1169     }
1170     else
1171         pixm = pixClone(pixt3);
1172     pixDestroy(&pixt1);
1173     pixDestroy(&pixt2);
1174     pixDestroy(&pixt3);
1175 
1176         /* Fill all the holes in the map. */
1177     nx = pixGetWidth(pixs) / reduction;
1178     ny = pixGetHeight(pixs) / reduction;
1179     if (pixFillMapHoles(pixm, nx, ny, L_FILL_BLACK)) {
1180         pixDestroy(&pixm);
1181         L_WARNING("can't make the map", procName);
1182         return 1;
1183     }
1184 
1185         /* Finally, for each connected region corresponding to the
1186          * fg mask, reset all pixels to their average value. */
1187     if (pixim && fgpixels) {
1188         pixSmoothConnectedRegions(pixm, pixims, 2);
1189         pixDestroy(&pixims);
1190     }
1191 
1192     *ppixm = pixm;
1193     return 0;
1194 }
1195 
1196 
1197 /*!
1198  *  pixGetBackgroundRGBMapMorph()
1199  *
1200  *      Input:  pixs (32 bpp rgb)
1201  *              pixim (<optional> 1 bpp 'image' mask; can be null; it
1202  *                     should not have all foreground pixels)
1203  *              reduction (factor at which closing is performed)
1204  *              size (of square Sel for the closing; use an odd number)
1205  *              &pixmr (<return> red component map)
1206  *              &pixmg (<return> green component map)
1207  *              &pixmb (<return> blue component map)
1208  *      Return: 0 if OK, 1 on error
1209  */
1210 l_int32
pixGetBackgroundRGBMapMorph(PIX * pixs,PIX * pixim,l_int32 reduction,l_int32 size,PIX ** ppixmr,PIX ** ppixmg,PIX ** ppixmb)1211 pixGetBackgroundRGBMapMorph(PIX     *pixs,
1212                             PIX     *pixim,
1213                             l_int32  reduction,
1214                             l_int32  size,
1215                             PIX    **ppixmr,
1216                             PIX    **ppixmg,
1217                             PIX    **ppixmb)
1218 {
1219 l_int32    nx, ny, empty, fgpixels;
1220 l_float32  scale;
1221 PIX       *pixm, *pixmr, *pixmg, *pixmb, *pixt1, *pixt2, *pixt3, *pixims;
1222 
1223     PROCNAME("pixGetBackgroundRGBMapMorph");
1224 
1225     if (!ppixmr || !ppixmg || !ppixmb)
1226         return ERROR_INT("&pixm* not all defined", procName, 1);
1227     *ppixmr = *ppixmg = *ppixmb = NULL;
1228     if (!pixs)
1229         return ERROR_INT("pixs not defined", procName, 1);
1230     if (pixGetDepth(pixs) != 32)
1231         return ERROR_INT("pixs not 32 bpp", procName, 1);
1232     if (pixim && pixGetDepth(pixim) != 1)
1233         return ERROR_INT("pixim not 1 bpp", procName, 1);
1234 
1235         /* Generate an 8 bpp version of the image mask, if it exists */
1236     scale = 1. / (l_float32)reduction;
1237     pixm = NULL;
1238     if (pixim) {
1239         pixims = pixScale(pixim, scale, scale);
1240         pixm = pixConvertTo8(pixims, FALSE);
1241     }
1242 
1243         /* Evaluate the mask pixim and make sure it is not all foreground. */
1244     fgpixels = 0;  /* boolean for existence of fg mask pixels */
1245     if (pixim) {
1246         pixInvert(pixim, pixim);  /* set background pixels to 1 */
1247         pixZero(pixim, &empty);
1248         if (empty)
1249             return ERROR_INT("pixim all fg; no background", procName, 1);
1250         pixInvert(pixim, pixim);  /* revert to original mask */
1251         pixZero(pixim, &empty);
1252         if (!empty)  /* there are fg pixels in pixim */
1253             fgpixels = 1;
1254     }
1255 
1256         /* Downscale as requested and do the closing to get the background.
1257          * Then remove the image mask pixels from the background.  They
1258          * will be filled in (twice) later.  Do this for all 3 components. */
1259     pixt1 = pixScaleRGBToGrayFast(pixs, reduction, COLOR_RED);
1260     pixt2 = pixCloseGray(pixt1, size, size);
1261     pixt3 = pixExtendByReplication(pixt2, 1, 1);
1262     if (pixim)
1263         pixmr = pixAnd(NULL, pixm, pixt3);
1264     else
1265         pixmr = pixClone(pixt3);
1266     pixDestroy(&pixt1);
1267     pixDestroy(&pixt2);
1268     pixDestroy(&pixt3);
1269 
1270     pixt1 = pixScaleRGBToGrayFast(pixs, reduction, COLOR_GREEN);
1271     pixt2 = pixCloseGray(pixt1, size, size);
1272     pixt3 = pixExtendByReplication(pixt2, 1, 1);
1273     if (pixim)
1274         pixmg = pixAnd(NULL, pixm, pixt3);
1275     else
1276         pixmg = pixClone(pixt3);
1277     pixDestroy(&pixt1);
1278     pixDestroy(&pixt2);
1279     pixDestroy(&pixt3);
1280 
1281     pixt1 = pixScaleRGBToGrayFast(pixs, reduction, COLOR_BLUE);
1282     pixt2 = pixCloseGray(pixt1, size, size);
1283     pixt3 = pixExtendByReplication(pixt2, 1, 1);
1284     if (pixim)
1285         pixmb = pixAnd(NULL, pixm, pixt3);
1286     else
1287         pixmb = pixClone(pixt3);
1288     pixDestroy(&pixm);
1289     pixDestroy(&pixt1);
1290     pixDestroy(&pixt2);
1291     pixDestroy(&pixt3);
1292 
1293         /* Fill all the holes in the three maps. */
1294     nx = pixGetWidth(pixs) / reduction;
1295     ny = pixGetHeight(pixs) / reduction;
1296     if (pixFillMapHoles(pixmr, nx, ny, L_FILL_BLACK) ||
1297         pixFillMapHoles(pixmg, nx, ny, L_FILL_BLACK) ||
1298         pixFillMapHoles(pixmb, nx, ny, L_FILL_BLACK)) {
1299         pixDestroy(&pixmr);
1300         pixDestroy(&pixmg);
1301         pixDestroy(&pixmb);
1302         L_WARNING("can't make the maps", procName);
1303         return 1;
1304     }
1305 
1306         /* Finally, for each connected region corresponding to the
1307          * fg mask in each component, reset all pixels to their
1308          * average value. */
1309     if (pixim && fgpixels) {
1310         pixSmoothConnectedRegions(pixmr, pixims, 2);
1311         pixSmoothConnectedRegions(pixmg, pixims, 2);
1312         pixSmoothConnectedRegions(pixmb, pixims, 2);
1313         pixDestroy(&pixims);
1314     }
1315 
1316     *ppixmr = pixmr;
1317     *ppixmg = pixmg;
1318     *ppixmb = pixmb;
1319     return 0;
1320 }
1321 
1322 
1323 /*!
1324  *  pixFillMapHoles()
1325  *
1326  *      Input:  pix (8 bpp; a map, with one pixel for each tile in
1327  *              a larger image)
1328  *              nx (number of horizontal pixel tiles that are entirely
1329  *                  covered with pixels in the original source image)
1330  *              ny (ditto for the number of vertical pixel tiles)
1331  *              filltype (L_FILL_WHITE or L_FILL_BLACK)
1332  *      Return: 0 if OK, 1 on error
1333  *
1334  *  Notes:
1335  *      (1) This is an in-place operation on pix (the map).  pix is
1336  *          typically a low-resolution version of some other image
1337  *          from which it was derived, where each pixel in pix
1338  *          corresponds to a rectangular tile (say, m x n) of pixels
1339  *          in the larger image.  All we need to know about the larger
1340  *          image is whether or not the rightmost column and bottommost
1341  *          row of pixels in pix correspond to tiles that are
1342  *          only partially covered by pixels in the larger image.
1343  *      (2) Typically, some number of pixels in the input map are
1344  *          not known, and their values must be determined by near
1345  *          pixels that are known.  These unknown pixels are the 'holes'.
1346  *          They can take on only two values, 0 and 255, and the
1347  *          instruction about which to fill is given by the filltype flag.
1348  *      (3) The "holes" can come from two sources.  The first is when there
1349  *          are not enough foreground or background pixels in a tile;
1350  *          the second is when a tile is at least partially covered
1351  *          by an image mask.  If we're filling holes in a fg mask,
1352  *          the holes are initialized to black (0) and use L_FILL_BLACK.
1353  *          For filling holes in a bg mask, initialize the holes to
1354  *          white (255) and use L_FILL_WHITE.
1355  *      (4) If w is the map width, nx = w or nx = w - 1; ditto for h and ny.
1356  */
1357 l_int32
pixFillMapHoles(PIX * pix,l_int32 nx,l_int32 ny,l_int32 filltype)1358 pixFillMapHoles(PIX     *pix,
1359                 l_int32  nx,
1360                 l_int32  ny,
1361                 l_int32  filltype)
1362 {
1363 l_int32   w, h, d, y, nmiss, goodcol, i, j, found, ival, valtest;
1364 l_uint32  val, lastval;
1365 NUMA     *na;  /* indicates if there is any data in the column */
1366 PIX      *pixt;
1367 
1368     PROCNAME("pixFillMapHoles");
1369 
1370     if (!pix)
1371         return ERROR_INT("pix not defined", procName, 1);
1372     pixGetDimensions(pix, &w, &h, &d);
1373     if (d != 8)
1374         return ERROR_INT("pix not 8 bpp", procName, 1);
1375 
1376     /* ------------- Fill holes in the mapping image columns ----------- */
1377     na = numaCreate(0);  /* holds flag for which columns have data */
1378     nmiss = 0;
1379     valtest = (filltype == L_FILL_WHITE) ? 255 : 0;
1380     for (j = 0; j < nx; j++) {  /* do it by columns */
1381         found = FALSE;
1382         for (i = 0; i < ny; i++) {
1383             pixGetPixel(pix, j, i, &val);
1384             if (val != valtest) {
1385                 y = i;
1386                 found = TRUE;
1387                 break;
1388             }
1389         }
1390         if (found == FALSE) {
1391             numaAddNumber(na, 0);  /* no data in the column */
1392             nmiss++;
1393         }
1394         else {
1395             numaAddNumber(na, 1);  /* data in the column */
1396             for (i = y - 1; i >= 0; i--)  /* replicate upwards to top */
1397                 pixSetPixel(pix, j, i, val);
1398             pixGetPixel(pix, j, 0, &lastval);
1399             for (i = 1; i < h; i++) {  /* set going down to bottom */
1400                 pixGetPixel(pix, j, i, &val);
1401                 if (val == valtest)
1402                     pixSetPixel(pix, j, i, lastval);
1403                 else
1404                     lastval = val;
1405             }
1406         }
1407     }
1408     numaAddNumber(na, 0);  /* last column */
1409 
1410     if (nmiss == nx) {  /* no data in any column! */
1411         numaDestroy(&na);
1412         L_WARNING("no bg found; no data in any column", procName);
1413         return 1;
1414     }
1415 
1416     /* ---------- Fill in missing columns by replication ----------- */
1417     if (nmiss > 0) {  /* replicate columns */
1418         pixt = pixCopy(NULL, pix);
1419             /* Find the first good column */
1420         goodcol = 0;
1421         for (j = 0; j < w; j++) {
1422             numaGetIValue(na, j, &ival);
1423             if (ival == 1) {
1424                 goodcol = j;
1425                 break;
1426             }
1427         }
1428         if (goodcol > 0) {  /* copy cols backward */
1429             for (j = goodcol - 1; j >= 0; j--) {
1430                 pixRasterop(pix, j, 0, 1, h, PIX_SRC, pixt, j + 1, 0);
1431                 pixRasterop(pixt, j, 0, 1, h, PIX_SRC, pix, j, 0);
1432             }
1433         }
1434         for (j = goodcol + 1; j < w; j++) {   /* copy cols forward */
1435             numaGetIValue(na, j, &ival);
1436             if (ival == 0) {
1437                     /* Copy the column to the left of j */
1438                 pixRasterop(pix, j, 0, 1, h, PIX_SRC, pixt, j - 1, 0);
1439                 pixRasterop(pixt, j, 0, 1, h, PIX_SRC, pix, j, 0);
1440             }
1441         }
1442         pixDestroy(&pixt);
1443     }
1444     if (w > nx) {  /* replicate the last column */
1445         for (i = 0; i < h; i++) {
1446             pixGetPixel(pix, w - 2, i, &val);
1447             pixSetPixel(pix, w - 1, i, val);
1448         }
1449     }
1450 
1451     numaDestroy(&na);
1452     return 0;
1453 }
1454 
1455 
1456 /*!
1457  *  pixExtendByReplication()
1458  *
1459  *      Input:  pixs (8 bpp)
1460  *              addw (number of extra pixels horizontally to add)
1461  *              addh (number of extra pixels vertically to add)
1462  *      Return: pixd (extended with replicated pixel values), or null on error
1463  *
1464  *  Notes:
1465  *      (1) The pixel values are extended to the left and down, as required.
1466  */
1467 PIX *
pixExtendByReplication(PIX * pixs,l_int32 addw,l_int32 addh)1468 pixExtendByReplication(PIX     *pixs,
1469                        l_int32  addw,
1470                        l_int32  addh)
1471 {
1472 l_int32   w, h, i, j;
1473 l_uint32  val;
1474 PIX      *pixd;
1475 
1476     PROCNAME("pixExtendByReplication");
1477 
1478     if (!pixs)
1479         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
1480     if (pixGetDepth(pixs) != 8)
1481         return (PIX *)ERROR_PTR("pixs not 8 bpp", procName, NULL);
1482 
1483     if (addw == 0 && addh == 0)
1484         return pixCopy(NULL, pixs);
1485 
1486     w = pixGetWidth(pixs);
1487     h = pixGetHeight(pixs);
1488     if ((pixd = pixCreate(w + addw, h + addh, 8)) == NULL)
1489         return (PIX *)ERROR_PTR("pixd not made", procName, NULL);
1490     pixRasterop(pixd, 0, 0, w, h, PIX_SRC, pixs, 0, 0);
1491 
1492     if (addw > 0) {
1493         for (i = 0; i < h; i++) {
1494             pixGetPixel(pixd, w - 1, i, &val);
1495             for (j = 0; j < addw; j++)
1496                 pixSetPixel(pixd, w + j, i, val);
1497         }
1498     }
1499 
1500     if (addh > 0) {
1501         for (j = 0; j < w + addw; j++) {
1502             pixGetPixel(pixd, j, h - 1, &val);
1503             for (i = 0; i < addh; i++)
1504                 pixSetPixel(pixd, j, h + i, val);
1505         }
1506     }
1507 
1508     return pixd;
1509 }
1510 
1511 
1512 /*!
1513  *  pixSmoothConnectedRegions()
1514  *
1515  *      Input:  pixs (8 bpp)
1516  *              pixm (<optional> 1 bpp; if null, this is a no-op)
1517  *              factor (subsampling factor for getting average; >= 1)
1518  *      Return: 0 if OK, 1 on error
1519  *
1520  *  Notes:
1521  *      (1) The pixels in pixs corresponding to those in each
1522  *          8-connected region in the mask are set to the average value.
1523  *      (2) This is required for adaptive mapping to avoid the
1524  *          generation of stripes in the background map, due to
1525  *          variations in the pixel values near the edges of mask regions.
1526  *      (3) This function is optimized for background smoothing, where
1527  *          there are a relatively small number of components.  It will
1528  *          be inefficient if used where there are many small components.
1529  */
1530 l_int32
pixSmoothConnectedRegions(PIX * pixs,PIX * pixm,l_int32 factor)1531 pixSmoothConnectedRegions(PIX     *pixs,
1532                           PIX     *pixm,
1533                           l_int32  factor)
1534 {
1535 l_int32    empty, i, n, x, y;
1536 l_float32  aveval;
1537 BOXA      *boxa;
1538 PIX       *pixmc;
1539 PIXA      *pixa;
1540 
1541     PROCNAME("pixSmoothConnectedRegions");
1542 
1543     if (!pixs)
1544         return ERROR_INT("pixs not defined", procName, 1);
1545     if (pixGetDepth(pixs) != 8)
1546         return ERROR_INT("pixs not 8 bpp", procName, 1);
1547     if (!pixm) {
1548         L_INFO("pixm not defined", procName);
1549         return 0;
1550     }
1551     if (pixGetDepth(pixm) != 1)
1552         return ERROR_INT("pixm not 1 bpp", procName, 1);
1553     pixZero(pixm, &empty);
1554     if (empty) {
1555         L_INFO("pixm has no fg pixels; nothing to do", procName);
1556         return 0;
1557     }
1558 
1559     boxa = pixConnComp(pixm, &pixa, 8);
1560     n = boxaGetCount(boxa);
1561     for (i = 0; i < n; i++) {
1562         if ((pixmc = pixaGetPix(pixa, i, L_CLONE)) == NULL) {
1563             L_WARNING("missing pixmc!", procName);
1564             continue;
1565         }
1566         boxaGetBoxGeometry(boxa, i, &x, &y, NULL, NULL);
1567         pixGetAverageMasked(pixs, pixmc, x, y, factor, L_MEAN_ABSVAL, &aveval);
1568         pixPaintThroughMask(pixs, pixmc, x, y, (l_int32)aveval);
1569         pixDestroy(&pixmc);
1570     }
1571 
1572     boxaDestroy(&boxa);
1573     pixaDestroy(&pixa);
1574     return 0;
1575 }
1576 
1577 
1578 /*------------------------------------------------------------------*
1579  *                 Measurement of local foreground                  *
1580  *------------------------------------------------------------------*/
1581 #if 0    /* Not working properly: do not use */
1582 
1583 /*!
1584  *  pixGetForegroundGrayMap()
1585  *
1586  *      Input:  pixs (8 bpp)
1587  *              pixim (<optional> 1 bpp 'image' mask; can be null)
1588  *              sx, sy (src tile size, in pixels)
1589  *              thresh (threshold for determining foreground)
1590  *              &pixd (<return> 8 bpp grayscale map)
1591  *      Return: 0 if OK, 1 on error
1592  *
1593  *  Notes:
1594  *      (1) Each (sx, sy) tile of pixs gets mapped to one pixel in pixd.
1595  *      (2) pixd is the estimate of the fg (darkest) value within each tile.
1596  *      (3) All pixels in pixd that are in 'image' regions, as specified
1597  *          by pixim, are given the background value 0.
1598  *      (4) For pixels in pixd that can't directly be given a fg value,
1599  *          the value is inferred by propagating from neighboring pixels.
1600  *      (5) In practice, pixd can be used to normalize the fg, and
1601  *          it can be done after background normalization.
1602  *      (6) The overall procedure is:
1603  *            - reduce 2x by sampling
1604  *            - paint all 'image' pixels white, so that they don't
1605  *              participate in the Min reduction
1606  *            - do a further (sx, sy) Min reduction -- think of
1607  *              it as a large opening followed by subsampling by the
1608  *              reduction factors
1609  *            - threshold the result to identify fg, and set the
1610  *              bg pixels to 255 (these are 'holes')
1611  *            - fill holes by propagation from fg values
1612  *            - replicatively expand by 2x, arriving at the final
1613  *              resolution of pixd
1614  *            - smooth with a 17x17 kernel
1615  *            - paint the 'image' regions black
1616  */
1617 l_int32
1618 pixGetForegroundGrayMap(PIX     *pixs,
1619                         PIX     *pixim,
1620                         l_int32  sx,
1621                         l_int32  sy,
1622                         l_int32  thresh,
1623                         PIX    **ppixd)
1624 {
1625 l_int32  w, h, d, wd, hd;
1626 l_int32  empty, fgpixels;
1627 PIX     *pixd, *piximi, *pixim2, *pixims, *pixs2, *pixb, *pixt1, *pixt2, *pixt3;
1628 
1629     PROCNAME("pixGetForegroundGrayMap");
1630 
1631     if (!ppixd)
1632         return ERROR_INT("&pixd not defined", procName, 1);
1633     *ppixd = NULL;
1634     if (!pixs)
1635         return ERROR_INT("pixs not defined", procName, 1);
1636     pixGetDimensions(pixs, &w, &h, &d);
1637     if (d != 8)
1638         return ERROR_INT("pixs not 8 bpp", procName, 1);
1639     if (pixim && pixGetDepth(pixim) != 1)
1640         return ERROR_INT("pixim not 1 bpp", procName, 1);
1641     if (sx < 2 || sy < 2)
1642         return ERROR_INT("sx and sy must be >= 2", procName, 1);
1643 
1644         /* Generate pixd, which is reduced by the factors (sx, sy). */
1645     wd = (w + sx - 1) / sx;
1646     hd = (h + sy - 1) / sy;
1647     pixd = pixCreate(wd, hd, 8);
1648     *ppixd = pixd;
1649 
1650         /* Evaluate the 'image' mask, pixim.  If it is all fg,
1651          * the output pixd has all pixels with value 0. */
1652     fgpixels = 0;  /* boolean for existence of fg pixels in the image mask. */
1653     if (pixim) {
1654         piximi = pixInvert(NULL, pixim);  /* set non-image pixels to 1 */
1655         pixZero(piximi, &empty);
1656         pixDestroy(&piximi);
1657         if (empty)  /* all 'image'; return with all pixels set to 0 */
1658             return 0;
1659         pixZero(pixim, &empty);
1660         if (!empty)  /* there are fg pixels in pixim */
1661             fgpixels = 1;
1662     }
1663 
1664         /* 2x subsampling; paint white through 'image' mask. */
1665     pixs2 = pixScaleBySampling(pixs, 0.5, 0.5);
1666     if (pixim && fgpixels) {
1667         pixim2 = pixReduceBinary2(pixim, NULL);
1668         pixPaintThroughMask(pixs2, pixim2, 0, 0, 255);
1669         pixDestroy(&pixim2);
1670     }
1671 
1672         /* Min (erosion) downscaling; total reduction (4 sx, 4 sy). */
1673     pixt1 = pixScaleGrayMinMax(pixs2, sx, sy, L_CHOOSE_MIN);
1674 
1675 /*    pixDisplay(pixt1, 300, 200); */
1676 
1677         /* Threshold to identify fg; paint bg pixels to white. */
1678     pixb = pixThresholdToBinary(pixt1, thresh);  /* fg pixels */
1679     pixInvert(pixb, pixb);
1680     pixPaintThroughMask(pixt1, pixb, 0, 0, 255);
1681     pixDestroy(&pixb);
1682 
1683         /* Replicative expansion by 2x to (sx, sy). */
1684     pixt2 = pixExpandReplicate(pixt1, 2);
1685 
1686 /*    pixDisplay(pixt2, 500, 200); */
1687 
1688         /* Fill holes in the fg by propagation */
1689     pixFillMapHoles(pixt2, w / sx, h / sy, L_FILL_WHITE);
1690 
1691 /*    pixDisplay(pixt2, 700, 200); */
1692 
1693         /* Smooth with 17x17 kernel. */
1694     pixt3 = pixBlockconv(pixt2, 8, 8);
1695     pixRasterop(pixd, 0, 0, wd, hd, PIX_SRC, pixt3, 0, 0);
1696 
1697         /* Paint the image parts black. */
1698     pixims = pixScaleBySampling(pixim, 1. / sx, 1. / sy);
1699     pixPaintThroughMask(pixd, pixims, 0, 0, 0);
1700 
1701     pixDestroy(&pixs2);
1702     pixDestroy(&pixt1);
1703     pixDestroy(&pixt2);
1704     pixDestroy(&pixt3);
1705     return 0;
1706 }
1707 #endif   /* Not working properly: do not use */
1708 
1709 
1710 /*------------------------------------------------------------------*
1711  *                  Generate inverted background map                *
1712  *------------------------------------------------------------------*/
1713 /*!
1714  *  pixGetInvBackgroundMap()
1715  *
1716  *      Input:  pixs (8 bpp)
1717  *              bgval (target bg val; typ. > 128)
1718  *              smoothx (half-width of block convolution kernel width)
1719  *              smoothy (half-width of block convolution kernel height)
1720  *      Return: pixd (16 bpp), or null on error
1721  *
1722  *  Note:
1723  *     - bgval should typically be > 120 and < 240
1724  *     - pixd is a normalization image; the original image is
1725  *       multiplied by pixd and the result is divided by 256.
1726  */
1727 PIX *
pixGetInvBackgroundMap(PIX * pixs,l_int32 bgval,l_int32 smoothx,l_int32 smoothy)1728 pixGetInvBackgroundMap(PIX     *pixs,
1729                        l_int32  bgval,
1730                        l_int32  smoothx,
1731                        l_int32  smoothy)
1732 {
1733 l_int32    w, h, wplsm, wpld, i, j;
1734 l_int32    val, val16;
1735 l_uint32  *datasm, *datad, *linesm, *lined;
1736 PIX       *pixsm, *pixd;
1737 
1738     PROCNAME("pixGetInvBackgroundMap");
1739 
1740     if (!pixs)
1741         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
1742     if (pixGetDepth(pixs) != 8)
1743         return (PIX *)ERROR_PTR("pixs not 8 bpp", procName, NULL);
1744     w = pixGetWidth(pixs);
1745     h = pixGetHeight(pixs);
1746     if (w < 5 || h < 5)
1747         return (PIX *)ERROR_PTR("w and h must be >= 5", procName, NULL);
1748 
1749         /* smooth the map image */
1750     pixsm = pixBlockconv(pixs, smoothx, smoothy);
1751     datasm = pixGetData(pixsm);
1752     wplsm = pixGetWpl(pixsm);
1753 
1754         /* invert the map image, scaling up to preserve dynamic range */
1755     pixd = pixCreate(w, h, 16);
1756     datad = pixGetData(pixd);
1757     wpld = pixGetWpl(pixd);
1758     for (i = 0; i < h; i++) {
1759         linesm = datasm + i * wplsm;
1760         lined = datad + i * wpld;
1761         for (j = 0; j < w; j++) {
1762             val = GET_DATA_BYTE(linesm, j);
1763             if (val > 0)
1764                 val16 = (256 * bgval) / val;
1765             else {  /* shouldn't happen */
1766                 L_WARNING("smoothed bg has 0 pixel!", procName);
1767                 val16 = bgval / 2;
1768             }
1769             SET_DATA_TWO_BYTES(lined, j, val16);
1770         }
1771     }
1772 
1773     pixDestroy(&pixsm);
1774     return pixd;
1775 }
1776 
1777 
1778 /*------------------------------------------------------------------*
1779  *                    Apply background map to image                 *
1780  *------------------------------------------------------------------*/
1781 /*!
1782  *  pixApplyInvBackgroundGrayMap()
1783  *
1784  *      Input:  pixs (8 bpp)
1785  *              pixm (16 bpp, inverse background map)
1786  *              sx (tile width in pixels)
1787  *              sy (tile height in pixels)
1788  *      Return: pixd (8 bpp), or null on error
1789  */
1790 PIX *
pixApplyInvBackgroundGrayMap(PIX * pixs,PIX * pixm,l_int32 sx,l_int32 sy)1791 pixApplyInvBackgroundGrayMap(PIX     *pixs,
1792                              PIX     *pixm,
1793                              l_int32  sx,
1794                              l_int32  sy)
1795 {
1796 l_int32    w, h, wm, hm, wpls, wpld, i, j, k, m, xoff, yoff;
1797 l_int32    vals, vald;
1798 l_uint32   val16;
1799 l_uint32  *datas, *datad, *lines, *lined, *flines, *flined;
1800 PIX       *pixd;
1801 
1802     PROCNAME("pixApplyInvBackgroundGrayMap");
1803 
1804     if (!pixs)
1805         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
1806     if (pixGetDepth(pixs) != 8)
1807         return (PIX *)ERROR_PTR("pixs not 8 bpp", procName, NULL);
1808     if (!pixm)
1809         return (PIX *)ERROR_PTR("pixm not defined", procName, NULL);
1810     if (pixGetDepth(pixm) != 16)
1811         return (PIX *)ERROR_PTR("pixm not 16 bpp", procName, NULL);
1812     if (sx == 0 || sy == 0)
1813         return (PIX *)ERROR_PTR("invalid sx and/or sy", procName, NULL);
1814 
1815     datas = pixGetData(pixs);
1816     wpls = pixGetWpl(pixs);
1817     w = pixGetWidth(pixs);
1818     h = pixGetHeight(pixs);
1819     wm = pixGetWidth(pixm);
1820     hm = pixGetHeight(pixm);
1821     pixd = pixCreateTemplate(pixs);
1822     datad = pixGetData(pixd);
1823     wpld = pixGetWpl(pixd);
1824     for (i = 0; i < hm; i++) {
1825         lines = datas + sy * i * wpls;
1826         lined = datad + sy * i * wpld;
1827         yoff = sy * i;
1828         for (j = 0; j < wm; j++) {
1829             pixGetPixel(pixm, j, i, &val16);
1830             xoff = sx * j;
1831             for (k = 0; k < sy && yoff + k < h; k++) {
1832                 flines = lines + k * wpls;
1833                 flined = lined + k * wpld;
1834                 for (m = 0; m < sx && xoff + m < w; m++) {
1835                     vals = GET_DATA_BYTE(flines, xoff + m);
1836                     vald = (vals * val16) / 256;
1837                     vald = L_MIN(vald, 255);
1838                     SET_DATA_BYTE(flined, xoff + m, vald);
1839                 }
1840             }
1841         }
1842     }
1843 
1844     return pixd;
1845 }
1846 
1847 
1848 /*!
1849  *  pixApplyInvBackgroundRGBMap()
1850  *
1851  *      Input:  pixs (32 bpp rbg)
1852  *              pixmr (16 bpp, red inverse background map)
1853  *              pixmg (16 bpp, green inverse background map)
1854  *              pixmb (16 bpp, blue inverse background map)
1855  *              sx (tile width in pixels)
1856  *              sy (tile height in pixels)
1857  *      Return: pixd (32 bpp rbg), or null on error
1858  */
1859 PIX *
pixApplyInvBackgroundRGBMap(PIX * pixs,PIX * pixmr,PIX * pixmg,PIX * pixmb,l_int32 sx,l_int32 sy)1860 pixApplyInvBackgroundRGBMap(PIX     *pixs,
1861                             PIX     *pixmr,
1862                             PIX     *pixmg,
1863                             PIX     *pixmb,
1864                             l_int32  sx,
1865                             l_int32  sy)
1866 {
1867 l_int32    w, h, wm, hm, wpls, wpld, i, j, k, m, xoff, yoff;
1868 l_int32    rvald, gvald, bvald;
1869 l_uint32   vals;
1870 l_uint32   rval16, gval16, bval16;
1871 l_uint32  *datas, *datad, *lines, *lined, *flines, *flined;
1872 PIX       *pixd;
1873 
1874     PROCNAME("pixApplyInvBackgroundRGBMap");
1875 
1876     if (!pixs)
1877         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
1878     if (pixGetDepth(pixs) != 32)
1879         return (PIX *)ERROR_PTR("pixs not 32 bpp", procName, NULL);
1880     if (!pixmr || !pixmg || !pixmb)
1881         return (PIX *)ERROR_PTR("pix maps not all defined", procName, NULL);
1882     if (pixGetDepth(pixmr) != 16 || pixGetDepth(pixmg) != 16 ||
1883         pixGetDepth(pixmb) != 16)
1884         return (PIX *)ERROR_PTR("pix maps not all 16 bpp", procName, NULL);
1885     if (sx == 0 || sy == 0)
1886         return (PIX *)ERROR_PTR("invalid sx and/or sy", procName, NULL);
1887 
1888     datas = pixGetData(pixs);
1889     wpls = pixGetWpl(pixs);
1890     w = pixGetWidth(pixs);
1891     h = pixGetHeight(pixs);
1892     wm = pixGetWidth(pixmr);
1893     hm = pixGetHeight(pixmr);
1894     pixd = pixCreateTemplate(pixs);
1895     datad = pixGetData(pixd);
1896     wpld = pixGetWpl(pixd);
1897     for (i = 0; i < hm; i++) {
1898         lines = datas + sy * i * wpls;
1899         lined = datad + sy * i * wpld;
1900         yoff = sy * i;
1901         for (j = 0; j < wm; j++) {
1902             pixGetPixel(pixmr, j, i, &rval16);
1903             pixGetPixel(pixmg, j, i, &gval16);
1904             pixGetPixel(pixmb, j, i, &bval16);
1905             xoff = sx * j;
1906             for (k = 0; k < sy && yoff + k < h; k++) {
1907                 flines = lines + k * wpls;
1908                 flined = lined + k * wpld;
1909                 for (m = 0; m < sx && xoff + m < w; m++) {
1910                     vals = *(flines + xoff + m);
1911                     rvald = ((vals >> 24) * rval16) / 256;
1912                     rvald = L_MIN(rvald, 255);
1913                     gvald = (((vals >> 16) & 0xff) * gval16) / 256;
1914                     gvald = L_MIN(gvald, 255);
1915                     bvald = (((vals >> 8) & 0xff) * bval16) / 256;
1916                     bvald = L_MIN(bvald, 255);
1917                     composeRGBPixel(rvald, gvald, bvald, flined + xoff + m);
1918                 }
1919             }
1920         }
1921     }
1922 
1923     return pixd;
1924 }
1925 
1926 
1927 /*------------------------------------------------------------------*
1928  *                         Apply variable map                       *
1929  *------------------------------------------------------------------*/
1930 /*!
1931  *  pixApplyVariableGrayMap()
1932  *
1933  *      Input:  pixs (8 bpp)
1934  *              pixg (8 bpp, variable map)
1935  *              target (typ. 128 for threshold)
1936  *      Return: pixd (8 bpp), or null on error
1937  *
1938  *  Notes:
1939  *      (1) Suppose you have an image that you want to transform based
1940  *          on some photometric measurement at each point, such as the
1941  *          threshold value for binarization.  Representing the photometric
1942  *          measurement as an image pixg, you can threshold in input image
1943  *          using pixVarThresholdToBinary().  Alternatively, you can map
1944  *          the input image pointwise so that the threshold over the
1945  *          entire image becomes a constant, such as 128.  For example,
1946  *          if a pixel in pixg is 150 and the target is 128, the
1947  *          corresponding pixel in pixs is mapped linearly to a value
1948  *          (128/150) of the input value.  If the resulting mapped image
1949  *          pixd were then thresholded at 128, you would obtain the
1950  *          same result as a direct binarization using pixg with
1951  *          pixVarThresholdToBinary().
1952  *      (2) The sizes of pixs and pixg must be equal.
1953  */
1954 PIX *
pixApplyVariableGrayMap(PIX * pixs,PIX * pixg,l_int32 target)1955 pixApplyVariableGrayMap(PIX     *pixs,
1956                         PIX     *pixg,
1957                         l_int32  target)
1958 {
1959 l_int32    i, j, w, h, d, wpls, wplg, wpld, vals, valg, vald;
1960 l_uint8   *lut;
1961 l_uint32  *datas, *datag, *datad, *lines, *lineg, *lined;
1962 l_float32  fval;
1963 PIX       *pixd;
1964 
1965     PROCNAME("pixApplyVariableGrayMap");
1966 
1967     if (!pixs)
1968         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
1969     if (!pixg)
1970         return (PIX *)ERROR_PTR("pixg not defined", procName, NULL);
1971     if (!pixSizesEqual(pixs, pixg))
1972         return (PIX *)ERROR_PTR("pix sizes not equal", procName, NULL);
1973     pixGetDimensions(pixs, &w, &h, &d);
1974     if (d != 8)
1975         return (PIX *)ERROR_PTR("depth not 8 bpp", procName, NULL);
1976 
1977         /* Generate a LUT for the mapping if the image is large enough
1978          * to warrant the overhead.  The LUT is of size 2^16.  For the
1979          * index to the table, get the MSB from pixs and the LSB from pixg.
1980          * Note: this LUT is bigger than the typical 32K L1 cache, so
1981          * we expect cache misses.  L2 latencies are about 5ns.  But
1982          * division is slooooow.  For large images, this function is about
1983          * 4x faster when using the LUT.  C'est la vie.  */
1984     lut = NULL;
1985     if (w * h > 100000) {  /* more pixels than 2^16 */
1986         if ((lut = (l_uint8 *)CALLOC(0x10000, sizeof(l_uint8))) == NULL)
1987             return (PIX *)ERROR_PTR("lut not made", procName, NULL);
1988         for (i = 0; i < 256; i++) {
1989             for (j = 0; j < 256; j++) {
1990                 fval = (l_float32)(i * target) / (j + 0.5);
1991                 lut[(i << 8) + j] = L_MIN(255, (l_int32)(fval + 0.5));
1992             }
1993         }
1994     }
1995 
1996     pixd = pixCreateNoInit(w, h, 8);
1997     datad = pixGetData(pixd);
1998     wpld = pixGetWpl(pixd);
1999     datas = pixGetData(pixs);
2000     wpls = pixGetWpl(pixs);
2001     datag = pixGetData(pixg);
2002     wplg = pixGetWpl(pixg);
2003     for (i = 0; i < h; i++) {
2004         lines = datas + i * wpls;
2005         lineg = datag + i * wplg;
2006         lined = datad + i * wpld;
2007         if (lut) {
2008             for (j = 0; j < w; j++) {
2009                 vals = GET_DATA_BYTE(lines, j);
2010                 valg = GET_DATA_BYTE(lineg, j);
2011                 vald = lut[(vals << 8) + valg];
2012                 SET_DATA_BYTE(lined, j, vald);
2013             }
2014         }
2015         else {
2016             for (j = 0; j < w; j++) {
2017                 vals = GET_DATA_BYTE(lines, j);
2018                 valg = GET_DATA_BYTE(lineg, j);
2019                 fval = (l_float32)(vals * target) / (valg + 0.5);
2020                 vald = L_MIN(255, (l_int32)(fval + 0.5));
2021                 SET_DATA_BYTE(lined, j, vald);
2022             }
2023         }
2024     }
2025 
2026     if (lut) FREE(lut);
2027     return pixd;
2028 }
2029 
2030 
2031 /*------------------------------------------------------------------*
2032  *                  Non-adaptive (global) mapping                   *
2033  *------------------------------------------------------------------*/
2034 /*!
2035  *  pixGlobalNormRGB()
2036  *
2037  *      Input:  pixd (<optional> null, existing or equal to pixs)
2038  *              pixs (32 bpp rgb, or colormapped)
2039  *              rval, gval, bval (pixel values in pixs that are
2040  *                                linearly mapped to mapval)
2041  *              mapval (use 255 for mapping to white)
2042  *      Return: pixd (32 bpp rgb or colormapped), or null on error
2043  *
2044  *  Notes:
2045  *    (1) The value of pixd determines if the results are written to a
2046  *        new pix (use NULL), in-place to pixs (use pixs), or to some
2047  *        other existing pix.
2048  *    (2) This does a global normalization of an image where the
2049  *        r,g,b color components are not balanced.  Thus, white in pixs is
2050  *        represented by a set of r,g,b values that are not all 255.
2051  *    (3) The input values (rval, gval, bval) should be chosen to
2052  *        represent the gray color (mapval, mapval, mapval) in src.
2053  *        Thus, this function will map (rval, gval, bval) to that gray color.
2054  *    (4) Typically, mapval = 255, so that (rval, gval, bval)
2055  *        corresponds to the white point of src.  In that case, these
2056  *        parameters should be chosen so that few pixels have higher values.
2057  *    (5) In all cases, we do a linear TRC separately on each of the
2058  *        components, saturating at 255.
2059  *    (6) If the input pix is 8 bpp without a colormap, you can get
2060  *        this functionality with mapval = 255 by calling:
2061  *            pixGammaTRC(pixd, pixs, 1.0, 0, bgval);
2062  *        where bgval is the value you want to be mapped to 255.
2063  *        Or more generally, if you want bgval to be mapped to mapval:
2064  *            pixGammaTRC(pixd, pixs, 1.0, 0, 255 * bgval / mapval);
2065  */
2066 PIX *
pixGlobalNormRGB(PIX * pixd,PIX * pixs,l_int32 rval,l_int32 gval,l_int32 bval,l_int32 mapval)2067 pixGlobalNormRGB(PIX     *pixd,
2068                  PIX     *pixs,
2069                  l_int32  rval,
2070                  l_int32  gval,
2071                  l_int32  bval,
2072                  l_int32  mapval)
2073 {
2074 l_int32    w, h, d, i, j, ncolors, rv, gv, bv, wpl;
2075 l_int32   *rarray, *garray, *barray;
2076 l_uint32  *data, *line;
2077 NUMA      *nar, *nag, *nab;
2078 PIXCMAP   *cmap;
2079 
2080     PROCNAME("pixGlobalNormRGB");
2081 
2082     if (!pixs)
2083         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
2084     cmap = pixGetColormap(pixs);
2085     pixGetDimensions(pixs, &w, &h, &d);
2086     if (!cmap && d != 32)
2087         return (PIX *)ERROR_PTR("pixs not cmapped or 32 bpp", procName, NULL);
2088     if (mapval <= 0) {
2089         L_WARNING("mapval must be > 0; setting to 255", procName);
2090         mapval = 255;
2091     }
2092 
2093         /* Prepare pixd to be a copy of pixs */
2094     if ((pixd = pixCopy(pixd, pixs)) == NULL)
2095         return (PIX *)ERROR_PTR("pixd not made", procName, NULL);
2096 
2097         /* Generate the TRC maps for each component.  Make sure the
2098          * upper range for each color is greater than zero. */
2099     nar = numaGammaTRC(1.0, 0, L_MAX(1, 255 * rval / mapval));
2100     nag = numaGammaTRC(1.0, 0, L_MAX(1, 255 * gval / mapval));
2101     nab = numaGammaTRC(1.0, 0, L_MAX(1, 255 * bval / mapval));
2102     if (!nar || !nag || !nab)
2103         return (PIX *)ERROR_PTR("trc maps not all made", procName, pixd);
2104 
2105         /* Extract copies of the internal arrays */
2106     rarray = numaGetIArray(nar);
2107     garray = numaGetIArray(nag);
2108     barray = numaGetIArray(nab);
2109     if (!rarray || !garray || !barray)
2110         return (PIX *)ERROR_PTR("*arrays not all made", procName, pixd);
2111 
2112     if (cmap) {
2113         ncolors = pixcmapGetCount(cmap);
2114         for (i = 0; i < ncolors; i++) {
2115             pixcmapGetColor(cmap, i, &rv, &gv, &bv);
2116             pixcmapResetColor(cmap, i, rarray[rv], garray[gv], barray[bv]);
2117         }
2118     }
2119     else {
2120         data = pixGetData(pixd);
2121         wpl = pixGetWpl(pixd);
2122         for (i = 0; i < h; i++) {
2123             line = data + i * wpl;
2124             for (j = 0; j < w; j++) {
2125                 extractRGBValues(line[j], &rv, &gv, &bv);
2126                 composeRGBPixel(rarray[rv], garray[gv], barray[bv], line + j);
2127             }
2128         }
2129     }
2130 
2131     numaDestroy(&nar);
2132     numaDestroy(&nag);
2133     numaDestroy(&nab);
2134     FREE(rarray);
2135     FREE(garray);
2136     FREE(barray);
2137     return pixd;
2138 }
2139 
2140 
2141 /*!
2142  *  pixGlobalNormNoSatRGB()
2143  *
2144  *      Input:  pixd (<optional> null, existing or equal to pixs)
2145  *              pixs (32 bpp rgb)
2146  *              rval, gval, bval (pixel values in pixs that are
2147  *                                linearly mapped to mapval; but see below)
2148  *              factor (subsampling factor; integer >= 1)
2149  *              rank (between 0.0 and 1.0; typ. use a value near 1.0)
2150  *      Return: pixd (32 bpp rgb), or null on error
2151  *
2152  *  Notes:
2153  *    (1) This is a version of pixGlobalNormRGB(), where the output
2154  *        intensity is scaled back so that a controlled fraction of
2155  *        pixel components is allowed to saturate.  See comments in
2156  *        pixGlobalNormRGB().
2157  *    (2) The value of pixd determines if the results are written to a
2158  *        new pix (use NULL), in-place to pixs (use pixs), or to some
2159  *        other existing pix.
2160  *    (3) This does a global normalization of an image where the
2161  *        r,g,b color components are not balanced.  Thus, white in pixs is
2162  *        represented by a set of r,g,b values that are not all 255.
2163  *    (4) The input values (rval, gval, bval) can be chosen to be the
2164  *        color that, after normalization, becomes white background.
2165  *        For images that are mostly background, the closer these values
2166  *        are to the median component values, the closer the resulting
2167  *        background will be to gray, becoming white at the brightest places.
2168  *    (5) The mapval used in pixGlobalNormRGB() is computed here to
2169  *        avoid saturation of any component in the image (save for a
2170  *        fraction of the pixels given by the input rank value).
2171  */
2172 PIX *
pixGlobalNormNoSatRGB(PIX * pixd,PIX * pixs,l_int32 rval,l_int32 gval,l_int32 bval,l_int32 factor,l_float32 rank)2173 pixGlobalNormNoSatRGB(PIX       *pixd,
2174                       PIX       *pixs,
2175                       l_int32    rval,
2176                       l_int32    gval,
2177                       l_int32    bval,
2178                       l_int32    factor,
2179                       l_float32  rank)
2180 {
2181 l_int32    mapval;
2182 l_float32  rankrval, rankgval, rankbval;
2183 l_float32  rfract, gfract, bfract, maxfract;
2184 
2185     PROCNAME("pixGlobalNormNoSatRGB");
2186 
2187     if (!pixs)
2188         return (PIX *)ERROR_PTR("pixs not defined", procName, NULL);
2189     if (pixGetDepth(pixs) != 32)
2190         return (PIX *)ERROR_PTR("pixs not 32 bpp", procName, NULL);
2191     if (factor < 1)
2192         return (PIX *)ERROR_PTR("sampling factor < 1", procName, NULL);
2193     if (rank < 0.0 || rank > 1.0)
2194         return (PIX *)ERROR_PTR("rank not in [0.0 ... 1.0]", procName, NULL);
2195     if (rval <= 0 || gval <= 0 || bval <= 0)
2196         return (PIX *)ERROR_PTR("invalid estim. color values", procName, NULL);
2197 
2198         /* The max value for each component may be larger than the
2199          * input estimated background value.  In that case, mapping
2200          * for those pixels would saturate.  To prevent saturation,
2201          * we compute the fraction for each component by which we
2202          * would oversaturate.  Then take the max of these, and
2203          * reduce, uniformly over all components, the output intensity
2204          * by this value.  Then no component will saturate.
2205          * In practice, if rank < 1.0, a fraction of pixels
2206          * may have a component saturate.  By keeping rank close to 1.0,
2207          * that fraction can be made arbitrarily small. */
2208     pixGetRankValueMaskedRGB(pixs, NULL, 0, 0, factor, rank, &rankrval,
2209                              &rankgval, &rankbval);
2210     rfract = rankrval / (l_float32)rval;
2211     gfract = rankgval / (l_float32)gval;
2212     bfract = rankbval / (l_float32)bval;
2213     maxfract = L_MAX(rfract, gfract);
2214     maxfract = L_MAX(maxfract, bfract);
2215 #if  DEBUG_GLOBAL
2216     fprintf(stderr, "rankrval = %7.2f, rankgval = %7.2f, rankbval = %7.2f\n",
2217             rankrval, rankgval, rankbval);
2218     fprintf(stderr, "rfract = %7.4f, gfract = %7.4f, bfract = %7.4f\n",
2219             rfract, gfract, bfract);
2220 #endif  /* DEBUG_GLOBAL */
2221 
2222     mapval = (l_int32)(255. / maxfract);
2223     pixd = pixGlobalNormRGB(pixd, pixs, rval, gval, bval, mapval);
2224     return pixd;
2225 }
2226 
2227 
2228 /*------------------------------------------------------------------*
2229  *              Adaptive threshold spread normalization             *
2230  *------------------------------------------------------------------*/
2231 /*!
2232  *  pixThresholdSpreadNorm()
2233  *
2234  *      Input:  pixs (8 bpp)
2235  *              filtertype (L_SOBEL_EDGE or L_TWO_SIDED_EDGE);
2236  *              edgethresh (threshold on magnitude of edge filter; typ 10-20)
2237  *              smoothx, smoothy (half-width of convolution kernel applied to
2238  *                                spread threshold: use 0 for no smoothing)
2239  *              gamma (gamma correction; typ. about 0.7)
2240  *              minval  (input value that gives 0 for output; typ. -25)
2241  *              maxval  (input value that gives 255 for output; typ. 255)
2242  *              targetthresh (target threshold for normalization)
2243  *              &pixth (<optional return> computed local threshold value)
2244  *              &pixb (<optional return> thresholded normalized image)
2245  *              &pixd (<optional return> normalized image)
2246  *      Return: 0 if OK, 1 on error
2247  *
2248  *  Notes:
2249  *      (1) The basis of this approach is the use of seed spreading
2250  *          on a (possibly) sparse set of estimates for the local threshold.
2251  *          The resulting dense estimates are smoothed by convolution
2252  *          and used to either threshold the input image or normalize it
2253  *          with a local transformation that linearly maps the pixels so
2254  *          that the local threshold estimate becomes constant over the
2255  *          resulting image.  This approach is one of several that
2256  *          have been suggested (and implemented) by Ray Smith.
2257  *      (2) You can use either the Sobel or TwoSided edge filters.
2258  *          The results appear to be similar, using typical values
2259  *          of edgethresh in the rang 10-20.
2260  *      (3) To skip the trc enhancement, use gamma = 1.0, minval = 0
2261  *          and maxval = 255.
2262  *      (4) For the normalized image pixd, each pixel is linearly mapped
2263  *          in such a way that the local threshold is equal to targetthresh.
2264  *      (5) The full width and height of the convolution kernel
2265  *          are (2 * smoothx + 1) and (2 * smoothy + 1).
2266  *      (6) This function can be used with the pixtiling utility if the
2267  *          images are too large.  See pixOtsuAdaptiveThreshold() for
2268  *          an example of this.
2269  */
2270 l_int32
pixThresholdSpreadNorm(PIX * pixs,l_int32 filtertype,l_int32 edgethresh,l_int32 smoothx,l_int32 smoothy,l_float32 gamma,l_int32 minval,l_int32 maxval,l_int32 targetthresh,PIX ** ppixth,PIX ** ppixb,PIX ** ppixd)2271 pixThresholdSpreadNorm(PIX       *pixs,
2272                        l_int32    filtertype,
2273                        l_int32    edgethresh,
2274                        l_int32    smoothx,
2275                        l_int32    smoothy,
2276                        l_float32  gamma,
2277                        l_int32    minval,
2278                        l_int32    maxval,
2279                        l_int32    targetthresh,
2280                        PIX      **ppixth,
2281                        PIX      **ppixb,
2282                        PIX      **ppixd)
2283 {
2284 l_int32  w, h, d;
2285 PIX     *pixe, *pixet, *pixsd, *pixg1, *pixg2, *pixth;
2286 
2287     PROCNAME("pixThresholdSpreadNorm");
2288 
2289     if (ppixth) *ppixth = NULL;
2290     if (ppixb) *ppixb = NULL;
2291     if (ppixd) *ppixd = NULL;
2292     if (!pixs)
2293         return ERROR_INT("pixs not defined", procName, 1);
2294     pixGetDimensions(pixs, &w, &h, &d);
2295     if (d != 8)
2296         return ERROR_INT("pixs not 8 bpp", procName, 1);
2297     if (!ppixth && !ppixb && !ppixd)
2298         return ERROR_INT("no output requested", procName, 1);
2299     if (filtertype != L_SOBEL_EDGE && filtertype != L_TWO_SIDED_EDGE)
2300         return ERROR_INT("invalid filter type", procName, 1);
2301 
2302         /* Get the thresholded edge pixels.  These are the ones
2303          * that have values in pixs near the local optimal fg/bg threshold. */
2304     if (filtertype == L_SOBEL_EDGE)
2305         pixe = pixSobelEdgeFilter(pixs, L_VERTICAL_EDGES);
2306     else  /* L_TWO_SIDED_EDGE */
2307         pixe = pixTwoSidedEdgeFilter(pixs, L_VERTICAL_EDGES);
2308     pixet = pixThresholdToBinary(pixe, edgethresh);
2309     pixInvert(pixet, pixet);
2310 
2311         /* Build a seed image whose only nonzero values are those
2312          * values of pixs corresponding to pixels in the fg of pixet. */
2313     pixsd = pixCreateTemplate(pixs);
2314     pixCombineMasked(pixsd, pixs, pixet);
2315 
2316         /* Spread the seed and optionally smooth to reduce noise */
2317     pixg1 = pixSeedspread(pixsd, 4);
2318     pixg2 = pixBlockconv(pixg1, smoothx, smoothy);
2319 
2320         /* Optionally do a gamma enhancement */
2321     pixth = pixGammaTRC(NULL, pixg2, gamma, minval, maxval);
2322 
2323         /* Do the mapping and thresholding */
2324     if (ppixd) {
2325         *ppixd = pixApplyVariableGrayMap(pixs, pixth, targetthresh);
2326         if (ppixb)
2327             *ppixb = pixThresholdToBinary(*ppixd, targetthresh);
2328     }
2329     else if (ppixb)
2330         *ppixb = pixVarThresholdToBinary(pixs, pixth);
2331 
2332     if (ppixth)
2333         *ppixth = pixth;
2334     else
2335         pixDestroy(&pixth);
2336 
2337     pixDestroy(&pixe);
2338     pixDestroy(&pixet);
2339     pixDestroy(&pixsd);
2340     pixDestroy(&pixg1);
2341     pixDestroy(&pixg2);
2342     return 0;
2343 }
2344 
2345 
2346 /*------------------------------------------------------------------*
2347  *      Adaptive background normalization (flexible adaptaption)    *
2348  *------------------------------------------------------------------*/
2349 /*!
2350  *  pixBackgroundNormFlex()
2351  *
2352  *      Input:  pixs (8 bpp)
2353  *              sx, sy (desired tile dimensions; actual size may vary; use
2354  *                      values between 3 and 10)
2355  *              smoothx, smoothy (half-width of convolution kernel applied to
2356  *                                threshold array: use values between 1 and 3)
2357  *              delta (difference parameter in basin filling; use 0
2358  *                     to skip)
2359  *      Return: pixd (8 bpp, background-normalized), or null on error)
2360  *
2361  *  Notes:
2362  *      (1) This does adaptation flexibly to a quickly varying background.
2363  *          For that reason, all input parameters should be small.
2364  *      (2) sx and sy give the tile size; they should be in [5 - 7].
2365  *      (3) The full width and height of the convolution kernel
2366  *          are (2 * smoothx + 1) and (2 * smoothy + 1).  They
2367  *          should be in [1 - 2].
2368  *      (4) Basin filling is used to fill the large fg regions.  The
2369  *          parameter @delta measures the height that the black
2370  *          background is raised from the local minima.  By raising
2371  *          the background, it is possible to threshold the large
2372  *          fg regions to foreground.  If @delta is too large,
2373  *          bg regions will be lifted, causing thickening of
2374  *          the fg regions.  Use 0 to skip.
2375  */
2376 PIX *
pixBackgroundNormFlex(PIX * pixs,l_int32 sx,l_int32 sy,l_int32 smoothx,l_int32 smoothy,l_int32 delta)2377 pixBackgroundNormFlex(PIX     *pixs,
2378                       l_int32  sx,
2379                       l_int32  sy,
2380                       l_int32  smoothx,
2381                       l_int32  smoothy,
2382                       l_int32  delta)
2383 {
2384 l_float32  scalex, scaley;
2385 PIX       *pixt, *pixsd, *pixmin, *pixbg, *pixbgi, *pixd;
2386 
2387     PROCNAME("pixBackgroundNormFlex");
2388 
2389     if (!pixs || pixGetDepth(pixs) != 8)
2390         return (PIX *)ERROR_PTR("pixs undefined or not 8 bpp", procName, NULL);
2391     if (sx < 3 || sy < 3)
2392         return (PIX *)ERROR_PTR("sx and/or sy less than 3", procName, NULL);
2393     if (sx > 10 || sy > 10)
2394         return (PIX *)ERROR_PTR("sx and/or sy exceed 10", procName, NULL);
2395     if (smoothx < 1 || smoothy < 1)
2396         return (PIX *)ERROR_PTR("smooth params less than 1", procName, NULL);
2397     if (smoothx > 3 || smoothy > 3)
2398         return (PIX *)ERROR_PTR("smooth params exceed 3", procName, NULL);
2399 
2400         /* Generate the bg estimate using smoothed average with subsampling */
2401     scalex = 1. / (l_float32)sx;
2402     scaley = 1. / (l_float32)sy;
2403     pixt = pixScaleSmooth(pixs, scalex, scaley);
2404 
2405         /* Do basin filling on the bg estimate if requested */
2406     if (delta <= 0)
2407         pixsd = pixClone(pixt);
2408     else {
2409         pixLocalExtrema(pixt, 0, 0, &pixmin, NULL);
2410         pixsd = pixSeedfillGrayBasin(pixmin, pixt, delta, 4);
2411         pixDestroy(&pixmin);
2412     }
2413     pixbg = pixExtendByReplication(pixsd, 1, 1);
2414 
2415         /* Map the bg to 200 */
2416     pixbgi = pixGetInvBackgroundMap(pixbg, 200, smoothx, smoothy);
2417     pixd = pixApplyInvBackgroundGrayMap(pixs, pixbgi, sx, sy);
2418 
2419     pixDestroy(&pixt);
2420     pixDestroy(&pixsd);
2421     pixDestroy(&pixbg);
2422     pixDestroy(&pixbgi);
2423     return pixd;
2424 }
2425 
2426 
2427 /*------------------------------------------------------------------*
2428  *                    Adaptive contrast normalization               *
2429  *------------------------------------------------------------------*/
2430 /*!
2431  *  pixContrastNorm()
2432  *
2433  *      Input:  pixd (<optional> 8 bpp; null or equal to pixs)
2434  *              pixs (8 bpp, not colormapped)
2435  *              sx, sy (tile dimensions)
2436  *              mindiff (minimum difference to accept as valid)
2437  *              smoothx, smoothy (half-width of convolution kernel applied to
2438  *                                min and max arrays: use 0 for no smoothing)
2439  *      Return: pixd always
2440  *
2441  *  Notes:
2442  *      (1) This function adaptively attempts to expand the contrast
2443  *          to the full dynamic range in each tile.  If the contrast in
2444  *          a tile is smaller than @mindiff, it uses the min and max
2445  *          pixel values from neighboring tiles.  It also can use
2446  *          convolution to smooth the min and max values from
2447  *          neighboring tiles.  After all that processing, it is
2448  *          possible that the actual pixel values in the tile are outside
2449  *          the computed [min ... max] range for local contrast
2450  *          normalization.  Such pixels are taken to be at either 0
2451  *          (if below the min) or 255 (if above the max).
2452  *      (2) pixd can be equal to pixs (in-place operation) or
2453  *          null (makes a new pixd).
2454  *      (3) sx and sy give the tile size; they are typically at least 20.
2455  *      (4) mindiff is used to eliminate results for tiles where it is
2456  *          likely that either fg or bg is missing.  A value around 50
2457  *          or more is reasonable.
2458  *      (5) The full width and height of the convolution kernel
2459  *          are (2 * smoothx + 1) and (2 * smoothy + 1).  Some smoothing
2460  *          is typically useful, and we limit the smoothing half-widths
2461  *          to the range from 0 to 8.
2462  *      (6) A linear TRC (gamma = 1.0) is applied to increase the contrast
2463  *          in each tile.  The result can subsequently be globally corrected,
2464  *          by applying pixGammaTRC() with arbitrary values of gamma
2465  *          and the 0 and 255 points of the mapping.
2466  */
2467 PIX *
pixContrastNorm(PIX * pixd,PIX * pixs,l_int32 sx,l_int32 sy,l_int32 mindiff,l_int32 smoothx,l_int32 smoothy)2468 pixContrastNorm(PIX       *pixd,
2469                 PIX       *pixs,
2470                 l_int32    sx,
2471                 l_int32    sy,
2472                 l_int32    mindiff,
2473                 l_int32    smoothx,
2474                 l_int32    smoothy)
2475 {
2476 PIX  *pixmin, *pixmax;
2477 
2478     PROCNAME("pixContrastNorm");
2479 
2480     if (!pixs || pixGetDepth(pixs) != 8)
2481         return (PIX *)ERROR_PTR("pixs undefined or not 8 bpp", procName, pixd);
2482     if (pixd && pixd != pixs)
2483         return (PIX *)ERROR_PTR("pixd not null or == pixs", procName, pixd);
2484     if (pixGetColormap(pixs))
2485         return (PIX *)ERROR_PTR("pixs is colormapped", procName, pixd);
2486     if (sx < 5 || sy < 5)
2487         return (PIX *)ERROR_PTR("sx and/or sy less than 5", procName, pixd);
2488     if (smoothx < 0 || smoothy < 0)
2489         return (PIX *)ERROR_PTR("smooth params less than 0", procName, pixd);
2490     if (smoothx > 8 || smoothy > 8)
2491         return (PIX *)ERROR_PTR("smooth params exceed 8", procName, pixd);
2492 
2493         /* Get the min and max pixel values in each tile, and represent
2494          * each value as a pixel in pixmin and pixmax, respectively. */
2495     pixMinMaxTiles(pixs, sx, sy, mindiff, smoothx, smoothy, &pixmin, &pixmax);
2496 
2497         /* For each tile, do a linear expansion of the dynamic range
2498          * of pixels so that the min value is mapped to 0 and the
2499          * max value is mapped to 255.  */
2500     pixd = pixLinearTRCTiled(pixd, pixs, sx, sy, pixmin, pixmax);
2501 
2502     pixDestroy(&pixmin);
2503     pixDestroy(&pixmax);
2504     return pixd;
2505 }
2506 
2507 
2508 /*!
2509  *  pixMinMaxTiles()
2510  *
2511  *      Input:  pixs (8 bpp, not colormapped)
2512  *              sx, sy (tile dimensions)
2513  *              mindiff (minimum difference to accept as valid)
2514  *              smoothx, smoothy (half-width of convolution kernel applied to
2515  *                                min and max arrays: use 0 for no smoothing)
2516  *              &pixmin (<return> tiled minima)
2517  *              &pixmax (<return> tiled maxima)
2518  *      Return: 0 if OK, 1 on error
2519  *
2520  *  Notes:
2521  *      (1) This computes filtered and smoothed values for the min and
2522  *          max pixel values in each tile of the image.
2523  *      (2) See pixContrastNorm() for usage.
2524  */
2525 l_int32
pixMinMaxTiles(PIX * pixs,l_int32 sx,l_int32 sy,l_int32 mindiff,l_int32 smoothx,l_int32 smoothy,PIX ** ppixmin,PIX ** ppixmax)2526 pixMinMaxTiles(PIX     *pixs,
2527                l_int32  sx,
2528                l_int32  sy,
2529                l_int32  mindiff,
2530                l_int32  smoothx,
2531                l_int32  smoothy,
2532                PIX    **ppixmin,
2533                PIX    **ppixmax)
2534 {
2535 l_int32  w, h;
2536 PIX     *pixmin1, *pixmax1, *pixmin2, *pixmax2;
2537 
2538     PROCNAME("pixMinMaxTiles");
2539 
2540     if (!ppixmin || !ppixmax)
2541         return ERROR_INT("&pixmin or &pixmax undefined", procName, 1);
2542     *ppixmin = *ppixmax = NULL;
2543     if (!pixs || pixGetDepth(pixs) != 8)
2544         return ERROR_INT("pixs undefined or not 8 bpp", procName, 1);
2545     if (pixGetColormap(pixs))
2546         return ERROR_INT("pixs is colormapped", procName, 1);
2547     if (sx < 5 || sy < 5)
2548         return ERROR_INT("sx and/or sy less than 3", procName, 1);
2549     if (smoothx < 0 || smoothy < 0)
2550         return ERROR_INT("smooth params less than 0", procName, 1);
2551     if (smoothx > 5 || smoothy > 5)
2552         return ERROR_INT("smooth params exceed 5", procName, 1);
2553 
2554         /* Get the min and max values in each tile */
2555     pixmin1 = pixScaleGrayMinMax(pixs, sx, sy, L_CHOOSE_MIN);
2556     pixmax1 = pixScaleGrayMinMax(pixs, sx, sy, L_CHOOSE_MAX);
2557 
2558     pixmin2 = pixExtendByReplication(pixmin1, 1, 1);
2559     pixmax2 = pixExtendByReplication(pixmax1, 1, 1);
2560     pixDestroy(&pixmin1);
2561     pixDestroy(&pixmax1);
2562 
2563         /* Make sure no value is 0 */
2564     pixAddConstantGray(pixmin2, 1);
2565     pixAddConstantGray(pixmax2, 1);
2566 
2567         /* Generate holes where the contrast is too small */
2568     pixSetLowContrast(pixmin2, pixmax2, mindiff);
2569 
2570         /* Fill the holes (0 values) */
2571     pixGetDimensions(pixmin2, &w, &h, NULL);
2572     pixFillMapHoles(pixmin2, w, h, L_FILL_BLACK);
2573     pixFillMapHoles(pixmax2, w, h, L_FILL_BLACK);
2574 
2575         /* Smooth if requested */
2576     if (smoothx > 0 || smoothy > 0) {
2577         smoothx = L_MIN(smoothx, (w - 1) / 2);
2578         smoothy = L_MIN(smoothy, (h - 1) / 2);
2579         *ppixmin = pixBlockconv(pixmin2, smoothx, smoothy);
2580         *ppixmax = pixBlockconv(pixmax2, smoothx, smoothy);
2581     }
2582     else {
2583         *ppixmin = pixClone(pixmin2);
2584         *ppixmax = pixClone(pixmax2);
2585     }
2586     pixDestroy(&pixmin2);
2587     pixDestroy(&pixmax2);
2588 
2589     return 0;
2590 }
2591 
2592 
2593 /*!
2594  *  pixSetLowContrast()
2595  *
2596  *      Input:  pixs1 (8 bpp)
2597  *              pixs2 (8 bpp)
2598  *              mindiff (minimum difference to accept as valid)
2599  *      Return: 0 if OK; 1 if no pixel diffs are large enough, or on error
2600  *
2601  *  Notes:
2602  *      (1) This compares corresponding pixels in pixs1 and pixs2.
2603  *          When they differ by less than @mindiff, set the pixel
2604  *          values to 0 in each.  Each pixel typically represents a tile
2605  *          in a larger image, and a very small difference between
2606  *          the min and max in the tile indicates that the min and max
2607  *          values are not to be trusted.
2608  *      (2) If contrast (pixel difference) detection is expected to fail,
2609  *          caller should check return value.
2610  */
2611 l_int32
pixSetLowContrast(PIX * pixs1,PIX * pixs2,l_int32 mindiff)2612 pixSetLowContrast(PIX     *pixs1,
2613                   PIX     *pixs2,
2614                   l_int32  mindiff)
2615 {
2616 l_int32    i, j, w, h, d, wpl, val1, val2, found;
2617 l_uint32  *data1, *data2, *line1, *line2;
2618 
2619     PROCNAME("pixSetLowContrast");
2620 
2621     if (!pixs1 || !pixs2)
2622         return ERROR_INT("pixs1 and pixs2 not both defined", procName, 1);
2623     if (pixSizesEqual(pixs1, pixs2) == 0)
2624         return ERROR_INT("pixs1 and pixs2 not equal size", procName, 1);
2625     pixGetDimensions(pixs1, &w, &h, &d);
2626     if (d != 8)
2627         return ERROR_INT("depth not 8 bpp", procName, 1);
2628     if (mindiff > 254) return 0;
2629 
2630     data1 = pixGetData(pixs1);
2631     data2 = pixGetData(pixs2);
2632     wpl = pixGetWpl(pixs1);
2633     found = 0;  /* init to not finding any diffs >= mindiff */
2634     for (i = 0; i < h; i++) {
2635         line1 = data1 + i * wpl;
2636         line2 = data2 + i * wpl;
2637         for (j = 0; j < w; j++) {
2638             val1 = GET_DATA_BYTE(line1, j);
2639             val2 = GET_DATA_BYTE(line2, j);
2640             if (L_ABS(val1 - val2) >= mindiff) {
2641                 found = 1;
2642                 break;
2643             }
2644         }
2645         if (found) break;
2646     }
2647     if (!found) {
2648         L_WARNING("no pixel pair diffs as large as mindiff", procName);
2649         pixClearAll(pixs1);
2650         pixClearAll(pixs2);
2651         return 1;
2652     }
2653 
2654     for (i = 0; i < h; i++) {
2655         line1 = data1 + i * wpl;
2656         line2 = data2 + i * wpl;
2657         for (j = 0; j < w; j++) {
2658             val1 = GET_DATA_BYTE(line1, j);
2659             val2 = GET_DATA_BYTE(line2, j);
2660             if (L_ABS(val1 - val2) < mindiff) {
2661                 SET_DATA_BYTE(line1, j, 0);
2662                 SET_DATA_BYTE(line2, j, 0);
2663             }
2664         }
2665     }
2666 
2667     return 0;
2668 }
2669 
2670 
2671 /*!
2672  *  pixLinearTRCTiled()
2673  *
2674  *      Input:  pixd (<optional> 8 bpp)
2675  *              pixs (8 bpp, not colormapped)
2676  *              sx, sy (tile dimensions)
2677  *              pixmin (pix of min values in tiles)
2678  *              pixmax (pix of max values in tiles)
2679  *      Return: pixd always
2680  *
2681  *  Notes:
2682  *      (1) pixd can be equal to pixs (in-place operation) or
2683  *          null (makes a new pixd).
2684  *      (2) sx and sy give the tile size; they are typically at least 20.
2685  *      (3) pixmin and pixmax are generated by pixMinMaxTiles()
2686  *      (4) For each tile, this does a linear expansion of the dynamic
2687  *          range so that the min value in the tile becomes 0 and the
2688  *          max value in the tile becomes 255.
2689  *      (5) The LUTs that do the mapping are generated as needed
2690  *          and stored for reuse in an integer array within the ptr array iaa[].
2691  */
2692 PIX *
pixLinearTRCTiled(PIX * pixd,PIX * pixs,l_int32 sx,l_int32 sy,PIX * pixmin,PIX * pixmax)2693 pixLinearTRCTiled(PIX       *pixd,
2694                   PIX       *pixs,
2695                   l_int32    sx,
2696                   l_int32    sy,
2697                   PIX       *pixmin,
2698                   PIX       *pixmax)
2699 {
2700 l_int32    i, j, k, m, w, h, wt, ht, wpl, wplt, xoff, yoff;
2701 l_int32    minval, maxval, val, sval;
2702 l_int32   *ia;
2703 l_int32  **iaa;
2704 l_uint32  *data, *datamin, *datamax, *line, *tline, *linemin, *linemax;
2705 
2706     PROCNAME("pixLinearTRCTiled");
2707 
2708     if (!pixs || pixGetDepth(pixs) != 8)
2709         return (PIX *)ERROR_PTR("pixs undefined or not 8 bpp", procName, pixd);
2710     if (pixd && pixd != pixs)
2711         return (PIX *)ERROR_PTR("pixd not null or == pixs", procName, pixd);
2712     if (pixGetColormap(pixs))
2713         return (PIX *)ERROR_PTR("pixs is colormapped", procName, pixd);
2714     if (!pixmin || !pixmax)
2715         return (PIX *)ERROR_PTR("pixmin & pixmax not defined", procName, pixd);
2716     if (sx < 5 || sy < 5)
2717         return (PIX *)ERROR_PTR("sx and/or sy less than 5", procName, pixd);
2718 
2719     pixd = pixCopy(pixd, pixs);
2720     iaa = (l_int32 **)CALLOC(256, sizeof(l_int32 *));
2721     pixGetDimensions(pixd, &w, &h, NULL);
2722 
2723     data = pixGetData(pixd);
2724     wpl = pixGetWpl(pixd);
2725     datamin = pixGetData(pixmin);
2726     datamax = pixGetData(pixmax);
2727     wplt = pixGetWpl(pixmin);
2728     pixGetDimensions(pixmin, &wt, &ht, NULL);
2729     for (i = 0; i < ht; i++) {
2730         line = data + sy * i * wpl;
2731         linemin = datamin + i * wplt;
2732         linemax = datamax + i * wplt;
2733         yoff = sy * i;
2734         for (j = 0; j < wt; j++) {
2735             xoff = sx * j;
2736             minval = GET_DATA_BYTE(linemin, j);
2737             maxval = GET_DATA_BYTE(linemax, j);
2738             if (maxval == minval) {  /* this is bad */
2739 /*                fprintf(stderr, "should't happen! i,j = %d,%d, minval = %d\n",
2740                         i, j, minval); */
2741                 continue;
2742             }
2743             ia = iaaGetLinearTRC(iaa, maxval - minval);
2744             for (k = 0; k < sy && yoff + k < h; k++) {
2745                 tline = line + k * wpl;
2746                 for (m = 0; m < sx && xoff + m < w; m++) {
2747                     val = GET_DATA_BYTE(tline, xoff + m);
2748                     sval = val - minval;
2749                     sval = L_MAX(0, sval);
2750                     SET_DATA_BYTE(tline, xoff + m, ia[sval]);
2751                 }
2752             }
2753         }
2754     }
2755 
2756     for (i = 0; i < 256; i++)
2757         if (iaa[i]) FREE(iaa[i]);
2758     FREE(iaa);
2759     return pixd;
2760 }
2761 
2762 
2763 /*!
2764  *  iaaGetLinearTRC()
2765  *
2766  *      Input:  iaa (bare array of ptrs to l_int32)
2767  *              diff (between min and max pixel values that are
2768  *                    to be mapped to 0 and 255)
2769  *      Return: ia (LUT with input (val - minval) and output a
2770  *                  value between 0 and 255)
2771  */
2772 static l_int32 *
iaaGetLinearTRC(l_int32 ** iaa,l_int32 diff)2773 iaaGetLinearTRC(l_int32  **iaa,
2774                 l_int32    diff)
2775 {
2776 l_int32    i;
2777 l_int32   *ia;
2778 l_float32  factor;
2779 
2780     PROCNAME("iaaGetLinearTRC");
2781 
2782     if (!iaa)
2783         return (l_int32 *)ERROR_PTR("iaa not defined", procName, NULL);
2784 
2785     if (iaa[diff] != NULL)  /* already have it */
2786        return iaa[diff];
2787 
2788     if ((ia = (l_int32 *)CALLOC(256, sizeof(l_int32))) == NULL)
2789         return (l_int32 *)ERROR_PTR("ia not made", procName, NULL);
2790     iaa[diff] = ia;
2791     if (diff == 0) {  /* shouldn't happen */
2792         for (i = 0; i < 256; i++)
2793             ia[i] = 128;
2794     }
2795     else {
2796         factor = 255. / (l_float32)diff;
2797         for (i = 0; i < diff + 1; i++)
2798             ia[i] = (l_int32)(factor * i + 0.5);
2799         for (i = diff + 1; i < 256; i++)
2800             ia[i] = 255;
2801     }
2802 
2803     return ia;
2804 }
2805 
2806 
2807