src/lib/openjp3d/t1_3d.c

   1 /*
   2  * Copyrigth (c) 2006, M�nica D�ez, LPI-UVA, Spain
   3  * All rights reserved.
   4  *
   5  * Redistribution and use in source and binary forms, with or without
   6  * modification, are permitted provided that the following conditions
   7  * are met:
   8  * 1. Redistributions of source code must retain the above copyright
   9  *    notice, this list of conditions and the following disclaimer.
  10  * 2. Redistributions in binary form must reproduce the above copyright
  11  *    notice, this list of conditions and the following disclaimer in the
  12  *    documentation and/or other materials provided with the distribution.
  13  *
  14  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS `AS IS'
  15  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  16  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
  17  * ARE DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
  18  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
  19  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
  20  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
  21  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
  22  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
  23  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
  24  * POSSIBILITY OF SUCH DAMAGE.
  25  */
  26
  27 #include "opj_includes.h"
  28
  29 /** @defgroup T1_3D T1_3D - Implementation of the tier-1 coding */
  30 /*@{*/
  31
  32 /** @name Local static functions */
  33 /*@{*/
  34
  35 static int t1_3d_getctxno_zc(unsigned int f, int orient);
  36 static int t1_3d_getctxno_sc(unsigned int f);
  37 static int t1_3d_getctxno_mag(unsigned int f, int fsvr);
  38 static int t1_3d_getspb(unsigned int f);
  39 static int t1_3d_getnmsedec_sig(opj_t1_3d_t *t1, int x, int bitpos);
  40 static int t1_3d_getnmsedec_ref(opj_t1_3d_t *t1, int x, int bitpos);
  41 static void t1_3d_updateflags(unsigned int *fp, int s);
  42 /**
  43 Encode significant pass
  44 */
  45 static void t1_3d_enc_sigpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
  46                                    int *dp, int orient, int bpno, int one, int *nmsedec, char type, int vsc);
  47 /**
  48 Decode significant pass
  49 */
  50 static void t1_3d_dec_sigpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
  51                                    int *dp, int orient, int oneplushalf, char type, int vsc);
  52 /**
  53 Encode significant pass
  54 */
  55 static void t1_3d_enc_sigpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
  56                               int orient, int *nmsedec, char type, int cblksty);
  57 /**
  58 Decode significant pass
  59 */
  60 static void t1_3d_dec_sigpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
  61                               int orient, char type, int cblksty);
  62 /**
  63 Encode refinement pass
  64 */
  65 static void t1_3d_enc_refpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
  66                                    int *dp, int bpno, int one, int *nmsedec, char type, int vsc);
  67 /**
  68 Decode refinement pass
  69 */
  70 static void t1_3d_dec_refpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
  71                                    int *dp, int poshalf, int neghalf, char type, int vsc);
  72 /**
  73 Encode refinement pass
  74 */
  75 static void t1_3d_enc_refpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
  76                               int *nmsedec, char type, int cblksty);
  77 /**
  78 Decode refinement pass
  79 */
  80 static void t1_3d_dec_refpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
  81                               char type, int cblksty);
  82 /**
  83 Encode clean-up pass
  84 */
  85 static void t1_3d_enc_clnpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
  86                                    int *dp, int orient, int bpno, int one, int *nmsedec, int partial, int vsc);
  87 /**
  88 Decode clean-up pass
  89 */
  90 static void t1_3d_dec_clnpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
  91                                    int *dp, int orient, int oneplushalf, int partial, int vsc);
  92 /**
  93 Encode clean-up pass
  94 */
  95 static void t1_3d_enc_clnpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
  96                               int orient, int *nmsedec, int cblksty);
  97 /**
  98 Decode clean-up pass
  99 */
 100 static void t1_3d_dec_clnpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
 101                               int orient, int cblksty);
 102 /**
 103 Encode 1 code-block
 104 @param t1 T1 handle
 105 @param cblk Code-block coding parameters
 106 @param orient
 107 @param compno Component number
 108 @param level[3]
 109 @param dwtid[3]
 110 @param stepsize
 111 @param cblksty Code-block style
 112 @param numcomps
 113 @param tile
 114 */
 115 static void t1_3d_encode_cblk(opj_t1_3d_t *t1, opj_tcd_cblk_t * cblk,
 116                               int orient, int compno,  int level[3], int dwtid[3], double stepsize,
 117                               int cblksty, int numcomps, opj_tcd_tile_t * tile);
 118 /**
 119 Decode 1 code-block
 120 @param t1 T1 handle
 121 @param cblk Code-block coding parameters
 122 @param orient
 123 @param roishift Region of interest shifting value
 124 @param cblksty Code-block style
 125 */
 126 static void t1_3d_decode_cblk(opj_t1_3d_t *t1, opj_tcd_cblk_t * cblk,
 127                               int orient, int roishift, int cblksty);
 128 static int t1_3d_init_ctxno_zc(unsigned int f, int orient);
 129 static int t1_3d_init_ctxno_sc(unsigned int f);
 130 static int t1_3d_init_ctxno_mag(unsigned int f, int f2);
 131 static int t1_3d_init_spb(unsigned int f);
 132 /**
 133 Initialize the look-up tables of the Tier-1 coder/decoder
 134 @param t1 T1 handle
 135 */
 136 static void t1_3d_init_luts(opj_t1_3d_t *t1);
 137
 138 /*@}*/
 139
 140 /*@}*/
 141
 142 /* ----------------------------------------------------------------------- */
 143
 144 static int t1_3d_getctxno_zc(unsigned int f, int orient)
 145 {
 146     return t1_3d_init_ctxno_zc((f & T1_3D_SIG_OTH), orient);
 147 }
 148
 149 static int t1_3d_getctxno_sc(unsigned int f)
 150 {
 151     return t1_3d_init_ctxno_sc((f & T1_3D_SIG_PRIM) | ((f >> 16) & T1_3D_SGN));
 152     /*return t1->lut_ctxno_sc[((f & T1_3D_SIG_PRIM) | ((f >> 16) & T1_3D_SGN)) >> 4];*/
 153 }
 154
 155 static int t1_3d_getctxno_mag(unsigned int f, int fsvr)
 156 {
 157     return t1_3d_init_ctxno_mag((f & T1_3D_SIG_OTH), fsvr);
 158 }
 159
 160 static int t1_3d_getspb(unsigned int f)
 161 {
 162     return t1_3d_init_spb((f & T1_3D_SIG_PRIM) | ((f >> 16) & T1_3D_SGN));
 163     /*return t1->lut_spb[((f & T1_3D_SIG_PRIM) | ((f >> 16) & T1_3D_SGN)) >> 4];*/
 164 }
 165
 166 static int t1_3d_getnmsedec_sig(opj_t1_3d_t *t1, int x, int bitpos)
 167 {
 168     if (bitpos > T1_NMSEDEC_FRACBITS) {
 169         return t1->lut_nmsedec_sig[(x >> (bitpos - T1_NMSEDEC_FRACBITS)) & ((
 170                                                                              1 << T1_NMSEDEC_BITS) - 1)];
 171     }
 172
 173     return t1->lut_nmsedec_sig0[x & ((1 << T1_NMSEDEC_BITS) - 1)];
 174 }
 175
 176 static int t1_3d_getnmsedec_ref(opj_t1_3d_t *t1, int x, int bitpos)
 177 {
 178     if (bitpos > T1_NMSEDEC_FRACBITS) {
 179         return t1->lut_nmsedec_ref[(x >> (bitpos - T1_NMSEDEC_FRACBITS)) & ((
 180                                                                              1 << T1_NMSEDEC_BITS) - 1)];
 181     }
 182
 183     return t1->lut_nmsedec_ref0[x & ((1 << T1_NMSEDEC_BITS) - 1)];
 184 }
 185
 186 static void t1_3d_updateflags(unsigned int *fp, int s)
 187 {
 188     unsigned int *np = fp - (T1_MAXCBLKW + 2);
 189     unsigned int *sp = fp + (T1_MAXCBLKW + 2);
 190
 191     unsigned int *bwp = fp + ((T1_MAXCBLKW + 2) * (T1_MAXCBLKH + 2));
 192     unsigned int *bnp = bwp - (T1_MAXCBLKW + 2);
 193     unsigned int *bsp = bwp + (T1_MAXCBLKW + 2);
 194
 195     unsigned int *fwp = fp - ((T1_MAXCBLKW + 2) * (T1_MAXCBLKH + 2));
 196     unsigned int *fnp = fwp - (T1_MAXCBLKW + 2);
 197     unsigned int *fsp = fwp + (T1_MAXCBLKW + 2);
 198
 199     np[-1] |= T1_3D_SIG_SE;
 200     np[1] |= T1_3D_SIG_SW;
 201     sp[-1] |= T1_3D_SIG_NE;
 202     sp[1] |= T1_3D_SIG_NW;
 203     *np |= T1_3D_SIG_S;
 204     *sp |= T1_3D_SIG_N;
 205     fp[-1] |= T1_3D_SIG_E;
 206     fp[1] |= T1_3D_SIG_W;
 207
 208     *fwp |= T1_3D_SIG_FC;
 209     *bwp |= T1_3D_SIG_BC;
 210
 211     fnp[-1] |= T1_3D_SIG_FSE;
 212     fnp[1] |= T1_3D_SIG_FSW;
 213     fsp[-1] |= T1_3D_SIG_FNE;
 214     fsp[1] |= T1_3D_SIG_FNW;
 215     *fnp |= T1_3D_SIG_FS;
 216     *fsp |= T1_3D_SIG_FN;
 217     fwp[-1] |= T1_3D_SIG_FE;
 218     fwp[1] |= T1_3D_SIG_FW;
 219
 220     bnp[-1] |= T1_3D_SIG_BSE;
 221     bnp[1] |= T1_3D_SIG_BSW;
 222     bsp[-1] |= T1_3D_SIG_BNE;
 223     bsp[1] |= T1_3D_SIG_BNW;
 224     *bnp |= T1_3D_SIG_BS;
 225     *bsp |= T1_3D_SIG_BN;
 226     bwp[-1] |= T1_3D_SIG_BE;
 227     bwp[1] |= T1_3D_SIG_BW;
 228
 229     if (s) {
 230         *np |= (T1_3D_SGN_S << 16);
 231         *sp |= (T1_3D_SGN_N << 16);
 232         fp[-1] |= (T1_3D_SGN_E << 16);
 233         fp[1] |= (T1_3D_SGN_W << 16);
 234         *fwp |= (T1_3D_SGN_F << 16);
 235         *bwp |= (T1_3D_SGN_B << 16);
 236     }
 237 }
 238
 239 static void t1_3d_enc_sigpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
 240                                    int *dp, int orient, int bpno, int one, int *nmsedec, char type, int vsc)
 241 {
 242     int v, flagsvr;
 243     unsigned int flag;
 244
 245     opj_mqc_t *mqc = t1->mqc;   /* MQC component */
 246
 247     flag = vsc ? ((*fp) & (~(T1_3D_SIG_S | T1_3D_SIG_SE | T1_3D_SIG_SW |
 248                              (T1_3D_SGN_S << 16)))) : (*fp);
 249     flagsvr = (*fsvr);
 250     if ((flag & T1_3D_SIG_OTH) && !(flagsvr & (T1_3D_SIG | T1_3D_VISIT))) {
 251         v = int_abs(*dp) & one ? 1 : 0;
 252         if (type == T1_TYPE_RAW) {  /* BYPASS/LAZY MODE */
 253             mqc_setcurctx(mqc, t1_3d_getctxno_zc(flag, orient));    /* ESSAI */
 254             mqc_bypass_enc(mqc, v);
 255         } else {
 256             mqc_setcurctx(mqc, t1_3d_getctxno_zc(flag, orient));
 257             mqc_encode(mqc, v);
 258         }
 259         if (v) {
 260             v = *dp < 0 ? 1 : 0;
 261             *nmsedec += t1_3d_getnmsedec_sig(t1, int_abs(*dp), bpno + T1_NMSEDEC_FRACBITS);
 262             if (type == T1_TYPE_RAW) {  /* BYPASS/LAZY MODE */
 263                 mqc_setcurctx(mqc, t1_3d_getctxno_sc(flag));    /* ESSAI */
 264                 mqc_bypass_enc(mqc, v);
 265             } else {
 266                 mqc_setcurctx(mqc, t1_3d_getctxno_sc(flag));
 267                 mqc_encode(mqc, v ^ t1_3d_getspb(flag));
 268             }
 269             t1_3d_updateflags(fp, v);
 270             *fsvr |= T1_3D_SIG;
 271         }
 272         *fsvr |= T1_3D_VISIT;
 273     }
 274 }
 275
 276 static void t1_3d_dec_sigpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
 277                                    int *dp, int orient, int oneplushalf, char type, int vsc)
 278 {
 279     int v, flagsvr;
 280     unsigned int flag;
 281
 282     opj_raw_t *raw = t1->raw;   /* RAW component */
 283     opj_mqc_t *mqc = t1->mqc;   /* MQC component */
 284
 285     flag = vsc ? ((*fp) & (~(T1_3D_SIG_S | T1_3D_SIG_SE | T1_3D_SIG_SW |
 286                              (T1_3D_SGN_S << 16)))) : (*fp);
 287     flagsvr = (*fsvr);
 288     if ((flag & T1_3D_SIG_OTH) && !(flagsvr & (T1_3D_SIG | T1_3D_VISIT))) {
 289         if (type == T1_TYPE_RAW) {
 290             if (raw_decode(raw)) {
 291                 v = raw_decode(raw);    /* ESSAI */
 292                 *dp = v ? -oneplushalf : oneplushalf;
 293                 t1_3d_updateflags(fp, v);
 294                 *fsvr |= T1_3D_SIG;
 295             }
 296         } else {
 297             mqc_setcurctx(mqc, t1_3d_getctxno_zc(flag, orient));
 298             if (mqc_decode(mqc)) {
 299                 mqc_setcurctx(mqc, t1_3d_getctxno_sc(flag));
 300                 v = mqc_decode(mqc) ^ t1_3d_getspb(flag);
 301                 *dp = v ? -oneplushalf : oneplushalf;
 302                 t1_3d_updateflags(fp, v);
 303                 *fsvr |= T1_3D_SIG;
 304             }
 305         }
 306         *fsvr |= T1_3D_VISIT;
 307     }
 308 }               /* VSC and  BYPASS by Antonin */
 309
 310 static void t1_3d_enc_sigpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
 311                               int orient, int *nmsedec, char type, int cblksty)
 312 {
 313     int i, j, k, m, one, vsc;
 314     *nmsedec = 0;
 315     one = 1 << (bpno + T1_NMSEDEC_FRACBITS);
 316     for (m = 0; m < l; m++) {
 317         for (k = 0; k < h; k += 4) {
 318             for (i = 0; i < w; i++) {
 319                 for (j = k; j < k + 4 && j < h; j++) {
 320                     vsc = ((cblksty & J3D_CCP_CBLKSTY_VSC) && (j == k + 3 || j == h - 1)) ? 1 : 0;
 321                     t1_3d_enc_sigpass_step(t1, &t1->flags[1 + m][1 + j][1 + i],
 322                                            &t1->flagSVR[1 + m][1 + j][1 + i], &t1->data[m][j][i], orient, bpno, one,
 323                                            nmsedec, type, vsc);
 324                 }
 325             }
 326         }
 327     }
 328 }
 329
 330 static void t1_3d_dec_sigpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
 331                               int orient, char type, int cblksty)
 332 {
 333     int i, j, k, m, one, half, oneplushalf, vsc;
 334     one = 1 << bpno;
 335     half = one >> 1;
 336     oneplushalf = one | half;
 337     for (m = 0; m < l; m++) {
 338         for (k = 0; k < h; k += 4) {
 339             for (i = 0; i < w; i++) {
 340                 for (j = k; j < k + 4 && j < h; j++) {
 341                     vsc = ((cblksty & J3D_CCP_CBLKSTY_VSC) && (j == k + 3 || j == h - 1)) ? 1 : 0;
 342                     t1_3d_dec_sigpass_step(t1, &t1->flags[1 + m][1 + j][1 + i],
 343                                            &t1->flagSVR[1 + m][1 + j][1 + i], &t1->data[m][j][i], orient, oneplushalf,
 344                                            type, vsc);
 345                 }
 346             }
 347         }
 348     }
 349 }               /* VSC and  BYPASS by Antonin */
 350
 351 static void t1_3d_enc_refpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
 352                                    int *dp, int bpno, int one, int *nmsedec, char type, int vsc)
 353 {
 354     int v, flagsvr;
 355     unsigned int flag;
 356
 357     opj_mqc_t *mqc = t1->mqc;   /* MQC component */
 358
 359     flag = vsc ? ((*fp) & (~(T1_3D_SIG_S | T1_3D_SIG_SE | T1_3D_SIG_SW |
 360                              (T1_3D_SGN_S << 16)))) : (*fp);
 361     flagsvr = (*fsvr);
 362     if ((flagsvr & (T1_3D_SIG | T1_3D_VISIT)) == T1_3D_SIG) {
 363         *nmsedec += t1_3d_getnmsedec_ref(t1, int_abs(*dp), bpno + T1_NMSEDEC_FRACBITS);
 364         v = int_abs(*dp) & one ? 1 : 0;
 365         if (type == T1_TYPE_RAW) {  /* BYPASS/LAZY MODE */
 366             mqc_setcurctx(mqc, t1_3d_getctxno_mag(flag, flagsvr));  /* ESSAI */
 367             mqc_bypass_enc(mqc, v);
 368         } else {
 369             mqc_setcurctx(mqc, t1_3d_getctxno_mag(flag, flagsvr));
 370             mqc_encode(mqc, v);
 371         }
 372         *fsvr |= T1_3D_REFINE;
 373     }
 374 }
 375
 376 static void t1_3d_dec_refpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
 377                                    int *dp, int poshalf, int neghalf, char type, int vsc)
 378 {
 379     int v, t, flagsvr;
 380     unsigned int flag;
 381
 382     opj_mqc_t *mqc = t1->mqc;   /* MQC component */
 383     opj_raw_t *raw = t1->raw;   /* RAW component */
 384
 385     flag = vsc ? ((*fp) & (~(T1_3D_SIG_S | T1_3D_SIG_SE | T1_3D_SIG_SW |
 386                              (T1_3D_SGN_S << 16)))) : (*fp);
 387     flagsvr = (*fsvr);
 388     if ((flagsvr & (T1_3D_SIG | T1_3D_VISIT)) == T1_3D_SIG) {
 389         if (type == T1_TYPE_RAW) {
 390             mqc_setcurctx(mqc, t1_3d_getctxno_mag(flag, flagsvr));  /* ESSAI */
 391             v = raw_decode(raw);
 392         } else {
 393             mqc_setcurctx(mqc, t1_3d_getctxno_mag(flag, flagsvr));
 394             v = mqc_decode(mqc);
 395         }
 396         t = v ? poshalf : neghalf;
 397         *dp += *dp < 0 ? -t : t;
 398         *fsvr |= T1_3D_REFINE;
 399     }
 400 }               /* VSC and  BYPASS by Antonin  */
 401
 402 static void t1_3d_enc_refpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
 403                               int *nmsedec, char type, int cblksty)
 404 {
 405     int i, j, k, m, one, vsc;
 406     *nmsedec = 0;
 407     one = 1 << (bpno + T1_NMSEDEC_FRACBITS);
 408     for (m = 0; m < l; m++) {
 409         for (k = 0; k < h; k += 4) {
 410             for (i = 0; i < w; i++) {
 411                 for (j = k; j < k + 4 && j < h; j++) {
 412                     vsc = ((cblksty & J3D_CCP_CBLKSTY_VSC) && (j == k + 3 || j == h - 1)) ? 1 : 0;
 413                     t1_3d_enc_refpass_step(t1, &t1->flags[1 + m][1 + j][1 + i],
 414                                            &t1->flagSVR[1 + m][1 + j][1 + i], &t1->data[m][j][i], bpno, one, nmsedec, type,
 415                                            vsc);
 416                 }
 417             }
 418         }
 419     }
 420 }
 421
 422 static void t1_3d_dec_refpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
 423                               char type, int cblksty)
 424 {
 425     int i, j, k, m, one, poshalf, neghalf;
 426     int vsc;
 427     one = 1 << bpno;
 428     poshalf = one >> 1;
 429     neghalf = bpno > 0 ? -poshalf : -1;
 430     for (m = 0; m < l; m++) {
 431         for (k = 0; k < h; k += 4) {
 432             for (i = 0; i < w; i++) {
 433                 for (j = k; j < k + 4 && j < h; j++) {
 434                     vsc = ((cblksty & J3D_CCP_CBLKSTY_VSC) && (j == k + 3 || j == h - 1)) ? 1 : 0;
 435                     t1_3d_dec_refpass_step(t1, &t1->flags[1 + m][1 + j][1 + i],
 436                                            &t1->flagSVR[1 + m][1 + j][1 + i], &t1->data[m][j][i], poshalf, neghalf, type,
 437                                            vsc);
 438                 }
 439             }
 440         }
 441     }
 442 }               /* VSC and  BYPASS by Antonin */
 443
 444 static void t1_3d_enc_clnpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
 445                                    int *dp, int orient, int bpno, int one, int *nmsedec, int partial, int vsc)
 446 {
 447     int v, flagsvr;
 448     unsigned int flag;
 449
 450     opj_mqc_t *mqc = t1->mqc;   /* MQC component */
 451
 452     flag = vsc ? ((*fp) & (~(T1_3D_SIG_S | T1_3D_SIG_SE | T1_3D_SIG_SW |
 453                              (T1_3D_SGN_S << 16)))) : (*fp);
 454     flagsvr = (*fsvr);
 455     if (partial) {
 456         goto LABEL_PARTIAL;
 457     }
 458     if (!(*fsvr & (T1_3D_SIG | T1_3D_VISIT))) {
 459         mqc_setcurctx(mqc, t1_3d_getctxno_zc(flag, orient));
 460         v = int_abs(*dp) & one ? 1 : 0;
 461         mqc_encode(mqc, v);
 462         if (v) {
 463 LABEL_PARTIAL:
 464             *nmsedec += t1_3d_getnmsedec_sig(t1, int_abs(*dp), bpno + T1_NMSEDEC_FRACBITS);
 465             mqc_setcurctx(mqc, t1_3d_getctxno_sc(flag));
 466             v = *dp < 0 ? 1 : 0;
 467             mqc_encode(mqc, v ^ t1_3d_getspb(flag));
 468             t1_3d_updateflags(fp, v);
 469             *fsvr |= T1_3D_SIG;
 470         }
 471     }
 472     *fsvr &= ~T1_3D_VISIT;
 473 }
 474
 475 static void t1_3d_dec_clnpass_step(opj_t1_3d_t *t1, unsigned int *fp, int *fsvr,
 476                                    int *dp, int orient, int oneplushalf, int partial, int vsc)
 477 {
 478     int v, flagsvr;
 479     unsigned int flag;
 480
 481     opj_mqc_t *mqc = t1->mqc;   /* MQC component */
 482
 483     flag = vsc ? ((*fp) & (~(T1_3D_SIG_S | T1_3D_SIG_SE | T1_3D_SIG_SW |
 484                              (T1_3D_SGN_S << 16)))) : (*fp);
 485     flagsvr = (*fsvr);
 486     if (partial) {
 487         goto LABEL_PARTIAL;
 488     }
 489     if (!(flagsvr & (T1_3D_SIG | T1_3D_VISIT))) {
 490         mqc_setcurctx(mqc, t1_3d_getctxno_zc(flag, orient));
 491         if (mqc_decode(mqc)) {
 492 LABEL_PARTIAL:
 493             mqc_setcurctx(mqc, t1_3d_getctxno_sc(flag));
 494             v = mqc_decode(mqc) ^ t1_3d_getspb(flag);
 495             *dp = v ? -oneplushalf : oneplushalf;
 496             t1_3d_updateflags(fp, v);
 497             *fsvr |= T1_3D_SIG;
 498         }
 499     }
 500     *fsvr &= ~T1_3D_VISIT;
 501 }               /* VSC and  BYPASS by Antonin */
 502
 503 static void t1_3d_enc_clnpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
 504                               int orient, int *nmsedec, int cblksty)
 505 {
 506     int i, j, k, m, one, agg, runlen, vsc;
 507
 508     opj_mqc_t *mqc = t1->mqc;   /* MQC component */
 509
 510     *nmsedec = 0;
 511     one = 1 << (bpno + T1_NMSEDEC_FRACBITS);
 512     for (m = 0; m < l; m++) {
 513         for (k = 0; k < h; k += 4) {
 514             for (i = 0; i < w; i++) {
 515                 if (k + 3 < h) {
 516                     if (cblksty & J3D_CCP_CBLKSTY_VSC) {
 517                         agg = !(((t1->flagSVR[1 + m][1 + k][1 + i] | (T1_3D_SIG | T1_3D_VISIT)) &
 518                                  (t1->flags[1 + m][1 + k][1 + i] & T1_3D_SIG_OTH))
 519                                 || ((t1->flagSVR[1 + m][1 + k + 1][1 + i] | (T1_3D_SIG | T1_3D_VISIT)) &
 520                                     (t1->flags[1 + m][1 + k + 1][1 + i] & T1_3D_SIG_OTH))
 521                                 || ((t1->flagSVR[1 + m][1 + k + 2][1 + i] | (T1_3D_SIG | T1_3D_VISIT)) &
 522                                     (t1->flags[1 + m][1 + k + 2][1 + i] & T1_3D_SIG_OTH))
 523                                 || ((t1->flagSVR[1 + m][1 + k + 3][1 + i] | (T1_3D_SIG | T1_3D_VISIT)) & ((
 524                                             t1->flags[1 + m][1 + k + 3][1 + i] & (~(T1_3D_SIG_S | T1_3D_SIG_SE |
 525                                                     T1_3D_SIG_SW | (T1_3D_SGN_S << 16)))) & (T1_3D_SIG_OTH)))
 526                                );
 527                     } else {
 528                         agg = !(
 529                                   ((t1->flagSVR[1 + m][1 + k][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) |
 530                                    (t1->flags[1 + m][1 + k][1 + i] & T1_3D_SIG_OTH))
 531                                   || ((t1->flagSVR[1 + m][1 + k + 1][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) |
 532                                       (t1->flags[1 + m][1 + k + 1][1 + i] & T1_3D_SIG_OTH))
 533                                   || ((t1->flagSVR[1 + m][1 + k + 2][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) |
 534                                       (t1->flags[1 + m][1 + k + 2][1 + i] & T1_3D_SIG_OTH))
 535                                   || ((t1->flagSVR[1 + m][1 + k + 3][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) |
 536                                       (t1->flags[1 + m][1 + k + 3][1 + i] & T1_3D_SIG_OTH))
 537                               );
 538                     }
 539                 } else {
 540                     agg = 0;
 541                 }
 542                 if (agg) {
 543                     for (runlen = 0; runlen < 4; runlen++) {
 544                         if (int_abs(t1->data[m][k + runlen][i]) & one) {
 545                             break;
 546                         }
 547                     }
 548                     mqc_setcurctx(mqc, T1_CTXNO_AGG);
 549                     mqc_encode(mqc, runlen != 4);
 550                     if (runlen == 4) {
 551                         continue;
 552                     }
 553                     mqc_setcurctx(mqc, T1_CTXNO_UNI);
 554                     mqc_encode(mqc, runlen >> 1);
 555                     mqc_encode(mqc, runlen & 1);
 556                 } else {
 557                     runlen = 0;
 558                 }
 559                 for (j = k + runlen; j < k + 4 && j < h; j++) {
 560                     vsc = ((cblksty & J3D_CCP_CBLKSTY_VSC) && (j == k + 3 || j == h - 1)) ? 1 : 0;
 561                     t1_3d_enc_clnpass_step(t1, &t1->flags[1 + m][1 + j][1 + i],
 562                                            &t1->flagSVR[1 + m][1 + j][1 + i], &t1->data[m][j][i], orient, bpno, one,
 563                                            nmsedec, agg && (j == k + runlen), vsc);
 564                 }
 565             }
 566         }
 567     }
 568 }
 569
 570 static void t1_3d_dec_clnpass(opj_t1_3d_t *t1, int w, int h, int l, int bpno,
 571                               int orient, int cblksty)
 572 {
 573     int i, j, k, m, one, half, oneplushalf, agg, runlen, vsc;
 574     int segsym = cblksty & J3D_CCP_CBLKSTY_SEGSYM;
 575
 576     opj_mqc_t *mqc = t1->mqc;   /* MQC component */
 577
 578     one = 1 << bpno;
 579     half = one >> 1;
 580     oneplushalf = one | half;
 581     for (m = 0; m < l; m++) {
 582         for (k = 0; k < h; k += 4) {
 583             for (i = 0; i < w; i++) {
 584                 if (k + 3 < h) {
 585                     if (cblksty & J3D_CCP_CBLKSTY_VSC) {
 586                         agg = !(
 587                                   ((t1->flagSVR[1 + m][1 + k][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) |
 588                                    (t1->flags[1 + m][1 + k][1 + i] & T1_3D_SIG_OTH))
 589                                   || ((t1->flagSVR[1 + m][1 + k + 1][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) |
 590                                       (t1->flags[1 + m][1 + k + 1][1 + i] & T1_3D_SIG_OTH))
 591                                   || ((t1->flagSVR[1 + m][1 + k + 2][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) |
 592                                       (t1->flags[1 + m][1 + k + 2][1 + i] & T1_3D_SIG_OTH))
 593                                   || ((t1->flagSVR[1 + m][1 + k + 3][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) | ((
 594                                               t1->flags[1 + m][1 + k + 3][1 + i] & (~(T1_3D_SIG_S | T1_3D_SIG_SE |
 595                                                       T1_3D_SIG_SW | (T1_3D_SGN_S << 16)))) & (T1_3D_SIG_OTH)))
 596                               );
 597                     } else {
 598                         agg = !(
 599                                   ((t1->flagSVR[1 + m][1 + k][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) |
 600                                    (t1->flags[1 + m][1 + k][1 + i] & T1_3D_SIG_OTH))
 601                                   || ((t1->flagSVR[1 + m][1 + k + 1][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) |
 602                                       (t1->flags[1 + m][1 + k + 1][1 + i] & T1_3D_SIG_OTH))
 603                                   || ((t1->flagSVR[1 + m][1 + k + 2][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) |
 604                                       (t1->flags[1 + m][1 + k + 2][1 + i] & T1_3D_SIG_OTH))
 605                                   || ((t1->flagSVR[1 + m][1 + k + 3][1 + i] & (T1_3D_SIG | T1_3D_VISIT)) |
 606                                       (t1->flags[1 + m][1 + k + 3][1 + i] & T1_3D_SIG_OTH))
 607                               );
 608                     }
 609                 } else {
 610                     agg = 0;
 611                 }
 612                 if (agg) {
 613                     mqc_setcurctx(mqc, T1_CTXNO_AGG);
 614                     if (!mqc_decode(mqc)) {
 615                         continue;
 616                     }
 617                     mqc_setcurctx(mqc, T1_CTXNO_UNI);
 618                     runlen = mqc_decode(mqc);
 619                     runlen = (runlen << 1) | mqc_decode(mqc);
 620                 } else {
 621                     runlen = 0;
 622                 }
 623                 for (j = k + runlen; j < k + 4 && j < h; j++) {
 624                     vsc = ((cblksty & J3D_CCP_CBLKSTY_VSC) && (j == k + 3 || j == h - 1)) ? 1 : 0;
 625                     t1_3d_dec_clnpass_step(t1, &t1->flags[1 + m][1 + j][1 + i],
 626                                            &t1->flagSVR[1 + m][1 + j][1 + i], &t1->data[m][j][i], orient, oneplushalf,
 627                                            agg && (j == k + runlen), vsc);
 628                 }
 629             }
 630         }
 631     }
 632     if (segsym) {
 633         int v = 0;
 634         mqc_setcurctx(mqc, T1_CTXNO_UNI);
 635         v = mqc_decode(mqc);
 636         v = (v << 1) | mqc_decode(mqc);
 637         v = (v << 1) | mqc_decode(mqc);
 638         v = (v << 1) | mqc_decode(mqc);
 639         /*
 640         if (v!=0xa) {
 641             opj_event_msg(t1->cinfo, EVT_WARNING, "Bad segmentation symbol %x\n", v);
 642         }
 643         */
 644     }
 645 }               /* VSC and  BYPASS by Antonin */
 646
 647
 648 static void t1_3d_encode_cblk(opj_t1_3d_t *t1, opj_tcd_cblk_t * cblk,
 649                               int orient, int compno, int level[3], int dwtid[3], double stepsize,
 650                               int cblksty, int numcomps, opj_tcd_tile_t * tile)
 651 {
 652     int i, j, k;
 653     int w, h, l;
 654     int passno;
 655     int bpno, passtype;
 656     int max;
 657     int nmsedec = 0;
 658     double cumwmsedec = 0;
 659     char type = T1_TYPE_MQ;
 660
 661     opj_mqc_t *mqc = t1->mqc;   /* MQC component */
 662
 663     w = cblk->x1 - cblk->x0;
 664     h = cblk->y1 - cblk->y0;
 665     l = cblk->z1 - cblk->z0;
 666
 667     max = 0;
 668     for (k = 0; k < l; k++) {
 669         for (j = 0; j < h; j++) {
 670             for (i = 0; i < w; i++) {
 671                 max = int_max(max, int_abs(t1->data[k][j][i]));
 672             }
 673         }
 674     }
 675     for (k = 0; k <= l; k++) {
 676         for (j = 0; j <= h; j++) {
 677             for (i = 0; i <= w; i++) {
 678                 t1->flags[k][j][i] = 0;
 679                 t1->flagSVR[k][j][i] = 0;
 680             }
 681         }
 682     }
 683
 684     cblk->numbps = max ? (int_floorlog2(max) + 1) - T1_NMSEDEC_FRACBITS : 0;
 685
 686     bpno = cblk->numbps - 1;
 687     passtype = 2;
 688
 689     mqc_reset_enc(mqc);
 690     mqc_init_enc(mqc, cblk->data);
 691
 692     for (passno = 0; bpno >= 0; passno++) {
 693         opj_tcd_pass_t *pass = &cblk->passes[passno];
 694         int correction = 3;
 695         double tmpwmsedec;
 696         type = ((bpno < (cblk->numbps - 4)) && (passtype < 2) &&
 697                 (cblksty & J3D_CCP_CBLKSTY_LAZY)) ? T1_TYPE_RAW : T1_TYPE_MQ;
 698
 699         switch (passtype) {
 700         case 0:
 701             t1_3d_enc_sigpass(t1, w, h, l, bpno, orient, &nmsedec, type, cblksty);
 702             break;
 703         case 1:
 704             t1_3d_enc_refpass(t1, w, h, l, bpno, &nmsedec, type, cblksty);
 705             break;
 706         case 2:
 707             t1_3d_enc_clnpass(t1, w, h, l, bpno, orient, &nmsedec, cblksty);
 708             /* code switch SEGMARK (i.e. SEGSYM) */
 709             if (cblksty & J3D_CCP_CBLKSTY_SEGSYM) {
 710                 mqc_segmark_enc(mqc);
 711             }
 712             break;
 713         }
 714
 715         /* fixed_quality */
 716         tmpwmsedec = t1_getwmsedec(nmsedec, compno, level, orient, bpno, stepsize,
 717                                    numcomps, dwtid);
 718         cumwmsedec += tmpwmsedec;
 719         tile->distotile += tmpwmsedec;
 720
 721         /* Code switch "RESTART" (i.e. TERMALL) */
 722         if ((cblksty & J3D_CCP_CBLKSTY_TERMALL) && !((passtype == 2) &&
 723                 (bpno - 1 < 0))) {
 724             if (type == T1_TYPE_RAW) {
 725                 mqc_flush(mqc);
 726                 correction = 1;
 727                 /* correction = mqc_bypass_flush_enc(); */
 728             } else {            /* correction = mqc_restart_enc(); */
 729                 mqc_flush(mqc);
 730                 correction = 1;
 731             }
 732             pass->term = 1;
 733         } else {
 734             if (((bpno < (cblk->numbps - 4) && (passtype > 0))
 735                     || ((bpno == (cblk->numbps - 4)) && (passtype == 2))) &&
 736                     (cblksty & J3D_CCP_CBLKSTY_LAZY)) {
 737                 if (type == T1_TYPE_RAW) {
 738                     mqc_flush(mqc);
 739                     correction = 1;
 740                 } else {
 741                     mqc_flush(mqc);
 742                     correction = 1;
 743                 }
 744                 pass->term = 1;
 745             } else {
 746                 pass->term = 0;
 747             }
 748         }
 749
 750         if (++passtype == 3) {
 751             passtype = 0;
 752             bpno--;
 753         }
 754
 755         if (pass->term && bpno > 0) {
 756             type = ((bpno < (cblk->numbps - 4)) && (passtype < 2) &&
 757                     (cblksty & J3D_CCP_CBLKSTY_LAZY)) ? T1_TYPE_RAW : T1_TYPE_MQ;
 758             if (type == T1_TYPE_RAW) {
 759                 mqc_bypass_init_enc(mqc);
 760             } else {
 761                 mqc_restart_init_enc(mqc);
 762             }
 763         }
 764
 765         pass->distortiondec = cumwmsedec;
 766         pass->rate = mqc_numbytes(mqc) + correction;    /* FIXME */
 767         pass->len = pass->rate - (passno == 0 ? 0 : cblk->passes[passno - 1].rate);
 768
 769         /* Code-switch "RESET" */
 770         if (cblksty & J3D_CCP_CBLKSTY_RESET) {
 771             mqc_reset_enc(mqc);
 772         }
 773     }
 774
 775     /* Code switch "ERTERM" (i.e. PTERM) */
 776     if (cblksty & J3D_CCP_CBLKSTY_PTERM) {
 777         mqc_erterm_enc(mqc);
 778     } else /* Default coding */ if (!(cblksty & J3D_CCP_CBLKSTY_LAZY)) {
 779         mqc_flush(mqc);
 780     }
 781
 782     cblk->totalpasses = passno;
 783 }
 784
 785 static void t1_3d_decode_cblk(opj_t1_3d_t *t1, opj_tcd_cblk_t * cblk,
 786                               int orient, int roishift, int cblksty)
 787 {
 788     int i, j, k;
 789     int w, h, l;
 790     int bpno, passtype;
 791     int segno, passno;
 792     char type = T1_TYPE_MQ; /* BYPASS mode */
 793     opj_raw_t *raw = t1->raw;   /* RAW component */
 794     opj_mqc_t *mqc = t1->mqc;   /* MQC component */
 795
 796     w = cblk->x1 - cblk->x0;
 797     h = cblk->y1 - cblk->y0;
 798     l = cblk->z1 - cblk->z0;
 799
 800     for (k = 0; k < l; k++) {
 801         for (j = 0; j < h; j++) {
 802             for (i = 0; i < w; i++) {
 803                 t1->data[k][j][i] = 0;
 804             }
 805         }
 806     }
 807
 808     for (k = 0; k <= l; k++) {
 809         for (j = 0; j <= h; j++) {
 810             for (i = 0; i <= w; i++) {
 811                 t1->flags[k][j][i] = 0;
 812                 t1->flagSVR[k][j][i] = 0;
 813             }
 814         }
 815     }
 816
 817
 818     bpno = roishift + cblk->numbps - 1;
 819     passtype = 2;
 820
 821     mqc_reset_enc(mqc);
 822
 823     for (segno = 0; segno < cblk->numsegs; segno++) {
 824         opj_tcd_seg_t *seg = &cblk->segs[segno];
 825
 826         /* BYPASS mode */
 827         type = ((bpno <= (cblk->numbps - 1) - 4) && (passtype < 2) &&
 828                 (cblksty & J3D_CCP_CBLKSTY_LAZY)) ? T1_TYPE_RAW : T1_TYPE_MQ;
 829         if (type == T1_TYPE_RAW) {
 830             raw_init_dec(raw, seg->data, seg->len);
 831         } else {
 832             mqc_init_dec(mqc, seg->data, seg->len);
 833         }
 834
 835         for (passno = 0; passno < seg->numpasses; passno++) {
 836             switch (passtype) {
 837             case 0:
 838                 t1_3d_dec_sigpass(t1, w, h, l, bpno + 1, orient, type, cblksty);
 839                 break;
 840             case 1:
 841                 t1_3d_dec_refpass(t1, w, h, l, bpno + 1, type, cblksty);
 842                 break;
 843             case 2:
 844                 t1_3d_dec_clnpass(t1, w, h, l, bpno + 1, orient, cblksty);
 845                 break;
 846             }
 847
 848             if ((cblksty & J3D_CCP_CBLKSTY_RESET) && type == T1_TYPE_MQ) {
 849                 mqc_reset_enc(mqc);
 850             }
 851             if (++passtype == 3) {
 852                 passtype = 0;
 853                 bpno--;
 854             }
 855
 856         }
 857     }
 858 }
 859
 860 static int t1_3d_init_ctxno_zc(unsigned int f, int orient)
 861 {
 862     unsigned int h, v, c;
 863     unsigned int d2xy, d2xz, d2yz, d3;
 864     int n;
 865     unsigned int hvc, hc, d2, d2xy2yz, d2xy2xz;
 866     n = 0;
 867     h = ((f & T1_3D_SIG_W) != 0) + ((f & T1_3D_SIG_E) != 0);
 868     v = ((f & T1_3D_SIG_N) != 0) + ((f & T1_3D_SIG_S) != 0);
 869     c = ((f & T1_3D_SIG_FC) != 0) + ((f & T1_3D_SIG_BC) != 0);
 870     d2xy = ((f & T1_3D_SIG_NW) != 0) + ((f & T1_3D_SIG_NE) != 0) + ((
 871                 f & T1_3D_SIG_SE) != 0) + ((f & T1_3D_SIG_SW) != 0);
 872     d2xz = ((f & T1_3D_SIG_FW) != 0) + ((f & T1_3D_SIG_BW) != 0) + ((
 873                 f & T1_3D_SIG_FE) != 0) + ((f & T1_3D_SIG_BE) != 0);
 874     d2yz = ((f & T1_3D_SIG_FN) != 0) + ((f & T1_3D_SIG_FS) != 0) + ((
 875                 f & T1_3D_SIG_BN) != 0) + ((f & T1_3D_SIG_BS) != 0);
 876     d3 = ((f & T1_3D_SIG_FNW) != 0) + ((f & T1_3D_SIG_FNE) != 0) + ((
 877                 f & T1_3D_SIG_FSE) != 0) + ((f & T1_3D_SIG_FSW) != 0)
 878          + ((f & T1_3D_SIG_BNW) != 0) + ((f & T1_3D_SIG_BNE) != 0) + ((
 879                      f & T1_3D_SIG_BSE) != 0) + ((f & T1_3D_SIG_BSW) != 0);
 880
 881     switch (orient) {
 882     case 0: /*LLL*/
 883     case 7: /*HHH*/
 884         hvc = h + v + c;
 885         d2 = d2xy + d2xz + d2yz;
 886         if (!hvc) {
 887             if (!d2) {
 888                 n = (!d3) ? 0 : 1;
 889             } else if (d2 == 1) {
 890                 n = (!d3) ? 2 : 3;
 891             } else {
 892                 n = (!d3) ? 4 : 5;
 893             }
 894         } else if (hvc == 1) {
 895             if (!d2) {
 896                 n = (!d3) ? 6 : 7;
 897             } else if (d2 == 1) {
 898                 n = (!d3) ? 8 : 9;
 899             } else {
 900                 n = 10;
 901             }
 902         } else if (hvc == 2) {
 903             if (!d2) {
 904                 n = (!d3) ? 11 : 12;
 905             } else {
 906                 n = 13;
 907             }
 908         } else if (hvc == 3) {
 909             n = 14;
 910         } else {
 911             n = 15;
 912         }
 913         break;
 914     /*LHL, HLL, LLH*/
 915     case 1:
 916     case 2:
 917     case 4:
 918         hc = h + c;
 919         d2xy2yz = d2xy + d2yz;
 920         if (!hc) {
 921             if (!v) {
 922                 if (!d2xy) {
 923                     n = (!d2xy2yz) ? ((!d3) ? 0 : 1) : ((!d3) ? 2 : 3);
 924                 } else if (d2xy == 1) {
 925                     n = (!d2xy2yz) ? ((!d3) ? 4 : 5) : 6;
 926                 } else { /*>=2*/
 927                     n = 7;
 928                 }
 929             } else {
 930                 n = (v == 1) ? 8 : 9; /* =1 or =2*/
 931             }
 932         } else if (hc == 1) {
 933             n = (!v) ? ((!d2xy) ? ((!d2xy2yz) ? ((!d3) ? 10 : 11) : (12)) : (13)) : (14);
 934         } else { /*if (hc >= 2)*/
 935             n = 15;
 936         }
 937         break;
 938     /*HLH, HHL, LHH*/
 939     case 3:
 940     case 5:
 941     case 6:
 942         hc = h + c;
 943         d2xy2xz = d2xy + d2xz;
 944         if (!v) {
 945             if (!d2xz) {
 946                 if (!hc && !d2xy2xz) {
 947                     n = (!d3) ? 0 : 1;
 948                 } else if (hc == 1) {
 949                     n = (!d2xy2xz) ?  2 : 3;
 950                 } else { /*if >= 2*/
 951                     n = 4;
 952                 }
 953             } else if (d2xz >= 1 && !hc) {
 954                 n = 5;
 955             } else if (hc >= 1) {
 956                 n = (d2xz == 1) ? 6 : 7;
 957             }
 958         } else if (v == 1) {
 959             if (!d2xz) {
 960                 n = (!hc) ? 8 : 9;
 961             } else if (d2xz == 1) {
 962                 n = (!hc) ? 10 : 11;
 963             } else if (d2xz == 2) {
 964                 n = (!hc) ? 12 : 13;
 965             } else { /* if (d2xz >= 3) {*/
 966                 n = 14;
 967             }
 968         } else if (v == 2) {
 969             n = 15;
 970         }
 971         break;
 972     }
 973
 974     return (T1_3D_CTXNO_ZC + n);
 975 }
 976
 977 static int t1_3d_init_ctxno_sc(unsigned int f)
 978 {
 979     int hc, vc, cc;
 980     int n = 0;
 981
 982     hc = int_min(((f & (T1_3D_SIG_E | T1_3D_SGN_E)) == T1_3D_SIG_E)
 983                  + ((f & (T1_3D_SIG_W | T1_3D_SGN_W)) == T1_3D_SIG_W), 1)
 984          - int_min(((f & (T1_3D_SIG_E | T1_3D_SGN_E)) == (T1_3D_SIG_E | T1_3D_SGN_E))
 985                    + ((f & (T1_3D_SIG_W | T1_3D_SGN_W)) == (T1_3D_SIG_W | T1_3D_SGN_W)), 1);
 986
 987     vc = int_min(((f & (T1_3D_SIG_N | T1_3D_SGN_N)) == T1_3D_SIG_N)
 988                  + ((f & (T1_3D_SIG_S | T1_3D_SGN_S)) == T1_3D_SIG_S), 1)
 989          - int_min(((f & (T1_3D_SIG_N | T1_3D_SGN_N)) == (T1_3D_SIG_N | T1_3D_SGN_N))
 990                    + ((f & (T1_3D_SIG_S | T1_3D_SGN_S)) == (T1_3D_SIG_S | T1_3D_SGN_S)), 1);
 991
 992     cc = int_min(((f & (T1_3D_SIG_FC | T1_3D_SGN_F)) == T1_3D_SIG_FC)
 993                  + ((f & (T1_3D_SIG_BC | T1_3D_SGN_B)) == T1_3D_SIG_BC), 1)
 994          - int_min(((f & (T1_3D_SIG_FC | T1_3D_SGN_F)) == (T1_3D_SIG_FC | T1_3D_SGN_F))
 995                    + ((f & (T1_3D_SIG_BC | T1_3D_SGN_B)) == (T1_3D_SIG_BC | T1_3D_SGN_B)), 1);
 996     if (hc < 0) {
 997         hc = -hc;
 998         vc = -vc;
 999         cc = -cc;
1000     }
1001
1002     if (!hc) {
1003         if (!vc) {
1004             n = (!cc) ? 0 : 1;
1005         } else if (vc == -1) {
1006             n = (!cc) ? 1 : ((cc > 0) ? 2 : 4);
1007         } else if (vc == 1) {
1008             n = (!cc) ? 1 : ((cc < 0) ? 2 : 4);
1009         }
1010     } else if (hc == 1) {
1011         if (!vc) {
1012             n = (!cc) ? 1 : ((cc < 0) ? 2 : 4);
1013         } else if (vc == 1) {
1014             n = (!cc) ? 4 : ((cc > 0) ? 5 : 3);
1015         } else if (vc == -1) {
1016             n = (!cc) ? 2 : 3;
1017         }
1018     } else if (hc == -1) {
1019         if (!vc) {
1020             n = (!cc) ? 1 : ((cc > 0) ? 2 : 4);
1021         } else if (vc == 1) {
1022             n = (!cc) ? 2 : 3;
1023         } else if (vc == -1) {
1024             n = (!cc) ? 4 : ((cc < 0) ? 5 : 3);
1025         }
1026     }
1027
1028     return (T1_3D_CTXNO_SC + n);
1029 }
1030
1031 static int t1_3d_init_ctxno_mag(unsigned int f, int f2)
1032 {
1033     int n;
1034     if (!(f2 &
1035             T1_3D_REFINE)) { /*First refinement for this coefficient (no previous refinement)*/
1036         n = (f & (T1_3D_SIG_PRIM)) ? 1 : 0;
1037     } else {
1038         n = 2;
1039     }
1040
1041     return (T1_3D_CTXNO_MAG + n);
1042 }
1043
1044 static int t1_3d_init_spb(unsigned int f)
1045 {
1046     int hc, vc, cc;
1047     int n = 0;
1048
1049     hc = int_min(((f & (T1_3D_SIG_E | T1_3D_SGN_E)) == T1_3D_SIG_E)
1050                  + ((f & (T1_3D_SIG_W | T1_3D_SGN_W)) == T1_3D_SIG_W), 1)
1051          - int_min(((f & (T1_3D_SIG_E | T1_3D_SGN_E)) == (T1_3D_SIG_E | T1_3D_SGN_E))
1052                    + ((f & (T1_3D_SIG_W | T1_3D_SGN_W)) == (T1_3D_SIG_W | T1_3D_SGN_W)), 1);
1053
1054     vc = int_min(((f & (T1_3D_SIG_N | T1_3D_SGN_N)) == T1_3D_SIG_N)
1055                  + ((f & (T1_3D_SIG_S | T1_3D_SGN_S)) == T1_3D_SIG_S), 1)
1056          - int_min(((f & (T1_3D_SIG_N | T1_3D_SGN_N)) == (T1_3D_SIG_N | T1_3D_SGN_N))
1057                    + ((f & (T1_3D_SIG_S | T1_3D_SGN_S)) == (T1_3D_SIG_S | T1_3D_SGN_S)), 1);
1058
1059     cc = int_min(((f & (T1_3D_SIG_FC | T1_3D_SGN_F)) == T1_3D_SIG_FC)
1060                  + ((f & (T1_3D_SIG_BC | T1_3D_SGN_B)) == T1_3D_SIG_BC), 1)
1061          - int_min(((f & (T1_3D_SIG_FC | T1_3D_SGN_F)) == (T1_3D_SIG_FC | T1_3D_SGN_F))
1062                    + ((f & (T1_3D_SIG_BC | T1_3D_SGN_B)) == (T1_3D_SIG_BC | T1_3D_SGN_B)), 1);
1063
1064     n = ((hc + vc + cc) < 0);
1065
1066     return n;
1067 }
1068
1069 static void t1_3d_init_luts(opj_t1_3d_t *t1)
1070 {
1071     int i;
1072     double u, v, t;
1073     /*for (j = 0; j < 4; j++) {
1074         for (i = 0; i < 256; ++i) {
1075             t1->lut_ctxno_zc[(j << 8) | i] = t1_3d_init_ctxno_zc(i, j);
1076         }
1077     }
1078     for (i = 0; i < 4096; i++) {
1079         t1->lut_ctxno_sc[i] = t1_3d_init_ctxno_sc(i << 4);
1080     }
1081     for (j = 0; j < 2; j++) {
1082         for (i = 0; i < 2048; ++i) {
1083             t1->lut_ctxno_mag[(j << 11) + i] = t1_3d_init_ctxno_mag((j ? T1_3D_REFINE : 0) | i);
1084         }
1085     }
1086     for (i = 0; i < 4096; ++i) {
1087         t1->lut_spb[i] = t1_3d_init_spb(i << 4);
1088     }*/
1089     /* FIXME FIXME FIXME */
1090     for (i = 0; i < (1 << T1_NMSEDEC_BITS); i++) {
1091         t = i / pow(2, T1_NMSEDEC_FRACBITS);
1092         u = t;
1093         v = t - 1.5;
1094         t1->lut_nmsedec_sig[i] =
1095             int_max(0,
1096                     (int)(floor((u * u - v * v) * pow(2, T1_NMSEDEC_FRACBITS) + 0.5) / pow(2,
1097                             T1_NMSEDEC_FRACBITS) * 8192.0));
1098         t1->lut_nmsedec_sig0[i] =
1099             int_max(0,
1100                     (int)(floor((u * u) * pow(2, T1_NMSEDEC_FRACBITS) + 0.5) / pow(2,
1101                             T1_NMSEDEC_FRACBITS) * 8192.0));
1102         u = t - 1.0;
1103         if (i & (1 << (T1_NMSEDEC_BITS - 1))) {
1104             v = t - 1.5;
1105         } else {
1106             v = t - 0.5;
1107         }
1108         t1->lut_nmsedec_ref[i] =
1109             int_max(0,
1110                     (int)(floor((u * u - v * v) * pow(2, T1_NMSEDEC_FRACBITS) + 0.5) / pow(2,
1111                             T1_NMSEDEC_FRACBITS) * 8192.0));
1112         t1->lut_nmsedec_ref0[i] =
1113             int_max(0,
1114                     (int)(floor((u * u) * pow(2, T1_NMSEDEC_FRACBITS) + 0.5) / pow(2,
1115                             T1_NMSEDEC_FRACBITS) * 8192.0));
1116     }
1117 }
1118
1119 /* ----------------------------------------------------------------------- */
1120
1121 opj_t1_3d_t* t1_3d_create(opj_common_ptr cinfo)
1122 {
1123     opj_t1_3d_t *t1 = (opj_t1_3d_t*)opj_malloc(sizeof(opj_t1_3d_t));
1124     if (t1) {
1125         t1->cinfo = cinfo;
1126         /* create MQC and RAW handles */
1127         t1->mqc = mqc_create();
1128         t1->raw = raw_create();
1129         /* initialize the look-up tables of the Tier-1 coder/decoder */
1130         t1_3d_init_luts(t1);
1131     }
1132     return t1;
1133 }
1134
1135 void t1_3d_destroy(opj_t1_3d_t *t1)
1136 {
1137     if (t1) {
1138         /* destroy MQC and RAW handles */
1139         mqc_destroy(t1->mqc);
1140         raw_destroy(t1->raw);
1141         opj_free(t1);
1142     }
1143 }
1144
1145 void t1_3d_encode_cblks(opj_t1_3d_t *t1, opj_tcd_tile_t *tile, opj_tcp_t *tcp)
1146 {
1147     int compno, resno, bandno, precno, cblkno;
1148     int x, y, z, i, j, k, orient;
1149     int level[3];
1150     tile->distotile = 0;        /* fixed_quality */
1151
1152     for (compno = 0; compno < tile->numcomps; compno++) {
1153         opj_tcd_tilecomp_t *tilec = &tile->comps[compno];
1154
1155         for (resno = 0; resno < tilec->numresolution[0]; resno++) {
1156             opj_tcd_resolution_t *res = &tilec->resolutions[resno];
1157
1158             for (bandno = 0; bandno < res->numbands; bandno++) {
1159                 opj_tcd_band_t *band = &res->bands[bandno];
1160
1161                 for (precno = 0; precno < res->prctno[0] * res->prctno[1] * res->prctno[2];
1162                         precno++) {
1163                     opj_tcd_precinct_t *prc = &band->precincts[precno];
1164
1165                     for (cblkno = 0; cblkno < prc->cblkno[0] * prc->cblkno[1] * prc->cblkno[2];
1166                             cblkno++) {
1167                         opj_tcd_cblk_t *cblk = &prc->cblks[cblkno];
1168
1169                         if (band->bandno == 0) {
1170                             x = cblk->x0 - band->x0;
1171                             y = cblk->y0 - band->y0;
1172                             z = cblk->z0 - band->z0;
1173                         } else if (band->bandno == 1) {
1174                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1175                             x = pres->x1 - pres->x0 + cblk->x0 - band->x0;
1176                             y = cblk->y0 - band->y0;
1177                             z = cblk->z0 - band->z0;
1178                         } else if (band->bandno == 2) {
1179                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1180                             x = cblk->x0 - band->x0;
1181                             y = pres->y1 - pres->y0 + cblk->y0 - band->y0;
1182                             z = cblk->z0 - band->z0;
1183                         } else if (band->bandno == 3) {
1184                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1185                             x = pres->x1 - pres->x0 + cblk->x0 - band->x0;
1186                             y = pres->y1 - pres->y0 + cblk->y0 - band->y0;
1187                             z = cblk->z0 - band->z0;
1188                         } else if (band->bandno == 4) {
1189                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1190                             x = cblk->x0 - band->x0;
1191                             y = cblk->y0 - band->y0;
1192                             z = pres->z1 - pres->z0 + cblk->z0 - band->z0;
1193                         } else if (band->bandno == 5) {
1194                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1195                             x = pres->x1 - pres->x0 + cblk->x0 - band->x0;
1196                             y = cblk->y0 - band->y0;
1197                             z = pres->z1 - pres->z0 + cblk->z0 - band->z0;
1198                         } else if (band->bandno == 6) {
1199                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1200                             x = cblk->x0 - band->x0;
1201                             y = pres->y1 - pres->y0 + cblk->y0 - band->y0;
1202                             z = pres->z1 - pres->z0 + cblk->z0 - band->z0;
1203                         } else if (band->bandno == 7) {
1204                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1205                             x = pres->x1 - pres->x0 + cblk->x0 - band->x0;
1206                             y = pres->y1 - pres->y0 + cblk->y0 - band->y0;
1207                             z = pres->z1 - pres->z0 + cblk->z0 - band->z0;
1208                         }
1209
1210                         if (tcp->tccps[compno].reversible == 1) {
1211                             for (k = 0; k < cblk->z1 - cblk->z0; k++) {
1212                                 for (j = 0; j < cblk->y1 - cblk->y0; j++) {
1213                                     for (i = 0; i < cblk->x1 - cblk->x0; i++) {
1214                                         t1->data[k][j][i] =
1215                                             tilec->data[(x + i) + (y + j) * (tilec->x1 - tilec->x0) + (z + k) *
1216                                                                 (tilec->x1 - tilec->x0) * (tilec->y1 - tilec->y0)] << T1_NMSEDEC_FRACBITS;
1217                                     }
1218                                 }
1219                             }
1220                         } else if (tcp->tccps[compno].reversible == 0) {
1221                             for (k = 0; k < cblk->z1 - cblk->z0; k++) {
1222                                 for (j = 0; j < cblk->y1 - cblk->y0; j++) {
1223                                     for (i = 0; i < cblk->x1 - cblk->x0; i++) {
1224                                         t1->data[k][j][i] = fix_mul(
1225                                                                 tilec->data[(x + i) + (y + j) * (tilec->x1 - tilec->x0) + (z + k) *
1226                                                                                     (tilec->x1 - tilec->x0) * (tilec->y1 - tilec->y0)],
1227                                                                 8192 * 8192 / ((int) floor(band->stepsize * 8192))) >> (13 -
1228                                                                         T1_NMSEDEC_FRACBITS);
1229                                     }
1230                                 }
1231                             }
1232                         }
1233                         orient = band->bandno;  /* FIXME */
1234                         for (i = 0; i < 3; i++) {
1235                             level[i] = tilec->numresolution[i] - 1 - resno;
1236                         }
1237
1238                         t1_3d_encode_cblk(t1, cblk, orient, compno, level, tcp->tccps[compno].dwtid,
1239                                           band->stepsize, tcp->tccps[compno].cblksty, tile->numcomps, tile);
1240
1241                     } /* cblkno */
1242                 } /* precno */
1243             } /* bandno */
1244         } /* resno  */
1245     } /* compno  */
1246 }
1247
1248 void t1_3d_decode_cblks(opj_t1_3d_t *t1, opj_tcd_tile_t *tile, opj_tcp_t *tcp)
1249 {
1250     int compno, resno, bandno, precno, cblkno;
1251
1252     for (compno = 0; compno < tile->numcomps; compno++) {
1253         opj_tcd_tilecomp_t *tilec = &tile->comps[compno];
1254
1255         for (resno = 0; resno < tilec->numresolution[0]; resno++) {
1256             opj_tcd_resolution_t *res = &tilec->resolutions[resno];
1257
1258             for (bandno = 0; bandno < res->numbands; bandno++) {
1259                 opj_tcd_band_t *band = &res->bands[bandno];
1260
1261                 for (precno = 0; precno < res->prctno[0] * res->prctno[1] * res->prctno[2];
1262                         precno++) {
1263                     opj_tcd_precinct_t *prc = &band->precincts[precno];
1264
1265                     for (cblkno = 0; cblkno < prc->cblkno[0] * prc->cblkno[1] * prc->cblkno[2];
1266                             cblkno++) {
1267                         int x, y, z, i, j, k, orient;
1268                         opj_tcd_cblk_t *cblk = &prc->cblks[cblkno];
1269
1270                         orient = band->bandno;  /* FIXME */
1271
1272                         /*fprintf(stdout,"[INFO] t1_3d_decode_cblk(t1, cblk, orient(%d), tcp->tccps[compno].roishift (%d), tcp->tccps[compno].cblksty (%d));\n",orient,tcp->tccps[compno].roishift, tcp->tccps[compno].cblksty);*/
1273                         t1_3d_decode_cblk(t1, cblk, orient, tcp->tccps[compno].roishift,
1274                                           tcp->tccps[compno].cblksty);
1275
1276                         if (band->bandno == 0) {
1277                             x = cblk->x0 - band->x0;
1278                             y = cblk->y0 - band->y0;
1279                             z = cblk->z0 - band->z0;
1280                         } else if (band->bandno == 1) {
1281                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1282                             x = pres->x1 - pres->x0 + cblk->x0 - band->x0;
1283                             y = cblk->y0 - band->y0;
1284                             z = cblk->z0 - band->z0;
1285                         } else if (band->bandno == 2) {
1286                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1287                             x = cblk->x0 - band->x0;
1288                             y = pres->y1 - pres->y0 + cblk->y0 - band->y0;
1289                             z = cblk->z0 - band->z0;
1290                         } else if (band->bandno == 3) {
1291                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1292                             x = pres->x1 - pres->x0 + cblk->x0 - band->x0;
1293                             y = pres->y1 - pres->y0 + cblk->y0 - band->y0;
1294                             z = cblk->z0 - band->z0;
1295                         } else if (band->bandno == 4) {
1296                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1297                             x = cblk->x0 - band->x0;
1298                             y = cblk->y0 - band->y0;
1299                             z = pres->z1 - pres->z0 + cblk->z0 - band->z0;
1300                         } else if (band->bandno == 5) {
1301                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1302                             x = pres->x1 - pres->x0 + cblk->x0 - band->x0;
1303                             y = cblk->y0 - band->y0;
1304                             z = pres->z1 - pres->z0 + cblk->z0 - band->z0;
1305                         } else if (band->bandno == 6) {
1306                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1307                             x = cblk->x0 - band->x0;
1308                             y = pres->y1 - pres->y0 + cblk->y0 - band->y0;
1309                             z = pres->z1 - pres->z0 + cblk->z0 - band->z0;
1310                         } else if (band->bandno == 7) {
1311                             opj_tcd_resolution_t *pres = &tilec->resolutions[resno - 1];
1312                             x = pres->x1 - pres->x0 + cblk->x0 - band->x0;
1313                             y = pres->y1 - pres->y0 + cblk->y0 - band->y0;
1314                             z = pres->z1 - pres->z0 + cblk->z0 - band->z0;
1315                         }
1316
1317                         if (tcp->tccps[compno].roishift) {
1318                             int thresh, val, mag;
1319                             thresh = 1 << tcp->tccps[compno].roishift;
1320                             for (k = 0; k < cblk->z1 - cblk->z0; k++) {
1321                                 for (j = 0; j < cblk->y1 - cblk->y0; j++) {
1322                                     for (i = 0; i < cblk->x1 - cblk->x0; i++) {
1323                                         val = t1->data[k][j][i];
1324                                         mag = int_abs(val);
1325                                         if (mag >= thresh) {
1326                                             mag >>= tcp->tccps[compno].roishift;
1327                                             t1->data[k][j][i] = val < 0 ? -mag : mag;
1328                                         }
1329                                     }
1330                                 }
1331                             }
1332                         }
1333
1334                         if (tcp->tccps[compno].reversible == 1) {
1335                             for (k = 0; k < cblk->z1 - cblk->z0; k++) {
1336                                 for (j = 0; j < cblk->y1 - cblk->y0; j++) {
1337                                     for (i = 0; i < cblk->x1 - cblk->x0; i++) {
1338                                         int tmp = t1->data[k][j][i];
1339                                         tilec->data[x + i + (y + j) * (tilec->x1 - tilec->x0) + (z + k) *
1340                                                       (tilec->x1 - tilec->x0) * (tilec->y1 - tilec->y0)] = tmp / 2;
1341                                     }
1342                                 }
1343                             }
1344                         } else {        /* if (tcp->tccps[compno].reversible == 0) */
1345                             for (k = 0; k < cblk->z1 - cblk->z0; k++) {
1346                                 for (j = 0; j < cblk->y1 - cblk->y0; j++) {
1347                                     for (i = 0; i < cblk->x1 - cblk->x0; i++) {
1348                                         double tmp = (double)(t1->data[k][j][i] * band->stepsize * 4096.0);
1349                                         if (t1->data[k][j][i] >> 1 == 0) {
1350                                             tilec->data[x + i + (y + j) * (tilec->x1 - tilec->x0) + (z + k) *
1351                                                           (tilec->x1 - tilec->x0) * (tilec->y1 - tilec->y0)] = 0;
1352                                         } else {
1353                                             int tmp2 = ((int)(floor(fabs(tmp)))) + ((int) floor(fabs(tmp * 2)) % 2);
1354                                             tilec->data[x + i + (y + j) * (tilec->x1 - tilec->x0) + (z + k) *
1355                                                           (tilec->x1 - tilec->x0) * (tilec->y1 - tilec->y0)] = ((tmp < 0) ? -tmp2 : tmp2);
1356                                         }
1357                                     }
1358                                 }
1359                             }
1360                         }
1361                     } /* cblkno */
1362                 } /* precno */
1363             } /* bandno */
1364         } /* resno */
1365     } /* compno */
1366 }