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