1 /* SPMG/JPEG-LS IMPLEMENTATION V.2.1
2 =====================================
3 These programs are Copyright (c) University of British Columbia. All rights reserved.
4 They may be freely redistributed in their entirety provided that this copyright
5 notice is not removed. THEY MAY NOT BE SOLD FOR PROFIT OR INCORPORATED IN
6 COMMERCIAL PROGRAMS WITHOUT THE WRITTEN PERMISSION OF THE COPYRIGHT HOLDER.
7 Each program is provided as is, without any express or implied warranty,
8 without even the warranty of fitness for a particular purpose.
10 =========================================================
11 THIS SOFTWARE IS BASED ON HP's implementation of jpeg-ls:
12 =========================================================
14 LOCO-I/JPEG-LS IMPLEMENTATION V.0.90
15 -------------------------------------------------------------------------------
16 (c) COPYRIGHT HEWLETT-PACKARD COMPANY, 1995-1999.
17 HEWLETT-PACKARD COMPANY ("HP") DOES NOT WARRANT THE ACCURACY OR
18 COMPLETENESS OF THE INFORMATION GIVEN HERE. ANY USE MADE OF, OR
19 RELIANCE ON, SUCH INFORMATION IS ENTIRELY AT USER'S OWN RISK.
20 BY DOWNLOADING THE LOCO-I/JPEG-LS COMPRESSORS/DECOMPRESSORS
21 ("THE SOFTWARE") YOU AGREE TO BE BOUND BY THE TERMS AND CONDITIONS
22 OF THIS LICENSING AGREEMENT.
23 YOU MAY DOWNLOAD AND USE THE SOFTWARE FOR NON-COMMERCIAL PURPOSES
24 FREE OF CHARGE OR FURTHER OBLIGATION. YOU MAY NOT, DIRECTLY OR
25 INDIRECTLY, DISTRIBUTE THE SOFTWARE FOR A FEE, INCORPORATE THIS
26 SOFTWARE INTO ANY PRODUCT OFFERED FOR SALE, OR USE THE SOFTWARE
27 TO PROVIDE A SERVICE FOR WHICH A FEE IS CHARGED.
28 YOU MAY MAKE COPIES OF THE SOFTWARE AND DISTRIBUTE SUCH COPIES TO
29 OTHER PERSONS PROVIDED THAT SUCH COPIES ARE ACCOMPANIED BY
30 HEWLETT-PACKARD'S COPYRIGHT NOTICE AND THIS AGREEMENT AND THAT
31 SUCH OTHER PERSONS AGREE TO BE BOUND BY THE TERMS OF THIS AGREEMENT.
32 THE SOFTWARE IS NOT OF PRODUCT QUALITY AND MAY HAVE ERRORS OR DEFECTS.
33 THE JPEG-LS STANDARD IS STILL UNDER DEVELOPMENT. THE SOFTWARE IS NOT A
34 FINAL OR FULL IMPLEMENTATION OF THE STANDARD. HP GIVES NO EXPRESS OR
35 IMPLIED WARRANTY OF ANY KIND AND ANY IMPLIED WARRANTIES OF
36 MERCHANTABILITY AND FITNESS FOR PURPOSE ARE DISCLAIMED.
37 HP SHALL NOT BE LIABLE FOR ANY DIRECT, INDIRECT, SPECIAL, INCIDENTAL,
38 OR CONSEQUENTIAL DAMAGES ARISING OUT OF ANY USE OF THE SOFTWARE.
39 -------------------------------------------------------------------------------
42 /* lossy_d.c --- the main pipeline which processes a scanline by doing
43 * prediction, context computation, context quantization,
44 * and statistics gathering. (for LOSSY images)
46 * Initial code by Alex Jakulin, Aug. 1995
48 * Modified and optimized: Gadiel Seroussi, October 1995
50 * Modified and added Restart marker and input tables by:
51 * David Cheng-Hsiu Chu, and Ismail R. Ismail march 1999
62 /* Do Golomb-Rice statistics and DECODING for LOSSY images*/
63 inline int lossy_regular_mode_d(int Q, int SIGN, int Px)
65 int At, Bt, Nt, Errval, absErrval;
68 /* This function is called only for regular contexts.
69 End_of_run context is treated separately */
75 register int nst = Nt;
76 for(k=0; nst < At; nst *=2, k++);
79 /* Get the number of leading zeros */
84 temp = zeroLUT[reg >> 24];
93 if ( absErrval < limit ) {
94 /* now add the binary part of the Rice code */
96 register unsigned long temp;
103 /* the original unary would have been too long:
104 (mapped value)-1 was sent verbatim */
105 GETBITS(absErrval, qbpp);
109 /* Do the Rice mapping */
110 if ( absErrval & 1 ) { /* negative */
111 absErrval = (absErrval + 1) / 2;
120 /* if ( k==0 && (2*Bt <= -qmul[Nt]) ) */
121 if ( k==0 && NEAR==0 && (2*Bt <= -Nt) )
123 /* special case: see encoder side */
124 Errval = -(Errval+1);
125 absErrval = (Errval<0)?(-Errval):Errval;
128 Errval = qmul[Errval]; /* dequantize prediction error */
130 /* center, clip if necessary, and mask final error */
134 current = Px - Errval;
139 current = Px + Errval;
142 /* first, we reduce mod beta in the range -NEAR <= x <= alpha-1+NEAR,
143 then we clip to [0,alpha] */
144 if (current < negNEAR)
146 else if (current > alpha1eps)
151 /* update bias stats */
152 B[Q] = (Bt += Errval);
154 /* update Golomb-Rice stats */
157 /* check reset (joint for Rice-Golomb and bias cancelation) */
164 /* Do bias estimation for NEXT pixel */
176 } else if ( Bt > 0 ) {
195 /* Do end of run DECODING for LOSSY images */
196 inline pixel lossy_end_of_run_d(pixel Ra, pixel Rb, int RItype)
216 for(k=0; Nt < At; Nt *=2, k++);
218 /* read and decode the Golomb code */
219 /* Get the number of leading zeros */
224 temp = zeroLUT[reg >> 24];
227 FILLBUFFER(temp + 1);
233 eor_limit = limit - limit_reduce;
235 if ( MErrval < eor_limit ) {
236 /* now add the binary part of the Golomb code */
238 register unsigned long temp;
245 /* the original unary would have been too long:
246 (mapped value)-1 was sent verbatim */
247 GETBITS(MErrval, qbpp);
251 oldmap = ( k==0 && (RItype||MErrval) && (2*B[Q]<Nt));
253 Note: the Boolean variable 'oldmap' is not
254 identical to the variable 'map' in the
255 JPEG-LS draft. We have
256 oldmap = (qdiff<0) ? (1-map) : map;
259 MErrval += ( RItype + oldmap );
261 if ( MErrval & 1 ) { /* negative */
262 Errval = oldmap - (MErrval+1)/2;
263 absErrval = -Errval-RItype;
266 else { /* nonnegative */
268 absErrval = Errval-RItype;
271 Errval = qmul[Errval]; /* de-quantize prediction error */
284 else if ( Ix > alpha1eps )
297 N[Q]++; /* for next pixel */
307 /* For DECODING line and plane interleaved mode for LOSSY images*/
308 int lossy_undoscanline( pixel *psl, /* previous scanline */
309 pixel *sl, /* current scanline */
310 int no, int color) /* number of values in it */
311 /*** watch it! actual pixels in the scan line are numbered 1 to no .
312 pixels with indices < 1 or > no are dummy "border" pixels */
315 pixel Ra, Rb, Rc, Rd;
320 /**********************************************/
321 /* Do for all pixels in the row in 8-bit mode */
322 /**********************************************/
335 /* Quantize the gradient */
336 cont = vLUT[0][Rd - Rb + LUTMAX8] +
337 vLUT[1][Rb - Rc + LUTMAX8] +
338 vLUT[2][Rc - Ra + LUTMAX8];
342 /********** RUN STATE *********/
346 /* get length of the run */
347 /* arg is # of pixels left */
348 m = n= process_run_dec(no-i+1, color);
350 if ( m > 0 ) { /* run of nonzero length, otherwise
351 we go directly to the end-of-run
361 /* update context pixels */
366 /* here we handle the "end-of-run" state */
367 run_int_type = ( (Rb-Ra) <= NEAR && (Rb-Ra) >= negNEAR);
368 Ra = lossy_end_of_run_d(Ra, Rb, run_int_type);
372 /****** REGULAR CONTEXT ******/
376 /* map symmetric contexts */
377 cont = classmap[cont];
387 /* decode a Rice code of a given context */
388 Ra = lossy_regular_mode_d(cont, SIGN, Px);
399 /***********************************************/
400 /* Do for all pixels in the row in 16-bit mode */
401 /***********************************************/
404 Rc = ENDIAN16(psl[0]);
405 Rb = ENDIAN16(psl[1]);
406 Ra = ENDIAN16(sl[0]);
412 Rd = ENDIAN16(psl[i + 1]);
414 /* Quantize the gradient */
418 /* Following segment assumes that T3 <= LUTMAX16 */
419 /* This condition should have been checked when the
420 lookup tables were built */
423 cont = (diff > -LUTMAX16) ? vLUT[0][diff + LUTMAX16] : 7*CREGIONS*CREGIONS;
425 cont = (diff < LUTMAX16) ? vLUT[0][diff + LUTMAX16] : 8*CREGIONS*CREGIONS;
429 cont += (diff > -LUTMAX16) ? vLUT[1][diff + LUTMAX16] : 7*CREGIONS;
431 cont += (diff < LUTMAX16) ? vLUT[1][diff + LUTMAX16] : 8*CREGIONS;
435 cont += (diff > -LUTMAX16) ? vLUT[2][diff + LUTMAX16] : 7;
437 cont += (diff < LUTMAX16) ? vLUT[2][diff + LUTMAX16] : 8;
442 /********* RUN STATE *********/
446 /* get length of the run */
447 /* arg is # of pixels left */
448 m = n = process_run_dec(no-i+1, color);
450 if ( m > 0 ) { /* run of nonzero length, otherwise
451 we go directly to the end-of-run
454 sl[i++] = ENDIAN16(Ra);
461 /* update context pixels */
462 Rb = ENDIAN16(psl[i]);
463 Rd = ENDIAN16(psl[i + 1]);
466 /* here we handle the "end-of-run" state, which is
467 treated separately from regular states */
468 run_int_type = ( (Rb-Ra) <= NEAR && (Rb-Ra) >= negNEAR);
469 Ra = lossy_end_of_run_d(Ra, Rb, run_int_type);
474 /******REGULAR CONTEXT ******/
478 /* map symmetric contexts */
479 cont = classmap[cont];
489 /* decode a Rice code of a given context */
490 Ra = lossy_regular_mode_d(cont, SIGN, Px);
493 sl[i] = ENDIAN16(Ra);
500 } /* ends "if 8/16 bit" */
511 /* For DECODING pixel interleaved mode in LOSSY mode */
512 int lossy_undoscanline_pixel( pixel *psl, /* previous scanline */
513 pixel *sl, /* current scanline */
514 int no) /* number of values in it */
515 /*** watch it! actual pixels in the scan line are numbered 1 to no .
516 pixels with indices < 1 or > no are dummy "border" pixels */
518 int i, n_c, color, was_in_run = 0,
520 pixel Ra, Rb, Rc, Rd;
521 pixel c_aa[MAX_COMPONENTS],
522 c_bb[MAX_COMPONENTS],
523 c_cc[MAX_COMPONENTS],
524 c_dd[MAX_COMPONENTS],
525 c_xx[MAX_COMPONENTS];
527 int cont,c_cont[MAX_COMPONENTS];
530 /**********************************************/
531 /* Do for all pixels in the row in 8-bit mode */
532 /**********************************************/
536 for (n_c=0; n_c<components; n_c++) {
537 c_cc[n_c] = psl[n_c];
538 c_bb[n_c] = psl[components+n_c];
548 if (!was_in_run) color = (color+1)%components;
552 for (n_c=0;n_c<components;n_c++) {
554 c_dd[n_c] = psl[i + components + n_c];
556 /* Quantize the gradient */
557 c_cont[n_c] = vLUT[0][c_dd[n_c] - c_bb[n_c] + LUTMAX8] +
558 vLUT[1][c_bb[n_c] - c_cc[n_c] + LUTMAX8] +
559 vLUT[2][c_cc[n_c] - c_aa[n_c] + LUTMAX8];
568 was_in_run = test_run = 0;
572 for (n_c=0;n_c<components;n_c++)
573 if (c_cont[n_c]!=0) {
581 /********* RUN STATE *********/
587 /* get length of the run */
588 /* arg is # of pixels left */
589 m = n = process_run_dec((no+components-1-i+1)/components, 0);
591 if ( m > 0 ) { /* run of nonzero length, otherwise
592 we go directly to the end-of-run
595 for (n_c=0;n_c<components;n_c++) {
600 if (i > no+components-1)
604 /* update context pixels */
605 for (n_c=0;n_c<components;n_c++) {
606 c_bb[n_c] = psl[i+n_c];
607 c_dd[n_c] = psl[i+components+n_c];
611 /* here we handle the "end-of-run" state */
612 for (n_c=0;n_c<components;n_c++) {
613 /* The end of run is processed for each component */
617 c_aa[n_c] = c_xx[n_c] = lossy_end_of_run_d(Ra, Rb, 0);
619 } /* Components loop */
623 /****** REGULAR CONTEXT *******/
627 cont = classmap[cont];
637 /* decode a Rice code of a given context */
638 c_aa[color] = Ra = lossy_regular_mode_d(cont, SIGN, Px);
649 for (n_c=0;n_c<components;n_c++) {
650 sl[i+n_c] = c_aa[n_c];
651 c_cc[n_c] = c_bb[n_c];
652 c_bb[n_c] = c_dd[n_c];
657 } while (i <= (no+components-1));
661 /***********************************************/
662 /* Do for all pixels in the row in 16-bit mode */
663 /***********************************************/
666 for (n_c=0; n_c<components; n_c++) {
667 c_cc[n_c] = ENDIAN16(psl[n_c]);
668 c_bb[n_c] = ENDIAN16(psl[components+n_c]);
669 c_aa[n_c] = ENDIAN16(sl[n_c]);
678 if (!was_in_run) color = (color+1)%components;
682 for (n_c=0;n_c<components;n_c++) {
684 c_dd[n_c] = ENDIAN16(psl[i + components + n_c]);
686 /* Quantize the gradient */
690 /* Following segment assumes that T3 <= LUTMAX16 */
691 /* This condition should have been checked when the
692 lookup tables were built */
693 diff = c_dd[n_c] - c_bb[n_c];
695 c_cont[n_c] = (diff > -LUTMAX16) ? vLUT[0][diff + LUTMAX16] : 7*CREGIONS*CREGIONS;
697 c_cont[n_c] = (diff < LUTMAX16) ? vLUT[0][diff + LUTMAX16] : 8*CREGIONS*CREGIONS;
699 diff = c_bb[n_c] - c_cc[n_c];
701 c_cont[n_c] += (diff > -LUTMAX16) ? vLUT[1][diff + LUTMAX16] : 7*CREGIONS;
703 c_cont[n_c] += (diff < LUTMAX16) ? vLUT[1][diff + LUTMAX16] : 8*CREGIONS;
705 diff = c_cc[n_c] - c_aa[n_c];
707 c_cont[n_c] += (diff > -LUTMAX16) ? vLUT[2][diff + LUTMAX16] : 7;
709 c_cont[n_c] += (diff < LUTMAX16) ? vLUT[2][diff + LUTMAX16] : 8;
719 was_in_run = test_run = 0;
723 for (n_c=0;n_c<components;n_c++)
724 if (c_cont[n_c]!=0) {
732 /********* RUN STATE *********/
738 /* get length of the run */
739 /* arg is # of pixels left */
740 m = n = process_run_dec((no+components-1-i+1)/components, 0);
742 if ( m > 0 ) { /* run of nonzero length, otherwise
743 we go directly to the end-of-run
746 for (n_c=0;n_c<components;n_c++) {
747 sl[i++] = ENDIAN16(c_aa[n_c]);
751 if (i > no+components-1)
755 /* update context pixels */
756 for (n_c=0;n_c<components;n_c++) {
757 c_bb[n_c] = ENDIAN16(psl[i+n_c]);
758 c_dd[n_c] = ENDIAN16(psl[i+components+n_c]);
762 /* here we handle the "end-of-run" state */
763 for (n_c=0;n_c<components;n_c++) {
764 /* The end of run is processed for each component */
767 c_aa[n_c] = c_xx[n_c] = lossy_end_of_run_d(Ra, Rb, 0);
768 } /* Components loop */
773 /******* REGULAR CONTEXT *******/
777 cont = classmap[cont];
786 /* decode a Rice code of a given context */
787 c_aa[color] = Ra = lossy_regular_mode_d(cont, SIGN, Px);
791 sl[i] = ENDIAN16(Ra);
797 for (n_c=0;n_c<components;n_c++) {
798 sl[i+n_c] = ENDIAN16(c_aa[n_c]);
799 c_cc[n_c] = c_bb[n_c];
800 c_bb[n_c] = c_dd[n_c];
805 } while (i <= (no+components-1));
807 } /* for "if 8/16 bit" mode */