]> Creatis software - gdcm.git/blob - src/gdcmjpegls/Encoder/initialize.c
* Fix compilation errors for the Python part
[gdcm.git] / 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 }