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;
322 /**********************************************/
323 /* Do for all pixels in the row in 8-bit mode */
324 /**********************************************/
337 /* Quantize the gradient */
338 cont = vLUT[0][Rd - Rb + LUTMAX8] +
339 vLUT[1][Rb - Rc + LUTMAX8] +
340 vLUT[2][Rc - Ra + LUTMAX8];
344 /********** RUN STATE *********/
348 /* get length of the run */
349 /* arg is # of pixels left */
350 m = n= process_run_dec(no-i+1, color);
352 if ( m > 0 ) { /* run of nonzero length, otherwise
353 we go directly to the end-of-run
363 /* update context pixels */
368 /* here we handle the "end-of-run" state */
369 run_int_type = ( (Rb-Ra) <= NEAR && (Rb-Ra) >= negNEAR);
370 Ra = lossy_end_of_run_d(Ra, Rb, run_int_type);
374 /****** REGULAR CONTEXT ******/
378 /* map symmetric contexts */
379 cont = classmap[cont];
389 /* decode a Rice code of a given context */
390 Ra = lossy_regular_mode_d(cont, SIGN, Px);
401 /***********************************************/
402 /* Do for all pixels in the row in 16-bit mode */
403 /***********************************************/
406 Rc = ENDIAN16(psl[0]);
407 Rb = ENDIAN16(psl[1]);
408 Ra = ENDIAN16(sl[0]);
414 Rd = ENDIAN16(psl[i + 1]);
416 /* Quantize the gradient */
420 /* Following segment assumes that T3 <= LUTMAX16 */
421 /* This condition should have been checked when the
422 lookup tables were built */
425 cont = (diff > -LUTMAX16) ? vLUT[0][diff + LUTMAX16] : 7*CREGIONS*CREGIONS;
427 cont = (diff < LUTMAX16) ? vLUT[0][diff + LUTMAX16] : 8*CREGIONS*CREGIONS;
431 cont += (diff > -LUTMAX16) ? vLUT[1][diff + LUTMAX16] : 7*CREGIONS;
433 cont += (diff < LUTMAX16) ? vLUT[1][diff + LUTMAX16] : 8*CREGIONS;
437 cont += (diff > -LUTMAX16) ? vLUT[2][diff + LUTMAX16] : 7;
439 cont += (diff < LUTMAX16) ? vLUT[2][diff + LUTMAX16] : 8;
444 /********* RUN STATE *********/
448 /* get length of the run */
449 /* arg is # of pixels left */
450 m = n = process_run_dec(no-i+1, color);
452 if ( m > 0 ) { /* run of nonzero length, otherwise
453 we go directly to the end-of-run
456 sl[i++] = ENDIAN16(Ra);
463 /* update context pixels */
464 Rb = ENDIAN16(psl[i]);
465 Rd = ENDIAN16(psl[i + 1]);
468 /* here we handle the "end-of-run" state, which is
469 treated separately from regular states */
470 run_int_type = ( (Rb-Ra) <= NEAR && (Rb-Ra) >= negNEAR);
471 Ra = lossy_end_of_run_d(Ra, Rb, run_int_type);
476 /******REGULAR CONTEXT ******/
480 /* map symmetric contexts */
481 cont = classmap[cont];
491 /* decode a Rice code of a given context */
492 Ra = lossy_regular_mode_d(cont, SIGN, Px);
495 sl[i] = ENDIAN16(Ra);
502 } /* ends "if 8/16 bit" */
513 /* For DECODING pixel interleaved mode in LOSSY mode */
514 int lossy_undoscanline_pixel( pixel *psl, /* previous scanline */
515 pixel *sl, /* current scanline */
516 int no) /* number of values in it */
517 /*** watch it! actual pixels in the scan line are numbered 1 to no .
518 pixels with indices < 1 or > no are dummy "border" pixels */
520 int i, psfix, n_c, color, enter_run=0, was_in_run = 0,
522 pixel Ra, Rb, Rc, Rd;
523 pixel c_aa[MAX_COMPONENTS],
524 c_bb[MAX_COMPONENTS],
525 c_cc[MAX_COMPONENTS],
526 c_dd[MAX_COMPONENTS],
527 c_xx[MAX_COMPONENTS];
529 int cont,c_cont[MAX_COMPONENTS];
534 /**********************************************/
535 /* Do for all pixels in the row in 8-bit mode */
536 /**********************************************/
540 for (n_c=0; n_c<components; n_c++) {
541 c_cc[n_c] = psl[n_c];
542 c_bb[n_c] = psl[components+n_c];
552 if (!was_in_run) color = (color+1)%components;
556 for (n_c=0;n_c<components;n_c++) {
558 c_dd[n_c] = psl[i + components + n_c];
560 /* Quantize the gradient */
561 c_cont[n_c] = vLUT[0][c_dd[n_c] - c_bb[n_c] + LUTMAX8] +
562 vLUT[1][c_bb[n_c] - c_cc[n_c] + LUTMAX8] +
563 vLUT[2][c_cc[n_c] - c_aa[n_c] + LUTMAX8];
572 enter_run = was_in_run = test_run = 0;
576 for (n_c=0;n_c<components;n_c++)
577 if (c_cont[n_c]!=0) {
585 /********* RUN STATE *********/
589 enter_run = was_in_run = 1;
591 /* get length of the run */
592 /* arg is # of pixels left */
593 m = n = process_run_dec((no+components-1-i+1)/components, 0);
595 if ( m > 0 ) { /* run of nonzero length, otherwise
596 we go directly to the end-of-run
599 for (n_c=0;n_c<components;n_c++) {
604 if (i > no+components-1)
608 /* update context pixels */
609 for (n_c=0;n_c<components;n_c++) {
610 c_bb[n_c] = psl[i+n_c];
611 c_dd[n_c] = psl[i+components+n_c];
615 /* here we handle the "end-of-run" state */
616 for (n_c=0;n_c<components;n_c++) {
617 /* The end of run is processed for each component */
621 c_aa[n_c] = c_xx[n_c] = lossy_end_of_run_d(Ra, Rb, 0);
623 } /* Components loop */
627 /****** REGULAR CONTEXT *******/
631 cont = classmap[cont];
641 /* decode a Rice code of a given context */
642 c_aa[color] = Ra = lossy_regular_mode_d(cont, SIGN, Px);
653 for (n_c=0;n_c<components;n_c++) {
654 sl[i+n_c] = c_aa[n_c];
655 c_cc[n_c] = c_bb[n_c];
656 c_bb[n_c] = c_dd[n_c];
661 } while (i <= (no+components-1));
665 /***********************************************/
666 /* Do for all pixels in the row in 16-bit mode */
667 /***********************************************/
670 for (n_c=0; n_c<components; n_c++) {
671 c_cc[n_c] = ENDIAN16(psl[n_c]);
672 c_bb[n_c] = ENDIAN16(psl[components+n_c]);
673 c_aa[n_c] = ENDIAN16(sl[n_c]);
682 if (!was_in_run) color = (color+1)%components;
686 for (n_c=0;n_c<components;n_c++) {
688 c_dd[n_c] = ENDIAN16(psl[i + components + n_c]);
690 /* Quantize the gradient */
694 /* Following segment assumes that T3 <= LUTMAX16 */
695 /* This condition should have been checked when the
696 lookup tables were built */
697 diff = c_dd[n_c] - c_bb[n_c];
699 c_cont[n_c] = (diff > -LUTMAX16) ? vLUT[0][diff + LUTMAX16] : 7*CREGIONS*CREGIONS;
701 c_cont[n_c] = (diff < LUTMAX16) ? vLUT[0][diff + LUTMAX16] : 8*CREGIONS*CREGIONS;
703 diff = c_bb[n_c] - c_cc[n_c];
705 c_cont[n_c] += (diff > -LUTMAX16) ? vLUT[1][diff + LUTMAX16] : 7*CREGIONS;
707 c_cont[n_c] += (diff < LUTMAX16) ? vLUT[1][diff + LUTMAX16] : 8*CREGIONS;
709 diff = c_cc[n_c] - c_aa[n_c];
711 c_cont[n_c] += (diff > -LUTMAX16) ? vLUT[2][diff + LUTMAX16] : 7;
713 c_cont[n_c] += (diff < LUTMAX16) ? vLUT[2][diff + LUTMAX16] : 8;
723 enter_run = was_in_run = test_run = 0;
727 for (n_c=0;n_c<components;n_c++)
728 if (c_cont[n_c]!=0) {
736 /********* RUN STATE *********/
740 enter_run = was_in_run = 1;
742 /* get length of the run */
743 /* arg is # of pixels left */
744 m = n = process_run_dec((no+components-1-i+1)/components, 0);
746 if ( m > 0 ) { /* run of nonzero length, otherwise
747 we go directly to the end-of-run
750 for (n_c=0;n_c<components;n_c++) {
751 sl[i++] = ENDIAN16(c_aa[n_c]);
755 if (i > no+components-1)
759 /* update context pixels */
760 for (n_c=0;n_c<components;n_c++) {
761 c_bb[n_c] = ENDIAN16(psl[i+n_c]);
762 c_dd[n_c] = ENDIAN16(psl[i+components+n_c]);
766 /* here we handle the "end-of-run" state */
767 for (n_c=0;n_c<components;n_c++) {
768 /* The end of run is processed for each component */
771 c_aa[n_c] = c_xx[n_c] = lossy_end_of_run_d(Ra, Rb, 0);
772 } /* Components loop */
777 /******* REGULAR CONTEXT *******/
781 cont = classmap[cont];
790 /* decode a Rice code of a given context */
791 c_aa[color] = Ra = lossy_regular_mode_d(cont, SIGN, Px);
795 sl[i] = ENDIAN16(Ra);
801 for (n_c=0;n_c<components;n_c++) {
802 sl[i+n_c] = ENDIAN16(c_aa[n_c]);
803 c_cc[n_c] = c_bb[n_c];
804 c_bb[n_c] = c_dd[n_c];
809 } while (i <= (no+components-1));
811 } /* for "if 8/16 bit" mode */