1 /*
2 * This file is part of the XForms library package.
3 *
4 * XForms is free software; you can redistribute it and/or modify it
5 * under the terms of the GNU Lesser General Public License as
6 * published by the Free Software Foundation; either version 2.1, or
7 * (at your option) any later version.
8 *
9 * XForms is distributed in the hope that it will be useful, but
10 * WITHOUT ANY WARRANTY; without even the implied warranty of
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 * Lesser General Public License for more details.
13 *
14 * You should have received a copy of the GNU Lesser General Public License
15 * along with XForms. If not, see <http://www.gnu.org/licenses/>.
16 */
17
18
19 /*
20 * Copyright (C) 1991-1996, Thomas G. Lane.
21 * This file is part of the XForms library package.
22 *
23 * The 2-pass quantizer from the JPEG distribution by the
24 * Independent JPEG group. Except for minor interface changes, the code
25 * here is almost verbatim copy of the IJG's code, which
26 * has the following copyright:
27 */
28
29 #ifdef HAVE_CONFIG_H
30 #include "config.h"
31 #endif
32
33 #include "include/forms.h"
34 #include "flimage.h"
35 #include "flimage_int.h"
36
37
38 /***************************************
39 ***************************************/
40
41 void
fl_select_mediancut_quantizer(void)42 fl_select_mediancut_quantizer( void )
43 {
44 flimage_quantize_rgb = j2pass_quantize_rgb;
45 flimage_quantize_packed = j2pass_quantize_packed;
46 }
47
48
49 #define R_SCALE 2 /* scale R distances by this much */
50 #define G_SCALE 3 /* scale G distances by this much */
51 #define B_SCALE 1 /* and B by this much */
52
53 #define C0_SCALE R_SCALE
54 #define C1_SCALE G_SCALE
55 #define C2_SCALE B_SCALE
56
57 #define BITS_IN_JSAMPLE FL_PCBITS
58 #define MAXJSAMPLE ( ( 1<<BITS_IN_JSAMPLE )-1 )
59
60 typedef unsigned char JSAMPLE;
61
62 #define GETJSAMPLE( a ) ( ( a ) &0xff )
63
64 #define MAXNUMCOLORS ( MAXJSAMPLE + 1 ) /* maximum size of colormap */
65
66 /* These will do the right thing for either R,G,B or B,G,R color order,
67 * but you may not like the results for other color orders.
68 */
69
70 #define HIST_C0_BITS 5 /* bits of precision in R/B histogram */
71 #define HIST_C1_BITS 6 /* bits of precision in G histogram */
72 #define HIST_C2_BITS 5 /* bits of precision in B/R histogram */
73
74 /* Number of elements along histogram axes. */
75
76 #define HIST_C0_ELEMS ( 1 << HIST_C0_BITS )
77 #define HIST_C1_ELEMS ( 1 << HIST_C1_BITS )
78 #define HIST_C2_ELEMS ( 1 << HIST_C2_BITS )
79
80 /* These are the amounts to shift an input value to get a histogram index. */
81
82 #define C0_SHIFT ( BITS_IN_JSAMPLE - HIST_C0_BITS )
83 #define C1_SHIFT ( BITS_IN_JSAMPLE - HIST_C1_BITS )
84 #define C2_SHIFT ( BITS_IN_JSAMPLE - HIST_C2_BITS )
85
86 typedef u_short histcell; /* histogram cell; prefer an unsigned int type */
87 typedef histcell *histptr; /* for pointers to histogram cells */
88 typedef histcell hist1d[HIST_C2_ELEMS]; /* typedefs for the array */
89 typedef hist1d *hist2d; /* type for the 2nd-level pointers */
90 typedef hist2d *hist3d; /* type for top-level pointer */
91
92 #if BITS_IN_JSAMPLE == 8
93 typedef short FSERROR; /* 16 bits should be enough */
94 typedef int LOCFSERROR; /* use 'int' for calculation temps */
95 #else
96 typedef int FSERROR; /* may need more than 16 bits */
97 typedef int LOCFSERROR; /* be sure calculation temps are big enough */
98 #endif
99
100 typedef FSERROR *FSERRPTR; /* pointer to error array (in FAR storage!) */
101
102
103 typedef struct
104 {
105 int c0min,
106 c0max;
107 int c1min,
108 c1max;
109 int c2min,
110 c2max;
111 int volume;
112 long colorcount;
113 } box;
114
115 typedef box *boxptr;
116
117
118 typedef struct
119 {
120 hist3d histogram; /* pointer to the 3D histogram array */
121 FSERRPTR fserrors; /* accumulated-errors array */
122 int * error_limiter; /* table for clamping applied error */
123 int on_odd_row; /* flag to remember which row we're on */
124 int * colormap[ 3 ]; /* selected colormap */
125 int actual_number_of_colors;/* number of selected colors */
126 FL_IMAGE * im; /* for progress monitor only */
127 }
128 SPEC;
129
130
131 static void init_error_limit( SPEC * );
132
133 static void prescan_quantize( SPEC *,
134 unsigned char **,
135 unsigned char **,
136 unsigned char **,
137 int,
138 int) ;
139
140 static void select_colors( SPEC *,
141 int );
142
143 static void pass2_fs_dither( SPEC *,
144 unsigned char **,
145 unsigned char **,
146 unsigned char **,
147 unsigned short **,
148 int,
149 int );
150
151
152 /***************************************
153 ***************************************/
154
155 static void
cleanup_spec(SPEC * sp)156 cleanup_spec( SPEC *sp )
157 {
158 if ( sp->fserrors )
159 fl_free( sp->fserrors );
160
161 if ( sp->error_limiter )
162 fl_free( sp->error_limiter - MAXJSAMPLE );
163
164 sp->error_limiter = NULL;
165 sp->fserrors = NULL;
166
167 if ( sp->histogram )
168 {
169 int i;
170
171 for ( i = 0; i < HIST_C0_ELEMS; i++ )
172 {
173 if ( sp->histogram[ i ] )
174 fl_free( sp->histogram[ i ] );
175 sp->histogram[ i ] = NULL;
176 }
177 }
178
179 fl_free( sp->histogram );
180 sp->histogram = NULL;
181
182 fl_free( sp );
183 }
184
185
186 /***************************************
187 ***************************************/
188
189 static SPEC *
alloc_spec(int w,int h FL_UNUSED_ARG,int * rlut,int * glut,int * blut)190 alloc_spec( int w,
191 int h FL_UNUSED_ARG,
192 int * rlut,
193 int * glut,
194 int * blut )
195 {
196 int fs_size = ( w + 2 ) * ( 3 * sizeof( FSERROR ) ),
197 i;
198 SPEC *sp = fl_calloc( 1, sizeof *sp );
199 int err = ! sp;
200
201 if ( ! err )
202 init_error_limit( sp );
203
204 err = err || ! ( sp->fserrors = fl_calloc( 1, fs_size ) );
205 err = err || ! ( sp->histogram = fl_calloc( 1,
206 HIST_C0_ELEMS * sizeof( hist2d ) ) );
207 for ( i = 0; ! err && i < HIST_C0_ELEMS; i++ )
208 err = ! ( sp->histogram[ i ] = fl_calloc( 1,
209 HIST_C1_ELEMS * HIST_C2_ELEMS * sizeof( histcell ) ) );
210
211 if ( err )
212 {
213 cleanup_spec( sp );
214 sp = NULL;
215 }
216 else
217 {
218 sp->colormap[ 0 ] = rlut;
219 sp->colormap[ 1 ] = glut;
220 sp->colormap[ 2 ] = blut;
221 }
222
223 return sp;
224 }
225
226
227 /***************************************
228 ***************************************/
229
230 int
j2pass_quantize_rgb(unsigned char ** red,unsigned char ** green,unsigned char ** blue,int w,int h,int max_color,unsigned short ** ci,int * actual_color,int * red_lut,int * green_lut,int * blue_lut,FL_IMAGE * im)231 j2pass_quantize_rgb( unsigned char ** red,
232 unsigned char ** green,
233 unsigned char ** blue,
234 int w,
235 int h,
236 int max_color,
237 unsigned short ** ci,
238 int * actual_color,
239 int * red_lut,
240 int * green_lut,
241 int * blue_lut,
242 FL_IMAGE * im )
243 {
244 SPEC *sp = alloc_spec( w, h, red_lut, green_lut, blue_lut );
245 int i;
246
247 if ( ! sp )
248 {
249 *actual_color = 0;
250 if ( im )
251 im->error_message( im, "Failed to allocate working memory" );
252 return -1;
253 }
254
255 if ( *actual_color > 256 )
256 *actual_color = 256;
257
258 sp->im = im;
259
260 /* get histogram */
261
262 prescan_quantize( sp, red, green, blue, w, h );
263
264 select_colors( sp, max_color );
265
266 /* re-init histogram for inverse lookup */
267
268 for ( i = 0; i < HIST_C0_ELEMS; i++ )
269 memset( sp->histogram[ i ], 0,
270 HIST_C1_ELEMS * HIST_C2_ELEMS * sizeof( histcell ) );
271
272 sp->on_odd_row = 0;
273 pass2_fs_dither( sp, red, green, blue, ci, w, h );
274 *actual_color = sp->actual_number_of_colors;
275 cleanup_spec( sp );
276
277 if ( im )
278 {
279 im->completed = im->h;
280 im->visual_cue( im, "Quantization Done" );
281 }
282
283 return 0;
284 }
285
286
287 /***************************************
288 ***************************************/
289
290 int
j2pass_quantize_packed(unsigned int ** packed,int w,int h,int max_color,unsigned short ** ci,int * actual_color,int * red_lut,int * green_lut,int * blue_lut,FL_IMAGE * im)291 j2pass_quantize_packed( unsigned int ** packed,
292 int w,
293 int h,
294 int max_color,
295 unsigned short ** ci,
296 int * actual_color,
297 int * red_lut,
298 int * green_lut,
299 int * blue_lut,
300 FL_IMAGE * im )
301 {
302 SPEC *sp = alloc_spec( w, h, red_lut, green_lut, blue_lut );
303 unsigned char **red = NULL,
304 **green = NULL,
305 **blue = NULL;
306 int i,
307 err;
308
309 if ( ! sp )
310 {
311 if ( im )
312 im->error_message( im, "Quantize: can't allocate memory" );
313 *actual_color = 0;
314 return -1;
315 }
316
317 sp->im = im;
318
319 /* we can process the image one piece a time to avoid the heavy memory
320 usage, but packed is not that common. For now, do it in one chunk */
321
322 err = ! ( red = fl_get_matrix( h, w, sizeof **red ) )
323 || ! ( green = fl_get_matrix( h, w, sizeof **green ) )
324 || ! ( blue = fl_get_matrix( h, w, sizeof **blue ) );
325
326 if ( err )
327 {
328 const char *s = "Quantize: can't allocate memory";
329 if ( im )
330 im->error_message( im, s );
331 else
332 fprintf( stderr, "%s\n", s );
333
334 fl_free_matrix( red );
335 fl_free_matrix( green );
336 fl_free_matrix( blue );
337
338 return -1;
339 }
340
341 for ( i = w * h; --i >= 0; )
342 FL_UNPACK( packed[ 0 ][ i ],
343 red[ 0 ][ i ], green[ 0 ][ i ], blue[ 0 ][ i ]);
344
345
346 /* get histogram */
347
348 prescan_quantize( sp, red, green, blue, w, h );
349
350 select_colors( sp, max_color );
351
352 /* re-init histogram for inverse lookup */
353
354 for ( i = 0; i < HIST_C0_ELEMS; i++ )
355 memset( sp->histogram[ i ], 0,
356 HIST_C1_ELEMS * HIST_C2_ELEMS * sizeof( histcell ) );
357
358 sp->on_odd_row = 0;
359 pass2_fs_dither( sp, red, green, blue, ci, w, h );
360 *actual_color = sp->actual_number_of_colors;
361
362 fl_free_matrix( red );
363 fl_free_matrix( green );
364 fl_free_matrix( blue );
365 cleanup_spec( sp );
366
367 if ( im )
368 {
369 im->completed = im->h;
370 im->visual_cue( im, "Quantization Done" );
371 }
372
373 return 0;
374 }
375
376
377 /* log2(histogram cells in update box) for each axis; this can be adjusted */
378
379 #define BOX_C0_LOG ( HIST_C0_BITS - 3 )
380 #define BOX_C1_LOG ( HIST_C1_BITS - 3 )
381 #define BOX_C2_LOG ( HIST_C2_BITS - 3 )
382
383 #define BOX_C0_ELEMS ( 1 << BOX_C0_LOG ) /* # of hist cells in update box */
384 #define BOX_C1_ELEMS ( 1 << BOX_C1_LOG )
385 #define BOX_C2_ELEMS ( 1 << BOX_C2_LOG )
386
387 #define BOX_C0_SHIFT ( C0_SHIFT + BOX_C0_LOG )
388 #define BOX_C1_SHIFT ( C1_SHIFT + BOX_C1_LOG )
389 #define BOX_C2_SHIFT ( C2_SHIFT + BOX_C2_LOG )
390
391
392 /***************************************
393 * Locate the colormap entries close enough to an update box to be candidates
394 * for the nearest entry to some cell(s) in the update box. The update box
395 * is specified by the center coordinates of its first cell. The number of
396 * candidate colormap entries is returned, and their colormap indexes are
397 * placed in colorlist[].
398 * This routine uses Heckbert's "locally sorted search" criterion to select
399 * the colors that need further consideration.
400 ***************************************/
401
402 static int
find_nearby_colors(SPEC * sp,int minc0,int minc1,int minc2,JSAMPLE colorlist[])403 find_nearby_colors( SPEC * sp,
404 int minc0,
405 int minc1,
406 int minc2,
407 JSAMPLE colorlist[ ] )
408 {
409 int numcolors = sp->actual_number_of_colors;
410 int maxc0,
411 maxc1,
412 maxc2;
413 int centerc0,
414 centerc1,
415 centerc2;
416 int i,
417 x,
418 ncolors;
419 int minmaxdist,
420 min_dist,
421 max_dist,
422 tdist;
423 int mindist[ MAXNUMCOLORS ]; /* min distance to colormap entry i */
424
425 /* Compute true coordinates of update box's upper corner and center. *
426 Actually we compute the coordinates of the center of the upper-corner
427 * histogram cell, which are the upper bounds of the volume we care
428 about. * Note that since ">>" rounds down, the "center" values may be
429 closer to * min than to max; hence comparisons to them must be "<=",
430 not "<". */
431
432 maxc0 = minc0 + ( ( 1 << BOX_C0_SHIFT ) - ( 1 << C0_SHIFT ) );
433 centerc0 = ( minc0 + maxc0 ) >> 1;
434 maxc1 = minc1 + ( ( 1 << BOX_C1_SHIFT ) - ( 1 << C1_SHIFT ) );
435 centerc1 = ( minc1 + maxc1 ) >> 1;
436 maxc2 = minc2 + ( ( 1 << BOX_C2_SHIFT ) - ( 1 << C2_SHIFT ) );
437 centerc2 = ( minc2 + maxc2 ) >> 1;
438
439 /* For each color in colormap, find:
440 1. its minimum squared-distance to any point in the update box
441 (zero if color is within update box)
442 2. its maximum squared-distance to any point in the update box.
443 Both of these can be found by considering only the corners of
444 the box. We save the minimum distance for each color in mindist[];
445 only the smallest maximum distance is of interest. */
446
447 minmaxdist = 0x7FFFFFFFL;
448
449 for ( i = 0; i < numcolors; i++ )
450 {
451 /* We compute the squared-c0-distance term, then add in the other
452 two. */
453
454 x = sp->colormap[ 0 ][ i ];
455 if ( x < minc0 )
456 {
457 tdist = ( x - minc0 ) * C0_SCALE;
458 min_dist = tdist * tdist;
459 tdist = ( x - maxc0 ) * C0_SCALE;
460 max_dist = tdist * tdist;
461 }
462 else if ( x > maxc0 )
463 {
464 tdist = ( x - maxc0 ) * C0_SCALE;
465 min_dist = tdist * tdist;
466 tdist = ( x - minc0 ) * C0_SCALE;
467 max_dist = tdist * tdist;
468 }
469 else
470 {
471 /* within cell range so no contribution to min_dist */
472
473 min_dist = 0;
474 if ( x <= centerc0 )
475 {
476 tdist = ( x - maxc0 ) * C0_SCALE;
477 max_dist = tdist * tdist;
478 }
479 else
480 {
481 tdist = ( x - minc0 ) * C0_SCALE;
482 max_dist = tdist * tdist;
483 }
484 }
485
486 x = sp->colormap[ 1 ][ i ];
487 if ( x < minc1 )
488 {
489 tdist = ( x - minc1 ) * C1_SCALE;
490 min_dist += tdist * tdist;
491 tdist = ( x - maxc1 ) * C1_SCALE;
492 max_dist += tdist * tdist;
493 }
494 else if ( x > maxc1 )
495 {
496 tdist = ( x - maxc1 ) * C1_SCALE;
497 min_dist += tdist * tdist;
498 tdist = (x - minc1 ) * C1_SCALE;
499 max_dist += tdist * tdist;
500 }
501 else
502 {
503 /* within cell range so no contribution to min_dist */
504
505 if ( x <= centerc1 )
506 {
507 tdist = ( x - maxc1 ) * C1_SCALE;
508 max_dist += tdist * tdist;
509 }
510 else
511 {
512 tdist = ( x - minc1 ) * C1_SCALE;
513 max_dist += tdist * tdist;
514 }
515 }
516
517 x = sp->colormap[ 2 ][ i ];
518 if ( x < minc2 )
519 {
520 tdist = ( x - minc2 ) * C2_SCALE;
521 min_dist += tdist * tdist;
522 tdist = ( x - maxc2 ) * C2_SCALE;
523 max_dist += tdist * tdist;
524 }
525 else if ( x > maxc2 )
526 {
527 tdist = ( x - maxc2 ) * C2_SCALE;
528 min_dist += tdist * tdist;
529 tdist = ( x - minc2 ) * C2_SCALE;
530 max_dist += tdist * tdist;
531 }
532 else
533 {
534 /* within cell range so no contribution to min_dist */
535
536 if ( x <= centerc2 )
537 {
538 tdist = ( x - maxc2 ) * C2_SCALE;
539 max_dist += tdist * tdist;
540 }
541 else
542 {
543 tdist = ( x - minc2 ) * C2_SCALE;
544 max_dist += tdist * tdist;
545 }
546 }
547
548 mindist[ i ] = min_dist; /* save away the results */
549 if ( max_dist < minmaxdist )
550 minmaxdist = max_dist;
551 }
552
553 /* Now we know that no cell in the update box is more than minmaxdist *
554 away from some colormap entry. Therefore, only colors that are *
555 within minmaxdist of some part of the box need be considered. */
556
557 ncolors = 0;
558
559 for ( i = 0; i < numcolors; i++ )
560 if ( mindist[ i ] <= minmaxdist )
561 colorlist[ ncolors++ ] = ( JSAMPLE ) i;
562
563 return ncolors;
564 }
565
566
567 /***************************************
568 * Find the closest colormap entry for each cell in the update box,
569 * given the list of candidate colors prepared by find_nearby_colors.
570 * Return the indexes of the closest entries in the bestcolor[] array.
571 * This routine uses Thomas' incremental distance calculation method to
572 * find the distance from a colormap entry to successive cells in the box.
573 ***************************************/
574
575 static void
find_best_colors(SPEC * sp,int minc0,int minc1,int minc2,int numcolors,JSAMPLE colorlist[],JSAMPLE bestcolor[])576 find_best_colors( SPEC * sp,
577 int minc0,
578 int minc1,
579 int minc2,
580 int numcolors,
581 JSAMPLE colorlist[ ],
582 JSAMPLE bestcolor[ ] )
583 {
584 int ic0,
585 ic1,
586 ic2;
587 int i,
588 icolor;
589 int *bptr; /* pointer into bestdist[] array */
590 JSAMPLE *cptr; /* pointer into bestcolor[] array */
591 int dist0,
592 dist1; /* initial distance values */
593 int dist2; /* current distance in inner loop */
594 int xx0,
595 xx1; /* distance increments */
596 int xx2;
597 int inc0,
598 inc1,
599 inc2; /* initial values for increments */
600 /* This array holds the distance to the nearest-so-far color for each
601 cell */
602 int bestdist[ BOX_C0_ELEMS * BOX_C1_ELEMS * BOX_C2_ELEMS ];
603
604 /* Initialize best-distance for each cell of the update box */
605
606 bptr = bestdist;
607 for ( i = BOX_C0_ELEMS * BOX_C1_ELEMS * BOX_C2_ELEMS - 1; i >= 0; i-- )
608 *bptr++ = 0x7FFFFFFFL;
609
610 /* For each color selected by find_nearby_colors, * compute its distance
611 to the center of each cell in the box. * If that's less than
612 best-so-far, update best distance and color number. */
613
614 /* Nominal steps between cell centers ("x" in Thomas article) */
615
616 #define STEP_C0 ( ( 1 << C0_SHIFT ) * C0_SCALE )
617 #define STEP_C1 ( ( 1 << C1_SHIFT ) * C1_SCALE )
618 #define STEP_C2 ( ( 1 << C2_SHIFT ) * C2_SCALE )
619
620 for ( i = 0; i < numcolors; i++ )
621 {
622 icolor = GETJSAMPLE( colorlist[ i ] );
623
624 /* Compute (square of) distance from minc0/c1/c2 to this color */
625
626 inc0 = ( minc0 - GETJSAMPLE( sp->colormap[ 0 ][ icolor] ) ) * C0_SCALE;
627 dist0 = inc0 * inc0;
628 inc1 = ( minc1 - GETJSAMPLE( sp->colormap[ 1 ][ icolor ] ) ) * C1_SCALE;
629 dist0 += inc1 * inc1;
630 inc2 = ( minc2 - GETJSAMPLE( sp->colormap[ 2 ][ icolor ] ) ) * C2_SCALE;
631 dist0 += inc2 * inc2;
632
633 /* Form the initial difference increments */
634
635 inc0 = inc0 * 2 * STEP_C0 + STEP_C0 * STEP_C0;
636 inc1 = inc1 * 2 * STEP_C1 + STEP_C1 * STEP_C1;
637 inc2 = inc2 * 2 * STEP_C2 + STEP_C2 * STEP_C2;
638
639 /* Now loop over all cells in box, updating distance per Thomas
640 method */
641
642 bptr = bestdist;
643 cptr = bestcolor;
644 xx0 = inc0;
645 for ( ic0 = BOX_C0_ELEMS - 1; ic0 >= 0; ic0-- )
646 {
647 dist1 = dist0;
648 xx1 = inc1;
649
650 for ( ic1 = BOX_C1_ELEMS - 1; ic1 >= 0; ic1-- )
651 {
652 dist2 = dist1;
653 xx2 = inc2;
654
655 for ( ic2 = BOX_C2_ELEMS - 1; ic2 >= 0; ic2-- )
656 {
657 if ( dist2 < *bptr )
658 {
659 *bptr = dist2;
660 *cptr = ( JSAMPLE ) icolor;
661 }
662
663 dist2 += xx2;
664 xx2 += 2 * STEP_C2 * STEP_C2;
665 bptr++;
666 cptr++;
667 }
668
669 dist1 += xx1;
670 xx1 += 2 * STEP_C1 * STEP_C1;
671 }
672
673 dist0 += xx0;
674 xx0 += 2 * STEP_C0 * STEP_C0;
675 }
676 }
677 }
678
679
680 /***************************************
681 * Fill the inverse-colormap entries in the update box that contains
682 * histogram cell c0/c1/c2. (Only that one cell MUST be filled, but
683 * we can fill as many others as we wish.)
684 ***************************************/
685
686 static void
fill_inverse_cmap(SPEC * cquantize,int c0,int c1,int c2)687 fill_inverse_cmap( SPEC * cquantize,
688 int c0,
689 int c1,
690 int c2 )
691 {
692 hist3d histogram = cquantize->histogram;
693 int minc0,
694 minc1,
695 minc2; /* lower left corner of update box */
696 int ic0,
697 ic1,
698 ic2;
699 JSAMPLE *cptr; /* pointer into bestcolor[] array */
700 histptr cachep; /* pointer into main cache array */
701 /* This array lists the candidate colormap indexes. */
702 JSAMPLE colorlist[ MAXNUMCOLORS ];
703 int numcolors; /* number of candidate colors */
704 /* This array holds the actually closest colormap index for each cell. */
705 JSAMPLE bestcolor[ BOX_C0_ELEMS * BOX_C1_ELEMS * BOX_C2_ELEMS ];
706
707 /* Convert cell coordinates to update box ID */
708
709 c0 >>= BOX_C0_LOG;
710 c1 >>= BOX_C1_LOG;
711 c2 >>= BOX_C2_LOG;
712
713 /* Compute true coordinates of update box's origin corner. * Actually we
714 compute the coordinates of the center of the corner * histogram cell,
715 which are the lower bounds of the volume we care about. */
716
717 minc0 = ( c0 << BOX_C0_SHIFT ) + ( ( 1 << C0_SHIFT ) >> 1 );
718 minc1 = ( c1 << BOX_C1_SHIFT ) + ( ( 1 << C1_SHIFT ) >> 1 );
719 minc2 = ( c2 << BOX_C2_SHIFT ) + ( ( 1 << C2_SHIFT ) >> 1 );
720
721 /* Determine which colormap entries are close enough to be candidates *
722 for the nearest entry to some cell in the update box. */
723
724 numcolors = find_nearby_colors( cquantize, minc0, minc1, minc2, colorlist );
725
726 /* Determine the actually nearest colors. */
727
728 find_best_colors( cquantize, minc0, minc1, minc2, numcolors, colorlist,
729 bestcolor );
730
731 /* Save the best color numbers (plus 1) in the main cache array */
732
733 c0 <<= BOX_C0_LOG; /* convert ID back to base cell indexes */
734 c1 <<= BOX_C1_LOG;
735 c2 <<= BOX_C2_LOG;
736 cptr = bestcolor;
737
738 for ( ic0 = 0; ic0 < BOX_C0_ELEMS; ic0++ )
739 {
740 for ( ic1 = 0; ic1 < BOX_C1_ELEMS; ic1++ )
741 {
742 cachep = &histogram[ c0 + ic0 ][ c1 + ic1 ][ c2 ];
743
744 for ( ic2 = 0; ic2 < BOX_C2_ELEMS; ic2++ )
745 *cachep++ = ( histcell ) ( GETJSAMPLE( *cptr++ ) + 1 );
746 }
747 }
748 }
749
750
751 #define RIGHT_SHIFT( x, shft) ( ( x ) >> ( shft ) )
752
753
754 /***************************************
755 * This version performs Floyd-Steinberg dithering
756 ***************************************/
757
758 static void
pass2_fs_dither(SPEC * sp,unsigned char ** red,unsigned char ** green,unsigned char ** blue,unsigned short ** output_buf,int width,int num_rows)759 pass2_fs_dither( SPEC * sp,
760 unsigned char ** red,
761 unsigned char ** green,
762 unsigned char ** blue,
763 unsigned short ** output_buf,
764 int width,
765 int num_rows )
766 {
767 hist3d histogram = sp->histogram;
768 LOCFSERROR cur0,
769 cur1,
770 cur2; /* current error or pixel value */
771 LOCFSERROR belowerr0,
772 belowerr1,
773 belowerr2; /* error for pixel below cur */
774 LOCFSERROR bpreverr0,
775 bpreverr1,
776 bpreverr2; /* error for below/prev col */
777 FSERRPTR errorptr; /* => fserrors[] at column before current */
778 unsigned short *outptr; /* => current output pixel */
779 histptr cachep;
780 int dir; /* +1 or -1 depending on direction */
781 int dir3; /* 3*dir, for advancing inptr & errorptr */
782 int row;
783 int col;
784 int *error_limit = sp->error_limiter;
785 int *colormap0 = sp->colormap[0];
786 int *colormap1 = sp->colormap[1];
787 int *colormap2 = sp->colormap[2];
788 unsigned char *r,
789 *g,
790 *b;
791
792 if ( sp->im )
793 {
794 sp->im->completed = -1;
795 sp->im->visual_cue( sp->im, "Dithering ..." );
796 }
797
798 for ( row = 0; row < num_rows; row++ )
799 {
800 r = red[ row ];
801 g = green[ row ];
802 b = blue[ row ];
803
804 outptr = output_buf[ row ];
805
806 if ( sp->on_odd_row )
807 {
808 /* work right to left in this row */
809
810 r += width - 1;
811 g += width - 1;
812 b += width - 1;
813 outptr += width - 1;
814 dir = -1;
815 dir3 = -3;
816 errorptr = sp->fserrors + (width + 1 ) * 3; /* => entry after
817 last column */
818 sp->on_odd_row = 0; /* flip for next time */
819 }
820 else
821 {
822 /* work left to right in this row */
823
824 dir = 1;
825 dir3 = 3;
826 errorptr = sp->fserrors; /* => entry before first real column */
827 sp->on_odd_row = 1; /* flip for next time */
828 }
829
830 /* Preset error values: no error propagated to first pixel from left */
831
832 cur0 = cur1 = cur2 = 0;
833
834 /* and no error propagated to row below yet */
835
836 belowerr0 = belowerr1 = belowerr2 = 0;
837 bpreverr0 = bpreverr1 = bpreverr2 = 0;
838
839 for ( col = 0; col < width; col++ )
840 {
841
842 /* curN holds the error propagated from the previous pixel on the
843 current line. Add the error propagated from the previous
844 line * to form the complete error correction term for this
845 pixel, and * round the error term (which is expressed * 16) to
846 an integer. * RIGHT_SHIFT rounds towards minus infinity, so
847 adding 8 is correct * for either sign of the error value. *
848 Note: errorptr points to *previous* column's array entry. */
849
850 cur0 = RIGHT_SHIFT( cur0 + errorptr[ dir3 + 0 ] + 8, 4 );
851 cur1 = RIGHT_SHIFT( cur1 + errorptr[ dir3 + 1 ] + 8, 4 );
852 cur2 = RIGHT_SHIFT( cur2 + errorptr[ dir3 + 2 ] + 8, 4 );
853
854 /* Limit the error using transfer function set by init_error_limit.
855 See comments with init_error_limit for rationale. */
856
857 cur0 = error_limit[ cur0 ];
858 cur1 = error_limit[ cur1 ];
859 cur2 = error_limit[ cur2 ];
860
861 /* Form pixel value + error, and range-limit to 0..MAXJSAMPLE.
862 The maximum error is +- MAXJSAMPLE (or less with error limiting);
863 this sets the required size of the range_limit array. */
864
865 cur0 += *r;
866 cur1 += *g;
867 cur2 += *b;
868
869 cur0 = FL_PCCLAMP( cur0 );
870 cur1 = FL_PCCLAMP( cur1 );
871 cur2 = FL_PCCLAMP( cur2 );
872
873 /* Index into the cache with adjusted pixel value */
874
875 cachep = &histogram[ cur0 >> C0_SHIFT ][ cur1 >> C1_SHIFT ][ cur2 >> C2_SHIFT ];
876
877 /* If we have not seen this color before, find nearest colormap */
878 /* entry and update the cache */
879
880 if ( *cachep == 0 )
881 fill_inverse_cmap( sp, cur0 >> C0_SHIFT, cur1 >> C1_SHIFT,
882 cur2 >> C2_SHIFT );
883
884 /* Now emit the colormap index for this cell */
885
886 {
887 int pixcode = *cachep - 1;
888 *outptr = ( JSAMPLE ) pixcode;
889
890 /* Compute representation error for this pixel */
891
892 cur0 -= colormap0[ pixcode ];
893 cur1 -= colormap1[ pixcode ];
894 cur2 -= colormap2[ pixcode ];
895 }
896
897 /* Compute error fractions to be propagated to adjacent pixels. *
898 Add these into the running sums, and simultaneously shift the
899 * next-line error sums left by 1 column. */
900
901 {
902 LOCFSERROR bnexterr,
903 delta;
904
905 bnexterr = cur0; /* Process component 0 */
906 delta = cur0 * 2;
907 cur0 += delta; /* form error * 3 */
908 errorptr[ 0 ] = ( FSERROR ) ( bpreverr0 + cur0 );
909 cur0 += delta; /* form error * 5 */
910 bpreverr0 = belowerr0 + cur0;
911 belowerr0 = bnexterr;
912 cur0 += delta; /* form error * 7 */
913 bnexterr = cur1; /* Process component 1 */
914
915 delta = cur1 * 2;
916 cur1 += delta; /* form error * 3 */
917 errorptr[ 1 ] = ( FSERROR ) ( bpreverr1 + cur1 );
918 cur1 += delta; /* form error * 5 */
919 bpreverr1 = belowerr1 + cur1;
920 belowerr1 = bnexterr;
921 cur1 += delta; /* form error * 7 */
922 bnexterr = cur2; /* Process component 2 */
923
924 delta = cur2 * 2;
925 cur2 += delta; /* form error * 3 */
926 errorptr[ 2 ] = ( FSERROR ) ( bpreverr2 + cur2 );
927 cur2 += delta; /* form error * 5 */
928 bpreverr2 = belowerr2 + cur2;
929 belowerr2 = bnexterr;
930 cur2 += delta; /* form error * 7 */
931 }
932
933 /* At this point curN contains the 7/16 error value to be propagated
934 to the next pixel on the current line, and all the errors for the
935 next line have been shifted over. We are therefore ready to move
936 on. */
937
938 r += dir;
939 g += dir;
940 b += dir;
941 outptr += dir;
942 errorptr += dir3; /* advance errorptr to current column */
943 }
944
945 /* Post-loop cleanup: we must unload the final error values into the
946 final fserrors[] entry. Note we need not unload belowerrN
947 because * it is for the dummy column before or after the actual
948 array. */
949
950 errorptr[ 0 ] = ( FSERROR ) bpreverr0; /* unload prev errs into
951 array */
952 errorptr[ 1 ] = ( FSERROR ) bpreverr1;
953 errorptr[ 2 ] = ( FSERROR ) bpreverr2;
954 }
955
956 if( sp->im )
957 {
958 sp->im->completed = sp->im->total = sp->im->h;
959 sp->im->visual_cue( sp->im, "Dithering done" );
960 }
961 }
962
963
964 /***************************************
965 * Shrink the min/max bounds of a box to enclose only nonzero elements,
966 * and recompute its volume and population
967 ***************************************/
968
969 static void
update_box(SPEC * sp,boxptr boxp)970 update_box( SPEC * sp,
971 boxptr boxp )
972 {
973 hist3d histogram = sp->histogram;
974 histptr histp;
975 int c0,
976 c1,
977 c2;
978 int c0min,
979 c0max,
980 c1min,
981 c1max,
982 c2min,
983 c2max;
984 int dist0,
985 dist1,
986 dist2;
987 long ccount;
988
989 c0min = boxp->c0min;
990 c0max = boxp->c0max;
991 c1min = boxp->c1min;
992 c1max = boxp->c1max;
993 c2min = boxp->c2min;
994 c2max = boxp->c2max;
995
996 if ( c0max > c0min )
997 for ( c0 = c0min; c0 <= c0max; c0++ )
998 for ( c1 = c1min; c1 <= c1max; c1++)
999 {
1000 histp = &histogram[ c0 ][ c1 ][ c2min ];
1001
1002 for ( c2 = c2min; c2 <= c2max; c2++ )
1003 if ( *histp++ != 0 )
1004 {
1005 boxp->c0min = c0min = c0;
1006 goto have_c0min;
1007 }
1008 }
1009
1010 have_c0min:
1011
1012 if ( c0max > c0min )
1013 for ( c0 = c0max; c0 >= c0min; c0-- )
1014 for ( c1 = c1min; c1 <= c1max; c1++ )
1015 {
1016 histp = &histogram[ c0 ][ c1 ][ c2min ];
1017
1018 for ( c2 = c2min; c2 <= c2max; c2++ )
1019 if ( *histp++ != 0 )
1020 {
1021 boxp->c0max = c0max = c0;
1022 goto have_c0max;
1023 }
1024 }
1025
1026 have_c0max:
1027
1028 if ( c1max > c1min )
1029 for ( c1 = c1min; c1 <= c1max; c1++ )
1030 for ( c0 = c0min; c0 <= c0max; c0++ )
1031 {
1032 histp = &histogram[ c0 ][ c1 ][ c2min ];
1033
1034 for ( c2 = c2min; c2 <= c2max; c2++ )
1035 if ( *histp++ != 0 )
1036 {
1037 boxp->c1min = c1min = c1;
1038 goto have_c1min;
1039 }
1040 }
1041
1042 have_c1min:
1043
1044 if ( c1max > c1min )
1045 for ( c1 = c1max; c1 >= c1min; c1-- )
1046 for ( c0 = c0min; c0 <= c0max; c0++ )
1047 {
1048 histp = &histogram[ c0 ][ c1 ][ c2min ];
1049
1050 for ( c2 = c2min; c2 <= c2max; c2++ )
1051 if ( *histp++ != 0 )
1052 {
1053 boxp->c1max = c1max = c1;
1054 goto have_c1max;
1055 }
1056 }
1057
1058 have_c1max:
1059
1060 if ( c2max > c2min )
1061 for ( c2 = c2min; c2 <= c2max; c2++ )
1062 for ( c0 = c0min; c0 <= c0max; c0++ )
1063 {
1064 histp = &histogram[ c0 ][ c1min ][ c2 ];
1065
1066 for ( c1 = c1min; c1 <= c1max; c1++, histp += HIST_C2_ELEMS )
1067 if ( *histp != 0 )
1068 {
1069 boxp->c2min = c2min = c2;
1070 goto have_c2min;
1071 }
1072 }
1073
1074 have_c2min:
1075
1076 if ( c2max > c2min )
1077 for ( c2 = c2max; c2 >= c2min; c2-- )
1078 for ( c0 = c0min; c0 <= c0max; c0++ )
1079 {
1080 histp = &histogram[ c0 ][ c1min ][ c2 ];
1081
1082 for ( c1 = c1min; c1 <= c1max; c1++, histp += HIST_C2_ELEMS )
1083 if ( *histp != 0 )
1084 {
1085 boxp->c2max = c2max = c2;
1086 goto have_c2max;
1087 }
1088 }
1089
1090 have_c2max:
1091
1092 /* Update box volume. * We use 2-norm rather than real volume here; this
1093 biases the method * against making long narrow boxes, and it has the
1094 side benefit that * a box is splittable iff norm > 0. * Since the
1095 differences are expressed in histogram-cell units, * we have to shift
1096 back to JSAMPLE units to get consistent distances; * after which, we
1097 scale according to the selected distance scale factors. */
1098
1099 dist0 = ( ( c0max - c0min ) << C0_SHIFT ) * C0_SCALE;
1100 dist1 = ( ( c1max - c1min ) << C1_SHIFT ) * C1_SCALE;
1101 dist2 = ( ( c2max - c2min ) << C2_SHIFT ) * C2_SCALE;
1102 boxp->volume = dist0 * dist0 + dist1 * dist1 + dist2 * dist2;
1103
1104 /* Now scan remaining volume of box and compute population */
1105
1106 ccount = 0;
1107 for ( c0 = c0min; c0 <= c0max; c0++ )
1108 for ( c1 = c1min; c1 <= c1max; c1++ )
1109 {
1110 histp = &histogram[ c0 ][ c1 ][ c2min ];
1111
1112 for ( c2 = c2min; c2 <= c2max; c2++, histp++ )
1113 if ( *histp != 0 )
1114 ccount++;
1115 }
1116
1117 boxp->colorcount = ccount;
1118 }
1119
1120
1121 /***************************************
1122 * Find the splittable box with the largest color population
1123 * Returns NULL if no splittable boxes remain
1124 ***************************************/
1125
1126 static boxptr
find_biggest_color_pop(boxptr boxlist,int numboxes)1127 find_biggest_color_pop( boxptr boxlist,
1128 int numboxes )
1129 {
1130 boxptr boxp;
1131 int i;
1132 long maxc = 0;
1133 boxptr which = NULL;
1134
1135 for ( i = 0, boxp = boxlist; i < numboxes; i++, boxp++ )
1136 if ( boxp->colorcount > maxc && boxp->volume > 0 )
1137 {
1138 which = boxp;
1139 maxc = boxp->colorcount;
1140 }
1141
1142 return which;
1143 }
1144
1145
1146 /***************************************
1147 * Find the splittable box with the largest (scaled) volume
1148 * Returns NULL if no splittable boxes remain
1149 ***************************************/
1150
1151 static boxptr
find_biggest_volume(boxptr boxlist,int numboxes)1152 find_biggest_volume( boxptr boxlist,
1153 int numboxes )
1154 {
1155 boxptr boxp;
1156 int i;
1157 int maxv = 0;
1158 boxptr which = NULL;
1159
1160 for ( i = 0, boxp = boxlist; i < numboxes; i++, boxp++ )
1161 if ( boxp->volume > maxv )
1162 {
1163 which = boxp;
1164 maxv = boxp->volume;
1165 }
1166
1167 return which;
1168 }
1169
1170
1171 /***************************************
1172 * Repeatedly select and split the largest box until we have enough boxes
1173 ***************************************/
1174
1175 static int
median_cut(SPEC * sp,boxptr boxlist,int numboxes,int desired_colors)1176 median_cut( SPEC * sp,
1177 boxptr boxlist,
1178 int numboxes,
1179 int desired_colors )
1180 {
1181 int n,
1182 lb;
1183 int c0,
1184 c1,
1185 c2,
1186 cmax;
1187 boxptr b1,
1188 b2;
1189
1190 while ( numboxes < desired_colors )
1191 {
1192 /* Select box to split. * Current algorithm: by population for first
1193 half, then by volume. */
1194
1195 if ( numboxes * 2 <= desired_colors )
1196 b1 = find_biggest_color_pop( boxlist, numboxes );
1197 else
1198 b1 = find_biggest_volume( boxlist, numboxes );
1199
1200 if ( b1 == NULL ) /* no splittable boxes left! */
1201 break;
1202
1203 b2 = &boxlist[ numboxes ]; /* where new box will go */
1204
1205 /* Copy the color bounds to the new box. */
1206
1207 b2->c0max = b1->c0max;
1208 b2->c1max = b1->c1max;
1209 b2->c2max = b1->c2max;
1210 b2->c0min = b1->c0min;
1211 b2->c1min = b1->c1min;
1212 b2->c2min = b1->c2min;
1213
1214 /* Choose which axis to split the box on. * Current algorithm:
1215 longest scaled axis. * See notes in update_box about scaling
1216 distances. */
1217
1218 c0 = ( ( b1->c0max - b1->c0min ) << C0_SHIFT ) * C0_SCALE;
1219 c1 = ( ( b1->c1max - b1->c1min ) << C1_SHIFT ) * C1_SCALE;
1220 c2 = ( ( b1->c2max - b1->c2min ) << C2_SHIFT ) * C2_SCALE;
1221
1222 /* We want to break any ties in favor of green, then red, blue last.
1223 * This code does the right thing for R,G,B or B,G,R color orders
1224 only. */
1225
1226 #if RGB_RED == 0
1227 cmax = c1;
1228 n = 1;
1229
1230 if ( c0 > cmax )
1231 {
1232 cmax = c0;
1233 n = 0;
1234 }
1235
1236 if ( c2 > cmax )
1237 n = 2;
1238 #else
1239 cmax = c1;
1240 n = 1;
1241
1242 if ( c2 > cmax )
1243 {
1244 cmax = c2;
1245 n = 2;
1246 }
1247
1248 if ( c0 > cmax )
1249 n = 0;
1250 #endif
1251
1252 /* Choose split point along selected axis, and update box bounds. *
1253 Current algorithm: split at halfway point. * (Since the box has
1254 been shrunk to minimum volume, * any split will produce two
1255 nonempty subboxes.) * Note that lb value is max for lower box, so
1256 must be < old max. */
1257
1258 switch ( n )
1259 {
1260 case 0 :
1261 lb = ( b1->c0max + b1->c0min ) / 2;
1262 b1->c0max = lb;
1263 b2->c0min = lb + 1;
1264 break;
1265
1266 case 1 :
1267 lb = ( b1->c1max + b1->c1min ) / 2;
1268 b1->c1max = lb;
1269 b2->c1min = lb + 1;
1270 break;
1271
1272 case 2 :
1273 lb = ( b1->c2max + b1->c2min ) / 2;
1274 b1->c2max = lb;
1275 b2->c2min = lb + 1;
1276 break;
1277 }
1278
1279 /* Update stats for boxes */
1280
1281 update_box( sp, b1 );
1282 update_box( sp, b2 );
1283 numboxes++;
1284 }
1285
1286 return numboxes;
1287 }
1288
1289
1290 /***************************************
1291 * Compute representative color for a box, put it in colormap[icolor]
1292 ***************************************/
1293
1294 static void
compute_color(SPEC * sp,boxptr boxp,int icolor)1295 compute_color( SPEC * sp,
1296 boxptr boxp,
1297 int icolor )
1298 {
1299 /* Current algorithm: mean weighted by pixels (not colors) */
1300 /* Note it is important to get the rounding correct! */
1301
1302 hist3d histogram = sp->histogram;
1303 histptr histp;
1304 int c0,
1305 c1,
1306 c2;
1307 int c0min,
1308 c0max,
1309 c1min,
1310 c1max,
1311 c2min, c2max;
1312 long count;
1313 long total = 0;
1314 long c0total = 0;
1315 long c1total = 0;
1316 long c2total = 0;
1317
1318 c0min = boxp->c0min;
1319 c0max = boxp->c0max;
1320 c1min = boxp->c1min;
1321 c1max = boxp->c1max;
1322 c2min = boxp->c2min;
1323 c2max = boxp->c2max;
1324
1325 for ( c0 = c0min; c0 <= c0max; c0++ )
1326 for ( c1 = c1min; c1 <= c1max; c1++ )
1327 {
1328 histp = &histogram[ c0][ c1 ][ c2min ];
1329
1330 for ( c2 = c2min; c2 <= c2max; c2++ )
1331 {
1332 if ( ( count = *histp++ ) != 0 )
1333 {
1334 total += count;
1335 c0total += ( ( c0 << C0_SHIFT )
1336 + ( ( 1 << C0_SHIFT ) >> 1 ) ) * count;
1337 c1total += ( ( c1 << C1_SHIFT )
1338 + ( ( 1 << C1_SHIFT ) >> 1 ) ) * count;
1339 c2total += ( ( c2 << C2_SHIFT )
1340 + ( ( 1 << C2_SHIFT ) >> 1 ) ) * count;
1341 }
1342 }
1343 }
1344
1345 sp->colormap[ 0 ][ icolor ] =
1346 ( JSAMPLE ) ( ( c0total + ( total >> 1 ) ) / total );
1347 sp->colormap[ 1 ][ icolor ] =
1348 ( JSAMPLE ) ( ( c1total + ( total >> 1 ) ) / total );
1349 sp->colormap[ 2 ][ icolor ] =
1350 ( JSAMPLE ) ( ( c2total + ( total >> 1 ) ) / total );
1351 }
1352
1353
1354 /***************************************
1355 * Master routine for color selection
1356 ***************************************/
1357
1358 static void
select_colors(SPEC * sp,int desired_colors)1359 select_colors( SPEC * sp,
1360 int desired_colors )
1361 {
1362 boxptr boxlist;
1363 int numboxes;
1364 int i;
1365
1366 if ( sp->im )
1367 sp->im->visual_cue( sp->im, "Selecting Colors ..." );
1368
1369 /* Allocate workspace for box list */
1370
1371 boxlist = fl_malloc( desired_colors * sizeof( box ) );
1372
1373 /* Initialize one box containing whole space */
1374
1375 numboxes = 1;
1376 boxlist[ 0 ].c0min = 0;
1377 boxlist[ 0 ].c0max = MAXJSAMPLE >> C0_SHIFT;
1378 boxlist[ 0 ].c1min = 0;
1379 boxlist[ 0 ].c1max = MAXJSAMPLE >> C1_SHIFT;
1380 boxlist[ 0 ].c2min = 0;
1381 boxlist[ 0 ].c2max = MAXJSAMPLE >> C2_SHIFT;
1382
1383 /* Shrink it to actually-used volume and set its statistics */
1384
1385 update_box( sp, &boxlist[ 0 ] );
1386
1387 /* Perform median-cut to produce final box list */
1388
1389 numboxes = median_cut( sp, boxlist, numboxes, desired_colors );
1390
1391 /* Compute the representative color for each box, fill colormap */
1392
1393 for ( i = 0; i < numboxes; i++ )
1394 compute_color( sp, boxlist+ i, i );
1395 sp->actual_number_of_colors = numboxes;
1396 fl_free( boxlist );
1397 }
1398
1399
1400 /***************************************
1401 * get histogram
1402 ***************************************/
1403
1404 static void
prescan_quantize(SPEC * sp,unsigned char ** r,unsigned char ** g,unsigned char ** b,int width,int num_rows)1405 prescan_quantize( SPEC * sp,
1406 unsigned char ** r,
1407 unsigned char ** g,
1408 unsigned char ** b,
1409 int width,
1410 int num_rows )
1411 {
1412 histptr histp;
1413 hist3d histogram = sp->histogram;
1414 int row, col;
1415
1416 if ( sp->im )
1417 {
1418 sp->im->completed = 0;
1419 sp->im->visual_cue( sp->im, "Getting Histogram ..." );
1420 }
1421
1422 for ( row = 0; row < num_rows; row++ )
1423 {
1424 for ( col = width; --col >= 0; )
1425 {
1426 /* get pixel value and index into the histogram */
1427
1428 histp = &histogram[ r[ row ][ col ] >> C0_SHIFT ]
1429 [ g[ row ][ col ] >> C1_SHIFT ]
1430 [ b[ row ][ col ] >> C2_SHIFT ];
1431
1432 /* increment, check for overflow and undo increment if so. */
1433
1434 if ( ++( *histp ) <= 0 )
1435 ( *histp )--;
1436 }
1437 }
1438 }
1439
1440
1441 /*
1442 * Initialize the error-limiting transfer function (lookup table).
1443 * The raw F-S error computation can potentially compute error values of up to
1444 * +- MAXJSAMPLE. But we want the maximum correction applied to a pixel to be
1445 * much less, otherwise obviously wrong pixels will be created. (Typical
1446 * effects include weird fringes at color-area boundaries, isolated bright
1447 * pixels in a dark area, etc.) The standard advice for avoiding this problem
1448 * is to ensure that the "corners" of the color cube are allocated as output
1449 * colors; then repeated errors in the same direction cannot cause cascading
1450 * error buildup. However, that only prevents the error from getting
1451 * completely out of hand; Aaron Giles reports that error limiting improves
1452 * the results even with corner colors allocated.
1453 * A simple clamping of the error values to about +- MAXJSAMPLE/8 works pretty
1454 * well, but the smoother transfer function used below is even better. Thanks
1455 * to Aaron Giles for this idea.
1456 */
1457
1458 /***************************************
1459 * Allocate and fill in the error_limiter table
1460 ***************************************/
1461
1462 static void
init_error_limit(SPEC * sp)1463 init_error_limit( SPEC *sp )
1464 {
1465 int *table;
1466 int in,
1467 out;
1468
1469 table = fl_malloc( ( MAXJSAMPLE * 2 + 1 ) * sizeof( int ) );
1470 table += MAXJSAMPLE; /* so can index -MAXJSAMPLE .. +MAXJSAMPLE */
1471 sp->error_limiter = table;
1472
1473 #define STEPSIZE ( ( MAXJSAMPLE+1 ) / 16 )
1474
1475 /* Map errors 1:1 up to +- MAXJSAMPLE/16 */
1476
1477 for ( out = in = 0; in < STEPSIZE; in++, out++ )
1478 {
1479 table[ in ] = out;
1480 table[ -in ] = -out;
1481 }
1482
1483 /* Map errors 1:2 up to +- 3*MAXJSAMPLE/16 */
1484
1485 for ( ; in < STEPSIZE * 3; in++, out += ( in & 1 ) ? 0 : 1 )
1486 {
1487 table[ in ] = out;
1488 table[ -in ] = -out;
1489 }
1490
1491 /* Clamp the rest to final out value (which is (MAXJSAMPLE+1)/8) */
1492
1493 for ( ; in <= MAXJSAMPLE; in++ )
1494 {
1495 table[ in ] = out;
1496 table[ -in ] = -out;
1497 }
1498 #undef STEPSIZE
1499 }
1500
1501
1502 /*
1503 * Local variables:
1504 * tab-width: 4
1505 * indent-tabs-mode: nil
1506 * End:
1507 */
1508