src/gdcmjpegls/Encoder/initialize.c

   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.
   9
  10    =========================================================
  11    THIS SOFTWARE IS BASED ON HP's implementation of jpeg-ls:
  12    =========================================================
  13
  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    -------------------------------------------------------------------------------
  40 */
  41
  42 /* initialize.c --- functions to initialize look up tables
  43  *                  and statistics tables
  44  *
  45  * Initial code by Alex Jakulin,  Aug. 1995
  46  *
  47  * Modified and optimized: Gadiel Seroussi, October 1995
  48  *
  49  * Modified and added Restart marker and input tables by:
  50  * David Cheng-Hsiu Chu, and Ismail R. Ismail march 1999
  51  */
  52
  53 #include <stdio.h>
  54 #include <math.h>
  55
  56 #include "global.h"
  57 #include "bitio.h"
  58
  59
  60 int vLUT[3][2 * LUTMAX16];
  61
  62 /*byte getk[65][3000];*/
  63 /*int clipPx[510];*/
  64
  65 int classmap[CONTEXTS1];
  66
  67 int  *qdiv0, *qdiv,        /* quantization table (division via look-up) */
  68   *qmul0, *qmul;        /* dequantization table */
  69
  70 int N[TOT_CONTEXTS],
  71     A[TOT_CONTEXTS],
  72     B[TOT_CONTEXTS],
  73   C[TOT_CONTEXTS];
  74
  75
  76
  77
  78
  79
  80
  81 /* Setup Look Up Tables for quantized gradient merging */
  82 void prepareLUTs()
  83 {
  84   int i, j, idx, lmax;
  85   byte k;
  86
  87   lmax = min(alpha,lutmax);
  88
  89   /* implementation limitation: */
  90   if ( T3 > lmax-1 ) {
  91     fprintf(stderr,"Sorry, current implementation does not support threshold T3 > %d, got %d\n",lmax-1,T3);
  92     exit(10);
  93   }
  94
  95
  96   /* Build classification tables (lossless or lossy) */
  97
  98   if (lossy==FALSE) {
  99
 100     for (i = -lmax + 1; i < lmax; i++) {
 101
 102       if  ( i <= -T3 )        /* ...... -T3  */
 103         idx = 7;
 104       else if ( i <= -T2 )    /* -(T3-1) ... -T2 */
 105         idx = 5;
 106       else if ( i <= -T1 )    /* -(T2-1) ... -T1 */
 107         idx = 3;
 108
 109       else if ( i <= -1 )     /* -(T1-1) ...  -1 */
 110         idx = 1;
 111       else if ( i == 0 )      /*  just 0 */
 112         idx = 0;
 113
 114       else if ( i < T1 )      /* 1 ... T1-1 */
 115         idx = 2;
 116       else if ( i < T2 )      /* T1 ... T2-1 */
 117         idx = 4;
 118       else if ( i < T3 )      /* T2 ... T3-1 */
 119         idx = 6;
 120       else                    /* T3 ... */
 121         idx = 8;
 122
 123       vLUT[0][i + lutmax] = CREGIONS * CREGIONS * idx;
 124       vLUT[1][i + lutmax] = CREGIONS * idx;
 125       vLUT[2][i + lutmax] = idx;
 126     }
 127
 128   } else {
 129
 130     for (i = -lmax + 1; i < lmax; i++) {
 131
 132       if ( NEAR >= (alpha-1) )
 133         idx = 0;   /* degenerate case, regardless of thresholds */
 134       else
 135
 136         if  ( i <= -T3 )        /* ...... -T3  */
 137           idx = 7;
 138         else if ( i <= -T2 )    /* -(T3-1) ... -T2 */
 139           idx = 5;
 140         else if ( i <= -T1 )    /* -(T2-1) ... -T1 */
 141           idx = 3;
 142
 143         else if ( i <= -NEAR-1 )     /* -(T1-1) ...  -NEAR-1 */
 144           idx = 1;
 145         else if ( i <= NEAR )      /*  within NEAR of 0 */
 146           idx = 0;
 147
 148         else if ( i < T1 )      /* 1 ... T1-1 */
 149           idx = 2;
 150         else if ( i < T2 )      /* T1 ... T2-1 */
 151           idx = 4;
 152         else if ( i < T3 )      /* T2 ... T3-1 */
 153           idx = 6;
 154         else                    /* T3 ... */
 155           idx = 8;
 156
 157       vLUT[0][i + lutmax] = CREGIONS * CREGIONS * idx;
 158       vLUT[1][i + lutmax] = CREGIONS * idx;
 159       vLUT[2][i + lutmax] = idx;
 160     }
 161
 162   }
 163
 164
 165   /*  prepare context mapping table (symmetric context merging) */
 166   classmap[0] = 0;
 167   for ( i=1, j=0; i<CONTEXTS1; i++) {
 168       int q1, q2, q3, n1=0, n2=0, n3=0, ineg, sgn;
 169
 170       if(classmap[i])
 171       continue;
 172
 173       q1 = i/(CREGIONS*CREGIONS);    /* first digit */
 174       q2 = (i/CREGIONS)%CREGIONS;    /* second digit */
 175       q3 = i%CREGIONS;          /* third digit */
 176
 177       if((q1%2)||(q1==0 && (q2%2))||(q1==0 && q2==0 && (q3%2)))
 178       sgn = -1;
 179       else
 180       sgn = 1;
 181
 182       /* compute negative context */
 183       if(q1) n1 = q1 + ((q1%2) ? 1 : -1);
 184       if(q2) n2 = q2 + ((q2%2) ? 1 : -1);
 185       if(q3) n3 = q3 + ((q3%2) ? 1 : -1);
 186
 187       ineg = (n1*CREGIONS+n2)*CREGIONS+n3;
 188       j++ ;    /* next class number */
 189       classmap[i] = sgn*j;
 190       classmap[ineg] = -sgn*j;
 191
 192   }
 193
 194   /* prepare table to find k */
 195
 196
 197 /*  for (i=1; i< 65; i++)  /* value of N[] */
 198 /*    for (j=1; j<2500; j++)  /* value of A[] */
 199 /*    {
 200       register nst = i;
 201       for(k=0; nst < j; nst<<=1, k++);
 202       getk[i][j]=k;
 203     }*/
 204
 205   /* prepare table to find clip of Px with 8-bit */
 206 /*  for (i=0; i<127; i++)
 207     clipPx[i] = 0;
 208   for (i=127; i<382; i++)
 209     clipPx[i] = i - 127;
 210   for (i=382; i<510; i++)
 211     clipPx[i] = 255;*/
 212
 213 }
 214
 215
 216
 217
 218
 219 /* prepare quantization tables for near-lossless quantization */
 220 void prepare_qtables(int absize, int NEAR)
 221 {
 222     int dif, qdiff;
 223     int beta, quant;
 224
 225     quant = 2*NEAR+1;
 226     beta = absize;
 227
 228     if ( (qdiv0 = (int *)calloc(2*absize-1,sizeof(int)))==NULL ) {
 229       perror("qdiv  table");
 230       exit(10);
 231     }
 232     qdiv = qdiv0+absize-1;
 233
 234     if ( (qmul0 = (int *)calloc(2*beta-1,sizeof(int)))==NULL ) {
 235       perror("qmul  table");
 236       exit(10);
 237     }
 238     qmul = qmul0+beta-1;
 239
 240     for ( dif = -(absize-1); dif<absize; dif++ ) {
 241       if ( dif<0 )
 242         qdiff = - ( (NEAR-dif)/quant );
 243       else
 244         qdiff = ( NEAR + dif )/quant;
 245       qdiv[dif] = qdiff;
 246     }
 247     for ( qdiff = -(beta-1); qdiff<beta; qdiff++ ) {
 248       dif = quant*qdiff;
 249       qmul[qdiff] = dif;
 250     }
 251 }
 252
 253
 254
 255 /* Initialize A[], B[], C[], and N[] arrays */
 256 void init_stats(int absize)
 257 {
 258   int i, initabstat, slack;
 259
 260   slack = 1<<INITABSLACK;
 261   initabstat = (absize + slack/2)/slack;
 262   if ( initabstat < MIN_INITABSTAT ) initabstat = MIN_INITABSTAT;
 263
 264   /* do the statistics initialization */
 265   for (i = 0; i < TOT_CONTEXTS; ++i) {
 266     C[i]= B[i] = 0;
 267     N[i] = INITNSTAT;
 268     A[i] = initabstat;
 269   }
 270 }