diff options
Diffstat (limited to 'leptonica/src/adaptmap.c')
-rw-r--r-- | leptonica/src/adaptmap.c | 2948 |
1 files changed, 2948 insertions, 0 deletions
diff --git a/leptonica/src/adaptmap.c b/leptonica/src/adaptmap.c new file mode 100644 index 00000000..c5c07162 --- /dev/null +++ b/leptonica/src/adaptmap.c @@ -0,0 +1,2948 @@ +/*====================================================================* + - Copyright (C) 2001 Leptonica. All rights reserved. + - + - Redistribution and use in source and binary forms, with or without + - modification, are permitted provided that the following conditions + - are met: + - 1. Redistributions of source code must retain the above copyright + - notice, this list of conditions and the following disclaimer. + - 2. Redistributions in binary form must reproduce the above + - copyright notice, this list of conditions and the following + - disclaimer in the documentation and/or other materials + - provided with the distribution. + - + - THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + - ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + - LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + - A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL ANY + - CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, + - EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, + - PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR + - PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY + - OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING + - NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS + - SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + *====================================================================*/ + +/*! + * \file adaptmap.c + * <pre> + * + * ------------------------------------------------------------------- + * + * Image binarization algorithms are found in: + * grayquant.c: standard, simple, general grayscale quantization + * adaptmap.c: local adaptive; mostly gray-to-gray in preparation + * for binarization + * binarize.c: special binarization methods, locally adaptive. + * + * ------------------------------------------------------------------- + * + * Clean background to white using background normalization + * PIX *pixCleanBackgroundToWhite() + * + * Adaptive background normalization (top-level functions) + * PIX *pixBackgroundNormSimple() 8 and 32 bpp + * PIX *pixBackgroundNorm() 8 and 32 bpp + * PIX *pixBackgroundNormMorph() 8 and 32 bpp + * + * Arrays of inverted background values for normalization (16 bpp) + * l_int32 pixBackgroundNormGrayArray() 8 bpp input + * l_int32 pixBackgroundNormRGBArrays() 32 bpp input + * l_int32 pixBackgroundNormGrayArrayMorph() 8 bpp input + * l_int32 pixBackgroundNormRGBArraysMorph() 32 bpp input + * + * Measurement of local background + * l_int32 pixGetBackgroundGrayMap() 8 bpp + * l_int32 pixGetBackgroundRGBMap() 32 bpp + * l_int32 pixGetBackgroundGrayMapMorph() 8 bpp + * l_int32 pixGetBackgroundRGBMapMorph() 32 bpp + * l_int32 pixFillMapHoles() + * PIX *pixExtendByReplication() 8 bpp + * l_int32 pixSmoothConnectedRegions() 8 bpp + * + * Measurement of local foreground + * l_int32 pixGetForegroundGrayMap() 8 bpp + * + * Generate inverted background map for each component + * PIX *pixGetInvBackgroundMap() 16 bpp + * + * Apply inverse background map to image + * PIX *pixApplyInvBackgroundGrayMap() 8 bpp + * PIX *pixApplyInvBackgroundRGBMap() 32 bpp + * + * Apply variable map + * PIX *pixApplyVariableGrayMap() 8 bpp + * + * Non-adaptive (global) mapping + * PIX *pixGlobalNormRGB() 32 bpp or cmapped + * PIX *pixGlobalNormNoSatRGB() 32 bpp + * + * Adaptive threshold spread normalization + * l_int32 pixThresholdSpreadNorm() 8 bpp + * + * Adaptive background normalization (flexible adaptaption) + * PIX *pixBackgroundNormFlex() 8 bpp + * + * Adaptive contrast normalization + * PIX *pixContrastNorm() 8 bpp + * static l_int32 pixMinMaxTiles() + * static l_int32 pixSetLowContrast() + * static PIX *pixLinearTRCTiled() + * static l_int32 *iaaGetLinearTRC() + * + * Background normalization is done by generating a reduced map (or set + * of maps) representing the estimated background value of the + * input image, and using this to shift the pixel values so that + * this background value is set to some constant value. + * + * Specifically, normalization has 3 steps: + * (1) Generate a background map at a reduced scale. + * (2) Make the array of inverted background values by inverting + * the map. The result is an array of local multiplicative factors. + * (3) Apply this inverse background map to the image + * + * The inverse background arrays can be generated in two different ways here: + * (1) Remove the 'foreground' pixels and average over the remaining + * pixels in each tile. Propagate values into tiles where + * values have not been assigned, either because there was not + * enough background in the tile or because the tile is covered + * by a foreground region described by an image mask. + * After the background map is made, the inverse map is generated by + * smoothing over some number of adjacent tiles + * (block convolution) and then inverting. + * (2) Remove the foreground pixels using a morphological closing + * on a subsampled version of the image. Propagate values + * into pixels covered by an optional image mask. Invert the + * background map without preconditioning by convolutional smoothing. + * + * Other methods for adaptively normalizing the image are also given here. + * + * (1) pixThresholdSpreadNorm() computes a local threshold over the image + * and normalizes the input pixel values so that this computed threshold + * is a constant across the entire image. + * + * (2) pixContrastNorm() computes and applies a local TRC so that the + * local dynamic range is expanded to the full 8 bits, where the + * darkest pixels are mapped to 0 and the lightest to 255. This is + * useful for improving the appearance of pages with very light + * foreground or very dark background, and where the local TRC + * function doesn't change rapidly with position. + * </pre> + */ + +#ifdef HAVE_CONFIG_H +#include <config_auto.h> +#endif /* HAVE_CONFIG_H */ + +#include "allheaders.h" + + /* Default input parameters for pixBackgroundNormSimple() + * Notes: + * (1) mincount must never exceed the tile area (width * height) + * (2) bgval must be sufficiently below 255 to avoid accidental + * saturation; otherwise it should be large to avoid + * shrinking the dynamic range + * (3) results should otherwise not be sensitive to these values + */ +static const l_int32 DefaultTileWidth = 10; /*!< default tile width */ +static const l_int32 DefaultTileHeight = 15; /*!< default tile height */ +static const l_int32 DefaultFgThreshold = 60; /*!< default fg threshold */ +static const l_int32 DefaultMinCount = 40; /*!< default minimum count */ +static const l_int32 DefaultBgVal = 200; /*!< default bg value */ +static const l_int32 DefaultXSmoothSize = 2; /*!< default x smooth size */ +static const l_int32 DefaultYSmoothSize = 1; /*!< default y smooth size */ + +static l_int32 pixMinMaxTiles(PIX *pixs, l_int32 sx, l_int32 sy, + l_int32 mindiff, l_int32 smoothx, l_int32 smoothy, + PIX **ppixmin, PIX **ppixmax); +static l_int32 pixSetLowContrast(PIX *pixs1, PIX *pixs2, l_int32 mindiff); +static PIX *pixLinearTRCTiled(PIX *pixd, PIX *pixs, l_int32 sx, l_int32 sy, + PIX *pixmin, PIX *pixmax); +static l_int32 *iaaGetLinearTRC(l_int32 **iaa, l_int32 diff); + +#ifndef NO_CONSOLE_IO +#define DEBUG_GLOBAL 0 /*!< set to 1 to debug pixGlobalNormNoSatRGB() */ +#endif /* ~NO_CONSOLE_IO */ + +/*------------------------------------------------------------------* + * Clean background to white using background normalization * + *------------------------------------------------------------------*/ +/*! + * \brief pixCleanBackgroundToWhite() + * + * \param[in] pixs 8 bpp grayscale or 32 bpp rgb + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null + * \param[in] pixg [optional] 8 bpp grayscale version; can be null + * \param[in] gamma gamma correction; must be > 0.0; typically ~1.0 + * \param[in] blackval dark value to set to black (0) + * \param[in] whiteval light value to set to white (255) + * \return pixd 8 bpp or 32 bpp rgb, or NULL on error + * + * <pre> + * Notes: + * (1) This is a simplified interface for cleaning an image. + * For comparison, see pixAdaptThresholdToBinaryGen(). + * (2) The suggested default values for the input parameters are: + * gamma: 1.0 (reduce this to increase the contrast; e.g., + * for light text) + * blackval 70 (a bit more than 60) + * whiteval 190 (a bit less than 200) + * </pre> + */ +PIX * +pixCleanBackgroundToWhite(PIX *pixs, + PIX *pixim, + PIX *pixg, + l_float32 gamma, + l_int32 blackval, + l_int32 whiteval) +{ +l_int32 d; +PIX *pixd; + + PROCNAME("pixCleanBackgroundToWhite"); + + if (!pixs) + return (PIX *)ERROR_PTR("pixs not defined", procName, NULL); + d = pixGetDepth(pixs); + if (d != 8 && d != 32) + return (PIX *)ERROR_PTR("depth not 8 or 32", procName, NULL); + + pixd = pixBackgroundNormSimple(pixs, pixim, pixg); + if (!pixd) + return (PIX *)ERROR_PTR("background norm failedd", procName, NULL); + pixGammaTRC(pixd, pixd, gamma, blackval, whiteval); + return pixd; +} + + +/*------------------------------------------------------------------* + * Adaptive background normalization * + *------------------------------------------------------------------*/ +/*! + * \brief pixBackgroundNormSimple() + * + * \param[in] pixs 8 bpp grayscale or 32 bpp rgb + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null + * \param[in] pixg [optional] 8 bpp grayscale version; can be null + * \return pixd 8 bpp or 32 bpp rgb, or NULL on error + * + * <pre> + * Notes: + * (1) This is a simplified interface to pixBackgroundNorm(), + * where seven parameters are defaulted. + * (2) The input image is either grayscale or rgb. + * (3) See pixBackgroundNorm() for usage and function. + * </pre> + */ +PIX * +pixBackgroundNormSimple(PIX *pixs, + PIX *pixim, + PIX *pixg) +{ + return pixBackgroundNorm(pixs, pixim, pixg, + DefaultTileWidth, DefaultTileHeight, + DefaultFgThreshold, DefaultMinCount, + DefaultBgVal, DefaultXSmoothSize, + DefaultYSmoothSize); +} + + +/*! + * \brief pixBackgroundNorm() + * + * \param[in] pixs 8 bpp grayscale or 32 bpp rgb + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null + * \param[in] pixg [optional] 8 bpp grayscale version; can be null + * \param[in] sx, sy tile size in pixels + * \param[in] thresh threshold for determining foreground + * \param[in] mincount min threshold on counts in a tile + * \param[in] bgval target bg val; typ. > 128 + * \param[in] smoothx half-width of block convolution kernel width + * \param[in] smoothy half-width of block convolution kernel height + * \return pixd 8 bpp or 32 bpp rgb, or NULL on error + * + * <pre> + * Notes: + * (1) This is a top-level interface for normalizing the image intensity + * by mapping the image so that the background is near the input + * value 'bgval'. + * (2) The input image is either grayscale or rgb. + * (3) For each component in the input image, the background value + * in each tile is estimated using the values in the tile that + * are not part of the foreground, where the foreground is + * determined by the input 'thresh' argument. + * (4) An optional binary mask can be specified, with the foreground + * pixels typically over image regions. The resulting background + * map values will be determined by surrounding pixels that are + * not under the mask foreground. The origin (0,0) of this mask + * is assumed to be aligned with the origin of the input image. + * This binary mask must not fully cover pixs, because then there + * will be no pixels in the input image available to compute + * the background. + * (5) An optional grayscale version of the input pixs can be supplied. + * The only reason to do this is if the input is RGB and this + * grayscale version can be used elsewhere. If the input is RGB + * and this is not supplied, it is made internally using only + * the green component, and destroyed after use. + * (6) The dimensions of the pixel tile (sx, sy) give the amount by + * by which the map is reduced in size from the input image. + * (7) The threshold is used to binarize the input image, in order to + * locate the foreground components. If this is set too low, + * some actual foreground may be used to determine the maps; + * if set too high, there may not be enough background + * to determine the map values accurately. Typically, it's + * better to err by setting the threshold too high. + * (8) A 'mincount' threshold is a minimum count of pixels in a + * tile for which a background reading is made, in order for that + * pixel in the map to be valid. This number should perhaps be + * at least 1/3 the size of the tile. + * (9) A 'bgval' target background value for the normalized image. This + * should be at least 128. If set too close to 255, some + * clipping will occur in the result. + * (10) Two factors, 'smoothx' and 'smoothy', are input for smoothing + * the map. Each low-pass filter kernel dimension is + * is 2 * (smoothing factor) + 1, so a + * value of 0 means no smoothing. A value of 1 or 2 is recommended. + * </pre> + */ +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) +{ +l_int32 d, allfg; +PIX *pixm, *pixmi, *pixd; +PIX *pixmr, *pixmg, *pixmb, *pixmri, *pixmgi, *pixmbi; + + PROCNAME("pixBackgroundNorm"); + + if (!pixs) + return (PIX *)ERROR_PTR("pixs not defined", procName, NULL); + d = pixGetDepth(pixs); + if (d != 8 && d != 32) + return (PIX *)ERROR_PTR("pixs not 8 or 32 bpp", procName, NULL); + if (sx < 4 || sy < 4) + return (PIX *)ERROR_PTR("sx and sy must be >= 4", procName, NULL); + if (mincount > sx * sy) { + L_WARNING("mincount too large for tile size\n", procName); + mincount = (sx * sy) / 3; + } + + /* If pixim exists, verify that it is not all foreground. */ + if (pixim) { + pixInvert(pixim, pixim); + pixZero(pixim, &allfg); + pixInvert(pixim, pixim); + if (allfg) + return (PIX *)ERROR_PTR("pixim all foreground", procName, NULL); + } + + pixd = NULL; + if (d == 8) { + pixm = NULL; + pixGetBackgroundGrayMap(pixs, pixim, sx, sy, thresh, mincount, &pixm); + if (!pixm) { + L_WARNING("map not made; return a copy of the source\n", procName); + return pixCopy(NULL, pixs); + } + + pixmi = pixGetInvBackgroundMap(pixm, bgval, smoothx, smoothy); + if (!pixmi) { + L_WARNING("pixmi not made; return a copy of source\n", procName); + pixDestroy(&pixm); + return pixCopy(NULL, pixs); + } else { + pixd = pixApplyInvBackgroundGrayMap(pixs, pixmi, sx, sy); + } + + pixDestroy(&pixm); + pixDestroy(&pixmi); + } + else { + pixmr = pixmg = pixmb = NULL; + pixGetBackgroundRGBMap(pixs, pixim, pixg, sx, sy, thresh, + mincount, &pixmr, &pixmg, &pixmb); + if (!pixmr || !pixmg || !pixmb) { + pixDestroy(&pixmr); + pixDestroy(&pixmg); + pixDestroy(&pixmb); + L_WARNING("map not made; return a copy of the source\n", procName); + return pixCopy(NULL, pixs); + } + + pixmri = pixGetInvBackgroundMap(pixmr, bgval, smoothx, smoothy); + pixmgi = pixGetInvBackgroundMap(pixmg, bgval, smoothx, smoothy); + pixmbi = pixGetInvBackgroundMap(pixmb, bgval, smoothx, smoothy); + if (!pixmri || !pixmgi || !pixmbi) { + L_WARNING("not all pixm*i are made; return src copy\n", procName); + pixd = pixCopy(NULL, pixs); + } else { + pixd = pixApplyInvBackgroundRGBMap(pixs, pixmri, pixmgi, pixmbi, + sx, sy); + } + + pixDestroy(&pixmr); + pixDestroy(&pixmg); + pixDestroy(&pixmb); + pixDestroy(&pixmri); + pixDestroy(&pixmgi); + pixDestroy(&pixmbi); + } + + if (!pixd) + ERROR_PTR("pixd not made", procName, NULL); + pixCopyResolution(pixd, pixs); + return pixd; +} + + +/*! + * \brief pixBackgroundNormMorph() + * + * \param[in] pixs 8 bpp grayscale or 32 bpp rgb + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null + * \param[in] reduction at which morph closings are done; between 2 and 16 + * \param[in] size of square Sel for the closing; use an odd number + * \param[in] bgval target bg val; typ. > 128 + * \return pixd 8 bpp, or NULL on error + * + * <pre> + * Notes: + * (1) This is a top-level interface for normalizing the image intensity + * by mapping the image so that the background is near the input + * value 'bgval'. + * (2) The input image is either grayscale or rgb. + * (3) For each component in the input image, the background value + * is estimated using a grayscale closing; hence the 'Morph' + * in the function name. + * (4) An optional binary mask can be specified, with the foreground + * pixels typically over image regions. The resulting background + * map values will be determined by surrounding pixels that are + * not under the mask foreground. The origin (0,0) of this mask + * is assumed to be aligned with the origin of the input image. + * This binary mask must not fully cover pixs, because then there + * will be no pixels in the input image available to compute + * the background. + * (5) The map is computed at reduced size (given by 'reduction') + * from the input pixs and optional pixim. At this scale, + * pixs is closed to remove the background, using a square Sel + * of odd dimension. The product of reduction * size should be + * large enough to remove most of the text foreground. + * (6) No convolutional smoothing needs to be done on the map before + * inverting it. + * (7) A 'bgval' target background value for the normalized image. This + * should be at least 128. If set too close to 255, some + * clipping will occur in the result. + * </pre> + */ +PIX * +pixBackgroundNormMorph(PIX *pixs, + PIX *pixim, + l_int32 reduction, + l_int32 size, + l_int32 bgval) +{ +l_int32 d, allfg; +PIX *pixm, *pixmi, *pixd; +PIX *pixmr, *pixmg, *pixmb, *pixmri, *pixmgi, *pixmbi; + + PROCNAME("pixBackgroundNormMorph"); + + if (!pixs) + return (PIX *)ERROR_PTR("pixs not defined", procName, NULL); + d = pixGetDepth(pixs); + if (d != 8 && d != 32) + return (PIX *)ERROR_PTR("pixs not 8 or 32 bpp", procName, NULL); + if (reduction < 2 || reduction > 16) + return (PIX *)ERROR_PTR("reduction must be between 2 and 16", + procName, NULL); + + /* If pixim exists, verify that it is not all foreground. */ + if (pixim) { + pixInvert(pixim, pixim); + pixZero(pixim, &allfg); + pixInvert(pixim, pixim); + if (allfg) + return (PIX *)ERROR_PTR("pixim all foreground", procName, NULL); + } + + pixd = NULL; + if (d == 8) { + pixGetBackgroundGrayMapMorph(pixs, pixim, reduction, size, &pixm); + if (!pixm) + return (PIX *)ERROR_PTR("pixm not made", procName, NULL); + pixmi = pixGetInvBackgroundMap(pixm, bgval, 0, 0); + if (!pixmi) + ERROR_PTR("pixmi not made", procName, NULL); + else + pixd = pixApplyInvBackgroundGrayMap(pixs, pixmi, + reduction, reduction); + pixDestroy(&pixm); + pixDestroy(&pixmi); + } + else { /* d == 32 */ + pixmr = pixmg = pixmb = NULL; + pixGetBackgroundRGBMapMorph(pixs, pixim, reduction, size, + &pixmr, &pixmg, &pixmb); + if (!pixmr || !pixmg || !pixmb) { + pixDestroy(&pixmr); + pixDestroy(&pixmg); + pixDestroy(&pixmb); + return (PIX *)ERROR_PTR("not all pixm*", procName, NULL); + } + + pixmri = pixGetInvBackgroundMap(pixmr, bgval, 0, 0); + pixmgi = pixGetInvBackgroundMap(pixmg, bgval, 0, 0); + pixmbi = pixGetInvBackgroundMap(pixmb, bgval, 0, 0); + if (!pixmri || !pixmgi || !pixmbi) + ERROR_PTR("not all pixm*i are made", procName, NULL); + else + pixd = pixApplyInvBackgroundRGBMap(pixs, pixmri, pixmgi, pixmbi, + reduction, reduction); + + pixDestroy(&pixmr); + pixDestroy(&pixmg); + pixDestroy(&pixmb); + pixDestroy(&pixmri); + pixDestroy(&pixmgi); + pixDestroy(&pixmbi); + } + + if (!pixd) + ERROR_PTR("pixd not made", procName, NULL); + pixCopyResolution(pixd, pixs); + return pixd; +} + + +/*-------------------------------------------------------------------------* + * Arrays of inverted background values for normalization * + *-------------------------------------------------------------------------* + * Notes for these four functions: * + * (1) They are useful if you need to save the actual mapping array. * + * (2) They could be used in the top-level functions but are * + * not because their use makes those functions less clear. * + * (3) Each component in the input pixs generates a 16 bpp pix array. * + *-------------------------------------------------------------------------*/ +/*! + * \brief pixBackgroundNormGrayArray() + * + * \param[in] pixs 8 bpp grayscale + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null + * \param[in] sx, sy tile size in pixels + * \param[in] thresh threshold for determining foreground + * \param[in] mincount min threshold on counts in a tile + * \param[in] bgval target bg val; typ. > 128 + * \param[in] smoothx half-width of block convolution kernel width + * \param[in] smoothy half-width of block convolution kernel height + * \param[out] ppixd 16 bpp array of inverted background value + * \return 0 if OK, 1 on error + * + * <pre> + * Notes: + * (1) See notes in pixBackgroundNorm(). + * (2) This returns a 16 bpp pix that can be used by + * pixApplyInvBackgroundGrayMap() to generate a normalized version + * of the input pixs. + * </pre> + */ +l_ok +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) +{ +l_int32 allfg; +PIX *pixm; + + PROCNAME("pixBackgroundNormGrayArray"); + + if (!ppixd) + return ERROR_INT("&pixd not defined", procName, 1); + *ppixd = NULL; + if (!pixs || pixGetDepth(pixs) != 8) + return ERROR_INT("pixs not defined or not 8 bpp", procName, 1); + if (pixGetColormap(pixs)) + return ERROR_INT("pixs is colormapped", procName, 1); + if (pixim && pixGetDepth(pixim) != 1) + return ERROR_INT("pixim not 1 bpp", procName, 1); + if (sx < 4 || sy < 4) + return ERROR_INT("sx and sy must be >= 4", procName, 1); + if (mincount > sx * sy) { + L_WARNING("mincount too large for tile size\n", procName); + mincount = (sx * sy) / 3; + } + + /* If pixim exists, verify that it is not all foreground. */ + if (pixim) { + pixInvert(pixim, pixim); + pixZero(pixim, &allfg); + pixInvert(pixim, pixim); + if (allfg) + return ERROR_INT("pixim all foreground", procName, 1); + } + + pixGetBackgroundGrayMap(pixs, pixim, sx, sy, thresh, mincount, &pixm); + if (!pixm) + return ERROR_INT("pixm not made", procName, 1); + *ppixd = pixGetInvBackgroundMap(pixm, bgval, smoothx, smoothy); + pixCopyResolution(*ppixd, pixs); + pixDestroy(&pixm); + return 0; +} + + +/*! + * \brief pixBackgroundNormRGBArrays() + * + * \param[in] pixs 32 bpp rgb + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null + * \param[in] pixg [optional] 8 bpp grayscale version; can be null + * \param[in] sx, sy tile size in pixels + * \param[in] thresh threshold for determining foreground + * \param[in] mincount min threshold on counts in a tile + * \param[in] bgval target bg val; typ. > 128 + * \param[in] smoothx half-width of block convolution kernel width + * \param[in] smoothy half-width of block convolution kernel height + * \param[out] ppixr 16 bpp array of inverted R background value + * \param[out] ppixg 16 bpp array of inverted G background value + * \param[out] ppixb 16 bpp array of inverted B background value + * \return 0 if OK, 1 on error + * + * <pre> + * Notes: + * (1) See notes in pixBackgroundNorm(). + * (2) This returns a set of three 16 bpp pix that can be used by + * pixApplyInvBackgroundGrayMap() to generate a normalized version + * of each component of the input pixs. + * </pre> + */ +l_ok +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) +{ +l_int32 allfg; +PIX *pixmr, *pixmg, *pixmb; + + PROCNAME("pixBackgroundNormRGBArrays"); + + if (!ppixr || !ppixg || !ppixb) + return ERROR_INT("&pixr, &pixg, &pixb not all defined", procName, 1); + *ppixr = *ppixg = *ppixb = NULL; + if (!pixs) + return ERROR_INT("pixs not defined", procName, 1); + if (pixGetDepth(pixs) != 32) + return ERROR_INT("pixs not 32 bpp", procName, 1); + if (pixim && pixGetDepth(pixim) != 1) + return ERROR_INT("pixim not 1 bpp", procName, 1); + if (sx < 4 || sy < 4) + return ERROR_INT("sx and sy must be >= 4", procName, 1); + if (mincount > sx * sy) { + L_WARNING("mincount too large for tile size\n", procName); + mincount = (sx * sy) / 3; + } + + /* If pixim exists, verify that it is not all foreground. */ + if (pixim) { + pixInvert(pixim, pixim); + pixZero(pixim, &allfg); + pixInvert(pixim, pixim); + if (allfg) + return ERROR_INT("pixim all foreground", procName, 1); + } + + pixGetBackgroundRGBMap(pixs, pixim, pixg, sx, sy, thresh, mincount, + &pixmr, &pixmg, &pixmb); + if (!pixmr || !pixmg || !pixmb) { + pixDestroy(&pixmr); + pixDestroy(&pixmg); + pixDestroy(&pixmb); + return ERROR_INT("not all pixm* made", procName, 1); + } + + *ppixr = pixGetInvBackgroundMap(pixmr, bgval, smoothx, smoothy); + *ppixg = pixGetInvBackgroundMap(pixmg, bgval, smoothx, smoothy); + *ppixb = pixGetInvBackgroundMap(pixmb, bgval, smoothx, smoothy); + pixDestroy(&pixmr); + pixDestroy(&pixmg); + pixDestroy(&pixmb); + return 0; +} + + +/*! + * \brief pixBackgroundNormGrayArrayMorph() + * + * \param[in] pixs 8 bpp grayscale + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null + * \param[in] reduction at which morph closings are done; between 2 and 16 + * \param[in] size of square Sel for the closing; use an odd number + * \param[in] bgval target bg val; typ. > 128 + * \param[out] ppixd 16 bpp array of inverted background value + * \return 0 if OK, 1 on error + * + * <pre> + * Notes: + * (1) See notes in pixBackgroundNormMorph(). + * (2) This returns a 16 bpp pix that can be used by + * pixApplyInvBackgroundGrayMap() to generate a normalized version + * of the input pixs. + * </pre> + */ +l_ok +pixBackgroundNormGrayArrayMorph(PIX *pixs, + PIX *pixim, + l_int32 reduction, + l_int32 size, + l_int32 bgval, + PIX **ppixd) +{ +l_int32 allfg; +PIX *pixm; + + PROCNAME("pixBackgroundNormGrayArrayMorph"); + + if (!ppixd) + return ERROR_INT("&pixd not defined", procName, 1); + *ppixd = NULL; + if (!pixs) + return ERROR_INT("pixs not defined", procName, 1); + if (pixGetDepth(pixs) != 8) + return ERROR_INT("pixs not 8 bpp", procName, 1); + if (pixim && pixGetDepth(pixim) != 1) + return ERROR_INT("pixim not 1 bpp", procName, 1); + if (reduction < 2 || reduction > 16) + return ERROR_INT("reduction must be between 2 and 16", procName, 1); + + /* If pixim exists, verify that it is not all foreground. */ + if (pixim) { + pixInvert(pixim, pixim); + pixZero(pixim, &allfg); + pixInvert(pixim, pixim); + if (allfg) + return ERROR_INT("pixim all foreground", procName, 1); + } + + pixGetBackgroundGrayMapMorph(pixs, pixim, reduction, size, &pixm); + if (!pixm) + return ERROR_INT("pixm not made", procName, 1); + *ppixd = pixGetInvBackgroundMap(pixm, bgval, 0, 0); + pixCopyResolution(*ppixd, pixs); + pixDestroy(&pixm); + return 0; +} + + +/*! + * \brief pixBackgroundNormRGBArraysMorph() + * + * \param[in] pixs 32 bpp rgb + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null + * \param[in] reduction at which morph closings are done; between 2 and 16 + * \param[in] size of square Sel for the closing; use an odd number + * \param[in] bgval target bg val; typ. > 128 + * \param[out] ppixr 16 bpp array of inverted R background value + * \param[out] ppixg 16 bpp array of inverted G background value + * \param[out] ppixb 16 bpp array of inverted B background value + * \return 0 if OK, 1 on error + * + * <pre> + * Notes: + * (1) See notes in pixBackgroundNormMorph(). + * (2) This returns a set of three 16 bpp pix that can be used by + * pixApplyInvBackgroundGrayMap() to generate a normalized version + * of each component of the input pixs. + * </pre> + */ +l_ok +pixBackgroundNormRGBArraysMorph(PIX *pixs, + PIX *pixim, + l_int32 reduction, + l_int32 size, + l_int32 bgval, + PIX **ppixr, + PIX **ppixg, + PIX **ppixb) +{ +l_int32 allfg; +PIX *pixmr, *pixmg, *pixmb; + + PROCNAME("pixBackgroundNormRGBArraysMorph"); + + if (!ppixr || !ppixg || !ppixb) + return ERROR_INT("&pixr, &pixg, &pixb not all defined", procName, 1); + *ppixr = *ppixg = *ppixb = NULL; + if (!pixs) + return ERROR_INT("pixs not defined", procName, 1); + if (pixGetDepth(pixs) != 32) + return ERROR_INT("pixs not 32 bpp", procName, 1); + if (pixim && pixGetDepth(pixim) != 1) + return ERROR_INT("pixim not 1 bpp", procName, 1); + if (reduction < 2 || reduction > 16) + return ERROR_INT("reduction must be between 2 and 16", procName, 1); + + /* If pixim exists, verify that it is not all foreground. */ + if (pixim) { + pixInvert(pixim, pixim); + pixZero(pixim, &allfg); + pixInvert(pixim, pixim); + if (allfg) + return ERROR_INT("pixim all foreground", procName, 1); + } + + pixGetBackgroundRGBMapMorph(pixs, pixim, reduction, size, + &pixmr, &pixmg, &pixmb); + if (!pixmr || !pixmg || !pixmb) { + pixDestroy(&pixmr); + pixDestroy(&pixmg); + pixDestroy(&pixmb); + return ERROR_INT("not all pixm* made", procName, 1); + } + + *ppixr = pixGetInvBackgroundMap(pixmr, bgval, 0, 0); + *ppixg = pixGetInvBackgroundMap(pixmg, bgval, 0, 0); + *ppixb = pixGetInvBackgroundMap(pixmb, bgval, 0, 0); + pixDestroy(&pixmr); + pixDestroy(&pixmg); + pixDestroy(&pixmb); + return 0; +} + + +/*------------------------------------------------------------------* + * Measurement of local background * + *------------------------------------------------------------------*/ +/*! + * \brief pixGetBackgroundGrayMap() + * + * \param[in] pixs 8 bpp grayscale; not cmapped + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null; + * it should not have only foreground pixels + * \param[in] sx, sy tile size in pixels + * \param[in] thresh threshold for determining foreground + * \param[in] mincount min threshold on counts in a tile + * \param[out] ppixd 8 bpp grayscale map + * \return 0 if OK, 1 on error + * + * <pre> + * Notes: + * (1) The background is measured in regions that don't have + * images. It is then propagated into the image regions, + * and finally smoothed in each image region. + * </pre> + */ +l_ok +pixGetBackgroundGrayMap(PIX *pixs, + PIX *pixim, + l_int32 sx, + l_int32 sy, + l_int32 thresh, + l_int32 mincount, + PIX **ppixd) +{ +l_int32 w, h, wd, hd, wim, him, wpls, wplim, wpld, wplf; +l_int32 xim, yim, delx, nx, ny, i, j, k, m; +l_int32 count, sum, val8; +l_int32 empty, fgpixels; +l_uint32 *datas, *dataim, *datad, *dataf, *lines, *lineim, *lined, *linef; +l_float32 scalex, scaley; +PIX *pixd, *piximi, *pixb, *pixf, *pixims; + + PROCNAME("pixGetBackgroundGrayMap"); + + if (!ppixd) + return ERROR_INT("&pixd not defined", procName, 1); + *ppixd = NULL; + if (!pixs || pixGetDepth(pixs) != 8) + return ERROR_INT("pixs not defined or not 8 bpp", procName, 1); + if (pixGetColormap(pixs)) + return ERROR_INT("pixs is colormapped", procName, 1); + if (pixim && pixGetDepth(pixim) != 1) + return ERROR_INT("pixim not 1 bpp", procName, 1); + if (sx < 4 || sy < 4) + return ERROR_INT("sx and sy must be >= 4", procName, 1); + if (mincount > sx * sy) { + L_WARNING("mincount too large for tile size\n", procName); + mincount = (sx * sy) / 3; + } + + /* Evaluate the 'image' mask, pixim, and make sure + * it is not all fg. */ + fgpixels = 0; /* boolean for existence of fg pixels in the image mask. */ + if (pixim) { + piximi = pixInvert(NULL, pixim); /* set non-'image' pixels to 1 */ + pixZero(piximi, &empty); + pixDestroy(&piximi); + if (empty) + return ERROR_INT("pixim all fg; no background", procName, 1); + pixZero(pixim, &empty); + if (!empty) /* there are fg pixels in pixim */ + fgpixels = 1; + } + + /* Generate the foreground mask, pixf, which is at + * full resolution. These pixels will be ignored when + * computing the background values. */ + pixb = pixThresholdToBinary(pixs, thresh); + pixf = pixMorphSequence(pixb, "d7.1 + d1.7", 0); + pixDestroy(&pixb); + if (!pixf) + return ERROR_INT("pixf not made", procName, 1); + + + /* ------------- Set up the output map pixd --------------- */ + /* Generate pixd, which is reduced by the factors (sx, sy). */ + w = pixGetWidth(pixs); + h = pixGetHeight(pixs); + wd = (w + sx - 1) / sx; + hd = (h + sy - 1) / sy; + pixd = pixCreate(wd, hd, 8); + + /* Note: we only compute map values in tiles that are complete. + * In general, tiles at right and bottom edges will not be + * complete, and we must fill them in later. */ + nx = w / sx; + ny = h / sy; + wpls = pixGetWpl(pixs); + datas = pixGetData(pixs); + wpld = pixGetWpl(pixd); + datad = pixGetData(pixd); + wplf = pixGetWpl(pixf); + dataf = pixGetData(pixf); + for (i = 0; i < ny; i++) { + lines = datas + sy * i * wpls; + linef = dataf + sy * i * wplf; + lined = datad + i * wpld; + for (j = 0; j < nx; j++) { + delx = j * sx; + sum = 0; + count = 0; + for (k = 0; k < sy; k++) { + for (m = 0; m < sx; m++) { + if (GET_DATA_BIT(linef + k * wplf, delx + m) == 0) { + sum += GET_DATA_BYTE(lines + k * wpls, delx + m); + count++; + } + } + } + if (count >= mincount) { + val8 = sum / count; + SET_DATA_BYTE(lined, j, val8); + } + } + } + pixDestroy(&pixf); + + /* If there is an optional mask with fg pixels, erase the previous + * calculation for the corresponding map pixels, setting the + * map values to 0. Then, when all the map holes are filled, + * these erased pixels will be set by the surrounding map values. + * + * The calculation here is relatively efficient: for each pixel + * in pixd (which corresponds to a tile of mask pixels in pixim) + * we look only at the pixel in pixim that is at the center + * of the tile. If the mask pixel is ON, we reset the map + * pixel in pixd to 0, so that it can later be filled in. */ + pixims = NULL; + if (pixim && fgpixels) { + wim = pixGetWidth(pixim); + him = pixGetHeight(pixim); + dataim = pixGetData(pixim); + wplim = pixGetWpl(pixim); + for (i = 0; i < ny; i++) { + yim = i * sy + sy / 2; + if (yim >= him) + break; + lineim = dataim + yim * wplim; + for (j = 0; j < nx; j++) { + xim = j * sx + sx / 2; + if (xim >= wim) + break; + if (GET_DATA_BIT(lineim, xim)) + pixSetPixel(pixd, j, i, 0); + } + } + } + + /* Fill all the holes in the map. */ + if (pixFillMapHoles(pixd, nx, ny, L_FILL_BLACK)) { + pixDestroy(&pixd); + L_WARNING("can't make the map\n", procName); + return 1; + } + + /* Finally, for each connected region corresponding to the + * 'image' mask, reset all pixels to their average value. + * Each of these components represents an image (or part of one) + * in the input, and this smooths the background values + * in each of these regions. */ + if (pixim && fgpixels) { + scalex = 1. / (l_float32)sx; + scaley = 1. / (l_float32)sy; + pixims = pixScaleBySampling(pixim, scalex, scaley); + pixSmoothConnectedRegions(pixd, pixims, 2); + pixDestroy(&pixims); + } + + *ppixd = pixd; + pixCopyResolution(*ppixd, pixs); + return 0; +} + + +/*! + * \brief pixGetBackgroundRGBMap() + * + * \param[in] pixs 32 bpp rgb + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null; it + * should not have all foreground pixels + * \param[in] pixg [optional] 8 bpp grayscale version; can be null + * \param[in] sx, sy tile size in pixels + * \param[in] thresh threshold for determining foreground + * \param[in] mincount min threshold on counts in a tile + * \param[out] ppixmr red component map + * \param[out] ppixmg green component map + * \param[out] ppixmb blue component map + * \return 0 if OK, 1 on error + * + * <pre> + * Notes: + * (1) If pixg, which is a grayscale version of pixs, is provided, + * use this internally to generate the foreground mask. + * Otherwise, a grayscale version of pixs will be generated + * from the green component only, used, and destroyed. + * </pre> + */ +l_ok +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) +{ +l_int32 w, h, wm, hm, wim, him, wpls, wplim, wplf; +l_int32 xim, yim, delx, nx, ny, i, j, k, m; +l_int32 count, rsum, gsum, bsum, rval, gval, bval; +l_int32 empty, fgpixels; +l_uint32 pixel; +l_uint32 *datas, *dataim, *dataf, *lines, *lineim, *linef; +l_float32 scalex, scaley; +PIX *piximi, *pixgc, *pixb, *pixf, *pixims; +PIX *pixmr, *pixmg, *pixmb; + + PROCNAME("pixGetBackgroundRGBMap"); + + if (!ppixmr || !ppixmg || !ppixmb) + return ERROR_INT("&pixm* not all defined", procName, 1); + *ppixmr = *ppixmg = *ppixmb = NULL; + if (!pixs) + return ERROR_INT("pixs not defined", procName, 1); + if (pixGetDepth(pixs) != 32) + return ERROR_INT("pixs not 32 bpp", procName, 1); + if (pixim && pixGetDepth(pixim) != 1) + return ERROR_INT("pixim not 1 bpp", procName, 1); + if (sx < 4 || sy < 4) + return ERROR_INT("sx and sy must be >= 4", procName, 1); + if (mincount > sx * sy) { + L_WARNING("mincount too large for tile size\n", procName); + mincount = (sx * sy) / 3; + } + + /* Evaluate the mask pixim and make sure it is not all foreground */ + fgpixels = 0; /* boolean for existence of fg mask pixels */ + if (pixim) { + piximi = pixInvert(NULL, pixim); /* set non-'image' pixels to 1 */ + pixZero(piximi, &empty); + pixDestroy(&piximi); + if (empty) + return ERROR_INT("pixim all fg; no background", procName, 1); + pixZero(pixim, &empty); + if (!empty) /* there are fg pixels in pixim */ + fgpixels = 1; + } + + /* Generate the foreground mask. These pixels will be + * ignored when computing the background values. */ + if (pixg) /* use the input grayscale version if it is provided */ + pixgc = pixClone(pixg); + else + pixgc = pixConvertRGBToGrayFast(pixs); + pixb = pixThresholdToBinary(pixgc, thresh); + pixf = pixMorphSequence(pixb, "d7.1 + d1.7", 0); + pixDestroy(&pixgc); + pixDestroy(&pixb); + + /* Generate the output mask images */ + w = pixGetWidth(pixs); + h = pixGetHeight(pixs); + wm = (w + sx - 1) / sx; + hm = (h + sy - 1) / sy; + pixmr = pixCreate(wm, hm, 8); + pixmg = pixCreate(wm, hm, 8); + pixmb = pixCreate(wm, hm, 8); + + /* ------------- Set up the mapping images --------------- */ + /* Note: we only compute map values in tiles that are complete. + * In general, tiles at right and bottom edges will not be + * complete, and we must fill them in later. */ + nx = w / sx; + ny = h / sy; + wpls = pixGetWpl(pixs); + datas = pixGetData(pixs); + wplf = pixGetWpl(pixf); + dataf = pixGetData(pixf); + for (i = 0; i < ny; i++) { + lines = datas + sy * i * wpls; + linef = dataf + sy * i * wplf; + for (j = 0; j < nx; j++) { + delx = j * sx; + rsum = gsum = bsum = 0; + count = 0; + for (k = 0; k < sy; k++) { + for (m = 0; m < sx; m++) { + if (GET_DATA_BIT(linef + k * wplf, delx + m) == 0) { + pixel = *(lines + k * wpls + delx + m); + rsum += (pixel >> 24); + gsum += ((pixel >> 16) & 0xff); + bsum += ((pixel >> 8) & 0xff); + count++; + } + } + } + if (count >= mincount) { + rval = rsum / count; + gval = gsum / count; + bval = bsum / count; + pixSetPixel(pixmr, j, i, rval); + pixSetPixel(pixmg, j, i, gval); + pixSetPixel(pixmb, j, i, bval); + } + } + } + pixDestroy(&pixf); + + /* If there is an optional mask with fg pixels, erase the previous + * calculation for the corresponding map pixels, setting the + * map values in each of the 3 color maps to 0. Then, when + * all the map holes are filled, these erased pixels will + * be set by the surrounding map values. */ + if (pixim) { + wim = pixGetWidth(pixim); + him = pixGetHeight(pixim); + dataim = pixGetData(pixim); + wplim = pixGetWpl(pixim); + for (i = 0; i < ny; i++) { + yim = i * sy + sy / 2; + if (yim >= him) + break; + lineim = dataim + yim * wplim; + for (j = 0; j < nx; j++) { + xim = j * sx + sx / 2; + if (xim >= wim) + break; + if (GET_DATA_BIT(lineim, xim)) { + pixSetPixel(pixmr, j, i, 0); + pixSetPixel(pixmg, j, i, 0); + pixSetPixel(pixmb, j, i, 0); + } + } + } + } + + /* ----------------- Now fill in the holes ----------------------- */ + if (pixFillMapHoles(pixmr, nx, ny, L_FILL_BLACK) || + pixFillMapHoles(pixmg, nx, ny, L_FILL_BLACK) || + pixFillMapHoles(pixmb, nx, ny, L_FILL_BLACK)) { + pixDestroy(&pixmr); + pixDestroy(&pixmg); + pixDestroy(&pixmb); + L_WARNING("can't make the maps\n", procName); + return 1; + } + + /* Finally, for each connected region corresponding to the + * fg mask, reset all pixels to their average value. */ + if (pixim && fgpixels) { + scalex = 1. / (l_float32)sx; + scaley = 1. / (l_float32)sy; + pixims = pixScaleBySampling(pixim, scalex, scaley); + pixSmoothConnectedRegions(pixmr, pixims, 2); + pixSmoothConnectedRegions(pixmg, pixims, 2); + pixSmoothConnectedRegions(pixmb, pixims, 2); + pixDestroy(&pixims); + } + + *ppixmr = pixmr; + *ppixmg = pixmg; + *ppixmb = pixmb; + pixCopyResolution(*ppixmr, pixs); + pixCopyResolution(*ppixmg, pixs); + pixCopyResolution(*ppixmb, pixs); + return 0; +} + + +/*! + * \brief pixGetBackgroundGrayMapMorph() + * + * \param[in] pixs 8 bpp grayscale; not cmapped + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null; it + * should not have all foreground pixels + * \param[in] reduction factor at which closing is performed + * \param[in] size of square Sel for the closing; use an odd number + * \param[out] ppixm grayscale map + * \return 0 if OK, 1 on error + */ +l_ok +pixGetBackgroundGrayMapMorph(PIX *pixs, + PIX *pixim, + l_int32 reduction, + l_int32 size, + PIX **ppixm) +{ +l_int32 nx, ny, empty, fgpixels; +l_float32 scale; +PIX *pixm, *pix1, *pix2, *pix3, *pixims; + + PROCNAME("pixGetBackgroundGrayMapMorph"); + + if (!ppixm) + return ERROR_INT("&pixm not defined", procName, 1); + *ppixm = NULL; + if (!pixs || pixGetDepth(pixs) != 8) + return ERROR_INT("pixs not defined or not 8 bpp", procName, 1); + if (pixGetColormap(pixs)) + return ERROR_INT("pixs is colormapped", procName, 1); + if (pixim && pixGetDepth(pixim) != 1) + return ERROR_INT("pixim not 1 bpp", procName, 1); + + /* Evaluate the mask pixim and make sure it is not all foreground. */ + fgpixels = 0; /* boolean for existence of fg mask pixels */ + if (pixim) { + pixInvert(pixim, pixim); /* set background pixels to 1 */ + pixZero(pixim, &empty); + if (empty) + return ERROR_INT("pixim all fg; no background", procName, 1); + pixInvert(pixim, pixim); /* revert to original mask */ + pixZero(pixim, &empty); + if (!empty) /* there are fg pixels in pixim */ + fgpixels = 1; + } + + /* Downscale as requested and do the closing to get the background. */ + scale = 1. / (l_float32)reduction; + pix1 = pixScaleBySampling(pixs, scale, scale); + pix2 = pixCloseGray(pix1, size, size); + pix3 = pixExtendByReplication(pix2, 1, 1); + pixDestroy(&pix1); + pixDestroy(&pix2); + + /* Downscale the image mask, if any, and remove it from the + * background. These pixels will be filled in (twice). */ + pixims = NULL; + if (pixim) { + pixims = pixScale(pixim, scale, scale); + pixm = pixConvertTo8(pixims, FALSE); + pixAnd(pixm, pixm, pix3); + } + else + pixm = pixClone(pix3); + pixDestroy(&pix3); + + /* Fill all the holes in the map. */ + nx = pixGetWidth(pixs) / reduction; + ny = pixGetHeight(pixs) / reduction; + if (pixFillMapHoles(pixm, nx, ny, L_FILL_BLACK)) { + pixDestroy(&pixm); + pixDestroy(&pixims); + L_WARNING("can't make the map\n", procName); + return 1; + } + + /* Finally, for each connected region corresponding to the + * fg mask, reset all pixels to their average value. */ + if (pixim && fgpixels) + pixSmoothConnectedRegions(pixm, pixims, 2); + pixDestroy(&pixims); + + *ppixm = pixm; + pixCopyResolution(*ppixm, pixs); + return 0; +} + + +/*! + * \brief pixGetBackgroundRGBMapMorph() + * + * \param[in] pixs 32 bpp rgb + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null; it + * should not have all foreground pixels + * \param[in] reduction factor at which closing is performed + * \param[in] size of square Sel for the closing; use an odd number + * \param[out] ppixmr red component map + * \param[out] ppixmg green component map + * \param[out] ppixmb blue component map + * \return 0 if OK, 1 on error + */ +l_ok +pixGetBackgroundRGBMapMorph(PIX *pixs, + PIX *pixim, + l_int32 reduction, + l_int32 size, + PIX **ppixmr, + PIX **ppixmg, + PIX **ppixmb) +{ +l_int32 nx, ny, empty, fgpixels; +l_float32 scale; +PIX *pixm, *pixmr, *pixmg, *pixmb, *pix1, *pix2, *pix3, *pixims; + + PROCNAME("pixGetBackgroundRGBMapMorph"); + + if (!ppixmr || !ppixmg || !ppixmb) + return ERROR_INT("&pixm* not all defined", procName, 1); + *ppixmr = *ppixmg = *ppixmb = NULL; + if (!pixs) + return ERROR_INT("pixs not defined", procName, 1); + if (pixGetDepth(pixs) != 32) + return ERROR_INT("pixs not 32 bpp", procName, 1); + if (pixim && pixGetDepth(pixim) != 1) + return ERROR_INT("pixim not 1 bpp", procName, 1); + + /* Evaluate the mask pixim and make sure it is not all foreground. */ + fgpixels = 0; /* boolean for existence of fg mask pixels */ + if (pixim) { + pixInvert(pixim, pixim); /* set background pixels to 1 */ + pixZero(pixim, &empty); + if (empty) + return ERROR_INT("pixim all fg; no background", procName, 1); + pixInvert(pixim, pixim); /* revert to original mask */ + pixZero(pixim, &empty); + if (!empty) /* there are fg pixels in pixim */ + fgpixels = 1; + } + + /* Generate an 8 bpp version of the image mask, if it exists */ + scale = 1. / (l_float32)reduction; + pixims = NULL; + pixm = NULL; + if (pixim) { + pixims = pixScale(pixim, scale, scale); + pixm = pixConvertTo8(pixims, FALSE); + } + + /* Downscale as requested and do the closing to get the background. + * Then remove the image mask pixels from the background. They + * will be filled in (twice) later. Do this for all 3 components. */ + pix1 = pixScaleRGBToGrayFast(pixs, reduction, COLOR_RED); + pix2 = pixCloseGray(pix1, size, size); + pix3 = pixExtendByReplication(pix2, 1, 1); + if (pixim) + pixmr = pixAnd(NULL, pixm, pix3); + else + pixmr = pixClone(pix3); + pixDestroy(&pix1); + pixDestroy(&pix2); + pixDestroy(&pix3); + + pix1 = pixScaleRGBToGrayFast(pixs, reduction, COLOR_GREEN); + pix2 = pixCloseGray(pix1, size, size); + pix3 = pixExtendByReplication(pix2, 1, 1); + if (pixim) + pixmg = pixAnd(NULL, pixm, pix3); + else + pixmg = pixClone(pix3); + pixDestroy(&pix1); + pixDestroy(&pix2); + pixDestroy(&pix3); + + pix1 = pixScaleRGBToGrayFast(pixs, reduction, COLOR_BLUE); + pix2 = pixCloseGray(pix1, size, size); + pix3 = pixExtendByReplication(pix2, 1, 1); + if (pixim) + pixmb = pixAnd(NULL, pixm, pix3); + else + pixmb = pixClone(pix3); + pixDestroy(&pixm); + pixDestroy(&pix1); + pixDestroy(&pix2); + pixDestroy(&pix3); + + /* Fill all the holes in the three maps. */ + nx = pixGetWidth(pixs) / reduction; + ny = pixGetHeight(pixs) / reduction; + if (pixFillMapHoles(pixmr, nx, ny, L_FILL_BLACK) || + pixFillMapHoles(pixmg, nx, ny, L_FILL_BLACK) || + pixFillMapHoles(pixmb, nx, ny, L_FILL_BLACK)) { + pixDestroy(&pixmr); + pixDestroy(&pixmg); + pixDestroy(&pixmb); + pixDestroy(&pixims); + L_WARNING("can't make the maps\n", procName); + return 1; + } + + /* Finally, for each connected region corresponding to the + * fg mask in each component, reset all pixels to their + * average value. */ + if (pixim && fgpixels) { + pixSmoothConnectedRegions(pixmr, pixims, 2); + pixSmoothConnectedRegions(pixmg, pixims, 2); + pixSmoothConnectedRegions(pixmb, pixims, 2); + pixDestroy(&pixims); + } + + *ppixmr = pixmr; + *ppixmg = pixmg; + *ppixmb = pixmb; + pixCopyResolution(*ppixmr, pixs); + pixCopyResolution(*ppixmg, pixs); + pixCopyResolution(*ppixmb, pixs); + return 0; +} + + +/*! + * \brief pixFillMapHoles() + * + * \param[in] pix 8 bpp; a map, with one pixel for each tile in + * a larger image + * \param[in] nx number of horizontal pixel tiles that are entirely + * covered with pixels in the original source image + * \param[in] ny ditto for the number of vertical pixel tiles + * \param[in] filltype L_FILL_WHITE or L_FILL_BLACK + * \return 0 if OK, 1 on error + * + * <pre> + * Notes: + * (1) This is an in-place operation on pix (the map). pix is + * typically a low-resolution version of some other image + * from which it was derived, where each pixel in pix + * corresponds to a rectangular tile (say, m x n) of pixels + * in the larger image. All we need to know about the larger + * image is whether or not the rightmost column and bottommost + * row of pixels in pix correspond to tiles that are + * only partially covered by pixels in the larger image. + * (2) Typically, some number of pixels in the input map are + * not known, and their values must be determined by near + * pixels that are known. These unknown pixels are the 'holes'. + * They can take on only two values, 0 and 255, and the + * instruction about which to fill is given by the filltype flag. + * (3) The "holes" can come from two sources. The first is when there + * are not enough foreground or background pixels in a tile; + * the second is when a tile is at least partially covered + * by an image mask. If we're filling holes in a fg mask, + * the holes are initialized to black (0) and use L_FILL_BLACK. + * For filling holes in a bg mask, initialize the holes to + * white (255) and use L_FILL_WHITE. + * (4) If w is the map width, nx = w or nx = w - 1; ditto for h and ny. + * </pre> + */ +l_ok +pixFillMapHoles(PIX *pix, + l_int32 nx, + l_int32 ny, + l_int32 filltype) +{ +l_int32 w, h, y, nmiss, goodcol, i, j, found, ival, valtest; +l_uint32 val, lastval; +NUMA *na; /* indicates if there is any data in the column */ + + PROCNAME("pixFillMapHoles"); + + if (!pix || pixGetDepth(pix) != 8) + return ERROR_INT("pix not defined or not 8 bpp", procName, 1); + if (pixGetColormap(pix)) + return ERROR_INT("pix is colormapped", procName, 1); + + /* ------------- Fill holes in the mapping image columns ----------- */ + pixGetDimensions(pix, &w, &h, NULL); + na = numaCreate(0); /* holds flag for which columns have data */ + nmiss = 0; + valtest = (filltype == L_FILL_WHITE) ? 255 : 0; + for (j = 0; j < nx; j++) { /* do it by columns */ + found = FALSE; + for (i = 0; i < ny; i++) { + pixGetPixel(pix, j, i, &val); + if (val != valtest) { + y = i; + found = TRUE; + break; + } + } + if (found == FALSE) { + numaAddNumber(na, 0); /* no data in the column */ + nmiss++; + } + else { + numaAddNumber(na, 1); /* data in the column */ + for (i = y - 1; i >= 0; i--) /* replicate upwards to top */ + pixSetPixel(pix, j, i, val); + pixGetPixel(pix, j, 0, &lastval); + for (i = 1; i < h; i++) { /* set going down to bottom */ + pixGetPixel(pix, j, i, &val); + if (val == valtest) + pixSetPixel(pix, j, i, lastval); + else + lastval = val; + } + } + } + + if (nmiss == nx) { /* no data in any column! */ + numaDestroy(&na); + L_WARNING("no bg found; no data in any column\n", procName); + return 1; + } + + /* ---------- Fill in missing columns by replication ----------- */ + if (nmiss > 0) { /* replicate columns */ + /* Find the first good column */ + goodcol = 0; + for (j = 0; j < w; j++) { + numaGetIValue(na, j, &ival); + if (ival == 1) { + goodcol = j; + break; + } + } + if (goodcol > 0) { /* copy cols backward */ + for (j = goodcol - 1; j >= 0; j--) + pixRasterop(pix, j, 0, 1, h, PIX_SRC, pix, j + 1, 0); + } + for (j = goodcol + 1; j < w; j++) { /* copy cols forward */ + numaGetIValue(na, j, &ival); + if (ival == 0) { + /* Copy the column to the left of j */ + pixRasterop(pix, j, 0, 1, h, PIX_SRC, pix, j - 1, 0); + } + } + } + if (w > nx) { /* replicate the last column */ + pixRasterop(pix, w - 1, 0, 1, h, PIX_SRC, pix, w - 2, 0); + } + + numaDestroy(&na); + return 0; +} + + +/*! + * \brief pixExtendByReplication() + * + * \param[in] pixs 8 bpp + * \param[in] addw number of extra pixels horizontally to add + * \param[in] addh number of extra pixels vertically to add + * \return pixd extended with replicated pixel values, or NULL on error + * + * <pre> + * Notes: + * (1) The pixel values are extended to the left and down, as required. + * </pre> + */ +PIX * +pixExtendByReplication(PIX *pixs, + l_int32 addw, + l_int32 addh) +{ +l_int32 w, h, i, j; +l_uint32 val; +PIX *pixd; + + PROCNAME("pixExtendByReplication"); + + if (!pixs || pixGetDepth(pixs) != 8) + return (PIX *)ERROR_PTR("pixs undefined or not 8 bpp", procName, NULL); + + if (addw == 0 && addh == 0) + return pixCopy(NULL, pixs); + + pixGetDimensions(pixs, &w, &h, NULL); + if ((pixd = pixCreate(w + addw, h + addh, 8)) == NULL) + return (PIX *)ERROR_PTR("pixd not made", procName, NULL); + pixRasterop(pixd, 0, 0, w, h, PIX_SRC, pixs, 0, 0); + + if (addw > 0) { + for (i = 0; i < h; i++) { + pixGetPixel(pixd, w - 1, i, &val); + for (j = 0; j < addw; j++) + pixSetPixel(pixd, w + j, i, val); + } + } + + if (addh > 0) { + for (j = 0; j < w + addw; j++) { + pixGetPixel(pixd, j, h - 1, &val); + for (i = 0; i < addh; i++) + pixSetPixel(pixd, j, h + i, val); + } + } + + pixCopyResolution(pixd, pixs); + return pixd; +} + + +/*! + * \brief pixSmoothConnectedRegions() + * + * \param[in] pixs 8 bpp grayscale; no colormap + * \param[in] pixm [optional] 1 bpp; if null, this is a no-op + * \param[in] factor subsampling factor for getting average; >= 1 + * \return 0 if OK, 1 on error + * + * <pre> + * Notes: + * (1) The pixels in pixs corresponding to those in each + * 8-connected region in the mask are set to the average value. + * (2) This is required for adaptive mapping to avoid the + * generation of stripes in the background map, due to + * variations in the pixel values near the edges of mask regions. + * (3) This function is optimized for background smoothing, where + * there are a relatively small number of components. It will + * be inefficient if used where there are many small components. + * </pre> + */ +l_ok +pixSmoothConnectedRegions(PIX *pixs, + PIX *pixm, + l_int32 factor) +{ +l_int32 empty, i, n, x, y; +l_float32 aveval; +BOXA *boxa; +PIX *pixmc; +PIXA *pixa; + + PROCNAME("pixSmoothConnectedRegions"); + + if (!pixs || pixGetDepth(pixs) != 8) + return ERROR_INT("pixs not defined or not 8 bpp", procName, 1); + if (pixGetColormap(pixs)) + return ERROR_INT("pixs has colormap", procName, 1); + if (!pixm) { + L_INFO("pixm not defined\n", procName); + return 0; + } + if (pixGetDepth(pixm) != 1) + return ERROR_INT("pixm not 1 bpp", procName, 1); + pixZero(pixm, &empty); + if (empty) { + L_INFO("pixm has no fg pixels; nothing to do\n", procName); + return 0; + } + + boxa = pixConnComp(pixm, &pixa, 8); + n = boxaGetCount(boxa); + for (i = 0; i < n; i++) { + if ((pixmc = pixaGetPix(pixa, i, L_CLONE)) == NULL) { + L_WARNING("missing pixmc!\n", procName); + continue; + } + boxaGetBoxGeometry(boxa, i, &x, &y, NULL, NULL); + pixGetAverageMasked(pixs, pixmc, x, y, factor, L_MEAN_ABSVAL, &aveval); + pixPaintThroughMask(pixs, pixmc, x, y, (l_int32)aveval); + pixDestroy(&pixmc); + } + + boxaDestroy(&boxa); + pixaDestroy(&pixa); + return 0; +} + + +/*------------------------------------------------------------------* + * Measurement of local foreground * + *------------------------------------------------------------------*/ +#if 0 /* Not working properly: do not use */ + +/*! + * \brief pixGetForegroundGrayMap() + * + * \param[in] pixs 8 bpp + * \param[in] pixim [optional] 1 bpp 'image' mask; can be null + * \param[in] sx, sy src tile size, in pixels + * \param[in] thresh threshold for determining foreground + * \param[out] ppixd 8 bpp grayscale map + * \return 0 if OK, 1 on error + * + * <pre> + * Notes: + * (1) Each (sx, sy) tile of pixs gets mapped to one pixel in pixd. + * (2) pixd is the estimate of the fg (darkest) value within each tile. + * (3) All pixels in pixd that are in 'image' regions, as specified + * by pixim, are given the background value 0. + * (4) For pixels in pixd that can't directly be given a fg value, + * the value is inferred by propagating from neighboring pixels. + * (5) In practice, pixd can be used to normalize the fg, and + * it can be done after background normalization. + * (6) The overall procedure is: + * ~ reduce 2x by sampling + * ~ paint all 'image' pixels white, so that they don't + * ~ participate in the Min reduction + * ~ do a further (sx, sy) Min reduction -- think of + * it as a large opening followed by subsampling by the + * reduction factors + * ~ threshold the result to identify fg, and set the + * bg pixels to 255 (these are 'holes') + * ~ fill holes by propagation from fg values + * ~ replicatively expand by 2x, arriving at the final + * resolution of pixd + * ~ smooth with a 17x17 kernel + * ~ paint the 'image' regions black + * </pre> + */ +l_ok +pixGetForegroundGrayMap(PIX *pixs, + PIX *pixim, + l_int32 sx, + l_int32 sy, + l_int32 thresh, + PIX **ppixd) +{ +l_int32 w, h, d, wd, hd; +l_int32 empty, fgpixels; +PIX *pixd, *piximi, *pixim2, *pixims, *pixs2, *pixb, *pixt1, *pixt2, *pixt3; + + PROCNAME("pixGetForegroundGrayMap"); + + if (!ppixd) + return ERROR_INT("&pixd not defined", procName, 1); + *ppixd = NULL; + if (!pixs) + return ERROR_INT("pixs not defined", procName, 1); + pixGetDimensions(pixs, &w, &h, &d); + if (d != 8) + return ERROR_INT("pixs not 8 bpp", procName, 1); + if (pixim && pixGetDepth(pixim) != 1) + return ERROR_INT("pixim not 1 bpp", procName, 1); + if (sx < 2 || sy < 2) + return ERROR_INT("sx and sy must be >= 2", procName, 1); + + /* Generate pixd, which is reduced by the factors (sx, sy). */ + wd = (w + sx - 1) / sx; + hd = (h + sy - 1) / sy; + pixd = pixCreate(wd, hd, 8); + *ppixd = pixd; + + /* Evaluate the 'image' mask, pixim. If it is all fg, + * the output pixd has all pixels with value 0. */ + fgpixels = 0; /* boolean for existence of fg pixels in the image mask. */ + if (pixim) { + piximi = pixInvert(NULL, pixim); /* set non-image pixels to 1 */ + pixZero(piximi, &empty); + pixDestroy(&piximi); + if (empty) /* all 'image'; return with all pixels set to 0 */ + return 0; + pixZero(pixim, &empty); + if (!empty) /* there are fg pixels in pixim */ + fgpixels = 1; + } + + /* 2x subsampling; paint white through 'image' mask. */ + pixs2 = pixScaleBySampling(pixs, 0.5, 0.5); + if (pixim && fgpixels) { + pixim2 = pixReduceBinary2(pixim, NULL); + pixPaintThroughMask(pixs2, pixim2, 0, 0, 255); + pixDestroy(&pixim2); + } + + /* Min (erosion) downscaling; total reduction (4 sx, 4 sy). */ + pixt1 = pixScaleGrayMinMax(pixs2, sx, sy, L_CHOOSE_MIN); + +/* pixDisplay(pixt1, 300, 200); */ + + /* Threshold to identify fg; paint bg pixels to white. */ + pixb = pixThresholdToBinary(pixt1, thresh); /* fg pixels */ + pixInvert(pixb, pixb); + pixPaintThroughMask(pixt1, pixb, 0, 0, 255); + pixDestroy(&pixb); + + /* Replicative expansion by 2x to (sx, sy). */ + pixt2 = pixExpandReplicate(pixt1, 2); + +/* pixDisplay(pixt2, 500, 200); */ + + /* Fill holes in the fg by propagation */ + pixFillMapHoles(pixt2, w / sx, h / sy, L_FILL_WHITE); + +/* pixDisplay(pixt2, 700, 200); */ + + /* Smooth with 17x17 kernel. */ + pixt3 = pixBlockconv(pixt2, 8, 8); + pixRasterop(pixd, 0, 0, wd, hd, PIX_SRC, pixt3, 0, 0); + + /* Paint the image parts black. */ + pixims = pixScaleBySampling(pixim, 1. / sx, 1. / sy); + pixPaintThroughMask(pixd, pixims, 0, 0, 0); + + pixDestroy(&pixs2); + pixDestroy(&pixt1); + pixDestroy(&pixt2); + pixDestroy(&pixt3); + return 0; +} +#endif /* Not working properly: do not use */ + + +/*------------------------------------------------------------------* + * Generate inverted background map * + *------------------------------------------------------------------*/ +/*! + * \brief pixGetInvBackgroundMap() + * + * \param[in] pixs 8 bpp grayscale; no colormap + * \param[in] bgval target bg val; typ. > 128 + * \param[in] smoothx half-width of block convolution kernel width + * \param[in] smoothy half-width of block convolution kernel height + * \return pixd 16 bpp, or NULL on error + * + * <pre> + * Notes: + * (1) bgval should typically be > 120 and < 240 + * (2) pixd is a normalization image; the original image is + * multiplied by pixd and the result is divided by 256. + * </pre> + */ +PIX * +pixGetInvBackgroundMap(PIX *pixs, + l_int32 bgval, + l_int32 smoothx, + l_int32 smoothy) +{ +l_int32 w, h, wplsm, wpld, i, j; +l_int32 val, val16; +l_uint32 *datasm, *datad, *linesm, *lined; +PIX *pixsm, *pixd; + + PROCNAME("pixGetInvBackgroundMap"); + + if (!pixs || pixGetDepth(pixs) != 8) + return (PIX *)ERROR_PTR("pixs undefined or not 8 bpp", procName, NULL); + if (pixGetColormap(pixs)) + return (PIX *)ERROR_PTR("pixs has colormap", procName, NULL); + pixGetDimensions(pixs, &w, &h, NULL); + if (w < 5 || h < 5) + return (PIX *)ERROR_PTR("w and h must be >= 5", procName, NULL); + + /* smooth the map image */ + pixsm = pixBlockconv(pixs, smoothx, smoothy); + datasm = pixGetData(pixsm); + wplsm = pixGetWpl(pixsm); + + /* invert the map image, scaling up to preserve dynamic range */ + pixd = pixCreate(w, h, 16); + datad = pixGetData(pixd); + wpld = pixGetWpl(pixd); + for (i = 0; i < h; i++) { + linesm = datasm + i * wplsm; + lined = datad + i * wpld; + for (j = 0; j < w; j++) { + val = GET_DATA_BYTE(linesm, j); + if (val > 0) + val16 = (256 * bgval) / val; + else { /* shouldn't happen */ + L_WARNING("smoothed bg has 0 pixel!\n", procName); + val16 = bgval / 2; + } + SET_DATA_TWO_BYTES(lined, j, val16); + } + } + + pixDestroy(&pixsm); + pixCopyResolution(pixd, pixs); + return pixd; +} + + +/*------------------------------------------------------------------* + * Apply background map to image * + *------------------------------------------------------------------*/ +/*! + * \brief pixApplyInvBackgroundGrayMap() + * + * \param[in] pixs 8 bpp grayscale; no colormap + * \param[in] pixm 16 bpp, inverse background map + * \param[in] sx tile width in pixels + * \param[in] sy tile height in pixels + * \return pixd 8 bpp, or NULL on error + */ +PIX * +pixApplyInvBackgroundGrayMap(PIX *pixs, + PIX *pixm, + l_int32 sx, + l_int32 sy) +{ +l_int32 w, h, wm, hm, wpls, wpld, i, j, k, m, xoff, yoff; +l_int32 vals, vald; +l_uint32 val16; +l_uint32 *datas, *datad, *lines, *lined, *flines, *flined; +PIX *pixd; + + PROCNAME("pixApplyInvBackgroundGrayMap"); + + if (!pixs || pixGetDepth(pixs) != 8) + return (PIX *)ERROR_PTR("pixs undefined or not 8 bpp", procName, NULL); + if (pixGetColormap(pixs)) + return (PIX *)ERROR_PTR("pixs has colormap", procName, NULL); + if (!pixm || pixGetDepth(pixm) != 16) + return (PIX *)ERROR_PTR("pixm undefined or not 16 bpp", procName, NULL); + if (sx == 0 || sy == 0) + return (PIX *)ERROR_PTR("invalid sx and/or sy", procName, NULL); + + datas = pixGetData(pixs); + wpls = pixGetWpl(pixs); + pixGetDimensions(pixs, &w, &h, NULL); + pixGetDimensions(pixm, &wm, &hm, NULL); + if ((pixd = pixCreateTemplate(pixs)) == NULL) + return (PIX *)ERROR_PTR("pixd not made", procName, NULL); + datad = pixGetData(pixd); + wpld = pixGetWpl(pixd); + for (i = 0; i < hm; i++) { + lines = datas + sy * i * wpls; + lined = datad + sy * i * wpld; + yoff = sy * i; + for (j = 0; j < wm; j++) { + pixGetPixel(pixm, j, i, &val16); + xoff = sx * j; + for (k = 0; k < sy && yoff + k < h; k++) { + flines = lines + k * wpls; + flined = lined + k * wpld; + for (m = 0; m < sx && xoff + m < w; m++) { + vals = GET_DATA_BYTE(flines, xoff + m); + vald = (vals * val16) / 256; + vald = L_MIN(vald, 255); + SET_DATA_BYTE(flined, xoff + m, vald); + } + } + } + } + + return pixd; +} + + +/*! + * \brief pixApplyInvBackgroundRGBMap() + * + * \param[in] pixs 32 bpp rbg + * \param[in] pixmr 16 bpp, red inverse background map + * \param[in] pixmg 16 bpp, green inverse background map + * \param[in] pixmb 16 bpp, blue inverse background map + * \param[in] sx tile width in pixels + * \param[in] sy tile height in pixels + * \return pixd 32 bpp rbg, or NULL on error + */ +PIX * +pixApplyInvBackgroundRGBMap(PIX *pixs, + PIX *pixmr, + PIX *pixmg, + PIX *pixmb, + l_int32 sx, + l_int32 sy) +{ +l_int32 w, h, wm, hm, wpls, wpld, i, j, k, m, xoff, yoff; +l_int32 rvald, gvald, bvald; +l_uint32 vals; +l_uint32 rval16, gval16, bval16; +l_uint32 *datas, *datad, *lines, *lined, *flines, *flined; +PIX *pixd; + + PROCNAME("pixApplyInvBackgroundRGBMap"); + + if (!pixs) + return (PIX *)ERROR_PTR("pixs not defined", procName, NULL); + if (pixGetDepth(pixs) != 32) + return (PIX *)ERROR_PTR("pixs not 32 bpp", procName, NULL); + if (!pixmr || !pixmg || !pixmb) + return (PIX *)ERROR_PTR("pix maps not all defined", procName, NULL); + if (pixGetDepth(pixmr) != 16 || pixGetDepth(pixmg) != 16 || + pixGetDepth(pixmb) != 16) + return (PIX *)ERROR_PTR("pix maps not all 16 bpp", procName, NULL); + if (sx == 0 || sy == 0) + return (PIX *)ERROR_PTR("invalid sx and/or sy", procName, NULL); + + datas = pixGetData(pixs); + wpls = pixGetWpl(pixs); + w = pixGetWidth(pixs); + h = pixGetHeight(pixs); + wm = pixGetWidth(pixmr); + hm = pixGetHeight(pixmr); + if ((pixd = pixCreateTemplate(pixs)) == NULL) + return (PIX *)ERROR_PTR("pixd not made", procName, NULL); + datad = pixGetData(pixd); + wpld = pixGetWpl(pixd); + for (i = 0; i < hm; i++) { + lines = datas + sy * i * wpls; + lined = datad + sy * i * wpld; + yoff = sy * i; + for (j = 0; j < wm; j++) { + pixGetPixel(pixmr, j, i, &rval16); + pixGetPixel(pixmg, j, i, &gval16); + pixGetPixel(pixmb, j, i, &bval16); + xoff = sx * j; + for (k = 0; k < sy && yoff + k < h; k++) { + flines = lines + k * wpls; + flined = lined + k * wpld; + for (m = 0; m < sx && xoff + m < w; m++) { + vals = *(flines + xoff + m); + rvald = ((vals >> 24) * rval16) / 256; + rvald = L_MIN(rvald, 255); + gvald = (((vals >> 16) & 0xff) * gval16) / 256; + gvald = L_MIN(gvald, 255); + bvald = (((vals >> 8) & 0xff) * bval16) / 256; + bvald = L_MIN(bvald, 255); + composeRGBPixel(rvald, gvald, bvald, flined + xoff + m); + } + } + } + } + + return pixd; +} + + +/*------------------------------------------------------------------* + * Apply variable map * + *------------------------------------------------------------------*/ +/*! + * \brief pixApplyVariableGrayMap() + * + * \param[in] pixs 8 bpp + * \param[in] pixg 8 bpp, variable map + * \param[in] target typ. 128 for threshold + * \return pixd 8 bpp, or NULL on error + * + * <pre> + * Notes: + * (1) Suppose you have an image that you want to transform based + * on some photometric measurement at each point, such as the + * threshold value for binarization. Representing the photometric + * measurement as an image pixg, you can threshold in input image + * using pixVarThresholdToBinary(). Alternatively, you can map + * the input image pointwise so that the threshold over the + * entire image becomes a constant, such as 128. For example, + * if a pixel in pixg is 150 and the target is 128, the + * corresponding pixel in pixs is mapped linearly to a value + * (128/150) of the input value. If the resulting mapped image + * pixd were then thresholded at 128, you would obtain the + * same result as a direct binarization using pixg with + * pixVarThresholdToBinary(). + * (2) The sizes of pixs and pixg must be equal. + * </pre> + */ +PIX * +pixApplyVariableGrayMap(PIX *pixs, + PIX *pixg, + l_int32 target) +{ +l_int32 i, j, w, h, d, wpls, wplg, wpld, vals, valg, vald; +l_uint8 *lut; +l_uint32 *datas, *datag, *datad, *lines, *lineg, *lined; +l_float32 fval; +PIX *pixd; + + PROCNAME("pixApplyVariableGrayMap"); + + if (!pixs) + return (PIX *)ERROR_PTR("pixs not defined", procName, NULL); + if (!pixg) + return (PIX *)ERROR_PTR("pixg not defined", procName, NULL); + if (!pixSizesEqual(pixs, pixg)) + return (PIX *)ERROR_PTR("pix sizes not equal", procName, NULL); + pixGetDimensions(pixs, &w, &h, &d); + if (d != 8) + return (PIX *)ERROR_PTR("depth not 8 bpp", procName, NULL); + + /* Generate a LUT for the mapping if the image is large enough + * to warrant the overhead. The LUT is of size 2^16. For the + * index to the table, get the MSB from pixs and the LSB from pixg. + * Note: this LUT is bigger than the typical 32K L1 cache, so + * we expect cache misses. L2 latencies are about 5ns. But + * division is slooooow. For large images, this function is about + * 4x faster when using the LUT. C'est la vie. */ + lut = NULL; + if (w * h > 100000) { /* more pixels than 2^16 */ + if ((lut = (l_uint8 *)LEPT_CALLOC(0x10000, sizeof(l_uint8))) == NULL) + return (PIX *)ERROR_PTR("lut not made", procName, NULL); + for (i = 0; i < 256; i++) { + for (j = 0; j < 256; j++) { + fval = (l_float32)(i * target) / (j + 0.5); + lut[(i << 8) + j] = L_MIN(255, (l_int32)(fval + 0.5)); + } + } + } + + if ((pixd = pixCreateNoInit(w, h, 8)) == NULL) { + LEPT_FREE(lut); + return (PIX *)ERROR_PTR("pixd not made", procName, NULL); + } + pixCopyResolution(pixd, pixs); + datad = pixGetData(pixd); + wpld = pixGetWpl(pixd); + datas = pixGetData(pixs); + wpls = pixGetWpl(pixs); + datag = pixGetData(pixg); + wplg = pixGetWpl(pixg); + for (i = 0; i < h; i++) { + lines = datas + i * wpls; + lineg = datag + i * wplg; + lined = datad + i * wpld; + if (lut) { + for (j = 0; j < w; j++) { + vals = GET_DATA_BYTE(lines, j); + valg = GET_DATA_BYTE(lineg, j); + vald = lut[(vals << 8) + valg]; + SET_DATA_BYTE(lined, j, vald); + } + } + else { + for (j = 0; j < w; j++) { + vals = GET_DATA_BYTE(lines, j); + valg = GET_DATA_BYTE(lineg, j); + fval = (l_float32)(vals * target) / (valg + 0.5); + vald = L_MIN(255, (l_int32)(fval + 0.5)); + SET_DATA_BYTE(lined, j, vald); + } + } + } + + LEPT_FREE(lut); + return pixd; +} + + +/*------------------------------------------------------------------* + * Non-adaptive (global) mapping * + *------------------------------------------------------------------*/ +/*! + * \brief pixGlobalNormRGB() + * + * \param[in] pixd [optional] null, existing or equal to pixs + * \param[in] pixs 32 bpp rgb, or colormapped + * \param[in] rval, gval, bval pixel values in pixs that are + * linearly mapped to mapval + * \param[in] mapval use 255 for mapping to white + * \return pixd 32 bpp rgb or colormapped, or NULL on error + * + * <pre> + * Notes: + * (1) The value of pixd determines if the results are written to a + * new pix (use NULL), in-place to pixs (use pixs), or to some + * other existing pix. + * (2) This does a global normalization of an image where the + * r,g,b color components are not balanced. Thus, white in pixs is + * represented by a set of r,g,b values that are not all 255. + * (3) The input values (rval, gval, bval) should be chosen to + * represent the gray color (mapval, mapval, mapval) in src. + * Thus, this function will map (rval, gval, bval) to that gray color. + * (4) Typically, mapval = 255, so that (rval, gval, bval) + * corresponds to the white point of src. In that case, these + * parameters should be chosen so that few pixels have higher values. + * (5) In all cases, we do a linear TRC separately on each of the + * components, saturating at 255. + * (6) If the input pix is 8 bpp without a colormap, you can get + * this functionality with mapval = 255 by calling: + * pixGammaTRC(pixd, pixs, 1.0, 0, bgval); + * where bgval is the value you want to be mapped to 255. + * Or more generally, if you want bgval to be mapped to mapval: + * pixGammaTRC(pixd, pixs, 1.0, 0, 255 * bgval / mapval); + * </pre> + */ +PIX * +pixGlobalNormRGB(PIX *pixd, + PIX *pixs, + l_int32 rval, + l_int32 gval, + l_int32 bval, + l_int32 mapval) +{ +l_int32 w, h, d, i, j, ncolors, rv, gv, bv, wpl; +l_int32 *rarray, *garray, *barray; +l_uint32 *data, *line; +NUMA *nar, *nag, *nab; +PIXCMAP *cmap; + + PROCNAME("pixGlobalNormRGB"); + + if (!pixs) + return (PIX *)ERROR_PTR("pixs not defined", procName, NULL); + cmap = pixGetColormap(pixs); + pixGetDimensions(pixs, &w, &h, &d); + if (!cmap && d != 32) + return (PIX *)ERROR_PTR("pixs not cmapped or 32 bpp", procName, NULL); + if (mapval <= 0) { + L_WARNING("mapval must be > 0; setting to 255\n", procName); + mapval = 255; + } + + /* Prepare pixd to be a copy of pixs */ + if ((pixd = pixCopy(pixd, pixs)) == NULL) + return (PIX *)ERROR_PTR("pixd not made", procName, NULL); + + /* Generate the TRC maps for each component. Make sure the + * upper range for each color is greater than zero. */ + nar = numaGammaTRC(1.0, 0, L_MAX(1, 255 * rval / mapval)); + nag = numaGammaTRC(1.0, 0, L_MAX(1, 255 * gval / mapval)); + nab = numaGammaTRC(1.0, 0, L_MAX(1, 255 * bval / mapval)); + + /* Extract copies of the internal arrays */ + rarray = numaGetIArray(nar); + garray = numaGetIArray(nag); + barray = numaGetIArray(nab); + if (!nar || !nag || !nab || !rarray || !garray || !barray) { + L_ERROR("allocation failure in arrays\n", procName); + goto cleanup_arrays; + } + + if (cmap) { + ncolors = pixcmapGetCount(cmap); + for (i = 0; i < ncolors; i++) { + pixcmapGetColor(cmap, i, &rv, &gv, &bv); + pixcmapResetColor(cmap, i, rarray[rv], garray[gv], barray[bv]); + } + } + else { + data = pixGetData(pixd); + wpl = pixGetWpl(pixd); + for (i = 0; i < h; i++) { + line = data + i * wpl; + for (j = 0; j < w; j++) { + extractRGBValues(line[j], &rv, &gv, &bv); + composeRGBPixel(rarray[rv], garray[gv], barray[bv], line + j); + } + } + } + +cleanup_arrays: + numaDestroy(&nar); + numaDestroy(&nag); + numaDestroy(&nab); + LEPT_FREE(rarray); + LEPT_FREE(garray); + LEPT_FREE(barray); + return pixd; +} + + +/*! + * \brief pixGlobalNormNoSatRGB() + * + * \param[in] pixd [optional] null, existing or equal to pixs + * \param[in] pixs 32 bpp rgb + * \param[in] rval, gval, bval pixel values in pixs that are + * linearly mapped to mapval; but see below + * \param[in] factor subsampling factor; integer >= 1 + * \param[in] rank between 0.0 and 1.0; typ. use a value near 1.0 + * \return pixd 32 bpp rgb, or NULL on error + * + * <pre> + * Notes: + * (1) This is a version of pixGlobalNormRGB(), where the output + * intensity is scaled back so that a controlled fraction of + * pixel components is allowed to saturate. See comments in + * pixGlobalNormRGB(). + * (2) The value of pixd determines if the results are written to a + * new pix (use NULL), in-place to pixs (use pixs), or to some + * other existing pix. + * (3) This does a global normalization of an image where the + * r,g,b color components are not balanced. Thus, white in pixs is + * represented by a set of r,g,b values that are not all 255. + * (4) The input values (rval, gval, bval) can be chosen to be the + * color that, after normalization, becomes white background. + * For images that are mostly background, the closer these values + * are to the median component values, the closer the resulting + * background will be to gray, becoming white at the brightest places. + * (5) The mapval used in pixGlobalNormRGB() is computed here to + * avoid saturation of any component in the image (save for a + * fraction of the pixels given by the input rank value). + * </pre> + */ +PIX * +pixGlobalNormNoSatRGB(PIX *pixd, + PIX *pixs, + l_int32 rval, + l_int32 gval, + l_int32 bval, + l_int32 factor, + l_float32 rank) +{ +l_int32 mapval; +l_float32 rankrval, rankgval, rankbval; +l_float32 rfract, gfract, bfract, maxfract; + + PROCNAME("pixGlobalNormNoSatRGB"); + + if (!pixs) + return (PIX *)ERROR_PTR("pixs not defined", procName, NULL); + if (pixGetDepth(pixs) != 32) + return (PIX *)ERROR_PTR("pixs not 32 bpp", procName, NULL); + if (factor < 1) + return (PIX *)ERROR_PTR("sampling factor < 1", procName, NULL); + if (rank < 0.0 || rank > 1.0) + return (PIX *)ERROR_PTR("rank not in [0.0 ... 1.0]", procName, NULL); + if (rval <= 0 || gval <= 0 || bval <= 0) + return (PIX *)ERROR_PTR("invalid estim. color values", procName, NULL); + + /* The max value for each component may be larger than the + * input estimated background value. In that case, mapping + * for those pixels would saturate. To prevent saturation, + * we compute the fraction for each component by which we + * would oversaturate. Then take the max of these, and + * reduce, uniformly over all components, the output intensity + * by this value. Then no component will saturate. + * In practice, if rank < 1.0, a fraction of pixels + * may have a component saturate. By keeping rank close to 1.0, + * that fraction can be made arbitrarily small. */ + pixGetRankValueMaskedRGB(pixs, NULL, 0, 0, factor, rank, &rankrval, + &rankgval, &rankbval); + rfract = rankrval / (l_float32)rval; + gfract = rankgval / (l_float32)gval; + bfract = rankbval / (l_float32)bval; + maxfract = L_MAX(rfract, gfract); + maxfract = L_MAX(maxfract, bfract); +#if DEBUG_GLOBAL + lept_stderr("rankrval = %7.2f, rankgval = %7.2f, rankbval = %7.2f\n", + rankrval, rankgval, rankbval); + lept_stderr("rfract = %7.4f, gfract = %7.4f, bfract = %7.4f\n", + rfract, gfract, bfract); +#endif /* DEBUG_GLOBAL */ + + mapval = (l_int32)(255. / maxfract); + pixd = pixGlobalNormRGB(pixd, pixs, rval, gval, bval, mapval); + return pixd; +} + + +/*------------------------------------------------------------------* + * Adaptive threshold spread normalization * + *------------------------------------------------------------------*/ +/*! + * \brief pixThresholdSpreadNorm() + * + * \param[in] pixs 8 bpp grayscale; not colormapped + * \param[in] filtertype L_SOBEL_EDGE or L_TWO_SIDED_EDGE; + * \param[in] edgethresh threshold on magnitude of edge filter; + * typ 10-20 + * \param[in] smoothx, smoothy half-width of convolution kernel applied to + * spread threshold: use 0 for no smoothing + * \param[in] gamma gamma correction; typ. about 0.7 + * \param[in] minval input value that gives 0 for output; typ. -25 + * \param[in] maxval input value that gives 255 for output; + * typ. 255 + * \param[in] targetthresh target threshold for normalization + * \param[out] ppixth [optional] computed local threshold value + * \param[out] ppixb [optional] thresholded normalized image + * \param[out] ppixd [optional] normalized image + * \return 0 if OK, 1 on error + * + * <pre> + * Notes: + * (1) The basis of this approach is the use of seed spreading + * on a (possibly) sparse set of estimates for the local threshold. + * The resulting dense estimates are smoothed by convolution + * and used to either threshold the input image or normalize it + * with a local transformation that linearly maps the pixels so + * that the local threshold estimate becomes constant over the + * resulting image. This approach is one of several that + * have been suggested (and implemented) by Ray Smith. + * (2) You can use either the Sobel or TwoSided edge filters. + * The results appear to be similar, using typical values + * of edgethresh in the rang 10-20. + * (3) To skip the trc enhancement, use gamma = 1.0, minval = 0 + * and maxval = 255. + * (4) For the normalized image pixd, each pixel is linearly mapped + * in such a way that the local threshold is equal to targetthresh. + * (5) The full width and height of the convolution kernel + * are (2 * smoothx + 1) and (2 * smoothy + 1). + * (6) This function can be used with the pixtiling utility if the + * images are too large. See pixOtsuAdaptiveThreshold() for + * an example of this. + * </pre> + */ +l_ok +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) +{ +PIX *pixe, *pixet, *pixsd, *pixg1, *pixg2, *pixth; + + PROCNAME("pixThresholdSpreadNorm"); + + if (ppixth) *ppixth = NULL; + if (ppixb) *ppixb = NULL; + if (ppixd) *ppixd = NULL; + if (!pixs || pixGetDepth(pixs) != 8) + return ERROR_INT("pixs not defined or not 8 bpp", procName, 1); + if (pixGetColormap(pixs)) + return ERROR_INT("pixs is colormapped", procName, 1); + if (!ppixth && !ppixb && !ppixd) + return ERROR_INT("no output requested", procName, 1); + if (filtertype != L_SOBEL_EDGE && filtertype != L_TWO_SIDED_EDGE) + return ERROR_INT("invalid filter type", procName, 1); + + /* Get the thresholded edge pixels. These are the ones + * that have values in pixs near the local optimal fg/bg threshold. */ + if (filtertype == L_SOBEL_EDGE) + pixe = pixSobelEdgeFilter(pixs, L_VERTICAL_EDGES); + else /* L_TWO_SIDED_EDGE */ + pixe = pixTwoSidedEdgeFilter(pixs, L_VERTICAL_EDGES); + pixet = pixThresholdToBinary(pixe, edgethresh); + pixInvert(pixet, pixet); + + /* Build a seed image whose only nonzero values are those + * values of pixs corresponding to pixels in the fg of pixet. */ + pixsd = pixCreateTemplate(pixs); + pixCombineMasked(pixsd, pixs, pixet); + + /* Spread the seed and optionally smooth to reduce noise */ + pixg1 = pixSeedspread(pixsd, 4); + pixg2 = pixBlockconv(pixg1, smoothx, smoothy); + + /* Optionally do a gamma enhancement */ + pixth = pixGammaTRC(NULL, pixg2, gamma, minval, maxval); + + /* Do the mapping and thresholding */ + if (ppixd) { + *ppixd = pixApplyVariableGrayMap(pixs, pixth, targetthresh); + if (ppixb) + *ppixb = pixThresholdToBinary(*ppixd, targetthresh); + } + else if (ppixb) + *ppixb = pixVarThresholdToBinary(pixs, pixth); + + if (ppixth) + *ppixth = pixth; + else + pixDestroy(&pixth); + + pixDestroy(&pixe); + pixDestroy(&pixet); + pixDestroy(&pixsd); + pixDestroy(&pixg1); + pixDestroy(&pixg2); + return 0; +} + + +/*------------------------------------------------------------------* + * Adaptive background normalization (flexible adaptaption) * + *------------------------------------------------------------------*/ +/*! + * \brief pixBackgroundNormFlex() + * + * \param[in] pixs 8 bpp grayscale; not colormapped + * \param[in] sx, sy desired tile dimensions; size may vary; + * use values between 3 and 10 + * \param[in] smoothx, smoothy half-width of convolution kernel applied to + * threshold array: use values between 1 and 3 + * \param[in] delta difference parameter in basin filling; + * use 0 to skip + * \return pixd 8 bpp, background-normalized), or NULL on error + * + * <pre> + * Notes: + * (1) This does adaptation flexibly to a quickly varying background. + * For that reason, all input parameters should be small. + * (2) sx and sy give the tile size; they should be in [5 - 7]. + * (3) The full width and height of the convolution kernel + * are (2 * smoothx + 1) and (2 * smoothy + 1). They + * should be in [1 - 2]. + * (4) Basin filling is used to fill the large fg regions. The + * parameter %delta measures the height that the black + * background is raised from the local minima. By raising + * the background, it is possible to threshold the large + * fg regions to foreground. If %delta is too large, + * bg regions will be lifted, causing thickening of + * the fg regions. Use 0 to skip. + * </pre> + */ +PIX * +pixBackgroundNormFlex(PIX *pixs, + l_int32 sx, + l_int32 sy, + l_int32 smoothx, + l_int32 smoothy, + l_int32 delta) +{ +l_float32 scalex, scaley; +PIX *pixt, *pixsd, *pixmin, *pixbg, *pixbgi, *pixd; + + PROCNAME("pixBackgroundNormFlex"); + + if (!pixs || pixGetDepth(pixs) != 8) + return (PIX *)ERROR_PTR("pixs undefined or not 8 bpp", procName, NULL); + if (pixGetColormap(pixs)) + return (PIX *)ERROR_PTR("pixs is colormapped", procName, NULL); + if (sx < 3 || sy < 3) + return (PIX *)ERROR_PTR("sx and/or sy less than 3", procName, NULL); + if (sx > 10 || sy > 10) + return (PIX *)ERROR_PTR("sx and/or sy exceed 10", procName, NULL); + if (smoothx < 1 || smoothy < 1) + return (PIX *)ERROR_PTR("smooth params less than 1", procName, NULL); + if (smoothx > 3 || smoothy > 3) + return (PIX *)ERROR_PTR("smooth params exceed 3", procName, NULL); + + /* Generate the bg estimate using smoothed average with subsampling */ + scalex = 1. / (l_float32)sx; + scaley = 1. / (l_float32)sy; + pixt = pixScaleSmooth(pixs, scalex, scaley); + + /* Do basin filling on the bg estimate if requested */ + if (delta <= 0) + pixsd = pixClone(pixt); + else { + pixLocalExtrema(pixt, 0, 0, &pixmin, NULL); + pixsd = pixSeedfillGrayBasin(pixmin, pixt, delta, 4); + pixDestroy(&pixmin); + } + pixbg = pixExtendByReplication(pixsd, 1, 1); + + /* Map the bg to 200 */ + pixbgi = pixGetInvBackgroundMap(pixbg, 200, smoothx, smoothy); + pixd = pixApplyInvBackgroundGrayMap(pixs, pixbgi, sx, sy); + + pixDestroy(&pixt); + pixDestroy(&pixsd); + pixDestroy(&pixbg); + pixDestroy(&pixbgi); + return pixd; +} + + +/*------------------------------------------------------------------* + * Adaptive contrast normalization * + *------------------------------------------------------------------*/ +/*! + * \brief pixContrastNorm() + * + * \param[in] pixd [optional] 8 bpp; null or equal to pixs + * \param[in] pixs 8 bpp grayscale; not colormapped + * \param[in] sx, sy tile dimensions + * \param[in] mindiff minimum difference to accept as valid + * \param[in] smoothx, smoothy half-width of convolution kernel applied to + * min and max arrays: use 0 for no smoothing + * \return pixd always + * + * <pre> + * Notes: + * (1) This function adaptively attempts to expand the contrast + * to the full dynamic range in each tile. If the contrast in + * a tile is smaller than %mindiff, it uses the min and max + * pixel values from neighboring tiles. It also can use + * convolution to smooth the min and max values from + * neighboring tiles. After all that processing, it is + * possible that the actual pixel values in the tile are outside + * the computed [min ... max] range for local contrast + * normalization. Such pixels are taken to be at either 0 + * (if below the min) or 255 (if above the max). + * (2) pixd can be equal to pixs (in-place operation) or + * null (makes a new pixd). + * (3) sx and sy give the tile size; they are typically at least 20. + * (4) mindiff is used to eliminate results for tiles where it is + * likely that either fg or bg is missing. A value around 50 + * or more is reasonable. + * (5) The full width and height of the convolution kernel + * are (2 * smoothx + 1) and (2 * smoothy + 1). Some smoothing + * is typically useful, and we limit the smoothing half-widths + * to the range from 0 to 8. + * (6) A linear TRC (gamma = 1.0) is applied to increase the contrast + * in each tile. The result can subsequently be globally corrected, + * by applying pixGammaTRC() with arbitrary values of gamma + * and the 0 and 255 points of the mapping. + * </pre> + */ +PIX * +pixContrastNorm(PIX *pixd, + PIX *pixs, + l_int32 sx, + l_int32 sy, + l_int32 mindiff, + l_int32 smoothx, + l_int32 smoothy) +{ +PIX *pixmin, *pixmax; + + PROCNAME("pixContrastNorm"); + + if (!pixs || pixGetDepth(pixs) != 8) + return (PIX *)ERROR_PTR("pixs undefined or not 8 bpp", procName, pixd); + if (pixd && pixd != pixs) + return (PIX *)ERROR_PTR("pixd not null or == pixs", procName, pixd); + if (pixGetColormap(pixs)) + return (PIX *)ERROR_PTR("pixs is colormapped", procName, pixd); + if (sx < 5 || sy < 5) + return (PIX *)ERROR_PTR("sx and/or sy less than 5", procName, pixd); + if (smoothx < 0 || smoothy < 0) + return (PIX *)ERROR_PTR("smooth params less than 0", procName, pixd); + if (smoothx > 8 || smoothy > 8) + return (PIX *)ERROR_PTR("smooth params exceed 8", procName, pixd); + + /* Get the min and max pixel values in each tile, and represent + * each value as a pixel in pixmin and pixmax, respectively. */ + pixMinMaxTiles(pixs, sx, sy, mindiff, smoothx, smoothy, &pixmin, &pixmax); + + /* For each tile, do a linear expansion of the dynamic range + * of pixels so that the min value is mapped to 0 and the + * max value is mapped to 255. */ + pixd = pixLinearTRCTiled(pixd, pixs, sx, sy, pixmin, pixmax); + + pixDestroy(&pixmin); + pixDestroy(&pixmax); + return pixd; +} + + +/*! + * \brief pixMinMaxTiles() + * + * \param[in] pixs 8 bpp grayscale; not colormapped + * \param[in] sx, sy tile dimensions + * \param[in] mindiff minimum difference to accept as valid + * \param[in] smoothx, smoothy half-width of convolution kernel applied to + * min and max arrays: use 0 for no smoothing + * \param[out] ppixmin tiled minima + * \param[out] ppixmax tiled maxima + * \return 0 if OK, 1 on error + * + * <pre> + * Notes: + * (1) This computes filtered and smoothed values for the min and + * max pixel values in each tile of the image. + * (2) See pixContrastNorm() for usage. + * </pre> + */ +static l_ok +pixMinMaxTiles(PIX *pixs, + l_int32 sx, + l_int32 sy, + l_int32 mindiff, + l_int32 smoothx, + l_int32 smoothy, + PIX **ppixmin, + PIX **ppixmax) +{ +l_int32 w, h; +PIX *pixmin1, *pixmax1, *pixmin2, *pixmax2; + + PROCNAME("pixMinMaxTiles"); + + if (ppixmin) *ppixmin = NULL; + if (ppixmax) *ppixmax = NULL; + if (!ppixmin || !ppixmax) + return ERROR_INT("&pixmin or &pixmax undefined", procName, 1); + if (!pixs || pixGetDepth(pixs) != 8) + return ERROR_INT("pixs undefined or not 8 bpp", procName, 1); + if (pixGetColormap(pixs)) + return ERROR_INT("pixs is colormapped", procName, 1); + if (sx < 5 || sy < 5) + return ERROR_INT("sx and/or sy less than 3", procName, 1); + if (smoothx < 0 || smoothy < 0) + return ERROR_INT("smooth params less than 0", procName, 1); + if (smoothx > 5 || smoothy > 5) + return ERROR_INT("smooth params exceed 5", procName, 1); + + /* Get the min and max values in each tile */ + pixmin1 = pixScaleGrayMinMax(pixs, sx, sy, L_CHOOSE_MIN); + pixmax1 = pixScaleGrayMinMax(pixs, sx, sy, L_CHOOSE_MAX); + + pixmin2 = pixExtendByReplication(pixmin1, 1, 1); + pixmax2 = pixExtendByReplication(pixmax1, 1, 1); + pixDestroy(&pixmin1); + pixDestroy(&pixmax1); + + /* Make sure no value is 0 */ + pixAddConstantGray(pixmin2, 1); + pixAddConstantGray(pixmax2, 1); + + /* Generate holes where the contrast is too small */ + pixSetLowContrast(pixmin2, pixmax2, mindiff); + + /* Fill the holes (0 values) */ + pixGetDimensions(pixmin2, &w, &h, NULL); + pixFillMapHoles(pixmin2, w, h, L_FILL_BLACK); + pixFillMapHoles(pixmax2, w, h, L_FILL_BLACK); + + /* Smooth if requested */ + if (smoothx > 0 || smoothy > 0) { + smoothx = L_MIN(smoothx, (w - 1) / 2); + smoothy = L_MIN(smoothy, (h - 1) / 2); + *ppixmin = pixBlockconv(pixmin2, smoothx, smoothy); + *ppixmax = pixBlockconv(pixmax2, smoothx, smoothy); + } + else { + *ppixmin = pixClone(pixmin2); + *ppixmax = pixClone(pixmax2); + } + pixCopyResolution(*ppixmin, pixs); + pixCopyResolution(*ppixmax, pixs); + pixDestroy(&pixmin2); + pixDestroy(&pixmax2); + + return 0; +} + + +/*! + * \brief pixSetLowContrast() + * + * \param[in] pixs1 8 bpp + * \param[in] pixs2 8 bpp + * \param[in] mindiff minimum difference to accept as valid + * \return 0 if OK; 1 if no pixel diffs are large enough, or on error + * + * <pre> + * Notes: + * (1) This compares corresponding pixels in pixs1 and pixs2. + * When they differ by less than %mindiff, set the pixel + * values to 0 in each. Each pixel typically represents a tile + * in a larger image, and a very small difference between + * the min and max in the tile indicates that the min and max + * values are not to be trusted. + * (2) If contrast (pixel difference) detection is expected to fail, + * caller should check return value. + * </pre> + */ +static l_ok +pixSetLowContrast(PIX *pixs1, + PIX *pixs2, + l_int32 mindiff) +{ +l_int32 i, j, w, h, d, wpl, val1, val2, found; +l_uint32 *data1, *data2, *line1, *line2; + + PROCNAME("pixSetLowContrast"); + + if (!pixs1 || !pixs2) + return ERROR_INT("pixs1 and pixs2 not both defined", procName, 1); + if (pixSizesEqual(pixs1, pixs2) == 0) + return ERROR_INT("pixs1 and pixs2 not equal size", procName, 1); + pixGetDimensions(pixs1, &w, &h, &d); + if (d != 8) + return ERROR_INT("depth not 8 bpp", procName, 1); + if (mindiff > 254) return 0; + + data1 = pixGetData(pixs1); + data2 = pixGetData(pixs2); + wpl = pixGetWpl(pixs1); + found = 0; /* init to not finding any diffs >= mindiff */ + for (i = 0; i < h; i++) { + line1 = data1 + i * wpl; + line2 = data2 + i * wpl; + for (j = 0; j < w; j++) { + val1 = GET_DATA_BYTE(line1, j); + val2 = GET_DATA_BYTE(line2, j); + if (L_ABS(val1 - val2) >= mindiff) { + found = 1; + break; + } + } + if (found) break; + } + if (!found) { + L_WARNING("no pixel pair diffs as large as mindiff\n", procName); + pixClearAll(pixs1); + pixClearAll(pixs2); + return 1; + } + + for (i = 0; i < h; i++) { + line1 = data1 + i * wpl; + line2 = data2 + i * wpl; + for (j = 0; j < w; j++) { + val1 = GET_DATA_BYTE(line1, j); + val2 = GET_DATA_BYTE(line2, j); + if (L_ABS(val1 - val2) < mindiff) { + SET_DATA_BYTE(line1, j, 0); + SET_DATA_BYTE(line2, j, 0); + } + } + } + + return 0; +} + + +/*! + * \brief pixLinearTRCTiled() + * + * \param[in] pixd [optional] 8 bpp + * \param[in] pixs 8 bpp, not colormapped + * \param[in] sx, sy tile dimensions + * \param[in] pixmin pix of min values in tiles + * \param[in] pixmax pix of max values in tiles + * \return pixd always + * + * <pre> + * Notes: + * (1) pixd can be equal to pixs (in-place operation) or + * null (makes a new pixd). + * (2) sx and sy give the tile size; they are typically at least 20. + * (3) pixmin and pixmax are generated by pixMinMaxTiles() + * (4) For each tile, this does a linear expansion of the dynamic + * range so that the min value in the tile becomes 0 and the + * max value in the tile becomes 255. + * (5) The LUTs that do the mapping are generated as needed + * and stored for reuse in an integer array within the ptr array iaa[]. + * </pre> + */ +static PIX * +pixLinearTRCTiled(PIX *pixd, + PIX *pixs, + l_int32 sx, + l_int32 sy, + PIX *pixmin, + PIX *pixmax) +{ +l_int32 i, j, k, m, w, h, wt, ht, wpl, wplt, xoff, yoff; +l_int32 minval, maxval, val, sval; +l_int32 *ia; +l_int32 **iaa; +l_uint32 *data, *datamin, *datamax, *line, *tline, *linemin, *linemax; + + PROCNAME("pixLinearTRCTiled"); + + if (!pixs || pixGetDepth(pixs) != 8) + return (PIX *)ERROR_PTR("pixs undefined or not 8 bpp", procName, pixd); + if (pixd && pixd != pixs) + return (PIX *)ERROR_PTR("pixd not null or == pixs", procName, pixd); + if (pixGetColormap(pixs)) + return (PIX *)ERROR_PTR("pixs is colormapped", procName, pixd); + if (!pixmin || !pixmax) + return (PIX *)ERROR_PTR("pixmin & pixmax not defined", procName, pixd); + if (sx < 5 || sy < 5) + return (PIX *)ERROR_PTR("sx and/or sy less than 5", procName, pixd); + + if ((iaa = (l_int32 **)LEPT_CALLOC(256, sizeof(l_int32 *))) == NULL) + return (PIX *)ERROR_PTR("iaa not made", procName, NULL); + if ((pixd = pixCopy(pixd, pixs)) == NULL) { + LEPT_FREE(iaa); + return (PIX *)ERROR_PTR("pixd not made", procName, NULL); + } + pixGetDimensions(pixd, &w, &h, NULL); + + data = pixGetData(pixd); + wpl = pixGetWpl(pixd); + datamin = pixGetData(pixmin); + datamax = pixGetData(pixmax); + wplt = pixGetWpl(pixmin); + pixGetDimensions(pixmin, &wt, &ht, NULL); + for (i = 0; i < ht; i++) { + line = data + sy * i * wpl; + linemin = datamin + i * wplt; + linemax = datamax + i * wplt; + yoff = sy * i; + for (j = 0; j < wt; j++) { + xoff = sx * j; + minval = GET_DATA_BYTE(linemin, j); + maxval = GET_DATA_BYTE(linemax, j); + if (maxval == minval) { + L_ERROR("shouldn't happen! i,j = %d,%d, minval = %d\n", + procName, i, j, minval); + continue; + } + if ((ia = iaaGetLinearTRC(iaa, maxval - minval)) == NULL) { + L_ERROR("failure to make ia for j = %d!\n", procName, j); + continue; + } + for (k = 0; k < sy && yoff + k < h; k++) { + tline = line + k * wpl; + for (m = 0; m < sx && xoff + m < w; m++) { + val = GET_DATA_BYTE(tline, xoff + m); + sval = val - minval; + sval = L_MAX(0, sval); + SET_DATA_BYTE(tline, xoff + m, ia[sval]); + } + } + } + } + + for (i = 0; i < 256; i++) + LEPT_FREE(iaa[i]); + LEPT_FREE(iaa); + return pixd; +} + + +/*! + * \brief iaaGetLinearTRC() + * + * \param[in] iaa bare array of ptrs to l_int32 + * \param[in] diff between min and max pixel values that are + * to be mapped to 0 and 255 + * \return ia LUT with input (val - minval) and output a + * value between 0 and 255) + */ +static l_int32 * +iaaGetLinearTRC(l_int32 **iaa, + l_int32 diff) +{ +l_int32 i; +l_int32 *ia; +l_float32 factor; + + PROCNAME("iaaGetLinearTRC"); + + if (!iaa) + return (l_int32 *)ERROR_PTR("iaa not defined", procName, NULL); + + if (iaa[diff] != NULL) /* already have it */ + return iaa[diff]; + + if ((ia = (l_int32 *)LEPT_CALLOC(256, sizeof(l_int32))) == NULL) + return (l_int32 *)ERROR_PTR("ia not made", procName, NULL); + iaa[diff] = ia; + if (diff == 0) { /* shouldn't happen */ + for (i = 0; i < 256; i++) + ia[i] = 128; + } + else { + factor = 255. / (l_float32)diff; + for (i = 0; i < diff + 1; i++) + ia[i] = (l_int32)(factor * i + 0.5); + for (i = diff + 1; i < 256; i++) + ia[i] = 255; + } + + return ia; +} |