Mercurial > hgrepos > Python2 > PyMuPDF

comparison mupdf-source/thirdparty/leptonica/src/pixarith.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 pixarith.c
29 * <pre>
30 *
31 * One-image grayscale arithmetic operations (8, 16, 32 bpp)
32 * l_int32 pixAddConstantGray()
33 * l_int32 pixMultConstantGray()
34 *
35 * Two-image grayscale arithmetic operations (8, 16, 32 bpp)
36 * PIX *pixAddGray()
37 * PIX *pixSubtractGray()
38 * PIX *pixMultiplyGray()
39 *
40 * Grayscale threshold operation (8, 16, 32 bpp)
41 * PIX *pixThresholdToValue()
42 *
43 * Image accumulator arithmetic operations
44 * PIX *pixInitAccumulate()
45 * PIX *pixFinalAccumulate()
46 * PIX *pixFinalAccumulateThreshold()
47 * l_int32 pixAccumulate()
48 * l_int32 pixMultConstAccumulate()
49 *
50 * Absolute value of difference
51 * PIX *pixAbsDifference()
52 *
53 * Sum of color images
54 * PIX *pixAddRGB()
55 *
56 * Two-image min and max operations (8 and 16 bpp)
57 * PIX *pixMinOrMax()
58 *
59 * Scale pix for maximum dynamic range
60 * PIX *pixMaxDynamicRange()
61 * PIX *pixMaxDynamicRangeRGB()
62 *
63 * RGB pixel value scaling
64 * l_uint32 linearScaleRGBVal()
65 * l_uint32 logScaleRGBVal()
66 *
67 * Log base2 lookup
68 * l_float32 *makeLogBase2Tab()
69 * l_float32 getLogBase2()
70 *
71 * The image accumulator operations are used when you expect
72 * overflow from 8 bits on intermediate results. For example,
73 * you might want a tophat contrast operator which is
74 * 3*I - opening(I,S) - closing(I,S)
75 * To use these operations, first use the init to generate
76 * a 16 bpp image, use the accumulate to add or subtract 8 bpp
77 * images from that, or the multiply constant to multiply
78 * by a small constant (much less than 256 -- we don't want
79 * overflow from the 16 bit images!), and when you're finished
80 * use final to bring the result back to 8 bpp, clipped
81 * if necessary. There is also a divide function, which
82 * can be used to divide one image by another, scaling the
83 * result for maximum dynamic range, and giving back the
84 * 8 bpp result.
85 *
86 * A simpler interface to the arithmetic operations is
87 * provided in pixacc.c.
88 * </pre>
89 */
90
91 #ifdef HAVE_CONFIG_H
92 #include <config_auto.h>
93 #endif /* HAVE_CONFIG_H */
94
95 #include <string.h>
96 #include <math.h>
97 #include "allheaders.h"
98
99 /*-------------------------------------------------------------*
100 * One-image grayscale arithmetic operations *
101 *-------------------------------------------------------------*/
102 /*!
103 * \brief pixAddConstantGray()
104 *
105 * \param[in] pixs 8, 16 or 32 bpp
106 * \param[in] val amount to add to each pixel
107 * \return 0 if OK, 1 on error
108 *
109 * <pre>
110 * Notes:
111 * (1) In-place operation.
112 * (2) No clipping for 32 bpp.
113 * (3) For 8 and 16 bpp, if val > 0 the result is clipped
114 * to 0xff and 0xffff, rsp.
115 * (4) For 8 and 16 bpp, if val < 0 the result is clipped to 0.
116 * </pre>
117 */
118 l_ok
119 pixAddConstantGray(PIX *pixs,
120 l_int32 val)
121 {
122 l_int32 i, j, w, h, d, wpl, pval;
123 l_uint32 *data, *line;
124
125 if (!pixs)
126 return ERROR_INT("pixs not defined", __func__, 1);
127 pixGetDimensions(pixs, &w, &h, &d);
128 if (d != 8 && d != 16 && d != 32)
129 return ERROR_INT("pixs not 8, 16 or 32 bpp", __func__, 1);
130
131 data = pixGetData(pixs);
132 wpl = pixGetWpl(pixs);
133 for (i = 0; i < h; i++) {
134 line = data + i * wpl;
135 if (d == 8) {
136 if (val < 0) {
137 for (j = 0; j < w; j++) {
138 pval = GET_DATA_BYTE(line, j);
139 pval = L_MAX(0, pval + val);
140 SET_DATA_BYTE(line, j, pval);
141 }
142 } else { /* val >= 0 */
143 for (j = 0; j < w; j++) {
144 pval = GET_DATA_BYTE(line, j);
145 pval = L_MIN(255, pval + val);
146 SET_DATA_BYTE(line, j, pval);
147 }
148 }
149 } else if (d == 16) {
150 if (val < 0) {
151 for (j = 0; j < w; j++) {
152 pval = GET_DATA_TWO_BYTES(line, j);
153 pval = L_MAX(0, pval + val);
154 SET_DATA_TWO_BYTES(line, j, pval);
155 }
156 } else { /* val >= 0 */
157 for (j = 0; j < w; j++) {
158 pval = GET_DATA_TWO_BYTES(line, j);
159 pval = L_MIN(0xffff, pval + val);
160 SET_DATA_TWO_BYTES(line, j, pval);
161 }
162 }
163 } else { /* d == 32; no check for overflow (< 0 or > 0xffffffff) */
164 for (j = 0; j < w; j++)
165 *(line + j) += val;
166 }
167 }
168
169 return 0;
170 }
171
172
173 /*!
174 * \brief pixMultConstantGray()
175 *
176 * \param[in] pixs 8, 16 or 32 bpp
177 * \param[in] val >= 0.0; amount to multiply by each pixel
178 * \return 0 if OK, 1 on error
179 *
180 * <pre>
181 * Notes:
182 * (1) In-place operation; val must be >= 0.
183 * (2) No clipping for 32 bpp.
184 * (3) For 8 and 16 bpp, the result is clipped to 0xff and 0xffff, rsp.
185 * </pre>
186 */
187 l_ok
188 pixMultConstantGray(PIX *pixs,
189 l_float32 val)
190 {
191 l_int32 i, j, w, h, d, wpl, pval;
192 l_uint32 upval;
193 l_uint32 *data, *line;
194
195 if (!pixs)
196 return ERROR_INT("pixs not defined", __func__, 1);
197 pixGetDimensions(pixs, &w, &h, &d);
198 if (d != 8 && d != 16 && d != 32)
199 return ERROR_INT("pixs not 8, 16 or 32 bpp", __func__, 1);
200 if (val < 0.0)
201 return ERROR_INT("val < 0.0", __func__, 1);
202
203 data = pixGetData(pixs);
204 wpl = pixGetWpl(pixs);
205 for (i = 0; i < h; i++) {
206 line = data + i * wpl;
207 if (d == 8) {
208 for (j = 0; j < w; j++) {
209 pval = GET_DATA_BYTE(line, j);
210 pval = (l_int32)(val * pval);
211 pval = L_MIN(255, pval);
212 SET_DATA_BYTE(line, j, pval);
213 }
214 } else if (d == 16) {
215 for (j = 0; j < w; j++) {
216 pval = GET_DATA_TWO_BYTES(line, j);
217 pval = (l_int32)(val * pval);
218 pval = L_MIN(0xffff, pval);
219 SET_DATA_TWO_BYTES(line, j, pval);
220 }
221 } else { /* d == 32; no clipping */
222 for (j = 0; j < w; j++) {
223 upval = *(line + j);
224 upval = (l_uint32)(val * upval);
225 *(line + j) = upval;
226 }
227 }
228 }
229
230 return 0;
231 }
232
233
234 /*-------------------------------------------------------------*
235 * Two-image grayscale arithmetic ops *
236 *-------------------------------------------------------------*/
237 /*!
238 * \brief pixAddGray()
239 *
240 * \param[in] pixd [optional]; this can be null, equal to pixs1, or
241 * different from pixs1
242 * \param[in] pixs1 can be equal to pixd
243 * \param[in] pixs2
244 * \return pixd always
245 *
246 * <pre>
247 * Notes:
248 * (1) Arithmetic addition of two 8, 16 or 32 bpp images.
249 * (2) For 8 and 16 bpp, we do explicit clipping to 0xff and 0xffff,
250 * respectively.
251 * (3) Alignment is to UL corner.
252 * (4) There are 3 cases. The result can go to a new dest,
253 * in-place to pixs1, or to an existing input dest:
254 * * pixd == null: (src1 + src2) --> new pixd
255 * * pixd == pixs1: (src1 + src2) --> src1 (in-place)
256 * * pixd != pixs1: (src1 + src2) --> input pixd
257 * (5) pixs2 must be different from both pixd and pixs1.
258 * </pre>
259 */
260 PIX *
261 pixAddGray(PIX *pixd,
262 PIX *pixs1,
263 PIX *pixs2)
264 {
265 l_int32 i, j, d, ws, hs, w, h, wpls, wpld, val, sum;
266 l_uint32 *datas, *datad, *lines, *lined;
267
268 if (!pixs1)
269 return (PIX *)ERROR_PTR("pixs1 not defined", __func__, pixd);
270 if (!pixs2)
271 return (PIX *)ERROR_PTR("pixs2 not defined", __func__, pixd);
272 if (pixs2 == pixs1)
273 return (PIX *)ERROR_PTR("pixs2 and pixs1 must differ", __func__, pixd);
274 if (pixs2 == pixd)
275 return (PIX *)ERROR_PTR("pixs2 and pixd must differ", __func__, pixd);
276 d = pixGetDepth(pixs1);
277 if (d != 8 && d != 16 && d != 32)
278 return (PIX *)ERROR_PTR("pix are not 8, 16 or 32 bpp", __func__, pixd);
279 if (pixGetDepth(pixs2) != d)
280 return (PIX *)ERROR_PTR("depths differ (pixs1, pixs2)", __func__, pixd);
281 if (pixd && (pixGetDepth(pixd) != d))
282 return (PIX *)ERROR_PTR("depths differ (pixs1, pixd)", __func__, pixd);
283
284 if (!pixSizesEqual(pixs1, pixs2))
285 L_WARNING("pixs1 and pixs2 not equal in size\n", __func__);
286 if (pixd && !pixSizesEqual(pixs1, pixd))
287 L_WARNING("pixs1 and pixd not equal in size\n", __func__);
288
289 if (pixs1 != pixd)
290 pixd = pixCopy(pixd, pixs1);
291
292 /* pixd + pixs2 ==> pixd */
293 datas = pixGetData(pixs2);
294 datad = pixGetData(pixd);
295 wpls = pixGetWpl(pixs2);
296 wpld = pixGetWpl(pixd);
297 pixGetDimensions(pixs2, &ws, &hs, NULL);
298 pixGetDimensions(pixd, &w, &h, NULL);
299 w = L_MIN(ws, w);
300 h = L_MIN(hs, h);
301 for (i = 0; i < h; i++) {
302 lined = datad + i * wpld;
303 lines = datas + i * wpls;
304 if (d == 8) {
305 for (j = 0; j < w; j++) {
306 sum = GET_DATA_BYTE(lines, j) + GET_DATA_BYTE(lined, j);
307 val = L_MIN(sum, 255);
308 SET_DATA_BYTE(lined, j, val);
309 }
310 } else if (d == 16) {
311 for (j = 0; j < w; j++) {
312 sum = GET_DATA_TWO_BYTES(lines, j)
313 + GET_DATA_TWO_BYTES(lined, j);
314 val = L_MIN(sum, 0xffff);
315 SET_DATA_TWO_BYTES(lined, j, val);
316 }
317 } else { /* d == 32; no clipping */
318 for (j = 0; j < w; j++)
319 *(lined + j) += *(lines + j);
320 }
321 }
322
323 return pixd;
324 }
325
326
327 /*!
328 * \brief pixSubtractGray()
329 *
330 * \param[in] pixd [optional]; this can be null, equal to pixs1, or
331 * different from pixs1
332 * \param[in] pixs1 can be equal to pixd
333 * \param[in] pixs2
334 * \return pixd always
335 *
336 * <pre>
337 * Notes:
338 * (1) Arithmetic subtraction of two 8, 16 or 32 bpp images.
339 * (2) Source pixs2 is always subtracted from source pixs1.
340 * (3) Do explicit clipping to 0.
341 * (4) Alignment is to UL corner.
342 * (5) There are 3 cases. The result can go to a new dest,
343 * in-place to pixs1, or to an existing input dest:
344 * (a) pixd == null (src1 - src2) --> new pixd
345 * (b) pixd == pixs1 (src1 - src2) --> src1 (in-place)
346 * (d) pixd != pixs1 (src1 - src2) --> input pixd
347 * (6) pixs2 must be different from both pixd and pixs1.
348 * </pre>
349 */
350 PIX *
351 pixSubtractGray(PIX *pixd,
352 PIX *pixs1,
353 PIX *pixs2)
354 {
355 l_int32 i, j, w, h, ws, hs, d, wpls, wpld, val, diff;
356 l_uint32 *datas, *datad, *lines, *lined;
357
358 if (!pixs1)
359 return (PIX *)ERROR_PTR("pixs1 not defined", __func__, pixd);
360 if (!pixs2)
361 return (PIX *)ERROR_PTR("pixs2 not defined", __func__, pixd);
362 if (pixs2 == pixs1)
363 return (PIX *)ERROR_PTR("pixs2 and pixs1 must differ", __func__, pixd);
364 if (pixs2 == pixd)
365 return (PIX *)ERROR_PTR("pixs2 and pixd must differ", __func__, pixd);
366 d = pixGetDepth(pixs1);
367 if (d != 8 && d != 16 && d != 32)
368 return (PIX *)ERROR_PTR("pix are not 8, 16 or 32 bpp", __func__, pixd);
369 if (pixGetDepth(pixs2) != d)
370 return (PIX *)ERROR_PTR("depths differ (pixs1, pixs2)", __func__, pixd);
371 if (pixd && (pixGetDepth(pixd) != d))
372 return (PIX *)ERROR_PTR("depths differ (pixs1, pixd)", __func__, pixd);
373
374 if (!pixSizesEqual(pixs1, pixs2))
375 L_WARNING("pixs1 and pixs2 not equal in size\n", __func__);
376 if (pixd && !pixSizesEqual(pixs1, pixd))
377 L_WARNING("pixs1 and pixd not equal in size\n", __func__);
378
379 if (pixs1 != pixd)
380 pixd = pixCopy(pixd, pixs1);
381
382 /* pixd - pixs2 ==> pixd */
383 datas = pixGetData(pixs2);
384 datad = pixGetData(pixd);
385 wpls = pixGetWpl(pixs2);
386 wpld = pixGetWpl(pixd);
387 pixGetDimensions(pixs2, &ws, &hs, NULL);
388 pixGetDimensions(pixd, &w, &h, NULL);
389 w = L_MIN(ws, w);
390 h = L_MIN(hs, h);
391 for (i = 0; i < h; i++) {
392 lined = datad + i * wpld;
393 lines = datas + i * wpls;
394 if (d == 8) {
395 for (j = 0; j < w; j++) {
396 diff = GET_DATA_BYTE(lined, j) - GET_DATA_BYTE(lines, j);
397 val = L_MAX(diff, 0);
398 SET_DATA_BYTE(lined, j, val);
399 }
400 } else if (d == 16) {
401 for (j = 0; j < w; j++) {
402 diff = GET_DATA_TWO_BYTES(lined, j)
403 - GET_DATA_TWO_BYTES(lines, j);
404 val = L_MAX(diff, 0);
405 SET_DATA_TWO_BYTES(lined, j, val);
406 }
407 } else { /* d == 32; no clipping */
408 for (j = 0; j < w; j++)
409 *(lined + j) -= *(lines + j);
410 }
411 }
412
413 return pixd;
414 }
415
416
417 /*!
418 * \brief pixMultiplyGray()
419 *
420 * \param[in] pixs 32 bpp rgb or 8 bpp gray
421 * \param[in] pixg 8 bpp gray
422 * \param[in] norm multiplicative factor to avoid overflow; 0 for default
423 * \return pixd, or null on error
424 *
425 * <pre>
426 * Notes:
427 * (1) This function can be used for correcting a scanned image
428 * under non-uniform illumination. For that application,
429 * %pixs is the scanned image, %pixg is an image whose values
430 * are inversely related to light from a uniform (say, white)
431 * target, and %norm is typically the inverse of the maximum
432 * pixel value in %pixg.
433 * (2) Set norm = 0 to get the default value, which is the inverse
434 * of the max value in %pixg. This avoids overflow in the product.
435 * (3) For 32 bpp %pixs, all 3 components are multiplied by the
436 * same number.
437 * (4) Alignment is to UL corner.
438 * </pre>
439 */
440 PIX *
441 pixMultiplyGray(PIX *pixs,
442 PIX *pixg,
443 l_float32 norm)
444 {
445 l_int32 i, j, w, h, d, ws, hs, ds, wpls, wplg, wpld;
446 l_int32 rval, gval, bval, rval2, gval2, bval2, vals, valg, val, maxgray;
447 l_uint32 val32;
448 l_uint32 *datas, *datag, *datad, *lines, *lineg, *lined;
449 PIX *pixd;
450
451 if (!pixs)
452 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
453 pixGetDimensions(pixs, &ws, &hs, &ds);
454 if (ds != 8 && ds != 32)
455 return (PIX *)ERROR_PTR("pixs not 8 or 32 bpp", __func__, NULL);
456 if (!pixg)
457 return (PIX *)ERROR_PTR("pixg not defined", __func__, NULL);
458 pixGetDimensions(pixg, &w, &h, &d);
459 if (d != 8)
460 return (PIX *)ERROR_PTR("pixg not 8 bpp", __func__, NULL);
461
462 if (norm <= 0.0) {
463 pixGetExtremeValue(pixg, 1, L_SELECT_MAX, NULL, NULL, NULL, &maxgray);
464 norm = (maxgray > 0) ? 1.0f / (l_float32)maxgray : 1.0f;
465 }
466
467 if ((pixd = pixCreateTemplate(pixs)) == NULL)
468 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
469 datas = pixGetData(pixs);
470 datag = pixGetData(pixg);
471 datad = pixGetData(pixd);
472 wpls = pixGetWpl(pixs);
473 wplg = pixGetWpl(pixg);
474 wpld = pixGetWpl(pixd);
475 w = L_MIN(ws, w);
476 h = L_MIN(hs, h);
477 for (i = 0; i < h; i++) {
478 lines = datas + i * wpls;
479 lineg = datag + i * wplg;
480 lined = datad + i * wpld;
481 if (ds == 8) {
482 for (j = 0; j < w; j++) {
483 vals = GET_DATA_BYTE(lines, j);
484 valg = GET_DATA_BYTE(lineg, j);
485 val = (l_int32)(vals * valg * norm + 0.5);
486 val = L_MIN(255, val);
487 SET_DATA_BYTE(lined, j, val);
488 }
489 } else { /* ds == 32 */
490 for (j = 0; j < w; j++) {
491 val32 = *(lines + j);
492 extractRGBValues(val32, &rval, &gval, &bval);
493 valg = GET_DATA_BYTE(lineg, j);
494 rval2 = (l_int32)(rval * valg * norm + 0.5);
495 rval2 = L_MIN(255, rval2);
496 gval2 = (l_int32)(gval * valg * norm + 0.5);
497 gval2 = L_MIN(255, gval2);
498 bval2 = (l_int32)(bval * valg * norm + 0.5);
499 bval2 = L_MIN(255, bval2);
500 composeRGBPixel(rval2, gval2, bval2, lined + j);
501 }
502 }
503 }
504
505 return pixd;
506 }
507
508
509 /*-------------------------------------------------------------*
510 * Grayscale threshold operation *
511 *-------------------------------------------------------------*/
512 /*!
513 * \brief pixThresholdToValue()
514 *
515 * \param[in] pixd [optional]; if not null, must be equal to pixs
516 * \param[in] pixs 8, 16, 32 bpp
517 * \param[in] threshval
518 * \param[in] setval
519 * \return pixd always
520 *
521 * <pre>
522 * Notes:
523 * ~ operation can be in-place (pixs == pixd) or to a new pixd
524 * ~ if %setval > %threshval, sets pixels with a value >= threshval to setval
525 * ~ if %setval < %threshval, sets pixels with a value <= threshval to setval
526 * ~ if %setval == %threshval, no-op
527 * </pre>
528 */
529 PIX *
530 pixThresholdToValue(PIX *pixd,
531 PIX *pixs,
532 l_int32 threshval,
533 l_int32 setval)
534 {
535 l_int32 i, j, w, h, d, wpld, setabove;
536 l_uint32 *datad, *lined;
537
538 if (!pixs)
539 return (PIX *)ERROR_PTR("pixs not defined", __func__, pixd);
540 d = pixGetDepth(pixs);
541 if (d != 8 && d != 16 && d != 32)
542 return (PIX *)ERROR_PTR("pixs not 8, 16 or 32 bpp", __func__, pixd);
543 if (pixd && (pixs != pixd))
544 return (PIX *)ERROR_PTR("pixd exists and is not pixs", __func__, pixd);
545 if (threshval < 0 || setval < 0)
546 return (PIX *)ERROR_PTR("threshval & setval not < 0", __func__, pixd);
547 if (d == 8 && setval > 255)
548 return (PIX *)ERROR_PTR("setval > 255 for 8 bpp", __func__, pixd);
549 if (d == 16 && setval > 0xffff)
550 return (PIX *)ERROR_PTR("setval > 0xffff for 16 bpp", __func__, pixd);
551
552 if (!pixd)
553 pixd = pixCopy(NULL, pixs);
554 if (setval == threshval) {
555 L_WARNING("setval == threshval; no operation\n", __func__);
556 return pixd;
557 }
558
559 datad = pixGetData(pixd);
560 pixGetDimensions(pixd, &w, &h, NULL);
561 wpld = pixGetWpl(pixd);
562 if (setval > threshval)
563 setabove = TRUE;
564 else
565 setabove = FALSE;
566
567 for (i = 0; i < h; i++) {
568 lined = datad + i * wpld;
569 if (setabove == TRUE) {
570 if (d == 8) {
571 for (j = 0; j < w; j++) {
572 if (GET_DATA_BYTE(lined, j) - threshval >= 0)
573 SET_DATA_BYTE(lined, j, setval);
574 }
575 } else if (d == 16) {
576 for (j = 0; j < w; j++) {
577 if (GET_DATA_TWO_BYTES(lined, j) - threshval >= 0)
578 SET_DATA_TWO_BYTES(lined, j, setval);
579 }
580 } else { /* d == 32 */
581 for (j = 0; j < w; j++) {
582 if (*(lined + j) >= threshval)
583 *(lined + j) = setval;
584 }
585 }
586 } else { /* set if below or at threshold */
587 if (d == 8) {
588 for (j = 0; j < w; j++) {
589 if (GET_DATA_BYTE(lined, j) - threshval <= 0)
590 SET_DATA_BYTE(lined, j, setval);
591 }
592 } else if (d == 16) {
593 for (j = 0; j < w; j++) {
594 if (GET_DATA_TWO_BYTES(lined, j) - threshval <= 0)
595 SET_DATA_TWO_BYTES(lined, j, setval);
596 }
597 } else { /* d == 32 */
598 for (j = 0; j < w; j++) {
599 if (*(lined + j) <= threshval)
600 *(lined + j) = setval;
601 }
602 }
603 }
604 }
605
606 return pixd;
607 }
608
609
610 /*-------------------------------------------------------------*
611 * Image accumulator arithmetic operations *
612 *-------------------------------------------------------------*/
613 /*!
614 * \brief pixInitAccumulate()
615 *
616 * \param[in] w, h of accumulate array
617 * \param[in] offset initialize the 32 bpp to have this
618 * value; not more than 0x40000000
619 * \return pixd 32 bpp, or NULL on error
620 *
621 * <pre>
622 * Notes:
623 * (1) %offset must be >= 0.
624 * (2) %offset is used so that we can do arithmetic
625 * with negative number results on l_uint32 data; it
626 * prevents the l_uint32 data from going negative.
627 * (3) Because we use l_int32 intermediate data results,
628 * these should never exceed the max of l_int32 (0x7fffffff).
629 * We do not permit the offset to be above 0x40000000,
630 * which is half way between 0 and the max of l_int32.
631 * (4) The same offset should be used for initialization,
632 * multiplication by a constant, and final extraction!
633 * (5) If you're only adding positive values, %offset can be 0.
634 * </pre>
635 */
636 PIX *
637 pixInitAccumulate(l_int32 w,
638 l_int32 h,
639 l_uint32 offset)
640 {
641 PIX *pixd;
642
643 if ((pixd = pixCreate(w, h, 32)) == NULL)
644 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
645 if (offset > 0x40000000)
646 offset = 0x40000000;
647 pixSetAllArbitrary(pixd, offset);
648 return pixd;
649 }
650
651
652 /*!
653 * \brief pixFinalAccumulate()
654 *
655 * \param[in] pixs 32 bpp
656 * \param[in] offset same as used for initialization
657 * \param[in] depth 8, 16 or 32 bpp, of destination
658 * \return pixd 8, 16 or 32 bpp, or NULL on error
659 *
660 * <pre>
661 * Notes:
662 * (1) %offset must be >= 0 and should not exceed 0x40000000.
663 * (2) %offset is subtracted from the src 32 bpp image
664 * (3) For 8 bpp dest, the result is clipped to [0, 0xff]
665 * (4) For 16 bpp dest, the result is clipped to [0, 0xffff]
666 * </pre>
667 */
668 PIX *
669 pixFinalAccumulate(PIX *pixs,
670 l_uint32 offset,
671 l_int32 depth)
672 {
673 l_int32 i, j, w, h, wpls, wpld, val;
674 l_uint32 *datas, *datad, *lines, *lined;
675 PIX *pixd;
676
677 if (!pixs)
678 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
679 if (pixGetDepth(pixs) != 32)
680 return (PIX *)ERROR_PTR("pixs not 32 bpp", __func__, NULL);
681 if (depth != 8 && depth != 16 && depth != 32)
682 return (PIX *)ERROR_PTR("dest depth not 8, 16, 32 bpp", __func__, NULL);
683 if (offset > 0x40000000)
684 offset = 0x40000000;
685
686 pixGetDimensions(pixs, &w, &h, NULL);
687 if ((pixd = pixCreate(w, h, depth)) == NULL)
688 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
689 pixCopyResolution(pixd, pixs); /* but how did pixs get it initially? */
690 datas = pixGetData(pixs);
691 datad = pixGetData(pixd);
692 wpls = pixGetWpl(pixs);
693 wpld = pixGetWpl(pixd);
694 if (depth == 8) {
695 for (i = 0; i < h; i++) {
696 lines = datas + i * wpls;
697 lined = datad + i * wpld;
698 for (j = 0; j < w; j++) {
699 val = lines[j] - offset;
700 val = L_MAX(0, val);
701 val = L_MIN(255, val);
702 SET_DATA_BYTE(lined, j, (l_uint8)val);
703 }
704 }
705 } else if (depth == 16) {
706 for (i = 0; i < h; i++) {
707 lines = datas + i * wpls;
708 lined = datad + i * wpld;
709 for (j = 0; j < w; j++) {
710 val = lines[j] - offset;
711 val = L_MAX(0, val);
712 val = L_MIN(0xffff, val);
713 SET_DATA_TWO_BYTES(lined, j, (l_uint16)val);
714 }
715 }
716 } else { /* depth == 32 */
717 for (i = 0; i < h; i++) {
718 lines = datas + i * wpls;
719 lined = datad + i * wpld;
720 for (j = 0; j < w; j++)
721 lined[j] = lines[j] - offset;
722 }
723 }
724
725 return pixd;
726 }
727
728
729 /*!
730 * \brief pixFinalAccumulateThreshold()
731 *
732 * \param[in] pixs 32 bpp
733 * \param[in] offset same as used for initialization
734 * \param[in] threshold values less than this are set in the destination
735 * \return pixd 1 bpp, or NULL on error
736 *
737 * <pre>
738 * Notes:
739 * (1) %offset must be >= 0 and should not exceed 0x40000000.
740 * (2) %offset is subtracted from the src 32 bpp image
741 * </pre>
742 */
743 PIX *
744 pixFinalAccumulateThreshold(PIX *pixs,
745 l_uint32 offset,
746 l_uint32 threshold)
747 {
748 l_int32 i, j, w, h, wpls, wpld, val;
749 l_uint32 *datas, *datad, *lines, *lined;
750 PIX *pixd;
751
752 if (!pixs)
753 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
754 if (pixGetDepth(pixs) != 32)
755 return (PIX *)ERROR_PTR("pixs not 32 bpp", __func__, NULL);
756 if (offset > 0x40000000)
757 offset = 0x40000000;
758
759 pixGetDimensions(pixs, &w, &h, NULL);
760 if ((pixd = pixCreate(w, h, 1)) == NULL)
761 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
762 pixCopyResolution(pixd, pixs); /* but how did pixs get it initially? */
763 datas = pixGetData(pixs);
764 datad = pixGetData(pixd);
765 wpls = pixGetWpl(pixs);
766 wpld = pixGetWpl(pixd);
767 for (i = 0; i < h; i++) {
768 lines = datas + i * wpls;
769 lined = datad + i * wpld;
770 for (j = 0; j < w; j++) {
771 val = lines[j] - offset;
772 if (val >= threshold) {
773 SET_DATA_BIT(lined, j);
774 }
775 }
776 }
777
778 return pixd;
779 }
780
781
782 /*!
783 * \brief pixAccumulate()
784 *
785 * \param[in] pixd 32 bpp
786 * \param[in] pixs 1, 8, 16 or 32 bpp
787 * \param[in] op L_ARITH_ADD or L_ARITH_SUBTRACT
788 * \return 0 if OK; 1 on error
789 *
790 * <pre>
791 * Notes:
792 * (1) This adds or subtracts each pixs value from pixd.
793 * (2) This clips to the minimum of pixs and pixd, so they
794 * do not need to be the same size.
795 * (3) The alignment is to the origin [UL corner] of pixs & pixd.
796 * </pre>
797 */
798 l_ok
799 pixAccumulate(PIX *pixd,
800 PIX *pixs,
801 l_int32 op)
802 {
803 l_int32 i, j, w, h, d, wd, hd, wpls, wpld;
804 l_uint32 *datas, *datad, *lines, *lined;
805
806
807 if (!pixd || (pixGetDepth(pixd) != 32))
808 return ERROR_INT("pixd not defined or not 32 bpp", __func__, 1);
809 if (!pixs)
810 return ERROR_INT("pixs not defined", __func__, 1);
811 d = pixGetDepth(pixs);
812 if (d != 1 && d != 8 && d != 16 && d != 32)
813 return ERROR_INT("pixs not 1, 8, 16 or 32 bpp", __func__, 1);
814 if (op != L_ARITH_ADD && op != L_ARITH_SUBTRACT)
815 return ERROR_INT("op must be in {L_ARITH_ADD, L_ARITH_SUBTRACT}",
816 __func__, 1);
817
818 datas = pixGetData(pixs);
819 datad = pixGetData(pixd);
820 wpls = pixGetWpl(pixs);
821 wpld = pixGetWpl(pixd);
822 pixGetDimensions(pixs, &w, &h, NULL);
823 pixGetDimensions(pixd, &wd, &hd, NULL);
824 w = L_MIN(w, wd);
825 h = L_MIN(h, hd);
826 if (d == 1) {
827 for (i = 0; i < h; i++) {
828 lines = datas + i * wpls;
829 lined = datad + i * wpld;
830 if (op == L_ARITH_ADD) {
831 for (j = 0; j < w; j++)
832 lined[j] += GET_DATA_BIT(lines, j);
833 } else { /* op == L_ARITH_SUBTRACT */
834 for (j = 0; j < w; j++)
835 lined[j] -= GET_DATA_BIT(lines, j);
836 }
837 }
838 } else if (d == 8) {
839 for (i = 0; i < h; i++) {
840 lines = datas + i * wpls;
841 lined = datad + i * wpld;
842 if (op == L_ARITH_ADD) {
843 for (j = 0; j < w; j++)
844 lined[j] += GET_DATA_BYTE(lines, j);
845 } else { /* op == L_ARITH_SUBTRACT */
846 for (j = 0; j < w; j++)
847 lined[j] -= GET_DATA_BYTE(lines, j);
848 }
849 }
850 } else if (d == 16) {
851 for (i = 0; i < h; i++) {
852 lines = datas + i * wpls;
853 lined = datad + i * wpld;
854 if (op == L_ARITH_ADD) {
855 for (j = 0; j < w; j++)
856 lined[j] += GET_DATA_TWO_BYTES(lines, j);
857 } else { /* op == L_ARITH_SUBTRACT */
858 for (j = 0; j < w; j++)
859 lined[j] -= GET_DATA_TWO_BYTES(lines, j);
860 }
861 }
862 } else { /* d == 32 */
863 for (i = 0; i < h; i++) {
864 lines = datas + i * wpls;
865 lined = datad + i * wpld;
866 if (op == L_ARITH_ADD) {
867 for (j = 0; j < w; j++)
868 lined[j] += lines[j];
869 } else { /* op == L_ARITH_SUBTRACT */
870 for (j = 0; j < w; j++)
871 lined[j] -= lines[j];
872 }
873 }
874 }
875
876 return 0;
877 }
878
879
880 /*!
881 * \brief pixMultConstAccumulate()
882 *
883 * \param[in] pixs 32 bpp
884 * \param[in] factor
885 * \param[in] offset same as used for initialization
886 * \return 0 if OK; 1 on error
887 *
888 * <pre>
889 * Notes:
890 * (1) %offset must be >= 0 and should not exceed 0x40000000.
891 * (2) This multiplies each pixel, relative to offset, by %factor.
892 * (3) The result is returned with %offset back in place.
893 * </pre>
894 */
895 l_ok
896 pixMultConstAccumulate(PIX *pixs,
897 l_float32 factor,
898 l_uint32 offset)
899 {
900 l_int32 i, j, w, h, wpl, val;
901 l_uint32 *data, *line;
902
903 if (!pixs)
904 return ERROR_INT("pixs not defined", __func__, 1);
905 if (pixGetDepth(pixs) != 32)
906 return ERROR_INT("pixs not 32 bpp", __func__, 1);
907 if (offset > 0x40000000)
908 offset = 0x40000000;
909
910 pixGetDimensions(pixs, &w, &h, NULL);
911 data = pixGetData(pixs);
912 wpl = pixGetWpl(pixs);
913 for (i = 0; i < h; i++) {
914 line = data + i * wpl;
915 for (j = 0; j < w; j++) {
916 val = line[j] - offset;
917 val = (l_int32)(val * factor);
918 val += offset;
919 line[j] = (l_uint32)val;
920 }
921 }
922
923 return 0;
924 }
925
926
927 /*-----------------------------------------------------------------------*
928 * Absolute value of difference *
929 *-----------------------------------------------------------------------*/
930 /*!
931 * \brief pixAbsDifference()
932 *
933 * \param[in] pixs1, pixs2 both either 8 or 16 bpp gray, or 32 bpp RGB
934 * \return pixd, or NULL on error
935 *
936 * <pre>
937 * Notes:
938 * (1) The depth of pixs1 and pixs2 must be equal.
939 * (2) Clips computation to the min size, aligning the UL corners
940 * (3) For 8 and 16 bpp, assumes one gray component.
941 * (4) For 32 bpp, assumes 3 color components, and ignores the
942 * LSB of each word (the alpha channel)
943 * (5) Computes the absolute value of the difference between
944 * each component value.
945 * </pre>
946 */
947 PIX *
948 pixAbsDifference(PIX *pixs1,
949 PIX *pixs2)
950 {
951 l_int32 i, j, w, h, w2, h2, d, wpls1, wpls2, wpld, val1, val2, diff;
952 l_int32 rval1, gval1, bval1, rval2, gval2, bval2, rdiff, gdiff, bdiff;
953 l_uint32 *datas1, *datas2, *datad, *lines1, *lines2, *lined;
954 PIX *pixd;
955
956 if (!pixs1)
957 return (PIX *)ERROR_PTR("pixs1 not defined", __func__, NULL);
958 if (!pixs2)
959 return (PIX *)ERROR_PTR("pixs2 not defined", __func__, NULL);
960 d = pixGetDepth(pixs1);
961 if (d != pixGetDepth(pixs2))
962 return (PIX *)ERROR_PTR("src1 and src2 depths unequal", __func__, NULL);
963 if (d != 8 && d != 16 && d != 32)
964 return (PIX *)ERROR_PTR("depths not in {8, 16, 32}", __func__, NULL);
965
966 pixGetDimensions(pixs1, &w, &h, NULL);
967 pixGetDimensions(pixs2, &w2, &h2, NULL);
968 w = L_MIN(w, w2);
969 h = L_MIN(h, h2);
970 if ((pixd = pixCreate(w, h, d)) == NULL)
971 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
972 pixCopyResolution(pixd, pixs1);
973 datas1 = pixGetData(pixs1);
974 datas2 = pixGetData(pixs2);
975 datad = pixGetData(pixd);
976 wpls1 = pixGetWpl(pixs1);
977 wpls2 = pixGetWpl(pixs2);
978 wpld = pixGetWpl(pixd);
979 if (d == 8) {
980 for (i = 0; i < h; i++) {
981 lines1 = datas1 + i * wpls1;
982 lines2 = datas2 + i * wpls2;
983 lined = datad + i * wpld;
984 for (j = 0; j < w; j++) {
985 val1 = GET_DATA_BYTE(lines1, j);
986 val2 = GET_DATA_BYTE(lines2, j);
987 diff = L_ABS(val1 - val2);
988 SET_DATA_BYTE(lined, j, diff);
989 }
990 }
991 } else if (d == 16) {
992 for (i = 0; i < h; i++) {
993 lines1 = datas1 + i * wpls1;
994 lines2 = datas2 + i * wpls2;
995 lined = datad + i * wpld;
996 for (j = 0; j < w; j++) {
997 val1 = GET_DATA_TWO_BYTES(lines1, j);
998 val2 = GET_DATA_TWO_BYTES(lines2, j);
999 diff = L_ABS(val1 - val2);
1000 SET_DATA_TWO_BYTES(lined, j, diff);
1001 }
1002 }
1003 } else { /* d == 32 */
1004 for (i = 0; i < h; i++) {
1005 lines1 = datas1 + i * wpls1;
1006 lines2 = datas2 + i * wpls2;
1007 lined = datad + i * wpld;
1008 for (j = 0; j < w; j++) {
1009 extractRGBValues(lines1[j], &rval1, &gval1, &bval1);
1010 extractRGBValues(lines2[j], &rval2, &gval2, &bval2);
1011 rdiff = L_ABS(rval1 - rval2);
1012 gdiff = L_ABS(gval1 - gval2);
1013 bdiff = L_ABS(bval1 - bval2);
1014 composeRGBPixel(rdiff, gdiff, bdiff, lined + j);
1015 }
1016 }
1017 }
1018
1019 return pixd;
1020 }
1021
1022
1023 /*-----------------------------------------------------------------------*
1024 * Sum of color images *
1025 *-----------------------------------------------------------------------*/
1026 /*!
1027 * \brief pixAddRGB()
1028 *
1029 * \param[in] pixs1, pixs2 32 bpp RGB, or colormapped
1030 * \return pixd, or NULL on error
1031 *
1032 * <pre>
1033 * Notes:
1034 * (1) Clips computation to the minimum size, aligning the UL corners.
1035 * (2) Removes any colormap to RGB, and ignores the LSB of each
1036 * pixel word (the alpha channel).
1037 * (3) Adds each component value, pixelwise, clipping to 255.
1038 * (4) This is useful to combine two images where most of the
1039 * pixels are essentially black, such as in pixPerceptualDiff().
1040 * </pre>
1041 */
1042 PIX *
1043 pixAddRGB(PIX *pixs1,
1044 PIX *pixs2)
1045 {
1046 l_int32 i, j, w, h, d, w2, h2, d2, wplc1, wplc2, wpld;
1047 l_int32 rval1, gval1, bval1, rval2, gval2, bval2, rval, gval, bval;
1048 l_uint32 *datac1, *datac2, *datad, *linec1, *linec2, *lined;
1049 PIX *pixc1, *pixc2, *pixd;
1050
1051 if (!pixs1)
1052 return (PIX *)ERROR_PTR("pixs1 not defined", __func__, NULL);
1053 if (!pixs2)
1054 return (PIX *)ERROR_PTR("pixs2 not defined", __func__, NULL);
1055 pixGetDimensions(pixs1, &w, &h, &d);
1056 pixGetDimensions(pixs2, &w2, &h2, &d2);
1057 if (!pixGetColormap(pixs1) && d != 32)
1058 return (PIX *)ERROR_PTR("pixs1 not cmapped or rgb", __func__, NULL);
1059 if (!pixGetColormap(pixs2) && d2 != 32)
1060 return (PIX *)ERROR_PTR("pixs2 not cmapped or rgb", __func__, NULL);
1061 if (pixGetColormap(pixs1))
1062 pixc1 = pixRemoveColormap(pixs1, REMOVE_CMAP_TO_FULL_COLOR);
1063 else
1064 pixc1 = pixClone(pixs1);
1065 if (pixGetColormap(pixs2))
1066 pixc2 = pixRemoveColormap(pixs2, REMOVE_CMAP_TO_FULL_COLOR);
1067 else
1068 pixc2 = pixClone(pixs2);
1069
1070 w = L_MIN(w, w2);
1071 h = L_MIN(h, h2);
1072 pixd = pixCreate(w, h, 32);
1073 pixCopyResolution(pixd, pixs1);
1074 datac1 = pixGetData(pixc1);
1075 datac2 = pixGetData(pixc2);
1076 datad = pixGetData(pixd);
1077 wplc1 = pixGetWpl(pixc1);
1078 wplc2 = pixGetWpl(pixc2);
1079 wpld = pixGetWpl(pixd);
1080 for (i = 0; i < h; i++) {
1081 linec1 = datac1 + i * wplc1;
1082 linec2 = datac2 + i * wplc2;
1083 lined = datad + i * wpld;
1084 for (j = 0; j < w; j++) {
1085 extractRGBValues(linec1[j], &rval1, &gval1, &bval1);
1086 extractRGBValues(linec2[j], &rval2, &gval2, &bval2);
1087 rval = L_MIN(255, rval1 + rval2);
1088 gval = L_MIN(255, gval1 + gval2);
1089 bval = L_MIN(255, bval1 + bval2);
1090 composeRGBPixel(rval, gval, bval, lined + j);
1091 }
1092 }
1093
1094 pixDestroy(&pixc1);
1095 pixDestroy(&pixc2);
1096 return pixd;
1097 }
1098
1099
1100 /*-----------------------------------------------------------------------*
1101 * Two-image min and max operations (8 and 16 bpp) *
1102 *-----------------------------------------------------------------------*/
1103 /*!
1104 * \brief pixMinOrMax()
1105 *
1106 * \param[in] pixd [optional] destination: this can be null,
1107 * equal to pixs1, or different from pixs1
1108 * \param[in] pixs1 can be equal to pixd
1109 * \param[in] pixs2
1110 * \param[in] type L_CHOOSE_MIN, L_CHOOSE_MAX
1111 * \return pixd always
1112 *
1113 * <pre>
1114 * Notes:
1115 * (1) This gives the min or max of two images, component-wise.
1116 * (2) The depth can be 8 or 16 bpp for 1 component, and 32 bpp
1117 * for a 3 component image. For 32 bpp, ignore the LSB
1118 * of each word (the alpha channel)
1119 * (3) There are 3 cases:
1120 * ~ if pixd == null, MinOrMax(src1, src2) --> new pixd
1121 * ~ if pixd == pixs1, MinOrMax(src1, src2) --> src1 (in-place)
1122 * ~ if pixd != pixs1, MinOrMax(src1, src2) --> input pixd
1123 * </pre>
1124 */
1125 PIX *
1126 pixMinOrMax(PIX *pixd,
1127 PIX *pixs1,
1128 PIX *pixs2,
1129 l_int32 type)
1130 {
1131 l_int32 d, ws, hs, w, h, wpls, wpld, i, j, vals, vald, val;
1132 l_int32 rval1, gval1, bval1, rval2, gval2, bval2, rval, gval, bval;
1133 l_uint32 *datas, *datad, *lines, *lined;
1134
1135 if (!pixs1)
1136 return (PIX *)ERROR_PTR("pixs1 not defined", __func__, pixd);
1137 if (!pixs2)
1138 return (PIX *)ERROR_PTR("pixs2 not defined", __func__, pixd);
1139 if (pixs1 == pixs2)
1140 return (PIX *)ERROR_PTR("pixs1 and pixs2 must differ", __func__, pixd);
1141 if (type != L_CHOOSE_MIN && type != L_CHOOSE_MAX)
1142 return (PIX *)ERROR_PTR("invalid type", __func__, pixd);
1143 d = pixGetDepth(pixs1);
1144 if (pixGetDepth(pixs2) != d)
1145 return (PIX *)ERROR_PTR("depths unequal", __func__, pixd);
1146 if (d != 8 && d != 16 && d != 32)
1147 return (PIX *)ERROR_PTR("depth not 8, 16 or 32 bpp", __func__, pixd);
1148
1149 if (pixs1 != pixd)
1150 pixd = pixCopy(pixd, pixs1);
1151
1152 pixGetDimensions(pixs2, &ws, &hs, NULL);
1153 pixGetDimensions(pixd, &w, &h, NULL);
1154 w = L_MIN(w, ws);
1155 h = L_MIN(h, hs);
1156 datas = pixGetData(pixs2);
1157 datad = pixGetData(pixd);
1158 wpls = pixGetWpl(pixs2);
1159 wpld = pixGetWpl(pixd);
1160 for (i = 0; i < h; i++) {
1161 lines = datas + i * wpls;
1162 lined = datad + i * wpld;
1163 if (d == 8) {
1164 for (j = 0; j < w; j++) {
1165 vals = GET_DATA_BYTE(lines, j);
1166 vald = GET_DATA_BYTE(lined, j);
1167 if (type == L_CHOOSE_MIN)
1168 val = L_MIN(vals, vald);
1169 else /* type == L_CHOOSE_MAX */
1170 val = L_MAX(vals, vald);
1171 SET_DATA_BYTE(lined, j, val);
1172 }
1173 } else if (d == 16) {
1174 for (j = 0; j < w; j++) {
1175 vals = GET_DATA_TWO_BYTES(lines, j);
1176 vald = GET_DATA_TWO_BYTES(lined, j);
1177 if (type == L_CHOOSE_MIN)
1178 val = L_MIN(vals, vald);
1179 else /* type == L_CHOOSE_MAX */
1180 val = L_MAX(vals, vald);
1181 SET_DATA_TWO_BYTES(lined, j, val);
1182 }
1183 } else { /* d == 32 */
1184 for (j = 0; j < w; j++) {
1185 extractRGBValues(lines[j], &rval1, &gval1, &bval1);
1186 extractRGBValues(lined[j], &rval2, &gval2, &bval2);
1187 if (type == L_CHOOSE_MIN) {
1188 rval = L_MIN(rval1, rval2);
1189 gval = L_MIN(gval1, gval2);
1190 bval = L_MIN(bval1, bval2);
1191 } else { /* type == L_CHOOSE_MAX */
1192 rval = L_MAX(rval1, rval2);
1193 gval = L_MAX(gval1, gval2);
1194 bval = L_MAX(bval1, bval2);
1195 }
1196 composeRGBPixel(rval, gval, bval, lined + j);
1197 }
1198 }
1199 }
1200
1201 return pixd;
1202 }
1203
1204
1205 /*-----------------------------------------------------------------------*
1206 * Scale for maximum dynamic range *
1207 *-----------------------------------------------------------------------*/
1208 /*!
1209 * \brief pixMaxDynamicRange()
1210 *
1211 * \param[in] pixs 4, 8, 16 or 32 bpp source
1212 * \param[in] type L_LINEAR_SCALE or L_LOG_SCALE
1213 * \return pixd 8 bpp, or NULL on error
1214 *
1215 * <pre>
1216 * Notes:
1217 * (1) Scales pixel values to fit maximally within the dest 8 bpp pixd
1218 * (2) Assumes the source 'pixels' are a 1-component scalar. For
1219 * a 32 bpp source, each pixel is treated as a single number --
1220 * not as a 3-component rgb pixel value.
1221 * (3) Uses a LUT for log scaling.
1222 * </pre>
1223 */
1224 PIX *
1225 pixMaxDynamicRange(PIX *pixs,
1226 l_int32 type)
1227 {
1228 l_uint8 dval;
1229 l_int32 i, j, w, h, d, wpls, wpld, max;
1230 l_uint32 *datas, *datad;
1231 l_uint32 word, sval;
1232 l_uint32 *lines, *lined;
1233 l_float32 factor;
1234 l_float32 *tab;
1235 PIX *pixd;
1236
1237 if (!pixs)
1238 return (PIX *)ERROR_PTR("pixs not defined", __func__, NULL);
1239 pixGetDimensions(pixs, &w, &h, &d);
1240 if (d != 4 && d != 8 && d != 16 && d != 32)
1241 return (PIX *)ERROR_PTR("pixs not in {4,8,16,32} bpp", __func__, NULL);
1242 if (type != L_LINEAR_SCALE && type != L_LOG_SCALE)
1243 return (PIX *)ERROR_PTR("invalid type", __func__, NULL);
1244
1245 if ((pixd = pixCreate(w, h, 8)) == NULL)
1246 return (PIX *)ERROR_PTR("pixd not made", __func__, NULL);
1247 pixCopyResolution(pixd, pixs);
1248 datas = pixGetData(pixs);
1249 datad = pixGetData(pixd);
1250 wpls = pixGetWpl(pixs);
1251 wpld = pixGetWpl(pixd);
1252
1253 /* Get max */
1254 max = 0;
1255 for (i = 0; i < h; i++) {
1256 lines = datas + i * wpls;
1257 for (j = 0; j < wpls; j++) {
1258 word = *(lines + j);
1259 if (d == 4) {
1260 max = L_MAX(max, word >> 28);
1261 max = L_MAX(max, (word >> 24) & 0xf);
1262 max = L_MAX(max, (word >> 20) & 0xf);
1263 max = L_MAX(max, (word >> 16) & 0xf);
1264 max = L_MAX(max, (word >> 12) & 0xf);
1265 max = L_MAX(max, (word >> 8) & 0xf);
1266 max = L_MAX(max, (word >> 4) & 0xf);
1267 max = L_MAX(max, word & 0xf);
1268 } else if (d == 8) {
1269 max = L_MAX(max, word >> 24);
1270 max = L_MAX(max, (word >> 16) & 0xff);
1271 max = L_MAX(max, (word >> 8) & 0xff);
1272 max = L_MAX(max, word & 0xff);
1273 } else if (d == 16) {
1274 max = L_MAX(max, word >> 16);
1275 max = L_MAX(max, word & 0xffff);
1276 } else { /* d == 32 (rgb) */
1277 max = L_MAX(max, word);
1278 }
1279 }
1280 }
1281
1282 /* Map to the full dynamic range */
1283 if (d == 4) {
1284 if (type == L_LINEAR_SCALE) {
1285 factor = 255.f / (l_float32)max;
1286 for (i = 0; i < h; i++) {
1287 lines = datas + i * wpls;
1288 lined = datad + i * wpld;
1289 for (j = 0; j < w; j++) {
1290 sval = GET_DATA_QBIT(lines, j);
1291 dval = (l_uint8)(factor * (l_float32)sval + 0.5);
1292 SET_DATA_QBIT(lined, j, dval);
1293 }
1294 }
1295 } else { /* type == L_LOG_SCALE) */
1296 tab = makeLogBase2Tab();
1297 factor = 255.f / getLogBase2(max, tab);
1298 for (i = 0; i < h; i++) {
1299 lines = datas + i * wpls;
1300 lined = datad + i * wpld;
1301 for (j = 0; j < w; j++) {
1302 sval = GET_DATA_QBIT(lines, j);
1303 dval = (l_uint8)(factor * getLogBase2(sval, tab) + 0.5);
1304 SET_DATA_BYTE(lined, j, dval);
1305 }
1306 }
1307 LEPT_FREE(tab);
1308 }
1309 } else if (d == 8) {
1310 if (type == L_LINEAR_SCALE) {
1311 factor = 255.f / (l_float32)max;
1312 for (i = 0; i < h; i++) {
1313 lines = datas + i * wpls;
1314 lined = datad + i * wpld;
1315 for (j = 0; j < w; j++) {
1316 sval = GET_DATA_BYTE(lines, j);
1317 dval = (l_uint8)(factor * (l_float32)sval + 0.5);
1318 SET_DATA_BYTE(lined, j, dval);
1319 }
1320 }
1321 } else { /* type == L_LOG_SCALE) */
1322 tab = makeLogBase2Tab();
1323 factor = 255.f / getLogBase2(max, tab);
1324 for (i = 0; i < h; i++) {
1325 lines = datas + i * wpls;
1326 lined = datad + i * wpld;
1327 for (j = 0; j < w; j++) {
1328 sval = GET_DATA_BYTE(lines, j);
1329 dval = (l_uint8)(factor * getLogBase2(sval, tab) + 0.5);
1330 SET_DATA_BYTE(lined, j, dval);
1331 }
1332 }
1333 LEPT_FREE(tab);
1334 }
1335 } else if (d == 16) {
1336 if (type == L_LINEAR_SCALE) {
1337 factor = 255.f / (l_float32)max;
1338 for (i = 0; i < h; i++) {
1339 lines = datas + i * wpls;
1340 lined = datad + i * wpld;
1341 for (j = 0; j < w; j++) {
1342 sval = GET_DATA_TWO_BYTES(lines, j);
1343 dval = (l_uint8)(factor * (l_float32)sval + 0.5);
1344 SET_DATA_BYTE(lined, j, dval);
1345 }
1346 }
1347 } else { /* type == L_LOG_SCALE) */
1348 tab = makeLogBase2Tab();
1349 factor = 255.f / getLogBase2(max, tab);
1350 for (i = 0; i < h; i++) {
1351 lines = datas + i * wpls;
1352 lined = datad + i * wpld;
1353 for (j = 0; j < w; j++) {
1354 sval = GET_DATA_TWO_BYTES(lines, j);
1355 dval = (l_uint8)(factor * getLogBase2(sval, tab) + 0.5);
1356 SET_DATA_BYTE(lined, j, dval);
1357 }
1358 }
1359 LEPT_FREE(tab);
1360 }
1361 } else { /* d == 32 */
1362 if (type == L_LINEAR_SCALE) {
1363 factor = 255.f / (l_float32)max;
1364 for (i = 0; i < h; i++) {
1365 lines = datas + i * wpls;
1366 lined = datad + i * wpld;
1367 for (j = 0; j < w; j++) {
1368 sval = lines[j];
1369 dval = (l_uint8)(factor * (l_float32)sval + 0.5);
1370 SET_DATA_BYTE(lined, j, dval);
1371 }
1372 }
1373 } else { /* type == L_LOG_SCALE) */
1374 tab = makeLogBase2Tab();
1375 factor = 255.f / getLogBase2(max, tab);
1376 for (i = 0; i < h; i++) {
1377 lines = datas + i * wpls;
1378 lined = datad + i * wpld;
1379 for (j = 0; j < w; j++) {
1380 sval = lines[j];
1381 dval = (l_uint8)(factor * getLogBase2(sval, tab) + 0.5);
1382 SET_DATA_BYTE(lined, j, dval);
1383 }
1384 }
1385 LEPT_FREE(tab);
1386 }
1387 }
1388
1389 return pixd;
1390 }
1391
1392
1393 /*!
1394 * \brief pixMaxDynamicRangeRGB()
1395 *
1396 * \param[in] pixs 32 bpp rgb source
1397 * \param[in] type L_LINEAR_SCALE or L_LOG_SCALE
1398 * \return pixd 32 bpp, or NULL on error
1399 *
1400 * <pre>
1401 * Notes:
1402 * (1) Scales pixel values to fit maximally within a 32 bpp dest pixd
1403 * (2) All color components are scaled with the same factor, based
1404 * on the maximum r, g or b component in the image. This should
1405 * not be used if the 32-bit value is a single number (e.g., a
1406 * count in a histogram generated by pixMakeHistoHS()).
1407 * (3) Uses a LUT for log scaling.
1408 * </pre>
1409 */
1410 PIX *
1411 pixMaxDynamicRangeRGB(PIX *pixs,
1412 l_int32 type)
1413 {
1414 l_int32 i, j, w, h, wpls, wpld, max;
1415 l_uint32 sval, dval, word;
1416 l_uint32 *datas, *datad;
1417 l_uint32 *lines, *lined;
1418 l_float32 factor;
1419 l_float32 *tab;
1420 PIX *pixd;
1421
1422 if (!pixs || pixGetDepth(pixs) != 32)
1423 return (PIX *)ERROR_PTR("pixs undefined or not 32 bpp", __func__, NULL);
1424 if (type != L_LINEAR_SCALE && type != L_LOG_SCALE)
1425 return (PIX *)ERROR_PTR("invalid type", __func__, NULL);
1426
1427 /* Get max */
1428 pixd = pixCreateTemplate(pixs);
1429 datas = pixGetData(pixs);
1430 datad = pixGetData(pixd);
1431 wpls = pixGetWpl(pixs);
1432 wpld = pixGetWpl(pixd);
1433 pixGetDimensions(pixs, &w, &h, NULL);
1434 max = 0;
1435 for (i = 0; i < h; i++) {
1436 lines = datas + i * wpls;
1437 for (j = 0; j < wpls; j++) {
1438 word = lines[j];
1439 max = L_MAX(max, word >> 24);
1440 max = L_MAX(max, (word >> 16) & 0xff);
1441 max = L_MAX(max, (word >> 8) & 0xff);
1442 }
1443 }
1444 if (max == 0) {
1445 L_WARNING("max = 0; setting to 1\n", __func__);
1446 max = 1;
1447 }
1448
1449 /* Map to the full dynamic range */
1450 if (type == L_LINEAR_SCALE) {
1451 factor = 255.f / (l_float32)max;
1452 for (i = 0; i < h; i++) {
1453 lines = datas + i * wpls;
1454 lined = datad + i * wpld;
1455 for (j = 0; j < w; j++) {
1456 sval = lines[j];
1457 dval = linearScaleRGBVal(sval, factor);
1458 lined[j] = dval;
1459 }
1460 }
1461 } else { /* type == L_LOG_SCALE) */
1462 tab = makeLogBase2Tab();
1463 factor = 255.f / getLogBase2(max, tab);
1464 for (i = 0; i < h; i++) {
1465 lines = datas + i * wpls;
1466 lined = datad + i * wpld;
1467 for (j = 0; j < w; j++) {
1468 sval = lines[j];
1469 dval = logScaleRGBVal(sval, tab, factor);
1470 lined[j] = dval;
1471 }
1472 }
1473 LEPT_FREE(tab);
1474 }
1475
1476 return pixd;
1477 }
1478
1479
1480 /*-----------------------------------------------------------------------*
1481 * RGB pixel value scaling *
1482 *-----------------------------------------------------------------------*/
1483 /*!
1484 * \brief linearScaleRGBVal()
1485 *
1486 * \param[in] sval 32-bit rgb pixel value
1487 * \param[in] factor multiplication factor on each component
1488 * \return dval linearly scaled version of %sval
1489 *
1490 * <pre>
1491 * Notes:
1492 * (1) %factor must be chosen to be not greater than (255 / maxcomp),
1493 * where maxcomp is the maximum value of the pixel components.
1494 * Otherwise, the product will overflow a uint8. In use, factor
1495 * is the same for all pixels in a pix.
1496 * (2) No scaling is performed on the transparency ("A") component.
1497 * </pre>
1498 */
1499 l_uint32
1500 linearScaleRGBVal(l_uint32 sval,
1501 l_float32 factor)
1502 {
1503 l_uint32 dval;
1504
1505 dval = ((l_uint8)(factor * (sval >> 24) + 0.5f) << 24) |
1506 ((l_uint8)(factor * ((sval >> 16) & 0xff) + 0.5f) << 16) |
1507 ((l_uint8)(factor * ((sval >> 8) & 0xff) + 0.5f) << 8) |
1508 (sval & 0xff);
1509 return dval;
1510 }
1511
1512
1513 /*!
1514 * \brief logScaleRGBVal()
1515 *
1516 * \param[in] sval 32-bit rgb pixel value
1517 * \param[in] tab 256 entry log-base-2 table
1518 * \param[in] factor multiplication factor on each component
1519 * \return dval log scaled version of %sval
1520 *
1521 * <pre>
1522 * Notes:
1523 * (1) %tab is made with makeLogBase2Tab().
1524 * (2) %factor must be chosen to be not greater than
1525 * 255.0 / log[base2](maxcomp), where maxcomp is the maximum
1526 * value of the pixel components. Otherwise, the product
1527 * will overflow a uint8. In use, factor is the same for
1528 * all pixels in a pix.
1529 * (3) No scaling is performed on the transparency ("A") component.
1530 * </pre>
1531 */
1532 l_uint32
1533 logScaleRGBVal(l_uint32 sval,
1534 l_float32 *tab,
1535 l_float32 factor)
1536 {
1537 l_uint32 dval;
1538
1539 dval = ((l_uint8)(factor * getLogBase2(sval >> 24, tab) + 0.5f) << 24) |
1540 ((l_uint8)(factor * getLogBase2(((sval >> 16) & 0xff), tab) + 0.5f)
1541 << 16) |
1542 ((l_uint8)(factor * getLogBase2(((sval >> 8) & 0xff), tab) + 0.5f)
1543 << 8) |
1544 (sval & 0xff);
1545 return dval;
1546 }
1547
1548
1549 /*-----------------------------------------------------------------------*
1550 * Log base2 lookup *
1551 *-----------------------------------------------------------------------*/
1552 /*
1553 * \brief makeLogBase2Tab()
1554 *
1555 * \return tab table giving the log[base2] of values from 1 to 255
1556 */
1557 l_float32 *
1558 makeLogBase2Tab(void)
1559 {
1560 l_int32 i;
1561 l_float32 log2;
1562 l_float32 *tab;
1563
1564 if ((tab = (l_float32 *)LEPT_CALLOC(256, sizeof(l_float32))) == NULL)
1565 return (l_float32 *)ERROR_PTR("tab not made", __func__, NULL);
1566
1567 log2 = (l_float32)log((l_float32)2);
1568 for (i = 0; i < 256; i++)
1569 tab[i] = (l_float32)log((l_float32)i) / log2;
1570
1571 return tab;
1572 }
1573
1574
1575 /*
1576 * \brief getLogBase2()
1577 *
1578 * \param[in] val in range [0 ... 255]
1579 * \param[in] logtab 256-entry table of logs
1580 * \return logval log[base2] of %val, or 0 on error
1581 */
1582 l_float32
1583 getLogBase2(l_int32 val,
1584 l_float32 *logtab)
1585 {
1586 if (!logtab)
1587 return ERROR_INT("logtab not defined", __func__, 0);
1588
1589 if (val < 0x100)
1590 return logtab[val];
1591 else if (val < 0x10000)
1592 return 8.0f + logtab[val >> 8];
1593 else if (val < 0x1000000)
1594 return 16.0f + logtab[val >> 16];
1595 else
1596 return 24.0f + logtab[val >> 24];
1597 }