Mercurial > hgrepos > Python2 > PyMuPDF

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

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
ADD: MuPDF v1.26.7: the MuPDF source as downloaded by a default build of PyMuPDF 1.26.4. The directory name has changed: no version number in the expanded directory now.
author Franz Glasner <fzglas.hg@dom66.de>
date Mon, 15 Sep 2025 11:43:07 +0200
parents
children
comparison
equal deleted inserted replaced
1:1d09e1dec1d9 2:b50eed0cc0ef
1 /*====================================================================*
2 - Copyright (C) 2001 Leptonica. All rights reserved.
3 -
4 - Redistribution and use in source and binary forms, with or without
5 - modification, are permitted provided that the following conditions
6 - are met:
7 - 1. Redistributions of source code must retain the above copyright
8 - notice, this list of conditions and the following disclaimer.
9 - 2. Redistributions in binary form must reproduce the above
10 - copyright notice, this list of conditions and the following
11 - disclaimer in the documentation and/or other materials
12 - provided with the distribution.
13 -
14 - THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
15 - ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
16 - LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
17 - A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL ANY
18 - CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
19 - EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
20 - PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
21 - PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY
22 - OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
23 - NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
24 - SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
25 *====================================================================*/
26
27 /*!
28 * \file enhance.c
29 * <pre>
30 *
31 * Gamma TRC (tone reproduction curve) mapping
32 * PIX *pixGammaTRC()
33 * PIX *pixGammaTRCMasked()
34 * PIX *pixGammaTRCWithAlpha()
35 * NUMA *numaGammaTRC()
36 *
37 * Contrast enhancement
38 * PIX *pixContrastTRC()
39 * PIX *pixContrastTRCMasked()
40 * NUMA *numaContrastTRC()
41 *
42 * Histogram equalization
43 * PIX *pixEqualizeTRC()
44 * NUMA *numaEqualizeTRC()
45 *
46 * Generic TRC mapper
47 * l_int32 pixTRCMap()
48 * l_int32 pixTRCMapGeneral()
49 *
50 * Unsharp-masking
51 * PIX *pixUnsharpMasking()
52 * PIX *pixUnsharpMaskingGray()
53 * PIX *pixUnsharpMaskingFast()
54 * PIX *pixUnsharpMaskingGrayFast()
55 * PIX *pixUnsharpMaskingGray1D()
56 * PIX *pixUnsharpMaskingGray2D()
57 *
58 * Hue and saturation modification
59 * PIX *pixModifyHue()
60 * PIX *pixModifySaturation()
61 * l_int32 pixMeasureSaturation()
62 * PIX *pixModifyBrightness()
63 *
64 * Color shifting
65 * PIX *pixMosaicColorShiftRGB()
66 * PIX *pixColorShiftRGB()
67 *
68 * Darken gray (unsaturated) pixels
69 * PIX *pixDarkenGray()
70 *
71 * General multiplicative constant color transform
72 * PIX *pixMultConstantColor()
73 * PIX *pixMultMatrixColor()
74 *
75 * Edge by bandpass
76 * PIX *pixHalfEdgeByBandpass()
77 *
78 * Gamma correction, contrast enhancement and histogram equalization
79 * apply a simple mapping function to each pixel (or, for color
80 * images, to each sample (i.e., r,g,b) of the pixel).
81 *
82 * ~ Gamma correction either lightens the image or darkens
83 * it, depending on whether the gamma factor is greater
84 * or less than 1.0, respectively.
85 *
86 * ~ Contrast enhancement darkens the pixels that are already
87 * darker than the middle of the dynamic range (128)
88 * and lightens pixels that are lighter than 128.
89 *
90 * ~ Histogram equalization remaps to have the same number
91 * of image pixels at each of 256 intensity values. This is
92 * a quick and dirty method of adjusting contrast and brightness
93 * to bring out details in both light and dark regions.
94 *
95 * Unsharp masking is a more complicated enhancement.
96 * A "high frequency" image, generated by subtracting
97 * the smoothed ("low frequency") part of the image from
98 * itself, has all the energy at the edges. This "edge image"
99 * has 0 average value. A fraction of the edge image is
100 * then added to the original, enhancing the differences
101 * between pixel values at edges. Because we represent
102 * images as l_uint8 arrays, we preserve dynamic range and
103 * handle negative values by doing all the arithmetic on
104 * shifted l_uint16 arrays; the l_uint8 values are recovered
105 * at the end.
106 *
107 * Hue and saturation modification work in HSV space. Because
108 * this is too large for efficient table lookup, each pixel value
109 * is transformed to HSV, modified, and transformed back.
110 * It's not the fastest way to do this, but the method is
111 * easily understood.
112 *
113 * Unsharp masking is never in-place, and returns a clone if no
114 * operation is to be performed.
115 * </pre>
116 */
117
118 #ifdef HAVE_CONFIG_H
119 #include <config_auto.h>
120 #endif /* HAVE_CONFIG_H */
121
122 #include <math.h>
123 #include "allheaders.h"
124
125 /* Scales contrast enhancement factor to have a useful range
126 * between 0.0 and 1.0 */
127 static const l_float32 EnhanceScaleFactor = 5.0;
128
129 /*-------------------------------------------------------------*
130 * Gamma TRC (tone reproduction curve) mapping *
131 *-------------------------------------------------------------*/
132 /*!
133 * \brief pixGammaTRC()
134 *
135 * \param[in] pixd [optional] null or equal to pixs
136 * \param[in] pixs 8 or 32 bpp; or 2, 4 or 8 bpp with colormap
137 * \param[in] gamma gamma correction; must be > 0.0
138 * \param[in] minval input value that gives 0 for output; can be < 0
139 * \param[in] maxval input value that gives 255 for output; can be > 255
140 * \return pixd always
141 *
142 * <pre>
143 * Notes:
144 * (1) pixd must either be null or equal to pixs.
145 * For in-place operation, set pixd == pixs:
146 * pixGammaTRC(pixs, pixs, ...);
147 * To get a new image, set pixd == null:
148 * pixd = pixGammaTRC(NULL, pixs, ...);
149 * (2) If pixs is colormapped, the colormap is transformed,
150 * either in-place or in a copy of pixs.
151 * (3) We use a gamma mapping between minval and maxval.
152 * (4) If gamma < 1.0, the image will appear darker;
153 * if gamma > 1.0, the image will appear lighter;
154 * (5) If gamma = 1.0 and minval = 0 and maxval = 255, no
155 * enhancement is performed; return a copy unless in-place,
156 * in which case this is a no-op.
157 * (6) For color images that are not colormapped, the mapping
158 * is applied to each component.
159 * (7) minval and maxval are not restricted to the interval [0, 255].
160 * If minval < 0, an input value of 0 is mapped to a
161 * nonzero output. This will turn black to gray.
162 * If maxval > 255, an input value of 255 is mapped to
163 * an output value less than 255. This will turn
164 * white (e.g., in the background) to gray.
165 * (8) Increasing minval darkens the image.
166 * (9) Decreasing maxval bleaches the image.
167 * (10) Simultaneously increasing minval and decreasing maxval
168 * will darken the image and make the colors more intense;
169 * e.g., minval = 50, maxval = 200.
170 * (11) See numaGammaTRC() for further examples of use.
171 * (12) Use pixTRCMapGeneral() if applying different mappings
172 * to each channel in an RGB image.
173 * </pre>
174 */
175 PIX *
176 pixGammaTRC(PIX *pixd,
177 PIX *pixs,
178 l_float32 gamma,
179 l_int32 minval,
180 l_int32 maxval)
181 {
182 l_int32 d;
183 NUMA *nag;
184 PIXCMAP *cmap;
185
186 if (!pixs)
187 return (PIX *)ERROR_PTR("pixs not defined", __func__, pixd);
188 if (pixd && (pixd != pixs))
189 return (PIX *)ERROR_PTR("pixd not null or pixs", __func__, pixd);
190 if (gamma <= 0.0) {
191 L_WARNING("gamma must be > 0.0; setting to 1.0\n", __func__);
192 gamma = 1.0;
193 }
194 if (minval >= maxval)
195 return (PIX *)ERROR_PTR("minval not < maxval", __func__, pixd);
196 cmap = pixGetColormap(pixs);
197 d = pixGetDepth(pixs);
198 if (!cmap && d != 8 && d != 32)
199 return (PIX *)ERROR_PTR("depth not 8 or 32 bpp", __func__, pixd);
200
201 if (gamma == 1.0 && minval == 0 && maxval == 255) /* no-op */
202 return pixCopy(pixd, pixs);
203
204 if (!pixd) /* start with a copy if not in-place */
205 pixd = pixCopy(NULL, pixs);
206
207 if (cmap) {
208 pixcmapGammaTRC(pixGetColormap(pixd), gamma, minval, maxval);
209 return pixd;
210 }
211
212 /* pixd is 8 or 32 bpp */
213 if ((nag = numaGammaTRC(gamma, minval, maxval)) == NULL)
214 return (PIX *)ERROR_PTR("nag not made", __func__, pixd);
215 pixTRCMap(pixd, NULL, nag);
216 numaDestroy(&nag);
217
218 return pixd;
219 }
220
221
222 /*!
223 * \brief pixGammaTRCMasked()
224 *
225 * \param[in] pixd [optional] null or equal to pixs
226 * \param[in] pixs 8 or 32 bpp; not colormapped
227 * \param[in] pixm [optional] null or 1 bpp
228 * \param[in] gamma gamma correction; must be > 0.0
229 * \param[in] minval input value that gives 0 for output; can be < 0
230 * \param[in] maxval input value that gives 255 for output; can be > 255
231 * \return pixd always
232 *
233 * <pre>
234 * Notes:
235 * (1) Same as pixGammaTRC() except mapping is optionally over
236 * a subset of pixels described by pixm.
237 * (2) Masking does not work for colormapped images.
238 * (3) See pixGammaTRC() for details on how to use the parameters.
239 * </pre>
240 */
241 PIX *
242 pixGammaTRCMasked(PIX *pixd,
243 PIX *pixs,
244 PIX *pixm,
245 l_float32 gamma,
246 l_int32 minval,
247 l_int32 maxval)
248 {
249 l_int32 d;
250 NUMA *nag;
251
252 if (!pixm)
253 return pixGammaTRC(pixd, pixs, gamma, minval, maxval);
254
255 if (!pixs)
256 return (PIX *)ERROR_PTR("pixs not defined", __func__, pixd);
257 if (pixGetColormap(pixs))
258 return (PIX *)ERROR_PTR("invalid: pixs has a colormap", __func__, pixd);
259 if (pixd && (pixd != pixs))
260 return (PIX *)ERROR_PTR("pixd not null or pixs", __func__, pixd);
261 d = pixGetDepth(pixs);
262 if (d != 8 && d != 32)
263 return (PIX *)ERROR_PTR("depth not 8 or 32 bpp", __func__, pixd);
264 if (minval >= maxval)
265 return (PIX *)ERROR_PTR("minval not < maxval", __func__, pixd);
266 if (gamma <= 0.0) {
267 L_WARNING("gamma must be > 0.0; setting to 1.0\n", __func__);
268 gamma = 1.0;
269 }
270
271 if (gamma == 1.0 && minval == 0 && maxval == 255)
272 return pixCopy(pixd, pixs);
273
274 if (!pixd) /* start with a copy if not in-place */
275 pixd = pixCopy(NULL, pixs);
276
277 if ((nag = numaGammaTRC(gamma, minval, maxval)) == NULL)
278 return (PIX *)ERROR_PTR("nag not made", __func__, pixd);
279 pixTRCMap(pixd, pixm, nag);
280 numaDestroy(&nag);
281
282 return pixd;
283 }
284
285
286 /*!
287 * \brief pixGammaTRCWithAlpha()
288 *
289 * \param[in] pixd [optional] null or equal to pixs
290 * \param[in] pixs 32 bpp
291 * \param[in] gamma gamma correction; must be > 0.0
292 * \param[in] minval input value that gives 0 for output; can be < 0
293 * \param[in] maxval input value that gives 255 for output; can be > 255
294 * \return pixd always
295 *
296 * <pre>
297 * Notes:
298 * (1) See usage notes in pixGammaTRC().
299 * (2) This version saves the alpha channel. It is only valid
300 * for 32 bpp (no colormap), and is a bit slower.
301 * </pre>
302 */
303 PIX *
304 pixGammaTRCWithAlpha(PIX *pixd,
305 PIX *pixs,
306 l_float32 gamma,
307 l_int32 minval,
308 l_int32 maxval)
309 {
310 NUMA *nag;
311 PIX *pixalpha;
312
313 if (!pixs || pixGetDepth(pixs) != 32)
314 return (PIX *)ERROR_PTR("pixs undefined or not 32 bpp", __func__, pixd);
315 if (pixd && (pixd != pixs))
316 return (PIX *)ERROR_PTR("pixd not null or pixs", __func__, pixd);
317 if (gamma <= 0.0) {
318 L_WARNING("gamma must be > 0.0; setting to 1.0\n", __func__);
319 gamma = 1.0;
320 }
321 if (minval >= maxval)
322 return (PIX *)ERROR_PTR("minval not < maxval", __func__, pixd);
323
324 if (gamma == 1.0 && minval == 0 && maxval == 255)
325 return pixCopy(pixd, pixs);
326 if (!pixd) /* start with a copy if not in-place */
327 pixd = pixCopy(NULL, pixs);
328
329 pixalpha = pixGetRGBComponent(pixs, L_ALPHA_CHANNEL); /* save */
330 if ((nag = numaGammaTRC(gamma, minval, maxval)) == NULL)
331 return (PIX *)ERROR_PTR("nag not made", __func__, pixd);
332 pixTRCMap(pixd, NULL, nag);
333 pixSetRGBComponent(pixd, pixalpha, L_ALPHA_CHANNEL); /* restore */
334 pixSetSpp(pixd, 4);
335
336 numaDestroy(&nag);
337 pixDestroy(&pixalpha);
338 return pixd;
339 }
340
341
342 /*!
343 * \brief numaGammaTRC()
344 *
345 * \param[in] gamma gamma factor; must be > 0.0
346 * \param[in] minval input value that gives 0 for output
347 * \param[in] maxval input value that gives 255 for output
348 * \return na, or NULL on error
349 *
350 * <pre>
351 * Notes:
352 * (1) The map is returned as a numa; values are clipped to [0, 255].
353 * (2) For a linear mapping, set gamma = 1.0.
354 * (3) To force all intensities into a range within fraction delta
355 * of white, use: minval = -256 * (1 - delta) / delta
356 * maxval = 255
357 * (4) To force all intensities into a range within fraction delta
358 * of black, use: minval = 0
359 * maxval = 256 * (1 - delta) / delta
360 * </pre>
361 */
362 NUMA *
363 numaGammaTRC(l_float32 gamma,
364 l_int32 minval,
365 l_int32 maxval)
366 {
367 l_int32 i, val;
368 l_float32 x, invgamma;
369 NUMA *na;
370
371 if (minval >= maxval)
372 return (NUMA *)ERROR_PTR("minval not < maxval", __func__, NULL);
373 if (gamma <= 0.0) {
374 L_WARNING("gamma must be > 0.0; setting to 1.0\n", __func__);
375 gamma = 1.0;
376 }
377
378 invgamma = 1. / gamma;
379 na = numaCreate(256);
380 for (i = 0; i < minval; i++)
381 numaAddNumber(na, 0);
382 for (i = minval; i <= maxval; i++) {
383 if (i < 0) continue;
384 if (i > 255) continue;
385 x = (l_float32)(i - minval) / (l_float32)(maxval - minval);
386 val = (l_int32)(255. * powf(x, invgamma) + 0.5);
387 val = L_MAX(val, 0);
388 val = L_MIN(val, 255);
389 numaAddNumber(na, val);
390 }
391 for (i = maxval + 1; i < 256; i++)
392 numaAddNumber(na, 255);
393
394 return na;
395 }
396
397
398 /*-------------------------------------------------------------*
399 * Contrast enhancement *
400 *-------------------------------------------------------------*/
401 /*!
402 * \brief pixContrastTRC()
403 *
404 * \param[in] pixd [optional] null or equal to pixs
405 * \param[in] pixs 8 or 32 bpp; or 2, 4 or 8 bpp with colormap
406 * \param[in] factor 0.0 is no enhancement
407 * \return pixd always
408 *
409 * <pre>
410 * Notes:
411 * (1) pixd must either be null or equal to pixs.
412 * For in-place operation, set pixd == pixs:
413 * pixContrastTRC(pixs, pixs, ...);
414 * To get a new image, set pixd == null:
415 * pixd = pixContrastTRC(NULL, pixs, ...);
416 * (2) If pixs is colormapped, the colormap is transformed,
417 * either in-place or in a copy of pixs.
418 * (3) Contrast is enhanced by mapping each color component
419 * using an atan function with maximum slope at 127.
420 * Pixels below 127 are lowered in intensity and pixels
421 * above 127 are increased.
422 * (4) The useful range for the contrast factor is scaled to
423 * be in (0.0 to 1.0), but larger values can also be used.
424 * (5) If factor == 0.0, no enhancement is performed; return a copy
425 * unless in-place, in which case this is a no-op.
426 * (6) For color images that are not colormapped, the mapping
427 * is applied to each component.
428 * </pre>
429 */
430 PIX *
431 pixContrastTRC(PIX *pixd,
432 PIX *pixs,
433 l_float32 factor)
434 {
435 l_int32 d;
436 NUMA *nac;
437 PIXCMAP *cmap;
438
439 if (!pixs)
440 return (PIX *)ERROR_PTR("pixs not defined", __func__, pixd);
441 if (pixd && (pixd != pixs))
442 return (PIX *)ERROR_PTR("pixd not null or pixs", __func__, pixd);
443 if (factor < 0.0) {
444 L_WARNING("factor must be >= 0.0; using 0.0\n", __func__);
445 factor = 0.0;
446 }
447 if (factor == 0.0)
448 return pixCopy(pixd, pixs);
449
450 cmap = pixGetColormap(pixs);
451 d = pixGetDepth(pixs);
452 if (!cmap && d != 8 && d != 32)
453 return (PIX *)ERROR_PTR("depth not 8 or 32 bpp", __func__, pixd);
454
455 if (!pixd) /* start with a copy if not in-place */
456 pixd = pixCopy(NULL, pixs);
457
458 if (cmap) {
459 pixcmapContrastTRC(pixGetColormap(pixd), factor);
460 return pixd;
461 }
462
463 /* pixd is 8 or 32 bpp */
464 if ((nac = numaContrastTRC(factor)) == NULL)
465 return (PIX *)ERROR_PTR("nac not made", __func__, pixd);
466 pixTRCMap(pixd, NULL, nac);
467 numaDestroy(&nac);
468
469 return pixd;
470 }
471
472
473 /*!
474 * \brief pixContrastTRCMasked()
475 *
476 * \param[in] pixd [optional] null or equal to pixs
477 * \param[in] pixs 8 or 32 bpp; or 2, 4 or 8 bpp with colormap
478 * \param[in] pixm [optional] null or 1 bpp
479 * \param[in] factor 0.0 is no enhancement
480 * \return pixd always
481 *
482 * <pre>
483 * Notes:
484 * (1) Same as pixContrastTRC() except mapping is optionally over
485 * a subset of pixels described by pixm.
486 * (2) Masking does not work for colormapped images.
487 * (3) See pixContrastTRC() for details on how to use the parameters.
488 * </pre>
489 */
490 PIX *
491 pixContrastTRCMasked(PIX *pixd,
492 PIX *pixs,
493 PIX *pixm,
494 l_float32 factor)
495 {
496 l_int32 d;
497 NUMA *nac;
498
499 if (!pixm)
500 return pixContrastTRC(pixd, pixs, factor);
501
502 if (!pixs)
503 return (PIX *)ERROR_PTR("pixs not defined", __func__, pixd);
504 if (pixGetColormap(pixs))
505 return (PIX *)ERROR_PTR("invalid: pixs has a colormap", __func__, pixd);
506 if (pixd && (pixd != pixs))
507 return (PIX *)ERROR_PTR("pixd not null or pixs", __func__, pixd);
508 d = pixGetDepth(pixs);
509 if (d != 8 && d != 32)
510 return (PIX *)ERROR_PTR("depth not 8 or 32 bpp", __func__, pixd);
511
512 if (factor < 0.0) {
513 L_WARNING("factor must be >= 0.0; using 0.0\n", __func__);
514 factor = 0.0;
515 }
516 if (factor == 0.0)
517 return pixCopy(pixd, pixs);
518
519 if (!pixd) /* start with a copy if not in-place */
520 pixd = pixCopy(NULL, pixs);
521
522 if ((nac = numaContrastTRC(factor)) == NULL)
523 return (PIX *)ERROR_PTR("nac not made", __func__, pixd);
524 pixTRCMap(pixd, pixm, nac);
525 numaDestroy(&nac);
526
527 return pixd;
528 }
529
530
531 /*!
532 * \brief numaContrastTRC()
533 *
534 * \param[in] factor generally between 0.0 [no enhancement]
535 * and 1.0, but can be larger than 1.0
536 * \return na, or NULL on error
537 *
538 * <pre>
539 * Notes:
540 * (1) The mapping is monotonic increasing, where 0 is mapped
541 * to 0 and 255 is mapped to 255.
542 * (2) As 'factor' is increased from 0.0 (where the mapping is linear),
543 * the map gets closer to its limit as a step function that
544 * jumps from 0 to 255 at the center (input value = 127).
545 * </pre>
546 */
547 NUMA *
548 numaContrastTRC(l_float32 factor)
549 {
550 l_int32 i, val;
551 l_float64 x, ymax, ymin, dely, scale;
552 NUMA *na;
553
554 if (factor < 0.0) {
555 L_WARNING("factor must be >= 0.0; using 0.0; no enhancement\n",
556 __func__);
557 factor = 0.0;
558 }
559 if (factor == 0.0)
560 return numaMakeSequence(0, 1, 256); /* linear map */
561
562 scale = EnhanceScaleFactor;
563 ymax = atan((l_float64)(1.0 * factor * scale));
564 ymin = atan((l_float64)(-127. * factor * scale / 128.));
565 dely = ymax - ymin;
566 na = numaCreate(256);
567 for (i = 0; i < 256; i++) {
568 x = (l_float64)i;
569 val = (l_int32)((255. / dely) *
570 (-ymin + atan((l_float64)(factor * scale * (x - 127.) / 128.))) +
571 0.5);
572 numaAddNumber(na, val);
573 }
574
575 return na;
576 }
577
578
579 /*-------------------------------------------------------------*
580 * Histogram equalization *
581 *-------------------------------------------------------------*/
582 /*!
583 * \brief pixEqualizeTRC()
584 *
585 * \param[in] pixd [optional] null or equal to pixs
586 * \param[in] pixs 8 bpp gray, 32 bpp rgb, or colormapped
587 * \param[in] fract fraction of equalization movement of pixel values
588 * \param[in] factor subsampling factor; integer >= 1
589 * \return pixd, or NULL on error
590 *
591 * <pre>
592 * Notes:
593 * (1) pixd must either be null or equal to pixs.
594 * For in-place operation, set pixd == pixs:
595 * pixEqualizeTRC(pixs, pixs, ...);
596 * To get a new image, set pixd == null:
597 * pixd = pixEqualizeTRC(NULL, pixs, ...);
598 * (2) In histogram equalization, a tone reproduction curve
599 * mapping is used to make the number of pixels at each
600 * intensity equal.
601 * (3) If fract == 0.0, no equalization is performed; return a copy
602 * unless in-place, in which case this is a no-op.
603 * If fract == 1.0, equalization is complete.
604 * (4) Set the subsampling factor > 1 to reduce the amount of computation.
605 * (5) If pixs is colormapped, the colormap is removed and
606 * converted to rgb or grayscale.
607 * (6) If pixs has color, equalization is done in each channel
608 * separately.
609 * (7) Note that even if there is a colormap, we can get an
610 * in-place operation because the intermediate image pixt
611 * is copied back to pixs (which for in-place is the same
612 * as pixd).
613 * </pre>
614 */
615 PIX *
616 pixEqualizeTRC(PIX *pixd,
617 PIX *pixs,
618 l_float32 fract,
619 l_int32 factor)
620 {
621 l_int32 d;
622 NUMA *na;
623 PIX *pixt, *pix8;
624 PIXCMAP *cmap;
625
626 if (!pixs)
627 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
628 if (pixd && (pixd != pixs))
629 return (PIX *)ERROR_PTR("pixd not null or pixs", __func__, pixd);
630 cmap = pixGetColormap(pixs);
631 d = pixGetDepth(pixs);
632 if (d != 8 && d != 32 && !cmap)
633 return (PIX *)ERROR_PTR("pixs not 8/32 bpp or cmapped", __func__, NULL);
634 if (fract < 0.0 || fract > 1.0)
635 return (PIX *)ERROR_PTR("fract not in [0.0 ... 1.0]", __func__, NULL);
636 if (factor < 1)
637 return (PIX *)ERROR_PTR("sampling factor < 1", __func__, NULL);
638
639 if (fract == 0.0)
640 return pixCopy(pixd, pixs);
641
642 /* If there is a colormap, remove it. */
643 if (cmap)
644 pixt = pixRemoveColormap(pixs, REMOVE_CMAP_BASED_ON_SRC);
645 else
646 pixt = pixClone(pixs);
647
648 /* Make a copy if necessary */
649 pixd = pixCopy(pixd, pixt);
650 pixDestroy(&pixt);
651
652 d = pixGetDepth(pixd);
653 if (d == 8) {
654 na = numaEqualizeTRC(pixd, fract, factor);
655 pixTRCMap(pixd, NULL, na);
656 numaDestroy(&na);
657 } else { /* 32 bpp */
658 pix8 = pixGetRGBComponent(pixd, COLOR_RED);
659 na = numaEqualizeTRC(pix8, fract, factor);
660 pixTRCMap(pix8, NULL, na);
661 pixSetRGBComponent(pixd, pix8, COLOR_RED);
662 numaDestroy(&na);
663 pixDestroy(&pix8);
664 pix8 = pixGetRGBComponent(pixd, COLOR_GREEN);
665 na = numaEqualizeTRC(pix8, fract, factor);
666 pixTRCMap(pix8, NULL, na);
667 pixSetRGBComponent(pixd, pix8, COLOR_GREEN);
668 numaDestroy(&na);
669 pixDestroy(&pix8);
670 pix8 = pixGetRGBComponent(pixd, COLOR_BLUE);
671 na = numaEqualizeTRC(pix8, fract, factor);
672 pixTRCMap(pix8, NULL, na);
673 pixSetRGBComponent(pixd, pix8, COLOR_BLUE);
674 numaDestroy(&na);
675 pixDestroy(&pix8);
676 }
677
678 return pixd;
679 }
680
681
682 /*!
683 * \brief numaEqualizeTRC()
684 *
685 * \param[in] pix 8 bpp, no colormap
686 * \param[in] fract fraction of equalization movement of pixel values
687 * \param[in] factor subsampling factor; integer >= 1
688 * \return nad, or NULL on error
689 *
690 * <pre>
691 * Notes:
692 * (1) If fract == 0.0, no equalization will be performed.
693 * If fract == 1.0, equalization is complete.
694 * (2) Set the subsampling factor > 1 to reduce the amount of computation.
695 * (3) The map is returned as a numa with 256 values, specifying
696 * the equalized value (array value) for every input value
697 * (the array index).
698 * </pre>
699 */
700 NUMA *
701 numaEqualizeTRC(PIX *pix,
702 l_float32 fract,
703 l_int32 factor)
704 {
705 l_int32 iin, iout, itarg;
706 l_float32 val, sum;
707 NUMA *nah, *nasum, *nad;
708
709 if (!pix)
710 return (NUMA *)ERROR_PTR("pix not defined", __func__, NULL);
711 if (pixGetDepth(pix) != 8)
712 return (NUMA *)ERROR_PTR("pix not 8 bpp", __func__, NULL);
713 if (fract < 0.0 || fract > 1.0)
714 return (NUMA *)ERROR_PTR("fract not in [0.0 ... 1.0]", __func__, NULL);
715 if (factor < 1)
716 return (NUMA *)ERROR_PTR("sampling factor < 1", __func__, NULL);
717
718 if (fract == 0.0)
719 L_WARNING("fract = 0.0; no equalization requested\n", __func__);
720
721 if ((nah = pixGetGrayHistogram(pix, factor)) == NULL)
722 return (NUMA *)ERROR_PTR("histogram not made", __func__, NULL);
723 numaGetSum(nah, &sum);
724 nasum = numaGetPartialSums(nah);
725
726 nad = numaCreate(256);
727 for (iin = 0; iin < 256; iin++) {
728 numaGetFValue(nasum, iin, &val);
729 itarg = (l_int32)(255. * val / sum + 0.5);
730 iout = iin + (l_int32)(fract * (itarg - iin));
731 iout = L_MIN(iout, 255); /* to be safe */
732 numaAddNumber(nad, iout);
733 }
734
735 numaDestroy(&nah);
736 numaDestroy(&nasum);
737 return nad;
738 }
739
740
741 /*-------------------------------------------------------------*
742 * Generic TRC mapping *
743 *-------------------------------------------------------------*/
744 /*!
745 * \brief pixTRCMap()
746 *
747 * \param[in] pixs 8 grayscale or 32 bpp rgb; not colormapped
748 * \param[in] pixm [optional] 1 bpp mask
749 * \param[in] na mapping array
750 * \return 0 if OK, 1 on error
751 *
752 * <pre>
753 * Notes:
754 * (1) This operation is in-place on pixs.
755 * (2) For 32 bpp, this applies the same map to each of the r,g,b
756 * components.
757 * (3) The mapping array is of size 256, and it maps the input
758 * index into values in the range [0, 255].
759 * (4) If defined, the optional 1 bpp mask pixm has its origin
760 * aligned with pixs, and the map function is applied only
761 * to pixels in pixs under the fg of pixm.
762 * (5) For 32 bpp, this does not save the alpha channel.
763 * </pre>
764 */
765 l_int32
766 pixTRCMap(PIX *pixs,
767 PIX *pixm,
768 NUMA *na)
769 {
770 l_int32 w, h, d, wm, hm, wpl, wplm, i, j, sval8, dval8;
771 l_uint32 sval32, dval32;
772 l_uint32 *data, *datam, *line, *linem, *tab;
773
774 if (!pixs)
775 return ERROR_INT("pixs not defined", __func__, 1);
776 if (pixGetColormap(pixs))
777 return ERROR_INT("pixs is colormapped", __func__, 1);
778 if (!na)
779 return ERROR_INT("na not defined", __func__, 1);
780 if (numaGetCount(na) != 256)
781 return ERROR_INT("na not of size 256", __func__, 1);
782 pixGetDimensions(pixs, &w, &h, &d);
783 if (d != 8 && d != 32)
784 return ERROR_INT("pixs not 8 or 32 bpp", __func__, 1);
785 if (pixm) {
786 if (pixGetDepth(pixm) != 1)
787 return ERROR_INT("pixm not 1 bpp", __func__, 1);
788 }
789
790 tab = (l_uint32 *)numaGetIArray(na); /* get the array for efficiency */
791 wpl = pixGetWpl(pixs);
792 data = pixGetData(pixs);
793 if (!pixm) {
794 if (d == 8) {
795 for (i = 0; i < h; i++) {
796 line = data + i * wpl;
797 for (j = 0; j < w; j++) {
798 sval8 = GET_DATA_BYTE(line, j);
799 dval8 = tab[sval8];
800 SET_DATA_BYTE(line, j, dval8);
801 }
802 }
803 } else { /* d == 32 */
804 for (i = 0; i < h; i++) {
805 line = data + i * wpl;
806 for (j = 0; j < w; j++) {
807 sval32 = *(line + j);
808 dval32 =
809 tab[(sval32 >> L_RED_SHIFT) & 0xff] << L_RED_SHIFT |
810 tab[(sval32 >> L_GREEN_SHIFT) & 0xff] << L_GREEN_SHIFT |
811 tab[(sval32 >> L_BLUE_SHIFT) & 0xff] << L_BLUE_SHIFT;
812 *(line + j) = dval32;
813 }
814 }
815 }
816 } else {
817 datam = pixGetData(pixm);
818 wplm = pixGetWpl(pixm);
819 pixGetDimensions(pixm, &wm, &hm, NULL);
820 if (d == 8) {
821 for (i = 0; i < h; i++) {
822 if (i >= hm)
823 break;
824 line = data + i * wpl;
825 linem = datam + i * wplm;
826 for (j = 0; j < w; j++) {
827 if (j >= wm)
828 break;
829 if (GET_DATA_BIT(linem, j) == 0)
830 continue;
831 sval8 = GET_DATA_BYTE(line, j);
832 dval8 = tab[sval8];
833 SET_DATA_BYTE(line, j, dval8);
834 }
835 }
836 } else { /* d == 32 */
837 for (i = 0; i < h; i++) {
838 if (i >= hm)
839 break;
840 line = data + i * wpl;
841 linem = datam + i * wplm;
842 for (j = 0; j < w; j++) {
843 if (j >= wm)
844 break;
845 if (GET_DATA_BIT(linem, j) == 0)
846 continue;
847 sval32 = *(line + j);
848 dval32 =
849 tab[(sval32 >> L_RED_SHIFT) & 0xff] << L_RED_SHIFT |
850 tab[(sval32 >> L_GREEN_SHIFT) & 0xff] << L_GREEN_SHIFT |
851 tab[(sval32 >> L_BLUE_SHIFT) & 0xff] << L_BLUE_SHIFT;
852 *(line + j) = dval32;
853 }
854 }
855 }
856 }
857
858 LEPT_FREE(tab);
859 return 0;
860 }
861
862
863 /*!
864 * \brief pixTRCMapGeneral()
865 *
866 * \param[in] pixs 32 bpp rgb; not colormapped
867 * \param[in] pixm [optional] 1 bpp mask
868 * \param[in] nar, nag, nab mapping arrays
869 * \return 0 if OK, 1 on error
870 *
871 * <pre>
872 * Notes:
873 * (1) This operation is in-place on %pixs.
874 * (2) Each of the r,g,b mapping arrays is of size 256. They map the
875 * input value for that color component into values in the
876 * range [0, 255].
877 * (3) In the special case where the r, g and b mapping arrays are
878 * all the same, call pixTRCMap() instead.
879 * (4) If defined, the optional 1 bpp mask %pixm has its origin
880 * aligned with %pixs, and the map function is applied only
881 * to pixels in %pixs under the fg of pixm.
882 * (5) The alpha channel is not saved.
883 * </pre>
884 */
885 l_int32
886 pixTRCMapGeneral(PIX *pixs,
887 PIX *pixm,
888 NUMA *nar,
889 NUMA *nag,
890 NUMA *nab)
891 {
892 l_int32 w, h, wm, hm, wpl, wplm, i, j;
893 l_uint32 sval32, dval32;
894 l_uint32 *data, *datam, *line, *linem, *tabr, *tabg, *tabb;
895
896 if (!pixs || pixGetDepth(pixs) != 32)
897 return ERROR_INT("pixs not defined or not 32 bpp", __func__, 1);
898 if (pixm && pixGetDepth(pixm) != 1)
899 return ERROR_INT("pixm defined and not 1 bpp", __func__, 1);
900 if (!nar || !nag || !nab)
901 return ERROR_INT("na{r,g,b} not all defined", __func__, 1);
902 if (numaGetCount(nar) != 256 || numaGetCount(nag) != 256 ||
903 numaGetCount(nab) != 256)
904 return ERROR_INT("na{r,g,b} not all of size 256", __func__, 1);
905
906 /* Get the arrays for efficiency */
907 tabr = (l_uint32 *)numaGetIArray(nar);
908 tabg = (l_uint32 *)numaGetIArray(nag);
909 tabb = (l_uint32 *)numaGetIArray(nab);
910 pixGetDimensions(pixs, &w, &h, NULL);
911 wpl = pixGetWpl(pixs);
912 data = pixGetData(pixs);
913 if (!pixm) {
914 for (i = 0; i < h; i++) {
915 line = data + i * wpl;
916 for (j = 0; j < w; j++) {
917 sval32 = *(line + j);
918 dval32 =
919 tabr[(sval32 >> L_RED_SHIFT) & 0xff] << L_RED_SHIFT |
920 tabg[(sval32 >> L_GREEN_SHIFT) & 0xff] << L_GREEN_SHIFT |
921 tabb[(sval32 >> L_BLUE_SHIFT) & 0xff] << L_BLUE_SHIFT;
922 *(line + j) = dval32;
923 }
924 }
925 } else {
926 datam = pixGetData(pixm);
927 wplm = pixGetWpl(pixm);
928 pixGetDimensions(pixm, &wm, &hm, NULL);
929 for (i = 0; i < h; i++) {
930 if (i >= hm)
931 break;
932 line = data + i * wpl;
933 linem = datam + i * wplm;
934 for (j = 0; j < w; j++) {
935 if (j >= wm)
936 break;
937 if (GET_DATA_BIT(linem, j) == 0)
938 continue;
939 sval32 = *(line + j);
940 dval32 =
941 tabr[(sval32 >> L_RED_SHIFT) & 0xff] << L_RED_SHIFT |
942 tabg[(sval32 >> L_GREEN_SHIFT) & 0xff] << L_GREEN_SHIFT |
943 tabb[(sval32 >> L_BLUE_SHIFT) & 0xff] << L_BLUE_SHIFT;
944 *(line + j) = dval32;
945 }
946 }
947 }
948
949 LEPT_FREE(tabr);
950 LEPT_FREE(tabg);
951 LEPT_FREE(tabb);
952 return 0;
953 }
954
955
956
957 /*-----------------------------------------------------------------------*
958 * Unsharp masking *
959 *-----------------------------------------------------------------------*/
960 /*!
961 * \brief pixUnsharpMasking()
962 *
963 * \param[in] pixs all depths except 1 bpp; with or without colormaps
964 * \param[in] halfwidth "half-width" of smoothing filter
965 * \param[in] fract fraction of edge added back into image
966 * \return pixd, or NULL on error
967 *
968 * <pre>
969 * Notes:
970 * (1) We use symmetric smoothing filters of odd dimension,
971 * typically use sizes of 3, 5, 7, etc. The %halfwidth parameter
972 * for these is (size - 1)/2; i.e., 1, 2, 3, etc.
973 * (2) The fract parameter is typically taken in the
974 * range: 0.2 < fract < 0.7
975 * (3) Returns a clone if no sharpening is requested.
976 * </pre>
977 */
978 PIX *
979 pixUnsharpMasking(PIX *pixs,
980 l_int32 halfwidth,
981 l_float32 fract)
982 {
983 l_int32 d;
984 PIX *pix1, *pixd, *pixr, *pixrs, *pixg, *pixgs, *pixb, *pixbs;
985
986 if (!pixs || (pixGetDepth(pixs) == 1))
987 return (PIX *)ERROR_PTR("pixs not defined or 1 bpp", __func__, NULL);
988 if (fract <= 0.0 || halfwidth <= 0) {
989 L_WARNING("no sharpening requested; clone returned\n", __func__);
990 return pixClone(pixs);
991 }
992
993 if (halfwidth == 1 || halfwidth == 2)
994 return pixUnsharpMaskingFast(pixs, halfwidth, fract, L_BOTH_DIRECTIONS);
995
996 /* Remove colormap; clone if possible; result is either 8 or 32 bpp */
997 if ((pix1 = pixConvertTo8Or32(pixs, L_CLONE, 0)) == NULL)
998 return (PIX *)ERROR_PTR("pix1 not made", __func__, NULL);
999
1000 /* Sharpen */
1001 d = pixGetDepth(pix1);
1002 if (d == 8) {
1003 pixd = pixUnsharpMaskingGray(pix1, halfwidth, fract);
1004 } else { /* d == 32 */
1005 pixr = pixGetRGBComponent(pix1, COLOR_RED);
1006 pixrs = pixUnsharpMaskingGray(pixr, halfwidth, fract);
1007 pixDestroy(&pixr);
1008 pixg = pixGetRGBComponent(pix1, COLOR_GREEN);
1009 pixgs = pixUnsharpMaskingGray(pixg, halfwidth, fract);
1010 pixDestroy(&pixg);
1011 pixb = pixGetRGBComponent(pix1, COLOR_BLUE);
1012 pixbs = pixUnsharpMaskingGray(pixb, halfwidth, fract);
1013 pixDestroy(&pixb);
1014 pixd = pixCreateRGBImage(pixrs, pixgs, pixbs);
1015 pixDestroy(&pixrs);
1016 pixDestroy(&pixgs);
1017 pixDestroy(&pixbs);
1018 if (pixGetSpp(pixs) == 4)
1019 pixCopyRGBComponent(pixd, pixs, L_ALPHA_CHANNEL);
1020 }
1021
1022 pixDestroy(&pix1);
1023 return pixd;
1024 }
1025
1026
1027 /*!
1028 * \brief pixUnsharpMaskingGray()
1029 *
1030 * \param[in] pixs 8 bpp; no colormap
1031 * \param[in] halfwidth "half-width" of smoothing filter
1032 * \param[in] fract fraction of edge added back into image
1033 * \return pixd, or NULL on error
1034 *
1035 * <pre>
1036 * Notes:
1037 * (1) We use symmetric smoothing filters of odd dimension,
1038 * typically use sizes of 3, 5, 7, etc. The %halfwidth parameter
1039 * for these is (size - 1)/2; i.e., 1, 2, 3, etc.
1040 * (2) The fract parameter is typically taken in the range:
1041 * 0.2 < fract < 0.7
1042 * (3) Returns a clone if no sharpening is requested.
1043 * </pre>
1044 */
1045 PIX *
1046 pixUnsharpMaskingGray(PIX *pixs,
1047 l_int32 halfwidth,
1048 l_float32 fract)
1049 {
1050 l_int32 w, h, d;
1051 PIX *pixc, *pixd;
1052 PIXACC *pixacc;
1053
1054 if (!pixs)
1055 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
1056 pixGetDimensions(pixs, &w, &h, &d);
1057 if (d != 8 || pixGetColormap(pixs) != NULL)
1058 return (PIX *)ERROR_PTR("pixs not 8 bpp or has cmap", __func__, NULL);
1059 if (fract <= 0.0 || halfwidth <= 0) {
1060 L_WARNING("no sharpening requested; clone returned\n", __func__);
1061 return pixClone(pixs);
1062 }
1063 if (halfwidth == 1 || halfwidth == 2)
1064 return pixUnsharpMaskingGrayFast(pixs, halfwidth, fract,
1065 L_BOTH_DIRECTIONS);
1066
1067 if ((pixc = pixBlockconvGray(pixs, NULL, halfwidth, halfwidth)) == NULL)
1068 return (PIX *)ERROR_PTR("pixc not made", __func__, NULL);
1069
1070 /* Steps:
1071 * (1) edge image is pixs - pixc (this is highpass part)
1072 * (2) multiply edge image by fract
1073 * (3) add fraction of edge to pixs
1074 *
1075 * To show how this is done with both interfaces to arithmetic
1076 * on integer Pix, here is the implementation in the lower-level
1077 * function calls:
1078 * pixt = pixInitAccumulate(w, h, 0x10000000)) == NULL)
1079 * pixAccumulate(pixt, pixs, L_ARITH_ADD);
1080 * pixAccumulate(pixt, pixc, L_ARITH_SUBTRACT);
1081 * pixMultConstAccumulate(pixt, fract, 0x10000000);
1082 * pixAccumulate(pixt, pixs, L_ARITH_ADD);
1083 * pixd = pixFinalAccumulate(pixt, 0x10000000, 8)) == NULL)
1084 * pixDestroy(&pixt);
1085 *
1086 * The code below does the same thing using the Pixacc accumulator,
1087 * hiding the details of the offset that is needed for subtraction.
1088 */
1089 pixacc = pixaccCreate(w, h, 1);
1090 pixaccAdd(pixacc, pixs);
1091 pixaccSubtract(pixacc, pixc);
1092 pixaccMultConst(pixacc, fract);
1093 pixaccAdd(pixacc, pixs);
1094 pixd = pixaccFinal(pixacc, 8);
1095 pixaccDestroy(&pixacc);
1096
1097 pixDestroy(&pixc);
1098 return pixd;
1099 }
1100
1101
1102 /*!
1103 * \brief pixUnsharpMaskingFast()
1104 *
1105 * \param[in] pixs all depths except 1 bpp; with or without colormaps
1106 * \param[in] halfwidth "half-width" of smoothing filter; 1 and 2 only
1107 * \param[in] fract fraction of high frequency added to image
1108 * \param[in] direction L_HORIZ, L_VERT, L_BOTH_DIRECTIONS
1109 * \return pixd, or NULL on error
1110 *
1111 * <pre>
1112 * Notes:
1113 * (1) The fast version uses separable 1-D filters directly on
1114 * the input image. The halfwidth is either 1 (full width = 3)
1115 * or 2 (full width = 5).
1116 * (2) The fract parameter is typically taken in the
1117 * range: 0.2 < fract < 0.7
1118 * (3) To skip horizontal sharpening, use %fracth = 0.0; ditto for %fractv
1119 * (4) For one dimensional filtering (as an example):
1120 * For %halfwidth = 1, the low-pass filter is
1121 * L: 1/3 1/3 1/3
1122 * and the high-pass filter is
1123 * H = I - L: -1/3 2/3 -1/3
1124 * For %halfwidth = 2, the low-pass filter is
1125 * L: 1/5 1/5 1/5 1/5 1/5
1126 * and the high-pass filter is
1127 * H = I - L: -1/5 -1/5 4/5 -1/5 -1/5
1128 * The new sharpened pixel value is found by adding some fraction
1129 * of the high-pass filter value (which sums to 0) to the
1130 * initial pixel value:
1131 * N = I + fract * H
1132 * (5) For 2D, the sharpening filter is not separable, because the
1133 * vertical filter depends on the horizontal location relative
1134 * to the filter origin, and v.v. So we either do the full
1135 * 2D filter (for %halfwidth == 1) or do the low-pass
1136 * convolution separably and then compose with the original pix.
1137 * (6) Returns a clone if no sharpening is requested.
1138 * </pre>
1139 */
1140 PIX *
1141 pixUnsharpMaskingFast(PIX *pixs,
1142 l_int32 halfwidth,
1143 l_float32 fract,
1144 l_int32 direction)
1145 {
1146 l_int32 d;
1147 PIX *pixt, *pixd, *pixr, *pixrs, *pixg, *pixgs, *pixb, *pixbs;
1148
1149 if (!pixs || (pixGetDepth(pixs) == 1))
1150 return (PIX *)ERROR_PTR("pixs not defined or 1 bpp", __func__, NULL);
1151 if (fract <= 0.0 || halfwidth <= 0) {
1152 L_WARNING("no sharpening requested; clone returned\n", __func__);
1153 return pixClone(pixs);
1154 }
1155 if (halfwidth != 1 && halfwidth != 2)
1156 return (PIX *)ERROR_PTR("halfwidth must be 1 or 2", __func__, NULL);
1157 if (direction != L_HORIZ && direction != L_VERT &&
1158 direction != L_BOTH_DIRECTIONS)
1159 return (PIX *)ERROR_PTR("invalid direction", __func__, NULL);
1160
1161 /* Remove colormap; clone if possible; result is either 8 or 32 bpp */
1162 if ((pixt = pixConvertTo8Or32(pixs, L_CLONE, 0)) == NULL)
1163 return (PIX *)ERROR_PTR("pixt not made", __func__, NULL);
1164
1165 /* Sharpen */
1166 d = pixGetDepth(pixt);
1167 if (d == 8) {
1168 pixd = pixUnsharpMaskingGrayFast(pixt, halfwidth, fract, direction);
1169 } else { /* d == 32 */
1170 pixr = pixGetRGBComponent(pixs, COLOR_RED);
1171 pixrs = pixUnsharpMaskingGrayFast(pixr, halfwidth, fract, direction);
1172 pixDestroy(&pixr);
1173 pixg = pixGetRGBComponent(pixs, COLOR_GREEN);
1174 pixgs = pixUnsharpMaskingGrayFast(pixg, halfwidth, fract, direction);
1175 pixDestroy(&pixg);
1176 pixb = pixGetRGBComponent(pixs, COLOR_BLUE);
1177 pixbs = pixUnsharpMaskingGrayFast(pixb, halfwidth, fract, direction);
1178 pixDestroy(&pixb);
1179 pixd = pixCreateRGBImage(pixrs, pixgs, pixbs);
1180 if (pixGetSpp(pixs) == 4)
1181 pixCopyRGBComponent(pixd, pixs, L_ALPHA_CHANNEL);
1182 pixDestroy(&pixrs);
1183 pixDestroy(&pixgs);
1184 pixDestroy(&pixbs);
1185 }
1186
1187 pixDestroy(&pixt);
1188 return pixd;
1189 }
1190
1191
1192
1193 /*!
1194 * \brief pixUnsharpMaskingGrayFast()
1195 *
1196 * \param[in] pixs 8 bpp; no colormap
1197 * \param[in] halfwidth "half-width" of smoothing filter: 1 or 2
1198 * \param[in] fract fraction of high frequency added to image
1199 * \param[in] direction L_HORIZ, L_VERT, L_BOTH_DIRECTIONS
1200 * \return pixd, or NULL on error
1201 *
1202 * <pre>
1203 * Notes:
1204 * (1) For usage and explanation of the algorithm, see notes
1205 * in pixUnsharpMaskingFast().
1206 * (2) Returns a clone if no sharpening is requested.
1207 * </pre>
1208 */
1209 PIX *
1210 pixUnsharpMaskingGrayFast(PIX *pixs,
1211 l_int32 halfwidth,
1212 l_float32 fract,
1213 l_int32 direction)
1214 {
1215 PIX *pixd;
1216
1217 if (!pixs)
1218 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
1219 if (pixGetDepth(pixs) != 8 || pixGetColormap(pixs) != NULL)
1220 return (PIX *)ERROR_PTR("pixs not 8 bpp or has cmap", __func__, NULL);
1221 if (fract <= 0.0 || halfwidth <= 0) {
1222 L_WARNING("no sharpening requested; clone returned\n", __func__);
1223 return pixClone(pixs);
1224 }
1225 if (halfwidth != 1 && halfwidth != 2)
1226 return (PIX *)ERROR_PTR("halfwidth must be 1 or 2", __func__, NULL);
1227 if (direction != L_HORIZ && direction != L_VERT &&
1228 direction != L_BOTH_DIRECTIONS)
1229 return (PIX *)ERROR_PTR("invalid direction", __func__, NULL);
1230
1231 if (direction != L_BOTH_DIRECTIONS)
1232 pixd = pixUnsharpMaskingGray1D(pixs, halfwidth, fract, direction);
1233 else /* 2D sharpening */
1234 pixd = pixUnsharpMaskingGray2D(pixs, halfwidth, fract);
1235
1236 return pixd;
1237 }
1238
1239
1240 /*!
1241 * \brief pixUnsharpMaskingGray1D()
1242 *
1243 * \param[in] pixs 8 bpp; no colormap
1244 * \param[in] halfwidth "half-width" of smoothing filter: 1 or 2
1245 * \param[in] fract fraction of high frequency added to image
1246 * \param[in] direction filtering direction; use L_HORIZ or L_VERT
1247 * \return pixd, or NULL on error
1248 *
1249 * <pre>
1250 * Notes:
1251 * (1) For usage and explanation of the algorithm, see notes
1252 * in pixUnsharpMaskingFast().
1253 * (2) Returns a clone if no sharpening is requested.
1254 * </pre>
1255 */
1256 PIX *
1257 pixUnsharpMaskingGray1D(PIX *pixs,
1258 l_int32 halfwidth,
1259 l_float32 fract,
1260 l_int32 direction)
1261 {
1262 l_int32 w, h, d, wpls, wpld, i, j, ival;
1263 l_uint32 *datas, *datad;
1264 l_uint32 *lines, *lines0, *lines1, *lines2, *lines3, *lines4, *lined;
1265 l_float32 val, a[5];
1266 PIX *pixd;
1267
1268 if (!pixs)
1269 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
1270 pixGetDimensions(pixs, &w, &h, &d);
1271 if (d != 8 || pixGetColormap(pixs) != NULL)
1272 return (PIX *)ERROR_PTR("pixs not 8 bpp or has cmap", __func__, NULL);
1273 if (fract <= 0.0 || halfwidth <= 0) {
1274 L_WARNING("no sharpening requested; clone returned\n", __func__);
1275 return pixClone(pixs);
1276 }
1277 if (halfwidth != 1 && halfwidth != 2)
1278 return (PIX *)ERROR_PTR("halfwidth must be 1 or 2", __func__, NULL);
1279
1280 /* Initialize pixd with pixels from pixs that will not be
1281 * set when computing the sharpened values. */
1282 pixd = pixCopyBorder(NULL, pixs, halfwidth, halfwidth,
1283 halfwidth, halfwidth);
1284 datas = pixGetData(pixs);
1285 datad = pixGetData(pixd);
1286 wpls = pixGetWpl(pixs);
1287 wpld = pixGetWpl(pixd);
1288
1289 if (halfwidth == 1) {
1290 a[0] = -fract / 3.0;
1291 a[1] = 1.0 + fract * 2.0 / 3.0;
1292 a[2] = a[0];
1293 } else { /* halfwidth == 2 */
1294 a[0] = -fract / 5.0;
1295 a[1] = a[0];
1296 a[2] = 1.0 + fract * 4.0 / 5.0;
1297 a[3] = a[0];
1298 a[4] = a[0];
1299 }
1300
1301 if (direction == L_HORIZ) {
1302 for (i = 0; i < h; i++) {
1303 lines = datas + i * wpls;
1304 lined = datad + i * wpld;
1305 if (halfwidth == 1) {
1306 for (j = 1; j < w - 1; j++) {
1307 val = a[0] * GET_DATA_BYTE(lines, j - 1) +
1308 a[1] * GET_DATA_BYTE(lines, j) +
1309 a[2] * GET_DATA_BYTE(lines, j + 1);
1310 ival = (l_int32)val;
1311 ival = L_MAX(0, ival);
1312 ival = L_MIN(255, ival);
1313 SET_DATA_BYTE(lined, j, ival);
1314 }
1315 } else { /* halfwidth == 2 */
1316 for (j = 2; j < w - 2; j++) {
1317 val = a[0] * GET_DATA_BYTE(lines, j - 2) +
1318 a[1] * GET_DATA_BYTE(lines, j - 1) +
1319 a[2] * GET_DATA_BYTE(lines, j) +
1320 a[3] * GET_DATA_BYTE(lines, j + 1) +
1321 a[4] * GET_DATA_BYTE(lines, j + 2);
1322 ival = (l_int32)val;
1323 ival = L_MAX(0, ival);
1324 ival = L_MIN(255, ival);
1325 SET_DATA_BYTE(lined, j, ival);
1326 }
1327 }
1328 }
1329 } else { /* direction == L_VERT */
1330 if (halfwidth == 1) {
1331 for (i = 1; i < h - 1; i++) {
1332 lines0 = datas + (i - 1) * wpls;
1333 lines1 = datas + i * wpls;
1334 lines2 = datas + (i + 1) * wpls;
1335 lined = datad + i * wpld;
1336 for (j = 0; j < w; j++) {
1337 val = a[0] * GET_DATA_BYTE(lines0, j) +
1338 a[1] * GET_DATA_BYTE(lines1, j) +
1339 a[2] * GET_DATA_BYTE(lines2, j);
1340 ival = (l_int32)val;
1341 ival = L_MAX(0, ival);
1342 ival = L_MIN(255, ival);
1343 SET_DATA_BYTE(lined, j, ival);
1344 }
1345 }
1346 } else { /* halfwidth == 2 */
1347 for (i = 2; i < h - 2; i++) {
1348 lines0 = datas + (i - 2) * wpls;
1349 lines1 = datas + (i - 1) * wpls;
1350 lines2 = datas + i * wpls;
1351 lines3 = datas + (i + 1) * wpls;
1352 lines4 = datas + (i + 2) * wpls;
1353 lined = datad + i * wpld;
1354 for (j = 0; j < w; j++) {
1355 val = a[0] * GET_DATA_BYTE(lines0, j) +
1356 a[1] * GET_DATA_BYTE(lines1, j) +
1357 a[2] * GET_DATA_BYTE(lines2, j) +
1358 a[3] * GET_DATA_BYTE(lines3, j) +
1359 a[4] * GET_DATA_BYTE(lines4, j);
1360 ival = (l_int32)val;
1361 ival = L_MAX(0, ival);
1362 ival = L_MIN(255, ival);
1363 SET_DATA_BYTE(lined, j, ival);
1364 }
1365 }
1366 }
1367 }
1368
1369 return pixd;
1370 }
1371
1372
1373 /*!
1374 * \brief pixUnsharpMaskingGray2D()
1375 *
1376 * \param[in] pixs 8 bpp; no colormap
1377 * \param[in] halfwidth "half-width" of smoothing filter: 1 or 2
1378 * \param[in] fract fraction of high frequency added to image
1379 * \return pixd, or NULL on error
1380 *
1381 * <pre>
1382 * Notes:
1383 * (1) This is for %halfwidth == 1, 2.
1384 * (2) The lowpass filter is implemented separably.
1385 * (3) Returns a clone if no sharpening is requested.
1386 * </pre>
1387 */
1388 PIX *
1389 pixUnsharpMaskingGray2D(PIX *pixs,
1390 l_int32 halfwidth,
1391 l_float32 fract)
1392 {
1393 l_int32 w, h, d, wpls, wpld, wplf, i, j, ival, sval;
1394 l_uint32 *datas, *datad, *lines, *lined;
1395 l_float32 val, norm;
1396 l_float32 *dataf, *linef, *linef0, *linef1, *linef2, *linef3, *linef4;
1397 PIX *pixd;
1398 FPIX *fpix;
1399
1400 if (!pixs)
1401 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
1402 pixGetDimensions(pixs, &w, &h, &d);
1403 if (d != 8 || pixGetColormap(pixs) != NULL)
1404 return (PIX *)ERROR_PTR("pixs not 8 bpp or has cmap", __func__, NULL);
1405 if (fract <= 0.0 || halfwidth <= 0) {
1406 L_WARNING("no sharpening requested; clone returned\n", __func__);
1407 return pixClone(pixs);
1408 }
1409 if (halfwidth != 1 && halfwidth != 2)
1410 return (PIX *)ERROR_PTR("halfwidth must be 1 or 2", __func__, NULL);
1411
1412 if ((pixd = pixCopyBorder(NULL, pixs, halfwidth, halfwidth,
1413 halfwidth, halfwidth)) == NULL)
1414 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
1415 datad = pixGetData(pixd);
1416 wpld = pixGetWpl(pixd);
1417 datas = pixGetData(pixs);
1418 wpls = pixGetWpl(pixs);
1419
1420 /* Do the low pass separably. Store the result of horizontal
1421 * smoothing in an intermediate fpix. */
1422 if ((fpix = fpixCreate(w, h)) == NULL) {
1423 pixDestroy(&pixd);
1424 return (PIX *)ERROR_PTR("fpix not made", __func__, NULL);
1425 }
1426 dataf = fpixGetData(fpix);
1427 wplf = fpixGetWpl(fpix);
1428 if (halfwidth == 1) {
1429 for (i = 0; i < h; i++) {
1430 lines = datas + i * wpls;
1431 linef = dataf + i * wplf;
1432 for (j = 1; j < w - 1; j++) {
1433 val = GET_DATA_BYTE(lines, j - 1) +
1434 GET_DATA_BYTE(lines, j) +
1435 GET_DATA_BYTE(lines, j + 1);
1436 linef[j] = val;
1437 }
1438 }
1439 } else {
1440 for (i = 0; i < h; i++) {
1441 lines = datas + i * wpls;
1442 linef = dataf + i * wplf;
1443 for (j = 2; j < w - 2; j++) {
1444 val = GET_DATA_BYTE(lines, j - 2) +
1445 GET_DATA_BYTE(lines, j - 1) +
1446 GET_DATA_BYTE(lines, j) +
1447 GET_DATA_BYTE(lines, j + 1) +
1448 GET_DATA_BYTE(lines, j + 2);
1449 linef[j] = val;
1450 }
1451 }
1452 }
1453
1454 /* Do vertical smoothing to finish the low-pass filter.
1455 * At each pixel, if L is the lowpass value, I is the
1456 * src pixel value and f is the fraction of highpass to
1457 * be added to I, then the highpass filter value is
1458 * H = I - L
1459 * and the new sharpened value is
1460 * N = I + f * H. */
1461 if (halfwidth == 1) {
1462 for (i = 1; i < h - 1; i++) {
1463 linef0 = dataf + (i - 1) * wplf;
1464 linef1 = dataf + i * wplf;
1465 linef2 = dataf + (i + 1) * wplf;
1466 lined = datad + i * wpld;
1467 lines = datas + i * wpls;
1468 norm = 1.0f / 9.0f;
1469 for (j = 1; j < w - 1; j++) {
1470 val = norm * (linef0[j] + linef1[j] +
1471 linef2[j]); /* L: lowpass filter value */
1472 sval = GET_DATA_BYTE(lines, j); /* I: source pixel */
1473 ival = (l_int32)(sval + fract * (sval - val) + 0.5);
1474 ival = L_MAX(0, ival);
1475 ival = L_MIN(255, ival);
1476 SET_DATA_BYTE(lined, j, ival);
1477 }
1478 }
1479 } else {
1480 for (i = 2; i < h - 2; i++) {
1481 linef0 = dataf + (i - 2) * wplf;
1482 linef1 = dataf + (i - 1) * wplf;
1483 linef2 = dataf + i * wplf;
1484 linef3 = dataf + (i + 1) * wplf;
1485 linef4 = dataf + (i + 2) * wplf;
1486 lined = datad + i * wpld;
1487 lines = datas + i * wpls;
1488 norm = 1.0f / 25.0f;
1489 for (j = 2; j < w - 2; j++) {
1490 val = norm * (linef0[j] + linef1[j] + linef2[j] + linef3[j] +
1491 linef4[j]); /* L: lowpass filter value */
1492 sval = GET_DATA_BYTE(lines, j); /* I: source pixel */
1493 ival = (l_int32)(sval + fract * (sval - val) + 0.5);
1494 ival = L_MAX(0, ival);
1495 ival = L_MIN(255, ival);
1496 SET_DATA_BYTE(lined, j, ival);
1497 }
1498 }
1499 }
1500
1501 fpixDestroy(&fpix);
1502 return pixd;
1503 }
1504
1505
1506 /*-----------------------------------------------------------------------*
1507 * Hue and saturation modification *
1508 *-----------------------------------------------------------------------*/
1509 /*!
1510 * \brief pixModifyHue()
1511 *
1512 * \param[in] pixd [optional] can be null or equal to pixs
1513 * \param[in] pixs 32 bpp rgb
1514 * \param[in] fract between -1.0 and 1.0
1515 * \return pixd, or NULL on error
1516 *
1517 * <pre>
1518 * Notes:
1519 * (1) pixd must either be null or equal to pixs.
1520 * For in-place operation, set pixd == pixs:
1521 * pixEqualizeTRC(pixs, pixs, ...);
1522 * To get a new image, set pixd == null:
1523 * pixd = pixEqualizeTRC(NULL, pixs, ...);
1524 * (2) Use fract > 0.0 to increase hue value; < 0.0 to decrease it.
1525 * 1.0 (or -1.0) represents a 360 degree rotation; i.e., no change.
1526 * (3) If no modification is requested (fract = -1.0 or 0 or 1.0),
1527 * return a copy unless in-place, in which case this is a no-op.
1528 * (4) This leaves saturation and intensity invariant.
1529 * (5) See discussion of color-modification methods, in coloring.c.
1530 * </pre>
1531 */
1532 PIX *
1533 pixModifyHue(PIX *pixd,
1534 PIX *pixs,
1535 l_float32 fract)
1536 {
1537 l_int32 w, h, d, i, j, wpl, delhue;
1538 l_int32 rval, gval, bval, hval, sval, vval;
1539 l_uint32 *data, *line;
1540
1541 if (!pixs)
1542 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
1543 if (pixGetColormap(pixs) != NULL)
1544 return (PIX *)ERROR_PTR("pixs colormapped", __func__, NULL);
1545 if (pixd && (pixd != pixs))
1546 return (PIX *)ERROR_PTR("pixd not null or pixs", __func__, pixd);
1547 pixGetDimensions(pixs, &w, &h, &d);
1548 if (d != 32)
1549 return (PIX *)ERROR_PTR("pixs not 32 bpp", __func__, NULL);
1550 if (L_ABS(fract) > 1.0)
1551 return (PIX *)ERROR_PTR("fract not in [-1.0 ... 1.0]", __func__, NULL);
1552
1553 pixd = pixCopy(pixd, pixs);
1554
1555 delhue = (l_int32)(240 * fract);
1556 if (delhue == 0 || delhue == 240 || delhue == -240) {
1557 L_WARNING("no change requested in hue\n", __func__);
1558 return pixd;
1559 }
1560 if (delhue < 0)
1561 delhue += 240;
1562
1563 data = pixGetData(pixd);
1564 wpl = pixGetWpl(pixd);
1565 for (i = 0; i < h; i++) {
1566 line = data + i * wpl;
1567 for (j = 0; j < w; j++) {
1568 extractRGBValues(line[j], &rval, &gval, &bval);
1569 convertRGBToHSV(rval, gval, bval, &hval, &sval, &vval);
1570 hval = (hval + delhue) % 240;
1571 convertHSVToRGB(hval, sval, vval, &rval, &gval, &bval);
1572 composeRGBPixel(rval, gval, bval, line + j);
1573 }
1574 }
1575 if (pixGetSpp(pixs) == 4)
1576 pixCopyRGBComponent(pixd, pixs, L_ALPHA_CHANNEL);
1577
1578 return pixd;
1579 }
1580
1581
1582 /*!
1583 * \brief pixModifySaturation()
1584 *
1585 * \param[in] pixd [optional] can be null, existing or equal to pixs
1586 * \param[in] pixs 32 bpp rgb
1587 * \param[in] fract between -1.0 and 1.0
1588 * \return pixd, or NULL on error
1589 *
1590 * <pre>
1591 * Notes:
1592 * (1) If fract > 0.0, it gives the fraction that the pixel
1593 * saturation is moved from its initial value toward 255.
1594 * If fract < 0.0, it gives the fraction that the pixel
1595 * saturation is moved from its initial value toward 0.
1596 * The limiting values for fract = -1.0 (1.0) thus set the
1597 * saturation to 0 (255).
1598 * (2) If fract = 0, no modification is requested; return a copy
1599 * unless in-place, in which case this is a no-op.
1600 * (3) This leaves hue and intensity invariant.
1601 * (4) See discussion of color-modification methods, in coloring.c.
1602 * </pre>
1603 */
1604 PIX *
1605 pixModifySaturation(PIX *pixd,
1606 PIX *pixs,
1607 l_float32 fract)
1608 {
1609 l_int32 w, h, d, i, j, wpl;
1610 l_int32 rval, gval, bval, hval, sval, vval;
1611 l_uint32 *data, *line;
1612
1613 if (!pixs)
1614 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
1615 pixGetDimensions(pixs, &w, &h, &d);
1616 if (d != 32)
1617 return (PIX *)ERROR_PTR("pixs not 32 bpp", __func__, NULL);
1618 if (L_ABS(fract) > 1.0)
1619 return (PIX *)ERROR_PTR("fract not in [-1.0 ... 1.0]", __func__, NULL);
1620
1621 pixd = pixCopy(pixd, pixs);
1622 if (fract == 0.0) {
1623 L_WARNING("no change requested in saturation\n", __func__);
1624 return pixd;
1625 }
1626
1627 data = pixGetData(pixd);
1628 wpl = pixGetWpl(pixd);
1629 for (i = 0; i < h; i++) {
1630 line = data + i * wpl;
1631 for (j = 0; j < w; j++) {
1632 extractRGBValues(line[j], &rval, &gval, &bval);
1633 convertRGBToHSV(rval, gval, bval, &hval, &sval, &vval);
1634 if (fract < 0.0)
1635 sval = (l_int32)(sval * (1.0 + fract));
1636 else
1637 sval = (l_int32)(sval + fract * (255 - sval));
1638 convertHSVToRGB(hval, sval, vval, &rval, &gval, &bval);
1639 composeRGBPixel(rval, gval, bval, line + j);
1640 }
1641 }
1642 if (pixGetSpp(pixs) == 4)
1643 pixCopyRGBComponent(pixd, pixs, L_ALPHA_CHANNEL);
1644
1645 return pixd;
1646 }
1647
1648
1649 /*!
1650 * \brief pixMeasureSaturation()
1651 *
1652 * \param[in] pixs 32 bpp rgb
1653 * \param[in] factor subsampling factor; integer >= 1
1654 * \param[out] psat average saturation
1655 * \return 0 if OK, 1 on error
1656 */
1657 l_int32
1658 pixMeasureSaturation(PIX *pixs,
1659 l_int32 factor,
1660 l_float32 *psat)
1661 {
1662 l_int32 w, h, d, i, j, wpl, sum, count;
1663 l_int32 rval, gval, bval, hval, sval, vval;
1664 l_uint32 *data, *line;
1665
1666 if (!psat)
1667 return ERROR_INT("pixs not defined", __func__, 1);
1668 *psat = 0.0;
1669 if (!pixs)
1670 return ERROR_INT("pixs not defined", __func__, 1);
1671 pixGetDimensions(pixs, &w, &h, &d);
1672 if (d != 32)
1673 return ERROR_INT("pixs not 32 bpp", __func__, 1);
1674 if (factor < 1)
1675 return ERROR_INT("subsampling factor < 1", __func__, 1);
1676
1677 data = pixGetData(pixs);
1678 wpl = pixGetWpl(pixs);
1679 for (i = 0, sum = 0, count = 0; i < h; i += factor) {
1680 line = data + i * wpl;
1681 for (j = 0; j < w; j += factor) {
1682 extractRGBValues(line[j], &rval, &gval, &bval);
1683 convertRGBToHSV(rval, gval, bval, &hval, &sval, &vval);
1684 sum += sval;
1685 count++;
1686 }
1687 }
1688
1689 if (count > 0)
1690 *psat = (l_float32)sum / (l_float32)count;
1691 return 0;
1692 }
1693
1694
1695 /*!
1696 * \brief pixModifyBrightness()
1697 *
1698 * \param[in] pixd [optional] can be null, existing or equal to pixs
1699 * \param[in] pixs 32 bpp rgb
1700 * \param[in] fract between -1.0 and 1.0
1701 * \return pixd, or NULL on error
1702 *
1703 * <pre>
1704 * Notes:
1705 * (1) If fract > 0.0, it gives the fraction that the v-parameter,
1706 * which is max(r,g,b), is moved from its initial value toward 255.
1707 * If fract < 0.0, it gives the fraction that the v-parameter
1708 * is moved from its initial value toward 0.
1709 * The limiting values for fract = -1.0 (1.0) thus set the
1710 * v-parameter to 0 (255).
1711 * (2) If fract = 0, no modification is requested; return a copy
1712 * unless in-place, in which case this is a no-op.
1713 * (3) This leaves hue and saturation invariant.
1714 * (4) See discussion of color-modification methods, in coloring.c.
1715 * </pre>
1716 */
1717 PIX *
1718 pixModifyBrightness(PIX *pixd,
1719 PIX *pixs,
1720 l_float32 fract)
1721 {
1722 l_int32 w, h, d, i, j, wpl;
1723 l_int32 rval, gval, bval, hval, sval, vval;
1724 l_uint32 *data, *line;
1725
1726 if (!pixs)
1727 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
1728 pixGetDimensions(pixs, &w, &h, &d);
1729 if (d != 32)
1730 return (PIX *)ERROR_PTR("pixs not 32 bpp", __func__, NULL);
1731 if (L_ABS(fract) > 1.0)
1732 return (PIX *)ERROR_PTR("fract not in [-1.0 ... 1.0]", __func__, NULL);
1733
1734 pixd = pixCopy(pixd, pixs);
1735 if (fract == 0.0) {
1736 L_WARNING("no change requested in brightness\n", __func__);
1737 return pixd;
1738 }
1739
1740 data = pixGetData(pixd);
1741 wpl = pixGetWpl(pixd);
1742 for (i = 0; i < h; i++) {
1743 line = data + i * wpl;
1744 for (j = 0; j < w; j++) {
1745 extractRGBValues(line[j], &rval, &gval, &bval);
1746 convertRGBToHSV(rval, gval, bval, &hval, &sval, &vval);
1747 if (fract > 0.0)
1748 vval = (l_int32)(vval + fract * (255.0 - vval));
1749 else
1750 vval = (l_int32)(vval * (1.0 + fract));
1751 convertHSVToRGB(hval, sval, vval, &rval, &gval, &bval);
1752 composeRGBPixel(rval, gval, bval, line + j);
1753 }
1754 }
1755 if (pixGetSpp(pixs) == 4)
1756 pixCopyRGBComponent(pixd, pixs, L_ALPHA_CHANNEL);
1757
1758 return pixd;
1759 }
1760
1761
1762 /*-----------------------------------------------------------------------*
1763 * Color shifting *
1764 *-----------------------------------------------------------------------*/
1765 /*!
1766 * \brief pixMosaicColorShiftRGB()
1767 *
1768 * \param[in] pixs 32 bpp rgb
1769 * \param[in] roff center offset of red component
1770 * \param[in] goff center offset of green component
1771 * \param[in] boff center offset of blue component
1772 * \param[in] delta increments from center offsets [0.0 - 0.1];
1773 * use 0.0 to get the default (0.04)
1774 * \param[in] nincr number of increments in each (positive and negative)
1775 * direction; use 0 to get the default (2).
1776 * \return pix, or NULL on error
1777 *
1778 * <pre>
1779 * Notes:
1780 * (1) This generates a mosaic view of the effect of shifting the RGB
1781 * components. See pixColorShiftRGB() for details on the shifting.
1782 * (2) The offsets (%roff, %goff, %boff) set the color center point,
1783 * and the deviations from this are shown separately for deltas
1784 * in r, g and b. For each component, we show 2 * %nincr + 1
1785 * images.
1786 * (3) The pix must have minimum dimensions of 100 and an aspect
1787 * ratio not exceeding 5.0.
1788 * (4) Usage: color prints differ from the original due to three factors:
1789 * illumination, calibration of the camera in acquisition,
1790 * and calibration of the printer. This function can be used
1791 * to iteratively match a color print to the original. On each
1792 * iteration, the center offsets are set to the best match so
1793 * far, and the %delta increments are typically reduced.
1794 * </pre>
1795 */
1796 PIX *
1797 pixMosaicColorShiftRGB(PIX *pixs,
1798 l_float32 roff,
1799 l_float32 goff,
1800 l_float32 boff,
1801 l_float32 delta,
1802 l_int32 nincr)
1803 {
1804 char buf[64];
1805 l_int32 i, w, h;
1806 l_float32 del, ratio;
1807 L_BMF *bmf;
1808 PIX *pix1, *pix2, *pix3;
1809 PIXA *pixa;
1810
1811 if (!pixs || pixGetDepth(pixs) != 32)
1812 return (PIX *)ERROR_PTR("pixs undefined or not rgb", __func__, NULL);
1813 if (roff < -1.0 || roff > 1.0)
1814 return (PIX *)ERROR_PTR("roff not in [-1.0, 1.0]", __func__, NULL);
1815 if (goff < -1.0 || goff > 1.0)
1816 return (PIX *)ERROR_PTR("goff not in [-1.0, 1.0]", __func__, NULL);
1817 if (boff < -1.0 || boff > 1.0)
1818 return (PIX *)ERROR_PTR("boff not in [-1.0, 1.0]", __func__, NULL);
1819 if (delta < 0.0 || delta > 0.1)
1820 return (PIX *)ERROR_PTR("delta not in [0.0, 0.1]", __func__, NULL);
1821 if (delta == 0.0) delta = 0.04f;
1822 if (nincr < 0 || nincr > 6)
1823 return (PIX *)ERROR_PTR("nincr not in [0, 6]", __func__, NULL);
1824 if (nincr == 0) nincr = 2;
1825
1826 /* Require width and height to be >= 100, and the aspect ratio <= 5.0 */
1827 pixGetDimensions(pixs, &w, &h, NULL);
1828 if (w < 100 || h < 100)
1829 return (PIX *)ERROR_PTR("w and h not both >= 100", __func__, NULL);
1830 pixMaxAspectRatio(pixs, &ratio);
1831 if (ratio < 1.0 || ratio > 5.0) {
1832 L_ERROR("invalid aspect ratio %5.1f\n", __func__, ratio);
1833 return NULL;
1834 }
1835
1836 pixa = pixaCreate(3 * (2 * nincr + 1));
1837 bmf = bmfCreate(NULL, 8);
1838 pix1 = pixScaleToSize(pixs, 400, 0);
1839 for (i = 0, del = - nincr * delta; i < 2 * nincr + 1; i++, del += delta) {
1840 pix2 = pixColorShiftRGB(pix1, roff + del, goff, boff);
1841 snprintf(buf, sizeof(buf), "%4.2f, %4.2f, %4.2f",
1842 roff + del, goff, boff);
1843 pix3 = pixAddSingleTextblock(pix2, bmf, buf, 0xff000000,
1844 L_ADD_BELOW, 0);
1845 pixaAddPix(pixa, pix3, L_INSERT);
1846 pixDestroy(&pix2);
1847 }
1848 for (i = 0, del = - nincr * delta; i < 2 * nincr + 1; i++, del += delta) {
1849 pix2 = pixColorShiftRGB(pix1, roff, goff + del, boff);
1850 snprintf(buf, sizeof(buf), "%4.2f, %4.2f, %4.2f",
1851 roff, goff + del, boff);
1852 pix3 = pixAddSingleTextblock(pix2, bmf, buf, 0xff000000,
1853 L_ADD_BELOW, 0);
1854 pixaAddPix(pixa, pix3, L_INSERT);
1855 pixDestroy(&pix2);
1856 }
1857 for (i = 0, del = - nincr * delta; i < 2 * nincr + 1; i++, del += delta) {
1858 pix2 = pixColorShiftRGB(pix1, roff, goff, boff + del);
1859 snprintf(buf, sizeof(buf), "%4.2f, %4.2f, %4.2f",
1860 roff, goff, boff + del);
1861 pix3 = pixAddSingleTextblock(pix2, bmf, buf, 0xff000000,
1862 L_ADD_BELOW, 0);
1863 pixaAddPix(pixa, pix3, L_INSERT);
1864 pixDestroy(&pix2);
1865 }
1866 pixDestroy(&pix1);
1867
1868 pix1 = pixaDisplayTiledAndScaled(pixa, 32, 300, 2 * nincr + 1, 0, 30, 2);
1869 pixaDestroy(&pixa);
1870 bmfDestroy(&bmf);
1871 return pix1;
1872 }
1873
1874
1875 /*!
1876 * \brief pixColorShiftRGB()
1877 *
1878 * \param[in] pixs 32 bpp rgb
1879 * \param[in] rfract fractional shift in red component
1880 * \param[in] gfract fractional shift in green component
1881 * \param[in] bfract fractional shift in blue component
1882 * \return pixd, or NULL on error
1883 *
1884 * <pre>
1885 * Notes:
1886 * (1) This allows independent fractional shifts of the r,g and b
1887 * components. A positive shift pushes to saturation (255);
1888 * a negative shift pushes toward 0 (black).
1889 * (2) The effect can be imagined using a color wheel that consists
1890 * (for our purposes) of these 6 colors, separated by 60 degrees:
1891 * red, magenta, blue, cyan, green, yellow
1892 * (3) So, for example, a negative shift of the blue component
1893 * (bfract < 0) could be accompanied by positive shifts
1894 * of red and green to make an image more yellow.
1895 * (4) Examples of limiting cases:
1896 * rfract = 1 ==> r = 255
1897 * rfract = -1 ==> r = 0
1898 * </pre>
1899 */
1900 PIX *
1901 pixColorShiftRGB(PIX *pixs,
1902 l_float32 rfract,
1903 l_float32 gfract,
1904 l_float32 bfract)
1905 {
1906 l_int32 w, h, i, j, wpls, wpld, rval, gval, bval;
1907 l_int32 *rlut, *glut, *blut;
1908 l_uint32 *datas, *datad, *lines, *lined;
1909 l_float32 fi;
1910 PIX *pixd;
1911
1912 if (!pixs)
1913 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
1914 if (pixGetDepth(pixs) != 32)
1915 return (PIX *)ERROR_PTR("pixs not 32 bpp", __func__, NULL);
1916 if (rfract < -1.0 || rfract > 1.0)
1917 return (PIX *)ERROR_PTR("rfract not in [-1.0, 1.0]", __func__, NULL);
1918 if (gfract < -1.0 || gfract > 1.0)
1919 return (PIX *)ERROR_PTR("gfract not in [-1.0, 1.0]", __func__, NULL);
1920 if (bfract < -1.0 || bfract > 1.0)
1921 return (PIX *)ERROR_PTR("bfract not in [-1.0, 1.0]", __func__, NULL);
1922 if (rfract == 0.0 && gfract == 0.0 && bfract == 0.0)
1923 return pixCopy(NULL, pixs);
1924
1925 rlut = (l_int32 *)LEPT_CALLOC(256, sizeof(l_int32));
1926 glut = (l_int32 *)LEPT_CALLOC(256, sizeof(l_int32));
1927 blut = (l_int32 *)LEPT_CALLOC(256, sizeof(l_int32));
1928 for (i = 0; i < 256; i++) {
1929 fi = i;
1930 if (rfract >= 0) {
1931 rlut[i] = (l_int32)(fi + (255.0 - fi) * rfract);
1932 } else {
1933 rlut[i] = (l_int32)(fi * (1.0 + rfract));
1934 }
1935 if (gfract >= 0) {
1936 glut[i] = (l_int32)(fi + (255.0 - fi) * gfract);
1937 } else {
1938 glut[i] = (l_int32)(fi * (1.0 + gfract));
1939 }
1940 if (bfract >= 0) {
1941 blut[i] = (l_int32)(fi + (255.0 - fi) * bfract);
1942 } else {
1943 blut[i] = (l_int32)(fi * (1.0 + bfract));
1944 }
1945 }
1946
1947 pixGetDimensions(pixs, &w, &h, NULL);
1948 datas = pixGetData(pixs);
1949 wpls = pixGetWpl(pixs);
1950 pixd = pixCreate(w, h, 32);
1951 datad = pixGetData(pixd);
1952 wpld = pixGetWpl(pixd);
1953 for (i = 0; i < h; i++) {
1954 lines = datas + i * wpls;
1955 lined = datad + i * wpld;
1956 for (j = 0; j < w; j++) {
1957 extractRGBValues(lines[j], &rval, &gval, &bval);
1958 composeRGBPixel(rlut[rval], glut[gval], blut[bval], lined + j);
1959 }
1960 }
1961
1962 LEPT_FREE(rlut);
1963 LEPT_FREE(glut);
1964 LEPT_FREE(blut);
1965 return pixd;
1966 }
1967
1968 /*-----------------------------------------------------------------------*
1969 * Darken gray (unsaturated) pixels
1970 *-----------------------------------------------------------------------*/
1971 /*!
1972 * \brief pixDarkenGray()
1973 *
1974 * \param[in] pixd [optional] can be null or equal to pixs
1975 * \param[in] pixs 32 bpp rgb
1976 * \param[in] thresh pixels with max component >= %thresh are unchanged
1977 * \param[in] satlimit pixels with saturation >= %satlimit are unchanged
1978 * \return pixd, or NULL on error
1979 *
1980 * <pre>
1981 * Notes:
1982 * (1) This darkens gray pixels, by a fraction (sat/%satlimit), where
1983 * the saturation, sat, is the component difference (max - min).
1984 * The pixel value is unchanged if sat >= %satlimit. A typical
1985 * value of %satlimit might be 40; the larger the value, the
1986 * more that pixels with a smaller saturation will be darkened.
1987 * (2) Pixels with max component >= %thresh are unchanged. This can be
1988 * used to prevent bright pixels with low saturation from being
1989 * darkened. Setting thresh == 0 is a no-op; setting %thresh == 255
1990 * causes the darkening to be applied to all pixels.
1991 * (3) This function is useful to enhance pixels relative to a
1992 * gray background.
1993 * (4) A related function that builds a 1 bpp mask over the gray
1994 * pixels is pixMaskOverGrayPixels().
1995 * </pre>
1996 */
1997 PIX *
1998 pixDarkenGray(PIX *pixd,
1999 PIX *pixs,
2000 l_int32 thresh,
2001 l_int32 satlimit)
2002 {
2003 l_int32 w, h, i, j, wpls, wpld;
2004 l_int32 rval, gval, bval, minrg, min, maxrg, max, sat;
2005 l_uint32 *datas, *datad, *lines, *lined;
2006 l_float32 ratio;
2007
2008 if (!pixs || pixGetDepth(pixs) != 32)
2009 return (PIX *)ERROR_PTR("pixs undefined or not 32 bpp", __func__, NULL);
2010 if (thresh < 0 || thresh > 255)
2011 return (PIX *)ERROR_PTR("invalid thresh", __func__, NULL);
2012 if (satlimit < 1)
2013 return (PIX *)ERROR_PTR("invalid satlimit", __func__, NULL);
2014 if (pixd && (pixs != pixd))
2015 return (PIX *)ERROR_PTR("not new or in-place", __func__, NULL);
2016
2017 pixGetDimensions(pixs, &w, &h, NULL);
2018 datas = pixGetData(pixs);
2019 wpls = pixGetWpl(pixs);
2020 if ((pixd = pixCopy(pixd, pixs)) == NULL)
2021 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
2022 datad = pixGetData(pixd);
2023 wpld = pixGetWpl(pixd);
2024
2025 for (i = 0; i < h; i++) {
2026 lines = datas + i * wpls;
2027 lined = datad + i * wpld;
2028 for (j = 0; j < w; j++) {
2029 extractRGBValues(lines[j], &rval, &gval, &bval);
2030 minrg = L_MIN(rval, gval);
2031 min = L_MIN(minrg, bval);
2032 maxrg = L_MAX(rval, gval);
2033 max = L_MAX(maxrg, bval);
2034 sat = max - min;
2035 if (max >= thresh || sat >= satlimit)
2036 continue;
2037 ratio = (l_float32)sat / (l_float32)satlimit;
2038 composeRGBPixel((l_int32)(ratio * rval), (l_int32)(ratio * gval),
2039 (l_int32)(ratio * bval), &lined[j]);
2040 }
2041 }
2042 return pixd;
2043 }
2044
2045
2046 /*-----------------------------------------------------------------------*
2047 * General multiplicative constant color transform *
2048 *-----------------------------------------------------------------------*/
2049 /*!
2050 * \brief pixMultConstantColor()
2051 *
2052 * \param[in] pixs colormapped or rgb
2053 * \param[in] rfact red multiplicative factor
2054 * \param[in] gfact green multiplicative factor
2055 * \param[in] bfact blue multiplicative factor
2056 * \return pixd colormapped or rgb, with colors scaled, or NULL on error
2057 *
2058 * <pre>
2059 * Notes:
2060 * (1) rfact, gfact and bfact can only have non-negative values.
2061 * They can be greater than 1.0. All transformed component
2062 * values are clipped to the interval [0, 255].
2063 * (2) For multiplication with a general 3x3 matrix of constants,
2064 * use pixMultMatrixColor().
2065 * </pre>
2066 */
2067 PIX *
2068 pixMultConstantColor(PIX *pixs,
2069 l_float32 rfact,
2070 l_float32 gfact,
2071 l_float32 bfact)
2072 {
2073 l_int32 i, j, w, h, d, wpls, wpld;
2074 l_int32 ncolors, rval, gval, bval, nrval, ngval, nbval;
2075 l_uint32 nval;
2076 l_uint32 *datas, *datad, *lines, *lined;
2077 PIX *pixd;
2078 PIXCMAP *cmap;
2079
2080 if (!pixs)
2081 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
2082 pixGetDimensions(pixs, &w, &h, &d);
2083 cmap = pixGetColormap(pixs);
2084 if (!cmap && d != 32)
2085 return (PIX *)ERROR_PTR("pixs not cmapped or 32 bpp", __func__, NULL);
2086 rfact = L_MAX(0.0, rfact);
2087 gfact = L_MAX(0.0, gfact);
2088 bfact = L_MAX(0.0, bfact);
2089
2090 if (cmap) {
2091 if ((pixd = pixCopy(NULL, pixs)) == NULL)
2092 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
2093 cmap = pixGetColormap(pixd);
2094 ncolors = pixcmapGetCount(cmap);
2095 for (i = 0; i < ncolors; i++) {
2096 pixcmapGetColor(cmap, i, &rval, &gval, &bval);
2097 nrval = (l_int32)(rfact * rval);
2098 ngval = (l_int32)(gfact * gval);
2099 nbval = (l_int32)(bfact * bval);
2100 nrval = L_MIN(255, nrval);
2101 ngval = L_MIN(255, ngval);
2102 nbval = L_MIN(255, nbval);
2103 pixcmapResetColor(cmap, i, nrval, ngval, nbval);
2104 }
2105 return pixd;
2106 }
2107
2108 if ((pixd = pixCreateTemplate(pixs)) == NULL)
2109 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
2110 datas = pixGetData(pixs);
2111 datad = pixGetData(pixd);
2112 wpls = pixGetWpl(pixs);
2113 wpld = pixGetWpl(pixd);
2114 for (i = 0; i < h; i++) {
2115 lines = datas + i * wpls;
2116 lined = datad + i * wpld;
2117 for (j = 0; j < w; j++) {
2118 extractRGBValues(lines[j], &rval, &gval, &bval);
2119 nrval = (l_int32)(rfact * rval);
2120 ngval = (l_int32)(gfact * gval);
2121 nbval = (l_int32)(bfact * bval);
2122 nrval = L_MIN(255, nrval);
2123 ngval = L_MIN(255, ngval);
2124 nbval = L_MIN(255, nbval);
2125 composeRGBPixel(nrval, ngval, nbval, &nval);
2126 *(lined + j) = nval;
2127 }
2128 }
2129
2130 return pixd;
2131 }
2132
2133
2134 /*!
2135 * \brief pixMultMatrixColor()
2136 *
2137 * \param[in] pixs colormapped or rgb
2138 * \param[in] kel kernel 3x3 matrix of floats
2139 * \return pixd colormapped or rgb, or NULL on error
2140 *
2141 * <pre>
2142 * Notes:
2143 * (1) The kernel is a data structure used mostly for floating point
2144 * convolution. Here it is a 3x3 matrix of floats that are used
2145 * to transform the pixel values by matrix multiplication:
2146 * nrval = a[0,0] * rval + a[0,1] * gval + a[0,2] * bval
2147 * ngval = a[1,0] * rval + a[1,1] * gval + a[1,2] * bval
2148 * nbval = a[2,0] * rval + a[2,1] * gval + a[2,2] * bval
2149 * (2) The matrix can be generated in several ways.
2150 * See kernel.c for details. Here are two of them:
2151 * (a) kel = kernelCreate(3, 3);
2152 * kernelSetElement(kel, 0, 0, val00);
2153 * kernelSetElement(kel, 0, 1, val01);
2154 * ...
2155 * (b) from a static string; e.g.,:
2156 * const char *kdata = " 0.6 0.3 -0.2 "
2157 * " 0.1 1.2 0.4 "
2158 * " -0.4 0.2 0.9 ";
2159 * kel = kernelCreateFromString(3, 3, 0, 0, kdata);
2160 * (3) For the special case where the matrix is diagonal, it is easier
2161 * to use pixMultConstantColor().
2162 * (4) Matrix entries can have positive and negative values, and can
2163 * be larger than 1.0. All transformed component values
2164 * are clipped to [0, 255].
2165 * </pre>
2166 */
2167 PIX *
2168 pixMultMatrixColor(PIX *pixs,
2169 L_KERNEL *kel)
2170 {
2171 l_int32 i, j, index, kw, kh, w, h, d, wpls, wpld;
2172 l_int32 ncolors, rval, gval, bval, nrval, ngval, nbval;
2173 l_uint32 nval;
2174 l_uint32 *datas, *datad, *lines, *lined;
2175 l_float32 v[9]; /* use linear array for convenience */
2176 PIX *pixd;
2177 PIXCMAP *cmap;
2178
2179 if (!pixs)
2180 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
2181 if (!kel)
2182 return (PIX *)ERROR_PTR("kel not defined", __func__, NULL);
2183 kernelGetParameters(kel, &kw, &kh, NULL, NULL);
2184 if (kw != 3 || kh != 3)
2185 return (PIX *)ERROR_PTR("matrix not 3x3", __func__, NULL);
2186 pixGetDimensions(pixs, &w, &h, &d);
2187 cmap = pixGetColormap(pixs);
2188 if (!cmap && d != 32)
2189 return (PIX *)ERROR_PTR("pixs not cmapped or 32 bpp", __func__, NULL);
2190
2191 for (i = 0, index = 0; i < 3; i++)
2192 for (j = 0; j < 3; j++, index++)
2193 kernelGetElement(kel, i, j, v + index);
2194
2195 if (cmap) {
2196 if ((pixd = pixCopy(NULL, pixs)) == NULL)
2197 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
2198 cmap = pixGetColormap(pixd);
2199 ncolors = pixcmapGetCount(cmap);
2200 for (i = 0; i < ncolors; i++) {
2201 pixcmapGetColor(cmap, i, &rval, &gval, &bval);
2202 nrval = (l_int32)(v[0] * rval + v[1] * gval + v[2] * bval);
2203 ngval = (l_int32)(v[3] * rval + v[4] * gval + v[5] * bval);
2204 nbval = (l_int32)(v[6] * rval + v[7] * gval + v[8] * bval);
2205 nrval = L_MAX(0, L_MIN(255, nrval));
2206 ngval = L_MAX(0, L_MIN(255, ngval));
2207 nbval = L_MAX(0, L_MIN(255, nbval));
2208 pixcmapResetColor(cmap, i, nrval, ngval, nbval);
2209 }
2210 return pixd;
2211 }
2212
2213 if ((pixd = pixCreateTemplate(pixs)) == NULL)
2214 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
2215 datas = pixGetData(pixs);
2216 datad = pixGetData(pixd);
2217 wpls = pixGetWpl(pixs);
2218 wpld = pixGetWpl(pixd);
2219 for (i = 0; i < h; i++) {
2220 lines = datas + i * wpls;
2221 lined = datad + i * wpld;
2222 for (j = 0; j < w; j++) {
2223 extractRGBValues(lines[j], &rval, &gval, &bval);
2224 nrval = (l_int32)(v[0] * rval + v[1] * gval + v[2] * bval);
2225 ngval = (l_int32)(v[3] * rval + v[4] * gval + v[5] * bval);
2226 nbval = (l_int32)(v[6] * rval + v[7] * gval + v[8] * bval);
2227 nrval = L_MAX(0, L_MIN(255, nrval));
2228 ngval = L_MAX(0, L_MIN(255, ngval));
2229 nbval = L_MAX(0, L_MIN(255, nbval));
2230 composeRGBPixel(nrval, ngval, nbval, &nval);
2231 *(lined + j) = nval;
2232 }
2233 }
2234
2235 return pixd;
2236 }
2237
2238
2239 /*-------------------------------------------------------------*
2240 * Half-edge by bandpass *
2241 *-------------------------------------------------------------*/
2242 /*!
2243 * \brief pixHalfEdgeByBandpass()
2244 *
2245 * \param[in] pixs 8 bpp gray or 32 bpp rgb
2246 * \param[in] sm1h, sm1v "half-widths" of smoothing filter sm1
2247 * \param[in] sm2h, sm2v "half-widths" of smoothing filter sm2;
2248 * require sm2 != sm1
2249 * \return pixd, or NULL on error
2250 *
2251 * <pre>
2252 * Notes:
2253 * (1) We use symmetric smoothing filters of odd dimension,
2254 * typically use 3, 5, 7, etc. The smoothing parameters
2255 * for these are 1, 2, 3, etc. The filter size is related
2256 * to the smoothing parameter by
2257 * size = 2 * smoothing + 1
2258 * (2) Because we take the difference of two lowpass filters,
2259 * this is actually a bandpass filter.
2260 * (3) We allow both filters to be anisotropic.
2261 * (4) Consider either the h or v component of the 2 filters.
2262 * Depending on whether sm1 > sm2 or sm2 > sm1, we get
2263 * different halves of the smoothed gradients (or "edges").
2264 * This difference of smoothed signals looks more like
2265 * a second derivative of a transition, which we rectify
2266 * by not allowing the signal to go below zero. If sm1 < sm2,
2267 * the sm2 transition is broader, so the difference between
2268 * sm1 and sm2 signals is positive on the upper half of
2269 * the transition. Likewise, if sm1 > sm2, the sm1 - sm2
2270 * signal difference is positive on the lower half of
2271 * the transition.
2272 * </pre>
2273 */
2274 PIX *
2275 pixHalfEdgeByBandpass(PIX *pixs,
2276 l_int32 sm1h,
2277 l_int32 sm1v,
2278 l_int32 sm2h,
2279 l_int32 sm2v)
2280 {
2281 l_int32 d;
2282 PIX *pixg, *pixacc, *pixc1, *pixc2;
2283
2284 if (!pixs)
2285 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
2286 if (sm1h == sm2h && sm1v == sm2v)
2287 return (PIX *)ERROR_PTR("sm2 = sm1", __func__, NULL);
2288 d = pixGetDepth(pixs);
2289 if (d != 8 && d != 32)
2290 return (PIX *)ERROR_PTR("pixs not 8 or 32 bpp", __func__, NULL);
2291 if (d == 32)
2292 pixg = pixConvertRGBToLuminance(pixs);
2293 else /* d == 8 */
2294 pixg = pixClone(pixs);
2295
2296 /* Make a convolution accumulator and use it twice */
2297 if ((pixacc = pixBlockconvAccum(pixg)) == NULL) {
2298 pixDestroy(&pixg);
2299 return (PIX *)ERROR_PTR("pixacc not made", __func__, NULL);
2300 }
2301 if ((pixc1 = pixBlockconvGray(pixg, pixacc, sm1h, sm1v)) == NULL) {
2302 pixDestroy(&pixg);
2303 pixDestroy(&pixacc);
2304 return (PIX *)ERROR_PTR("pixc1 not made", __func__, NULL);
2305 }
2306 pixc2 = pixBlockconvGray(pixg, pixacc, sm2h, sm2v);
2307 pixDestroy(&pixg);
2308 pixDestroy(&pixacc);
2309 if (!pixc2) {
2310 pixDestroy(&pixc1);
2311 return (PIX *)ERROR_PTR("pixc2 not made", __func__, NULL);
2312 }
2313
2314 /* Compute the half-edge using pixc1 - pixc2. */
2315 pixSubtractGray(pixc1, pixc1, pixc2);
2316 pixDestroy(&pixc2);
2317 return pixc1;
2318 }