libopenjpeg/tcd.c

   1 /*
   2  * Copyright (c) 2001-2002, David Janssens
   3  * Copyright (c) 2002-2004, Yannick Verschueren
   4  * Copyright (c) 2002-2004, Communications and remote sensing Laboratory, Universite catholique de Louvain, Belgium
   5  * All rights reserved.
   6  *
   7  * Redistribution and use in source and binary forms, with or without
   8  * modification, are permitted provided that the following conditions
   9  * are met:
  10  * 1. Redistributions of source code must retain the above copyright
  11  *    notice, this list of conditions and the following disclaimer.
  12  * 2. Redistributions in binary form must reproduce the above copyright
  13  *    notice, this list of conditions and the following disclaimer in the
  14  *    documentation and/or other materials provided with the distribution.
  15  *
  16  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS `AS IS'
  17  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  18  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
  19  * ARE DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
  20  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
  21  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
  22  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
  23  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
  24  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
  25  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
  26  * POSSIBILITY OF SUCH DAMAGE.
  27  */
  28
  29 #include "tcd.h"
  30 #include "int.h"
  31 #include "t1.h"
  32 #include "t2.h"
  33 #include "dwt.h"
  34 #include "mct.h"
  35 #include <setjmp.h>
  36 #include <float.h>
  37 #include <stdio.h>
  38 #include <time.h>
  39 #include <math.h>
  40 #include <stdlib.h>
  41 #include <string.h>
  42
  43 static tcd_image_t tcd_image;
  44
  45 static j2k_image_t *tcd_img;
  46 static j2k_cp_t *tcd_cp;
  47
  48 static tcd_tile_t *tcd_tile;
  49 static j2k_tcp_t *tcd_tcp;
  50 static int tcd_tileno;
  51
  52 static tcd_tile_t *tile;
  53 static tcd_tilecomp_t *tilec;
  54 static tcd_resolution_t *res;
  55 static tcd_band_t *band;
  56 static tcd_precinct_t *prc;
  57 static tcd_cblk_t *cblk;
  58
  59 extern jmp_buf j2k_error;
  60
  61 void tcd_dump(tcd_image_t * img, int curtileno)
  62 {
  63   int tileno, compno, resno, bandno, precno, cblkno;
  64   fprintf(stdout, "image {\n");
  65   fprintf(stdout, "  tw=%d, th=%d x0=%d x1=%d y0=%d y1=%d\n", img->tw,
  66           img->th, tcd_img->x0, tcd_img->x1, tcd_img->y0, tcd_img->y1);
  67   for (tileno = 0; tileno < img->th * img->tw; tileno++) {
  68     tcd_tile_t *tile = &tcd_image.tiles[tileno];
  69     fprintf(stdout, "  tile {\n");
  70     fprintf(stdout, "    x0=%d, y0=%d, x1=%d, y1=%d, numcomps=%d\n",
  71             tile->x0, tile->y0, tile->x1, tile->y1, tile->numcomps);
  72     for (compno = 0; compno < tile->numcomps; compno++) {
  73       tcd_tilecomp_t *tilec = &tile->comps[compno];
  74       fprintf(stdout, "    tilec {\n");
  75       fprintf(stdout,
  76               "      x0=%d, y0=%d, x1=%d, y1=%d, numresolutions=%d\n",
  77               tilec->x0, tilec->y0, tilec->x1, tilec->y1,
  78               tilec->numresolutions);
  79       for (resno = 0; resno < tilec->numresolutions; resno++) {
  80         tcd_resolution_t *res = &tilec->resolutions[resno];
  81         fprintf(stdout, "\n   res {\n");
  82         fprintf(stdout,
  83                 "          x0=%d, y0=%d, x1=%d, y1=%d, pw=%d, ph=%d, numbands=%d\n",
  84                 res->x0, res->y0, res->x1, res->y1, res->pw, res->ph,
  85                 res->numbands);
  86         for (bandno = 0; bandno < res->numbands; bandno++) {
  87           tcd_band_t *band = &res->bands[bandno];
  88           fprintf(stdout, "        band {\n");
  89           fprintf(stdout,
  90                   "          x0=%d, y0=%d, x1=%d, y1=%d, stepsize=%d, numbps=%d\n",
  91                   band->x0, band->y0, band->x1, band->y1,
  92                   band->stepsize, band->numbps);
  93           for (precno = 0; precno < res->pw * res->ph; precno++) {
  94             tcd_precinct_t *prec = &band->precincts[precno];
  95             fprintf(stdout, "          prec {\n");
  96             fprintf(stdout,
  97                     "            x0=%d, y0=%d, x1=%d, y1=%d, cw=%d, ch=%d\n",
  98                     prec->x0, prec->y0, prec->x1, prec->y1,
  99                     prec->cw, prec->ch);
 100             for (cblkno = 0; cblkno < prec->cw * prec->ch; cblkno++) {
 101               tcd_cblk_t *cblk = &prec->cblks[cblkno];
 102               fprintf(stdout, "            cblk {\n");
 103               fprintf(stdout,
 104                       "              x0=%d, y0=%d, x1=%d, y1=%d\n",
 105                       cblk->x0, cblk->y0, cblk->x1, cblk->y1);
 106               fprintf(stdout, "            }\n");
 107             }
 108             fprintf(stdout, "          }\n");
 109           }
 110           fprintf(stdout, "        }\n");
 111         }
 112         fprintf(stdout, "      }\n");
 113       }
 114       fprintf(stdout, "    }\n");
 115     }
 116     fprintf(stdout, "  }\n");
 117   }
 118   fprintf(stdout, "}\n");
 119 }
 120
 121 void tcd_malloc_encode(j2k_image_t * img, j2k_cp_t * cp, int curtileno)
 122 {
 123   int tileno, compno, resno, bandno, precno, cblkno;
 124   tcd_img = img;
 125   tcd_cp = cp;
 126   tcd_image.tw = cp->tw;
 127   tcd_image.th = cp->th;
 128   tcd_image.tiles = (tcd_tile_t *) malloc(sizeof(tcd_tile_t));
 129
 130   for (tileno = 0; tileno < 1; tileno++) {
 131     j2k_tcp_t *tcp = &cp->tcps[curtileno];
 132     int j;
 133     /* cfr p59 ISO/IEC FDIS15444-1 : 2000 (18 august 2000) */
 134     int p = curtileno % cp->tw; /* si numerotation matricielle .. */
 135     int q = curtileno / cp->tw; /* .. coordonnees de la tile (q,p) q pour ligne et p pour colonne */
 136     /* tcd_tile_t *tile=&tcd_image.tiles[tileno]; */
 137     tile = tcd_image.tiles;
 138     /* 4 borders of the tile rescale on the image if necessary */
 139     tile->x0 = int_max(cp->tx0 + p * cp->tdx, img->x0);
 140     tile->y0 = int_max(cp->ty0 + q * cp->tdy, img->y0);
 141     tile->x1 = int_min(cp->tx0 + (p + 1) * cp->tdx, img->x1);
 142     tile->y1 = int_min(cp->ty0 + (q + 1) * cp->tdy, img->y1);
 143     tile->numcomps = img->numcomps;
 144     /* tile->PPT=img->PPT;  */
 145     /* Modification of the RATE >> */
 146     for (j = 0; j < tcp->numlayers; j++) {
 147       tcp->rates[j] = tcp->rates[j] ? int_ceildiv(tile->numcomps * (tile->x1 - tile->x0) * (tile->y1 - tile->y0) * img->comps[0].prec, (tcp->rates[j] * 8 * img->comps[0].dx * img->comps[0].dy)) : 0;  //Mod antonin losslessbug
 148       if (tcp->rates[j]) {
 149         if (j && tcp->rates[j] < tcp->rates[j - 1] + 10) {
 150           tcp->rates[j] = tcp->rates[j - 1] + 20;
 151         } else {
 152           if (!j && tcp->rates[j] < 30)
 153             tcp->rates[j] = 30;
 154         }
 155       }
 156     }
 157     /* << Modification of the RATE */
 158
 159     tile->comps =
 160       (tcd_tilecomp_t *) malloc(img->numcomps * sizeof(tcd_tilecomp_t));
 161     for (compno = 0; compno < tile->numcomps; compno++) {
 162       j2k_tccp_t *tccp = &tcp->tccps[compno];
 163       /* tcd_tilecomp_t *tilec=&tile->comps[compno]; */
 164       tilec = &tile->comps[compno];
 165       /* border of each tile component (global) */
 166       tilec->x0 = int_ceildiv(tile->x0, img->comps[compno].dx);
 167
 168       tilec->y0 = int_ceildiv(tile->y0, img->comps[compno].dy);
 169       tilec->x1 = int_ceildiv(tile->x1, img->comps[compno].dx);
 170       tilec->y1 = int_ceildiv(tile->y1, img->comps[compno].dy);
 171
 172       tilec->data =
 173         (int *) malloc((tilec->x1 - tilec->x0) *
 174                        (tilec->y1 - tilec->y0) * sizeof(int));
 175       tilec->numresolutions = tccp->numresolutions;
 176
 177       tilec->resolutions =
 178         (tcd_resolution_t *) malloc(tilec->numresolutions *
 179                                     sizeof(tcd_resolution_t));
 180
 181       for (resno = 0; resno < tilec->numresolutions; resno++) {
 182         int pdx, pdy;
 183         int levelno = tilec->numresolutions - 1 - resno;
 184         int tlprcxstart, tlprcystart, brprcxend, brprcyend;
 185         int tlcbgxstart, tlcbgystart, brcbgxend, brcbgyend;
 186         int cbgwidthexpn, cbgheightexpn;
 187         int cblkwidthexpn, cblkheightexpn;
 188         /* tcd_resolution_t *res=&tilec->resolutions[resno]; */
 189
 190         res = &tilec->resolutions[resno];
 191
 192         /* border for each resolution level (global) */
 193         res->x0 = int_ceildivpow2(tilec->x0, levelno);
 194         res->y0 = int_ceildivpow2(tilec->y0, levelno);
 195         res->x1 = int_ceildivpow2(tilec->x1, levelno);
 196         res->y1 = int_ceildivpow2(tilec->y1, levelno);
 197
 198         res->numbands = resno == 0 ? 1 : 3;
 199         /* p. 35, table A-23, ISO/IEC FDIS154444-1 : 2000 (18 august 2000) */
 200         if (tccp->csty & J2K_CCP_CSTY_PRT) {
 201           pdx = tccp->prcw[resno];
 202           pdy = tccp->prch[resno];
 203         } else {
 204           pdx = 15;
 205           pdy = 15;
 206         }
 207         /* p. 64, B.6, ISO/IEC FDIS15444-1 : 2000 (18 august 2000)  */
 208         tlprcxstart = int_floordivpow2(res->x0, pdx) << pdx;
 209         tlprcystart = int_floordivpow2(res->y0, pdy) << pdy;
 210         brprcxend = int_ceildivpow2(res->x1, pdx) << pdx;
 211         brprcyend = int_ceildivpow2(res->y1, pdy) << pdy;
 212
 213         res->pw = (brprcxend - tlprcxstart) >> pdx;
 214         res->ph = (brprcyend - tlprcystart) >> pdy;
 215
 216         if (resno == 0) {
 217           tlcbgxstart = tlprcxstart;
 218           tlcbgystart = tlprcystart;
 219           brcbgxend = brprcxend;
 220           brcbgyend = brprcyend;
 221           cbgwidthexpn = pdx;
 222           cbgheightexpn = pdy;
 223         } else {
 224           tlcbgxstart = int_ceildivpow2(tlprcxstart, 1);
 225           tlcbgystart = int_ceildivpow2(tlprcystart, 1);
 226           brcbgxend = int_ceildivpow2(brprcxend, 1);
 227           brcbgyend = int_ceildivpow2(brprcyend, 1);
 228           cbgwidthexpn = pdx - 1;
 229           cbgheightexpn = pdy - 1;
 230         }
 231
 232         cblkwidthexpn = int_min(tccp->cblkw, cbgwidthexpn);
 233         cblkheightexpn = int_min(tccp->cblkh, cbgheightexpn);
 234
 235         for (bandno = 0; bandno < res->numbands; bandno++) {
 236           int x0b, y0b, i;
 237           int gain, numbps;
 238           j2k_stepsize_t *ss;
 239           band = &res->bands[bandno];
 240           band->bandno = resno == 0 ? 0 : bandno + 1;
 241           x0b = (band->bandno == 1) || (band->bandno == 3) ? 1 : 0;
 242           y0b = (band->bandno == 2) || (band->bandno == 3) ? 1 : 0;
 243
 244           if (band->bandno == 0) {
 245             /* band border (global) */
 246             band->x0 = int_ceildivpow2(tilec->x0, levelno);
 247             band->y0 = int_ceildivpow2(tilec->y0, levelno);
 248             band->x1 = int_ceildivpow2(tilec->x1, levelno);
 249             band->y1 = int_ceildivpow2(tilec->y1, levelno);
 250           } else {
 251             /* band border (global) */
 252             band->x0 =
 253               int_ceildivpow2(tilec->x0 -
 254                               (1 << levelno) * x0b, levelno + 1);
 255             band->y0 =
 256               int_ceildivpow2(tilec->y0 -
 257                               (1 << levelno) * y0b, levelno + 1);
 258             band->x1 =
 259               int_ceildivpow2(tilec->x1 -
 260                               (1 << levelno) * x0b, levelno + 1);
 261             band->y1 =
 262               int_ceildivpow2(tilec->y1 -
 263                               (1 << levelno) * y0b, levelno + 1);
 264
 265           }
 266
 267           ss = &tccp->stepsizes[resno ==
 268                                 0 ? 0 : 3 * (resno - 1) + bandno + 1];
 269           gain =
 270             tccp->qmfbid ==
 271             0 ? dwt_getgain_real(band->bandno) : dwt_getgain(band->bandno);
 272           numbps = img->comps[compno].prec + gain;
 273           band->stepsize =
 274             (int) floor((1.0 + ss->mant / 2048.0) *
 275                         pow(2.0, numbps - ss->expn) * 8192.0);
 276           band->numbps = ss->expn + tccp->numgbits - 1; /* WHY -1 ? */
 277
 278           band->precincts =
 279             (tcd_precinct_t *) malloc(3 * res->pw * res->ph *
 280                                       sizeof(tcd_precinct_t));
 281
 282           for (i = 0; i < res->pw * res->ph * 3; i++) {
 283             band->precincts[i].imsbtree = NULL;
 284             band->precincts[i].incltree = NULL;
 285           }
 286
 287           for (precno = 0; precno < res->pw * res->ph; precno++) {
 288             int tlcblkxstart, tlcblkystart, brcblkxend, brcblkyend;
 289             int cbgxstart =
 290               tlcbgxstart + (precno % res->pw) * (1 << cbgwidthexpn);
 291             int cbgystart =
 292               tlcbgystart + (precno / res->pw) * (1 << cbgheightexpn);
 293             int cbgxend = cbgxstart + (1 << cbgwidthexpn);
 294             int cbgyend = cbgystart + (1 << cbgheightexpn);
 295             /* tcd_precinct_t *prc=&band->precincts[precno]; */
 296             prc = &band->precincts[precno];
 297             /* precinct size (global) */
 298             prc->x0 = int_max(cbgxstart, band->x0);
 299             prc->y0 = int_max(cbgystart, band->y0);
 300             prc->x1 = int_min(cbgxend, band->x1);
 301             prc->y1 = int_min(cbgyend, band->y1);
 302
 303             tlcblkxstart =
 304               int_floordivpow2(prc->x0, cblkwidthexpn) << cblkwidthexpn;
 305             tlcblkystart =
 306               int_floordivpow2(prc->y0, cblkheightexpn) << cblkheightexpn;
 307             brcblkxend =
 308               int_ceildivpow2(prc->x1, cblkwidthexpn) << cblkwidthexpn;
 309             brcblkyend =
 310               int_ceildivpow2(prc->y1, cblkheightexpn) << cblkheightexpn;
 311             prc->cw = (brcblkxend - tlcblkxstart) >> cblkwidthexpn;
 312             prc->ch = (brcblkyend - tlcblkystart) >> cblkheightexpn;
 313
 314             prc->cblks =
 315               (tcd_cblk_t *) malloc((prc->cw * prc->ch) *
 316                                     sizeof(tcd_cblk_t));
 317             prc->incltree = tgt_create(prc->cw, prc->ch);
 318             prc->imsbtree = tgt_create(prc->cw, prc->ch);
 319
 320             for (cblkno = 0; cblkno < prc->cw * prc->ch; cblkno++) {
 321               int cblkxstart =
 322                 tlcblkxstart + (cblkno % prc->cw) * (1 << cblkwidthexpn);
 323               int cblkystart =
 324                 tlcblkystart + (cblkno / prc->cw) * (1 << cblkheightexpn);
 325               int cblkxend = cblkxstart + (1 << cblkwidthexpn);
 326               int cblkyend = cblkystart + (1 << cblkheightexpn);
 327
 328               cblk = &prc->cblks[cblkno];
 329               /* code-block size (global) */
 330               cblk->x0 = int_max(cblkxstart, prc->x0);
 331               cblk->y0 = int_max(cblkystart, prc->y0);
 332               cblk->x1 = int_min(cblkxend, prc->x1);
 333               cblk->y1 = int_min(cblkyend, prc->y1);
 334             }
 335           }
 336         }
 337       }
 338     }
 339   }
 340   /* tcd_dump(&tcd_image,curtileno); */
 341 }
 342
 343 void tcd_free_encode(j2k_image_t * img, j2k_cp_t * cp, int curtileno)
 344 {
 345   int tileno, compno, resno, bandno, precno;
 346   tcd_img = img;
 347   tcd_cp = cp;
 348   tcd_image.tw = cp->tw;
 349   tcd_image.th = cp->th;
 350   for (tileno = 0; tileno < 1; tileno++) {
 351     /* j2k_tcp_t *tcp=&cp->tcps[curtileno]; */
 352     tile = tcd_image.tiles;
 353     for (compno = 0; compno < tile->numcomps; compno++) {
 354       tilec = &tile->comps[compno];
 355       for (resno = 0; resno < tilec->numresolutions; resno++) {
 356         res = &tilec->resolutions[resno];
 357         for (bandno = 0; bandno < res->numbands; bandno++) {
 358           band = &res->bands[bandno];
 359           for (precno = 0; precno < res->pw * res->ph; precno++) {
 360             prc = &band->precincts[precno];
 361
 362             if (prc->incltree != NULL)
 363               tgt_destroy(prc->incltree);
 364             if (prc->imsbtree != NULL)
 365               tgt_destroy(prc->imsbtree);
 366             free(prc->cblks);
 367           }                     /* for (precno */
 368           free(band->precincts);
 369         }                       /* for (bandno */
 370       }                         /* for (resno */
 371       free(tilec->resolutions);
 372     }                           /* for (compno */
 373     free(tile->comps);
 374   }                             /* for (tileno */
 375   free(tcd_image.tiles);
 376 }
 377
 378 void tcd_init_encode(j2k_image_t * img, j2k_cp_t * cp, int curtileno)
 379 {
 380   int tileno, compno, resno, bandno, precno, cblkno;
 381
 382   for (tileno = 0; tileno < 1; tileno++) {
 383     j2k_tcp_t *tcp = &cp->tcps[curtileno];
 384     int j;
 385     //              int previous_x0, previous_x1, previous_y0, previous_y1;
 386     /* cfr p59 ISO/IEC FDIS15444-1 : 2000 (18 august 2000) */
 387     int p = curtileno % cp->tw;
 388     int q = curtileno / cp->tw;
 389     tile = tcd_image.tiles;
 390
 391     /* 4 borders of the tile rescale on the image if necessary */
 392     tile->x0 = int_max(cp->tx0 + p * cp->tdx, img->x0);
 393     tile->y0 = int_max(cp->ty0 + q * cp->tdy, img->y0);
 394     tile->x1 = int_min(cp->tx0 + (p + 1) * cp->tdx, img->x1);
 395     tile->y1 = int_min(cp->ty0 + (q + 1) * cp->tdy, img->y1);
 396
 397     tile->numcomps = img->numcomps;
 398     /* tile->PPT=img->PPT; */
 399     /* Modification of the RATE >> */
 400     for (j = 0; j < tcp->numlayers; j++) {
 401       tcp->rates[j] = tcp->rates[j] ? int_ceildiv(tile->numcomps * (tile->x1 - tile->x0) * (tile->y1 - tile->y0) * img->comps[0].prec, (tcp->rates[j] * 8 * img->comps[0].dx * img->comps[0].dy)) : 0;  //Mod antonin losslessbug
 402       if (tcp->rates[j]) {
 403         if (j && tcp->rates[j] < tcp->rates[j - 1] + 10) {
 404           tcp->rates[j] = tcp->rates[j - 1] + 20;
 405         } else {
 406           if (!j && tcp->rates[j] < 30)
 407             tcp->rates[j] = 30;
 408         }
 409       }
 410     }
 411     /* << Modification of the RATE */
 412     /* tile->comps=(tcd_tilecomp_t*)realloc(tile->comps,img->numcomps*sizeof(tcd_tilecomp_t)); */
 413     for (compno = 0; compno < tile->numcomps; compno++) {
 414       j2k_tccp_t *tccp = &tcp->tccps[compno];
 415       /* int realloc_op; */
 416
 417       tilec = &tile->comps[compno];
 418       /* border of each tile component (global) */
 419       tilec->x0 = int_ceildiv(tile->x0, img->comps[compno].dx);
 420       tilec->y0 = int_ceildiv(tile->y0, img->comps[compno].dy);
 421       tilec->x1 = int_ceildiv(tile->x1, img->comps[compno].dx);
 422       tilec->y1 = int_ceildiv(tile->y1, img->comps[compno].dy);
 423
 424       tilec->data =
 425         (int *) malloc((tilec->x1 - tilec->x0) *
 426                        (tilec->y1 - tilec->y0) * sizeof(int));
 427       tilec->numresolutions = tccp->numresolutions;
 428       /* tilec->resolutions=(tcd_resolution_t*)realloc(tilec->resolutions,tilec->numresolutions*sizeof(tcd_resolution_t)); */
 429       for (resno = 0; resno < tilec->numresolutions; resno++) {
 430         int pdx, pdy;
 431         int levelno = tilec->numresolutions - 1 - resno;
 432         int tlprcxstart, tlprcystart, brprcxend, brprcyend;
 433         int tlcbgxstart, tlcbgystart, brcbgxend, brcbgyend;
 434         int cbgwidthexpn, cbgheightexpn;
 435         int cblkwidthexpn, cblkheightexpn;
 436
 437         res = &tilec->resolutions[resno];
 438         /* border for each resolution level (global) */
 439         res->x0 = int_ceildivpow2(tilec->x0, levelno);
 440         res->y0 = int_ceildivpow2(tilec->y0, levelno);
 441         res->x1 = int_ceildivpow2(tilec->x1, levelno);
 442         res->y1 = int_ceildivpow2(tilec->y1, levelno);
 443
 444         res->numbands = resno == 0 ? 1 : 3;
 445         /* p. 35, table A-23, ISO/IEC FDIS154444-1 : 2000 (18 august 2000) */
 446         if (tccp->csty & J2K_CCP_CSTY_PRT) {
 447           pdx = tccp->prcw[resno];
 448           pdy = tccp->prch[resno];
 449         } else {
 450           pdx = 15;
 451           pdy = 15;
 452         }
 453         /* p. 64, B.6, ISO/IEC FDIS15444-1 : 2000 (18 august 2000)  */
 454         tlprcxstart = int_floordivpow2(res->x0, pdx) << pdx;
 455         tlprcystart = int_floordivpow2(res->y0, pdy) << pdy;
 456         brprcxend = int_ceildivpow2(res->x1, pdx) << pdx;
 457         brprcyend = int_ceildivpow2(res->y1, pdy) << pdy;
 458
 459         res->pw = (brprcxend - tlprcxstart) >> pdx;
 460         res->ph = (brprcyend - tlprcystart) >> pdy;
 461
 462         if (resno == 0) {
 463           tlcbgxstart = tlprcxstart;
 464           tlcbgystart = tlprcystart;
 465           brcbgxend = brprcxend;
 466           brcbgyend = brprcyend;
 467           cbgwidthexpn = pdx;
 468           cbgheightexpn = pdy;
 469         } else {
 470           tlcbgxstart = int_ceildivpow2(tlprcxstart, 1);
 471           tlcbgystart = int_ceildivpow2(tlprcystart, 1);
 472           brcbgxend = int_ceildivpow2(brprcxend, 1);
 473           brcbgyend = int_ceildivpow2(brprcyend, 1);
 474           cbgwidthexpn = pdx - 1;
 475           cbgheightexpn = pdy - 1;
 476         }
 477
 478         cblkwidthexpn = int_min(tccp->cblkw, cbgwidthexpn);
 479         cblkheightexpn = int_min(tccp->cblkh, cbgheightexpn);
 480
 481         for (bandno = 0; bandno < res->numbands; bandno++) {
 482           int x0b, y0b;
 483           int gain, numbps;
 484           j2k_stepsize_t *ss;
 485           band = &res->bands[bandno];
 486           band->bandno = resno == 0 ? 0 : bandno + 1;
 487           x0b = (band->bandno == 1) || (band->bandno == 3) ? 1 : 0;
 488           y0b = (band->bandno == 2) || (band->bandno == 3) ? 1 : 0;
 489
 490           if (band->bandno == 0) {
 491             /* band border */
 492             band->x0 = int_ceildivpow2(tilec->x0, levelno);
 493             band->y0 = int_ceildivpow2(tilec->y0, levelno);
 494             band->x1 = int_ceildivpow2(tilec->x1, levelno);
 495             band->y1 = int_ceildivpow2(tilec->y1, levelno);
 496           } else {
 497             band->x0 =
 498               int_ceildivpow2(tilec->x0 -
 499                               (1 << levelno) * x0b, levelno + 1);
 500             band->y0 =
 501               int_ceildivpow2(tilec->y0 -
 502                               (1 << levelno) * y0b, levelno + 1);
 503             band->x1 =
 504               int_ceildivpow2(tilec->x1 -
 505                               (1 << levelno) * x0b, levelno + 1);
 506             band->y1 =
 507               int_ceildivpow2(tilec->y1 -
 508                               (1 << levelno) * y0b, levelno + 1);
 509           }
 510
 511           ss = &tccp->stepsizes[resno ==
 512                                 0 ? 0 : 3 * (resno - 1) + bandno + 1];
 513           gain =
 514             tccp->qmfbid ==
 515             0 ? dwt_getgain_real(band->bandno) : dwt_getgain(band->bandno);
 516           numbps = img->comps[compno].prec + gain;
 517           band->stepsize =
 518             (int) floor((1.0 + ss->mant / 2048.0) *
 519                         pow(2.0, numbps - ss->expn) * 8192.0);
 520           band->numbps = ss->expn + tccp->numgbits - 1; /* WHY -1 ? */
 521
 522           for (precno = 0; precno < res->pw * res->ph; precno++) {
 523             int tlcblkxstart, tlcblkystart, brcblkxend, brcblkyend;
 524             int cbgxstart =
 525               tlcbgxstart + (precno % res->pw) * (1 << cbgwidthexpn);
 526             int cbgystart =
 527               tlcbgystart + (precno / res->pw) * (1 << cbgheightexpn);
 528             int cbgxend = cbgxstart + (1 << cbgwidthexpn);
 529             int cbgyend = cbgystart + (1 << cbgheightexpn);
 530
 531             prc = &band->precincts[precno];
 532             /* precinct size (global) */
 533             prc->x0 = int_max(cbgxstart, band->x0);
 534             prc->y0 = int_max(cbgystart, band->y0);
 535             prc->x1 = int_min(cbgxend, band->x1);
 536             prc->y1 = int_min(cbgyend, band->y1);
 537
 538             tlcblkxstart =
 539               int_floordivpow2(prc->x0, cblkwidthexpn) << cblkwidthexpn;
 540             tlcblkystart =
 541               int_floordivpow2(prc->y0, cblkheightexpn) << cblkheightexpn;
 542             brcblkxend =
 543               int_ceildivpow2(prc->x1, cblkwidthexpn) << cblkwidthexpn;
 544             brcblkyend =
 545               int_ceildivpow2(prc->y1, cblkheightexpn) << cblkheightexpn;
 546             prc->cw = (brcblkxend - tlcblkxstart) >> cblkwidthexpn;
 547             prc->ch = (brcblkyend - tlcblkystart) >> cblkheightexpn;
 548
 549             free(prc->cblks);
 550             prc->cblks =
 551               (tcd_cblk_t *) malloc(prc->cw * prc->ch *
 552                                     sizeof(tcd_cblk_t));
 553
 554             if (prc->incltree != NULL)
 555               tgt_destroy(prc->incltree);
 556             if (prc->imsbtree != NULL)
 557               tgt_destroy(prc->imsbtree);
 558
 559             prc->incltree = tgt_create(prc->cw, prc->ch);
 560             prc->imsbtree = tgt_create(prc->cw, prc->ch);
 561
 562             for (cblkno = 0; cblkno < prc->cw * prc->ch; cblkno++) {
 563               int cblkxstart =
 564                 tlcblkxstart + (cblkno % prc->cw) * (1 << cblkwidthexpn);
 565               int cblkystart =
 566                 tlcblkystart + (cblkno / prc->cw) * (1 << cblkheightexpn);
 567               int cblkxend = cblkxstart + (1 << cblkwidthexpn);
 568               int cblkyend = cblkystart + (1 << cblkheightexpn);
 569               cblk = &prc->cblks[cblkno];
 570
 571               /* code-block size (global) */
 572               cblk->x0 = int_max(cblkxstart, prc->x0);
 573               cblk->y0 = int_max(cblkystart, prc->y0);
 574               cblk->x1 = int_min(cblkxend, prc->x1);
 575               cblk->y1 = int_min(cblkyend, prc->y1);
 576
 577             }
 578           }
 579         }
 580       }
 581     }
 582   }
 583   /* tcd_dump(&tcd_image,0); */
 584 }
 585
 586 void tcd_init(j2k_image_t * img, j2k_cp_t * cp)
 587 {
 588   int tileno, compno, resno, bandno, precno, cblkno, i, j;
 589   unsigned int x0 = 0, y0 = 0, x1 = 0, y1 = 0, w, h, p, q;
 590   tcd_img = img;
 591   tcd_cp = cp;
 592   tcd_image.tw = cp->tw;
 593   tcd_image.th = cp->th;
 594   tcd_image.tiles =
 595     (tcd_tile_t *) malloc(cp->tw * cp->th * sizeof(tcd_tile_t));
 596
 597   /*for (tileno = 0; tileno < cp->tw * cp->th; tileno++) {
 598      j2k_tcp_t *tcp = &cp->tcps[tileno];
 599      tcd_tile_t *tile = &tcd_image.tiles[tileno]; */
 600
 601   for (i = 0; i < cp->tileno_size; i++) {
 602     j2k_tcp_t *tcp = &cp->tcps[cp->tileno[i]];
 603     tcd_tile_t *tile = &tcd_image.tiles[cp->tileno[i]];
 604     tileno = cp->tileno[i];
 605
 606
 607     //              int previous_x0, previous_x1, previous_y0, previous_y1;
 608     /* cfr p59 ISO/IEC FDIS15444-1 : 2000 (18 august 2000) */
 609     p = tileno % cp->tw;        /* si numerotation matricielle .. */
 610     q = tileno / cp->tw;        /* .. coordonnees de la tile (q,p) q pour ligne et p pour colonne */
 611
 612     /* 4 borders of the tile rescale on the image if necessary */
 613     tile->x0 = int_max(cp->tx0 + p * cp->tdx, img->x0);
 614     tile->y0 = int_max(cp->ty0 + q * cp->tdy, img->y0);
 615     tile->x1 = int_min(cp->tx0 + (p + 1) * cp->tdx, img->x1);
 616     tile->y1 = int_min(cp->ty0 + (q + 1) * cp->tdy, img->y1);
 617
 618     tile->numcomps = img->numcomps;
 619     tile->comps =
 620       (tcd_tilecomp_t *) malloc(img->numcomps * sizeof(tcd_tilecomp_t));
 621     for (compno = 0; compno < tile->numcomps; compno++) {
 622       j2k_tccp_t *tccp = &tcp->tccps[compno];
 623       tcd_tilecomp_t *tilec = &tile->comps[compno];
 624       /* border of each tile component (global) */
 625       tilec->x0 = int_ceildiv(tile->x0, img->comps[compno].dx);
 626       tilec->y0 = int_ceildiv(tile->y0, img->comps[compno].dy);
 627       tilec->x1 = int_ceildiv(tile->x1, img->comps[compno].dx);
 628       tilec->y1 = int_ceildiv(tile->y1, img->comps[compno].dy);
 629
 630       tilec->data =
 631         (int *) malloc((tilec->x1 - tilec->x0) *
 632                        (tilec->y1 - tilec->y0) * sizeof(int));
 633       tilec->numresolutions = tccp->numresolutions;
 634       tilec->resolutions =
 635         (tcd_resolution_t *) malloc(tilec->numresolutions *
 636                                     sizeof(tcd_resolution_t));
 637       for (resno = 0; resno < tilec->numresolutions; resno++) {
 638         int pdx, pdy;
 639         int levelno = tilec->numresolutions - 1 - resno;
 640         int tlprcxstart, tlprcystart, brprcxend, brprcyend;
 641         int tlcbgxstart, tlcbgystart, brcbgxend, brcbgyend;
 642         int cbgwidthexpn, cbgheightexpn;
 643         int cblkwidthexpn, cblkheightexpn;
 644         tcd_resolution_t *res = &tilec->resolutions[resno];
 645
 646         /* border for each resolution level (global) */
 647         res->x0 = int_ceildivpow2(tilec->x0, levelno);
 648         res->y0 = int_ceildivpow2(tilec->y0, levelno);
 649         res->x1 = int_ceildivpow2(tilec->x1, levelno);
 650         res->y1 = int_ceildivpow2(tilec->y1, levelno);
 651
 652         res->numbands = resno == 0 ? 1 : 3;
 653         /* p. 35, table A-23, ISO/IEC FDIS154444-1 : 2000 (18 august 2000) */
 654         if (tccp->csty & J2K_CCP_CSTY_PRT) {
 655           pdx = tccp->prcw[resno];
 656           pdy = tccp->prch[resno];
 657         } else {
 658           pdx = 15;
 659           pdy = 15;
 660         }
 661         /* p. 64, B.6, ISO/IEC FDIS15444-1 : 2000 (18 august 2000)  */
 662         tlprcxstart = int_floordivpow2(res->x0, pdx) << pdx;
 663         tlprcystart = int_floordivpow2(res->y0, pdy) << pdy;
 664         brprcxend = int_ceildivpow2(res->x1, pdx) << pdx;
 665         brprcyend = int_ceildivpow2(res->y1, pdy) << pdy;
 666         res->pw = (res->x0 == res->x1) ? 0 : ((brprcxend - tlprcxstart) >> pdx);        // Mod Antonin : sizebug1
 667         res->ph = (res->y0 == res->y1) ? 0 : ((brprcyend - tlprcystart) >> pdy);        // Mod Antonin : sizebug1
 668
 669         if (resno == 0) {
 670           tlcbgxstart = tlprcxstart;
 671           tlcbgystart = tlprcystart;
 672           brcbgxend = brprcxend;
 673           brcbgyend = brprcyend;
 674           cbgwidthexpn = pdx;
 675           cbgheightexpn = pdy;
 676         } else {
 677           tlcbgxstart = int_ceildivpow2(tlprcxstart, 1);
 678           tlcbgystart = int_ceildivpow2(tlprcystart, 1);
 679           brcbgxend = int_ceildivpow2(brprcxend, 1);
 680           brcbgyend = int_ceildivpow2(brprcyend, 1);
 681           cbgwidthexpn = pdx - 1;
 682           cbgheightexpn = pdy - 1;
 683         }
 684
 685         cblkwidthexpn = int_min(tccp->cblkw, cbgwidthexpn);
 686         cblkheightexpn = int_min(tccp->cblkh, cbgheightexpn);
 687
 688         for (bandno = 0; bandno < res->numbands; bandno++) {
 689           int x0b, y0b;
 690           int gain, numbps;
 691           j2k_stepsize_t *ss;
 692           tcd_band_t *band = &res->bands[bandno];
 693           band->bandno = resno == 0 ? 0 : bandno + 1;
 694           x0b = (band->bandno == 1) || (band->bandno == 3) ? 1 : 0;
 695           y0b = (band->bandno == 2) || (band->bandno == 3) ? 1 : 0;
 696
 697           if (band->bandno == 0) {
 698             /* band border (global) */
 699             band->x0 = int_ceildivpow2(tilec->x0, levelno);
 700             band->y0 = int_ceildivpow2(tilec->y0, levelno);
 701             band->x1 = int_ceildivpow2(tilec->x1, levelno);
 702             band->y1 = int_ceildivpow2(tilec->y1, levelno);
 703           } else {
 704             /* band border (global) */
 705             band->x0 =
 706               int_ceildivpow2(tilec->x0 -
 707                               (1 << levelno) * x0b, levelno + 1);
 708             band->y0 =
 709               int_ceildivpow2(tilec->y0 -
 710                               (1 << levelno) * y0b, levelno + 1);
 711             band->x1 =
 712               int_ceildivpow2(tilec->x1 -
 713                               (1 << levelno) * x0b, levelno + 1);
 714             band->y1 =
 715               int_ceildivpow2(tilec->y1 -
 716                               (1 << levelno) * y0b, levelno + 1);
 717           }
 718
 719           ss = &tccp->stepsizes[resno ==
 720                                 0 ? 0 : 3 * (resno - 1) + bandno + 1];
 721           gain =
 722             tccp->qmfbid ==
 723             0 ? dwt_getgain_real(band->bandno) : dwt_getgain(band->bandno);
 724           numbps = img->comps[compno].prec + gain;
 725           band->stepsize =
 726             (int) floor((1.0 + ss->mant / 2048.0) *
 727                         pow(2.0, numbps - ss->expn) * 8192.0);
 728           band->numbps = ss->expn + tccp->numgbits - 1; /* WHY -1 ? */
 729
 730           band->precincts =
 731             (tcd_precinct_t *) malloc(res->pw * res->ph *
 732                                       sizeof(tcd_precinct_t));
 733
 734           for (precno = 0; precno < res->pw * res->ph; precno++) {
 735             int tlcblkxstart, tlcblkystart, brcblkxend, brcblkyend;
 736             int cbgxstart =
 737               tlcbgxstart + (precno % res->pw) * (1 << cbgwidthexpn);
 738             int cbgystart =
 739               tlcbgystart + (precno / res->pw) * (1 << cbgheightexpn);
 740             int cbgxend = cbgxstart + (1 << cbgwidthexpn);
 741             int cbgyend = cbgystart + (1 << cbgheightexpn);
 742             tcd_precinct_t *prc = &band->precincts[precno];
 743             /* precinct size (global) */
 744             prc->x0 = int_max(cbgxstart, band->x0);
 745             prc->y0 = int_max(cbgystart, band->y0);
 746             prc->x1 = int_min(cbgxend, band->x1);
 747             prc->y1 = int_min(cbgyend, band->y1);
 748
 749             tlcblkxstart =
 750               int_floordivpow2(prc->x0, cblkwidthexpn) << cblkwidthexpn;
 751             tlcblkystart =
 752               int_floordivpow2(prc->y0, cblkheightexpn) << cblkheightexpn;
 753             brcblkxend =
 754               int_ceildivpow2(prc->x1, cblkwidthexpn) << cblkwidthexpn;
 755             brcblkyend =
 756               int_ceildivpow2(prc->y1, cblkheightexpn) << cblkheightexpn;
 757             prc->cw = (brcblkxend - tlcblkxstart) >> cblkwidthexpn;
 758             prc->ch = (brcblkyend - tlcblkystart) >> cblkheightexpn;
 759
 760             prc->cblks =
 761               (tcd_cblk_t *) malloc(prc->cw * prc->ch *
 762                                     sizeof(tcd_cblk_t));
 763
 764             prc->incltree = tgt_create(prc->cw, prc->ch);
 765             prc->imsbtree = tgt_create(prc->cw, prc->ch);
 766
 767             for (cblkno = 0; cblkno < prc->cw * prc->ch; cblkno++) {
 768               int cblkxstart =
 769                 tlcblkxstart + (cblkno % prc->cw) * (1 << cblkwidthexpn);
 770               int cblkystart =
 771                 tlcblkystart + (cblkno / prc->cw) * (1 << cblkheightexpn);
 772               int cblkxend = cblkxstart + (1 << cblkwidthexpn);
 773               int cblkyend = cblkystart + (1 << cblkheightexpn);
 774               tcd_cblk_t *cblk = &prc->cblks[cblkno];
 775               /* code-block size (global) */
 776               cblk->x0 = int_max(cblkxstart, prc->x0);
 777               cblk->y0 = int_max(cblkystart, prc->y0);
 778               cblk->x1 = int_min(cblkxend, prc->x1);
 779               cblk->y1 = int_min(cblkyend, prc->y1);
 780
 781
 782
 783               cblk->lastbp = 0; // Add Antonin : quantizbug1
 784             }
 785           }
 786         }
 787       }
 788     }
 789   }
 790   //tcd_dump(&tcd_image,0);
 791
 792
 793   /* Allocate place to store the data decoded = final image */
 794   /* Place limited by the tile really present in the codestream */
 795
 796
 797   for (i = 0; i < img->numcomps; i++) {
 798     for (j = 0; j < cp->tileno_size; j++) {
 799       tileno = cp->tileno[j];
 800       x0 = j == 0 ? tcd_image.tiles[tileno].comps[i].x0 : int_min(x0,
 801                                                                   tcd_image.
 802                                                                   tiles
 803                                                                   [tileno].
 804                                                                   comps
 805                                                                   [i].x0);
 806       y0 =
 807         j == 0 ? tcd_image.tiles[tileno].comps[i].y0 : int_min(y0,
 808                                                                tcd_image.
 809                                                                tiles
 810                                                                [tileno].
 811                                                                comps[i].
 812                                                                y0);
 813       x1 =
 814         j == 0 ? tcd_image.tiles[tileno].comps[i].x1 : int_max(x1,
 815                                                                tcd_image.
 816                                                                tiles
 817                                                                [tileno].
 818                                                                comps[i].
 819                                                                x1);
 820       y1 =
 821         j == 0 ? tcd_image.tiles[tileno].comps[i].y1 : int_max(y1,
 822                                                                tcd_image.
 823                                                                tiles
 824                                                                [tileno].
 825                                                                comps[i].
 826                                                                y1);
 827     }
 828     //w = int_ceildiv(x1 - x0, img->comps[i].dx);
 829     //h = int_ceildiv(y1 - y0, img->comps[i].dy);
 830
 831     w = x1 - x0;
 832
 833     h = y1 - y0;
 834     img->comps[i].data = (int *) calloc(w * h, sizeof(int));
 835     img->comps[i].w = w;
 836     img->comps[i].h = h;
 837     img->comps[i].x0 = x0;
 838     img->comps[i].y0 = y0;
 839   }
 840 }
 841
 842 void tcd_makelayer_fixed(int layno, int final)
 843 {
 844   int compno, resno, bandno, precno, cblkno;
 845   int value;                    //, matrice[tcd_tcp->numlayers][tcd_tile->comps[0].numresolutions][3];
 846   int matrice[10][10][3];
 847   int i, j, k;
 848
 849   /*matrice=(int*)malloc(tcd_tcp->numlayers*tcd_tile->comps[0].numresolutions*3*sizeof(int)); */
 850
 851   for (compno = 0; compno < tcd_tile->numcomps; compno++) {
 852     tcd_tilecomp_t *tilec = &tcd_tile->comps[compno];
 853     for (i = 0; i < tcd_tcp->numlayers; i++) {
 854       for (j = 0; j < tilec->numresolutions; j++) {
 855         for (k = 0; k < 3; k++) {
 856           matrice[i][j][k] =
 857             (int) (tcd_cp->
 858                    matrice[i * tilec->numresolutions * 3 +
 859                            j * 3 +
 860                            k] *
 861                    (float) (tcd_img->comps[compno].prec / 16.0));
 862     }}}
 863
 864     for (resno = 0; resno < tilec->numresolutions; resno++) {
 865       tcd_resolution_t *res = &tilec->resolutions[resno];
 866       for (bandno = 0; bandno < res->numbands; bandno++) {
 867         tcd_band_t *band = &res->bands[bandno];
 868         for (precno = 0; precno < res->pw * res->ph; precno++) {
 869           tcd_precinct_t *prc = &band->precincts[precno];
 870           for (cblkno = 0; cblkno < prc->cw * prc->ch; cblkno++) {
 871             tcd_cblk_t *cblk = &prc->cblks[cblkno];
 872             tcd_layer_t *layer = &cblk->layers[layno];
 873             int n;
 874             int imsb = tcd_img->comps[compno].prec - cblk->numbps;      /* number of bit-plan equal to zero */
 875             /* Correction of the matrix of coefficient to include the IMSB information */
 876
 877             if (layno == 0) {
 878               value = matrice[layno][resno][bandno];
 879               if (imsb >= value)
 880                 value = 0;
 881               else
 882                 value -= imsb;
 883             } else {
 884               value =
 885                 matrice[layno][resno][bandno] -
 886                 matrice[layno - 1][resno][bandno];
 887               if (imsb >= matrice[layno - 1][resno][bandno]) {
 888                 value -= (imsb - matrice[layno - 1][resno][bandno]);
 889                 if (value < 0)
 890                   value = 0;
 891               }
 892             }
 893
 894             if (layno == 0)
 895               cblk->numpassesinlayers = 0;
 896
 897             n = cblk->numpassesinlayers;
 898             if (cblk->numpassesinlayers == 0) {
 899               if (value != 0)
 900                 n = 3 * value - 2 + cblk->numpassesinlayers;
 901               else
 902                 n = cblk->numpassesinlayers;
 903             } else
 904               n = 3 * value + cblk->numpassesinlayers;
 905
 906             layer->numpasses = n - cblk->numpassesinlayers;
 907
 908             if (!layer->numpasses)
 909               continue;
 910
 911             if (cblk->numpassesinlayers == 0) {
 912               layer->len = cblk->passes[n - 1].rate;
 913               layer->data = cblk->data;
 914             } else {
 915               layer->len =
 916                 cblk->passes[n - 1].rate -
 917                 cblk->passes[cblk->numpassesinlayers - 1].rate;
 918               layer->data =
 919                 cblk->data +
 920                 cblk->passes[cblk->numpassesinlayers - 1].rate;
 921             }
 922             if (final)
 923               cblk->numpassesinlayers = n;
 924           }
 925         }
 926       }
 927     }
 928   }
 929 }
 930
 931 void tcd_rateallocate_fixed()
 932 {
 933   int layno;
 934
 935   for (layno = 0; layno < tcd_tcp->numlayers; layno++) {
 936     tcd_makelayer_fixed(layno, 1);
 937   }
 938 }
 939
 940 void tcd_makelayer(int layno, double thresh, int final)
 941 {
 942   int compno, resno, bandno, precno, cblkno, passno;
 943
 944   tcd_tile->distolayer[layno] = 0;      //add fixed_quality
 945
 946   for (compno = 0; compno < tcd_tile->numcomps; compno++) {
 947     tcd_tilecomp_t *tilec = &tcd_tile->comps[compno];
 948     for (resno = 0; resno < tilec->numresolutions; resno++) {
 949       tcd_resolution_t *res = &tilec->resolutions[resno];
 950       for (bandno = 0; bandno < res->numbands; bandno++) {
 951         tcd_band_t *band = &res->bands[bandno];
 952         for (precno = 0; precno < res->pw * res->ph; precno++) {
 953           tcd_precinct_t *prc = &band->precincts[precno];
 954           for (cblkno = 0; cblkno < prc->cw * prc->ch; cblkno++) {
 955             tcd_cblk_t *cblk = &prc->cblks[cblkno];
 956             tcd_layer_t *layer = &cblk->layers[layno];
 957             int n;
 958
 959             if (layno == 0) {
 960               cblk->numpassesinlayers = 0;
 961             }
 962             n = cblk->numpassesinlayers;
 963             for (passno = cblk->numpassesinlayers;
 964                  passno < cblk->totalpasses; passno++) {
 965               int dr;
 966               double dd;
 967               tcd_pass_t *pass = &cblk->passes[passno];
 968               if (n == 0) {
 969                 dr = pass->rate;
 970                 dd = pass->distortiondec;
 971               } else {
 972                 dr = pass->rate - cblk->passes[n - 1].rate;
 973                 dd = pass->distortiondec - cblk->passes[n -
 974                                                         1].distortiondec;
 975               }
 976               if (!dr) {
 977                 if (dd)
 978                   n = passno + 1;
 979                 continue;
 980               }
 981
 982               if (dd / dr >= thresh)
 983                 n = passno + 1;
 984             }
 985             layer->numpasses = n - cblk->numpassesinlayers;
 986
 987             if (!layer->numpasses) {
 988               layer->disto = 0;
 989               continue;
 990             }
 991
 992             if (cblk->numpassesinlayers == 0) {
 993               layer->len = cblk->passes[n - 1].rate;
 994               layer->data = cblk->data;
 995               layer->disto = cblk->passes[n - 1].distortiondec;
 996             } else {
 997               layer->len = cblk->passes[n - 1].rate -
 998                 cblk->passes[cblk->numpassesinlayers - 1].rate;
 999               layer->data =
1000                 cblk->data +
1001                 cblk->passes[cblk->numpassesinlayers - 1].rate;
1002               layer->disto =
1003                 cblk->passes[n - 1].distortiondec -
1004                 cblk->passes[cblk->numpassesinlayers - 1].distortiondec;
1005             }
1006
1007             tcd_tile->distolayer[layno] += layer->disto;        //add fixed_quality
1008
1009             if (final)
1010               cblk->numpassesinlayers = n;
1011           }
1012         }
1013       }
1014     }
1015   }
1016 }
1017
1018 void tcd_rateallocate(unsigned char *dest, int len, info_image * info_IM)
1019 {
1020   int compno, resno, bandno, precno, cblkno, passno, layno;
1021   double min, max;
1022   double cumdisto[100];         //add fixed_quality
1023   const double K = 1;           // 1.1; //add fixed_quality
1024
1025   double maxSE = 0;
1026   min = DBL_MAX;
1027   max = 0;
1028
1029   tcd_tile->nbpix = 0;          //add fixed_quality
1030
1031   for (compno = 0; compno < tcd_tile->numcomps; compno++) {
1032     tcd_tilecomp_t *tilec = &tcd_tile->comps[compno];
1033
1034     tilec->nbpix = 0;
1035     for (resno = 0; resno < tilec->numresolutions; resno++) {
1036       tcd_resolution_t *res = &tilec->resolutions[resno];
1037       for (bandno = 0; bandno < res->numbands; bandno++) {
1038         tcd_band_t *band = &res->bands[bandno];
1039         for (precno = 0; precno < res->pw * res->ph; precno++) {
1040           tcd_precinct_t *prc = &band->precincts[precno];
1041           for (cblkno = 0; cblkno < prc->cw * prc->ch; cblkno++) {
1042             tcd_cblk_t *cblk = &prc->cblks[cblkno];
1043             for (passno = 0; passno < cblk->totalpasses; passno++) {
1044               tcd_pass_t *pass = &cblk->passes[passno];
1045               int dr;
1046               double dd, rdslope;
1047               if (passno == 0) {
1048                 dr = pass->rate;
1049                 dd = pass->distortiondec;
1050               } else {
1051                 dr = pass->rate - cblk->passes[passno - 1].rate;
1052                 dd = pass->distortiondec -
1053                   cblk->passes[passno - 1].distortiondec;
1054               }
1055               if (dr == 0) {
1056                 continue;
1057               }
1058
1059               rdslope = dd / dr;
1060
1061               if (rdslope < min) {
1062                 min = rdslope;
1063               }
1064               if (rdslope > max) {
1065                 max = rdslope;
1066               }
1067             }                   /* passno */
1068
1069             tcd_tile->nbpix += ((cblk->x1 - cblk->x0) * (cblk->y1 - cblk->y0)); //add fixed_quality
1070
1071             tilec->nbpix += ((cblk->x1 - cblk->x0) * (cblk->y1 - cblk->y0));    //add fixed_quality
1072
1073           }                     /* cbklno */
1074         }                       /* precno */
1075       }                         /* bandno */
1076     }                           /* resno */
1077
1078     maxSE +=
1079       (double) (((1 << tcd_img->comps[compno].prec) -
1080                  1) * ((1 << tcd_img->comps[compno].prec) -
1081                        1)) * (tilec->nbpix);
1082   }                             /* compno */
1083
1084   /* add antonin index */
1085   if (info_IM->index_on) {
1086     info_tile *info_TL = &info_IM->tile[tcd_tileno];
1087     info_TL->nbpix = tcd_tile->nbpix;
1088     info_TL->distotile = tcd_tile->distotile;
1089     info_TL->thresh =
1090       (double *) malloc(tcd_tcp->numlayers * sizeof(double));
1091   }
1092   /* dda */
1093
1094   for (layno = 0; layno < tcd_tcp->numlayers; layno++) {
1095     volatile double lo = min;
1096     volatile double hi = max;
1097     volatile int success = 0;
1098     volatile int maxlen = tcd_tcp->rates[layno] ? int_min(tcd_tcp->rates[layno], len) : len;    //Mod antonin losslessbug
1099     volatile double goodthresh;
1100     volatile int i;
1101     double distotarget;         //add fixed_quality
1102
1103     distotarget = tcd_tile->distotile - ((K * maxSE) / pow(10, tcd_tcp->distoratio[layno] / 10));       // add fixed_quality
1104
1105     if (tcd_tcp->rates[layno]) {
1106       for (i = 0; i < 32; i++) {
1107         volatile double thresh = (lo + hi) / 2;
1108         int l = 0;
1109         double distoachieved = 0;       // add fixed_quality
1110
1111         tcd_makelayer(layno, thresh, 0);
1112
1113         if (tcd_cp->fixed_quality) {    // add fixed_quality
1114           distoachieved =
1115             layno ==
1116             0 ? tcd_tile->distolayer[0] : cumdisto[layno - 1] +
1117             tcd_tile->distolayer[layno];
1118           if (distoachieved < distotarget) {
1119             hi = thresh;
1120             continue;
1121           }
1122           lo = thresh;
1123         } else {
1124           l =
1125             t2_encode_packets(tcd_img, tcd_cp, tcd_tileno, tcd_tile,
1126                               layno + 1, dest, maxlen, info_IM);
1127           /* fprintf(stderr, "rate alloc: len=%d, max=%d\n", l, maxlen); */
1128           if (l == -999) {
1129             lo = thresh;
1130             continue;
1131           }
1132           hi = thresh;
1133         }
1134
1135         success = 1;
1136         goodthresh = thresh;
1137       }
1138     } else {
1139       success = 1;
1140       goodthresh = min;
1141     }
1142
1143     if (!success) {
1144       longjmp(j2k_error, 1);
1145     }
1146
1147     if (info_IM->index_on) {    /* Threshold for Marcela Index */
1148       info_IM->tile[tcd_tileno].thresh[layno] = goodthresh;
1149     }
1150     tcd_makelayer(layno, goodthresh, 1);
1151
1152     cumdisto[layno] = layno == 0 ? tcd_tile->distolayer[0] : cumdisto[layno - 1] + tcd_tile->distolayer[layno]; // add fixed_quality
1153   }
1154 }
1155
1156 int
1157 tcd_encode_tile_pxm(int tileno, unsigned char *dest, int len,
1158                     info_image * info_IM)
1159 {
1160   int compno;
1161   int l, i;
1162   clock_t time7;
1163   tcd_tile_t *tile;
1164   j2k_tcp_t *tcp = &tcd_cp->tcps[0];
1165   j2k_tccp_t *tccp = &tcp->tccps[0];
1166
1167   tcd_tileno = tileno;
1168   tcd_tile = tcd_image.tiles;
1169   tcd_tcp = &tcd_cp->tcps[tileno];
1170   tile = tcd_tile;
1171   /* INDEX >> "Precinct_nb_X et Precinct_nb_Y" */
1172   if (info_IM->index_on) {
1173     tcd_tilecomp_t *tilec_idx = &tile->comps[0];        //Based on Component 0
1174
1175     for (i = 0; i < tilec_idx->numresolutions; i++) {
1176
1177       tcd_resolution_t *res_idx = &tilec_idx->resolutions[i];
1178
1179       info_IM->tile[tileno].pw[i] = res_idx->pw;
1180       info_IM->tile[tileno].ph[i] = res_idx->ph;
1181
1182       info_IM->tile[tileno].pdx[i] = tccp->prcw[i];
1183       info_IM->tile[tileno].pdy[i] = tccp->prch[i];
1184
1185     }
1186   }
1187   /* << INDEX */
1188
1189 /*---------------TILE-------------------*/
1190
1191   time7 = clock();
1192
1193   for (compno = 0; compno < tile->numcomps; compno++) {
1194     FILE *src;
1195     char tmp[256];
1196     int k;
1197     unsigned char elmt;
1198     int i, j;
1199     int tw, w;
1200     tcd_tilecomp_t *tilec = &tile->comps[compno];
1201     int adjust =
1202       tcd_img->comps[compno].sgnd ? 0 : 1 << (tcd_img->comps[compno].
1203                                               prec - 1);
1204     int offset_x, offset_y;
1205
1206     offset_x = int_ceildiv(tcd_img->x0, tcd_img->comps[compno].dx);
1207     offset_y = int_ceildiv(tcd_img->y0, tcd_img->comps[compno].dy);
1208     tw = tilec->x1 - tilec->x0;
1209     w = int_ceildiv(tcd_img->x1 - tcd_img->x0, tcd_img->comps[compno].dx);
1210     sprintf(tmp, "Compo%d", compno);    /* component file */
1211     src = fopen(tmp, "rb");
1212     if (!src) {
1213       fprintf(stderr, "failed to open %s for reading\n", tmp);
1214       return 1;
1215     }
1216
1217     /* read the Compo file to extract data of the tile */
1218     k = 0;
1219     fseek(src, (tilec->x0 - offset_x) + (tilec->y0 - offset_y) * w,
1220           SEEK_SET);
1221     k = (tilec->x0 - offset_x) + (tilec->y0 - offset_y) * w;
1222     for (j = tilec->y0; j < tilec->y1; j++) {
1223       for (i = tilec->x0; i < tilec->x1; i++) {
1224         if (tcd_tcp->tccps[compno].qmfbid == 1) {
1225           elmt = fgetc(src);
1226           tilec->data[i - tilec->x0 + (j - tilec->y0) * tw] =
1227             elmt - adjust;
1228           k++;
1229         } else if (tcd_tcp->tccps[compno].qmfbid == 0) {
1230           elmt = fgetc(src);
1231           tilec->data[i - tilec->x0 + (j - tilec->y0) * tw] =
1232             (elmt - adjust) << 13;
1233           k++;
1234         }
1235       }
1236       fseek(src, (tilec->x0 - offset_x) + (j + 1 - offset_y) * w - k,
1237             SEEK_CUR);
1238       k = tilec->x0 - offset_x + (j + 1 - offset_y) * w;
1239
1240     }
1241     fclose(src);
1242   }
1243
1244 /*----------------MCT-------------------*/
1245
1246   if (tcd_tcp->mct) {
1247     if (tcd_tcp->tccps[0].qmfbid == 0) {
1248       mct_encode_real(tile->comps[0].data, tile->comps[1].data,
1249                       tile->comps[2].data,
1250                       (tile->comps[0].x1 -
1251                        tile->comps[0].x0) * (tile->comps[0].y1 -
1252                                              tile->comps[0].y0));
1253     } else {
1254       mct_encode(tile->comps[0].data, tile->comps[1].data,
1255                  tile->comps[2].data,
1256                  (tile->comps[0].x1 -
1257                   tile->comps[0].x0) * (tile->comps[0].y1 -
1258                                         tile->comps[0].y0));
1259     }
1260   }
1261 /*----------------DWT---------------------*/
1262
1263   /* time3=clock(); */
1264   for (compno = 0; compno < tile->numcomps; compno++) {
1265     tcd_tilecomp_t *tilec = &tile->comps[compno];
1266     if (tcd_tcp->tccps[compno].qmfbid == 1) {
1267       dwt_encode(tilec->data, tilec->x1 - tilec->x0,
1268                  tilec->y1 - tilec->y0, tilec, tilec->numresolutions - 1);
1269     } else if (tcd_tcp->tccps[compno].qmfbid == 0) {
1270       dwt_encode_real(tilec->data, tilec->x1 - tilec->x0,
1271                       tilec->y1 - tilec->y0, tilec,
1272                       tilec->numresolutions - 1);
1273     }
1274   }
1275 /*------------------TIER1-----------------*/
1276
1277   t1_init_luts();
1278   t1_encode_cblks(tile, tcd_tcp);
1279
1280 /*-----------RATE-ALLOCATE------------------*/
1281   info_IM->index_write = 0;     /* INDEX     */
1282
1283   if (tcd_cp->disto_alloc || tcd_cp->fixed_quality)     // mod fixed_quality
1284     /* Normal Rate/distortion allocation */
1285     tcd_rateallocate(dest, len, info_IM);
1286   else
1287     /* Fixed layer allocation */
1288     tcd_rateallocate_fixed();
1289
1290 /*--------------TIER2------------------*/
1291   info_IM->index_write = 1;     /* INDEX     */
1292   l = t2_encode_packets(tcd_img, tcd_cp, tileno, tile,
1293                         tcd_tcp->numlayers, dest, len, info_IM);
1294 /*---------------CLEAN-------------------*/
1295
1296   time7 = clock() - time7;
1297   fprintf(stdout,"total:     %ld.%.3ld s\n", time7 / CLOCKS_PER_SEC,
1298          (time7 % CLOCKS_PER_SEC) * 1000 / CLOCKS_PER_SEC);
1299
1300   /* cleaning memory */
1301   for (compno = 0; compno < tile->numcomps; compno++) {
1302     tilec = &tile->comps[compno];
1303     free(tilec->data);
1304   }
1305
1306   return l;
1307 }
1308
1309 int
1310 tcd_encode_tile_pgx(int tileno, unsigned char *dest, int len,
1311                     info_image * info_IM)
1312 {
1313   int compno;
1314   int l, i;
1315   clock_t time;
1316   tcd_tile_t *tile;
1317   j2k_tcp_t *tcp = &tcd_cp->tcps[0];
1318   j2k_tccp_t *tccp = &tcp->tccps[0];
1319
1320   tcd_tileno = tileno;
1321   tcd_tile = tcd_image.tiles;
1322   tcd_tcp = &tcd_cp->tcps[tileno];
1323   tile = tcd_tile;
1324   /* INDEX >> "Precinct_nb_X et Precinct_nb_Y" */
1325
1326   if (info_IM->index_on) {
1327
1328     tcd_tilecomp_t *tilec_idx = &tile->comps[0];        //Based on Component 0
1329
1330     for (i = 0; i < tilec_idx->numresolutions; i++) {
1331
1332       tcd_resolution_t *res_idx = &tilec_idx->resolutions[i];
1333
1334
1335
1336       info_IM->tile[tileno].pw[i] = res_idx->pw;
1337
1338       info_IM->tile[tileno].ph[i] = res_idx->ph;
1339
1340
1341
1342       info_IM->tile[tileno].pdx[i] = tccp->prcw[i];
1343
1344       info_IM->tile[tileno].pdy[i] = tccp->prch[i];
1345
1346     }
1347
1348   }
1349
1350   /* << INDEX */
1351 /*---------------TILE-------------------*/
1352   time = clock();
1353
1354   for (compno = 0; compno < tile->numcomps; compno++) {
1355     FILE *src;
1356     char tmp[256];
1357     int k;
1358     int elmt;
1359     int i, j;
1360     int tw, w;
1361     tcd_tilecomp_t *tilec = &tile->comps[compno];
1362     int adjust =
1363       tcd_img->comps[compno].sgnd ? 0 : 1 << (tcd_img->comps[compno].
1364                                               prec - 1);
1365     int offset_x, offset_y;
1366
1367     offset_x = int_ceildiv(tcd_img->x0, tcd_img->comps[compno].dx);
1368     offset_y = int_ceildiv(tcd_img->y0, tcd_img->comps[compno].dy);
1369     tw = tilec->x1 - tilec->x0;
1370     w = int_ceildiv(tcd_img->x1 - tcd_img->x0, tcd_img->comps[compno].dx);
1371     sprintf(tmp, "bandtile%d", tileno / tcd_cp->tw + 1);        /* bandtile file opening */
1372     src = fopen(tmp, "rb");
1373     if (!src) {
1374       fprintf(stderr, "failed to open %s for reading\n", tmp);
1375       return 1;
1376     }
1377     /* Extract data from bandtile file limited to the current tile */
1378     k = 0;
1379     while (k < tilec->x0 - offset_x) {
1380       k++;
1381       fscanf(src, "%d", &elmt);
1382     }
1383
1384     for (j = 0; j < tilec->y1 - tilec->y0; j++) {
1385       for (i = tilec->x0; i < tilec->x1; i++) {
1386         if (tcd_tcp->tccps[compno].qmfbid == 1) {
1387           fscanf(src, "%d", &elmt);
1388           tilec->data[i - tilec->x0 + (j) * tw] = elmt - adjust;
1389           k++;
1390         } else if (tcd_tcp->tccps[compno].qmfbid == 0) {
1391           fscanf(src, "%d", &elmt);
1392           tilec->data[i - tilec->x0 + (j) * tw] = (elmt - adjust) << 13;
1393           k++;
1394         }
1395       }
1396       while (k < tilec->x0 - offset_x + (j + 1) * w) {
1397         k++;
1398         fscanf(src, "%d", &elmt);
1399       }
1400     }
1401     fclose(src);
1402   }
1403
1404 /*----------------MCT-------------------*/
1405
1406   if (tcd_tcp->mct) {
1407     if (tcd_tcp->tccps[0].qmfbid == 0) {
1408       mct_encode_real(tile->comps[0].data, tile->comps[1].data,
1409                       tile->comps[2].data,
1410                       (tile->comps[0].x1 -
1411                        tile->comps[0].x0) * (tile->comps[0].y1 -
1412                                              tile->comps[0].y0));
1413     } else {
1414       mct_encode(tile->comps[0].data, tile->comps[1].data,
1415                  tile->comps[2].data,
1416                  (tile->comps[0].x1 -
1417                   tile->comps[0].x0) * (tile->comps[0].y1 -
1418                                         tile->comps[0].y0));
1419     }
1420   }
1421
1422 /*----------------DWT---------------------*/
1423
1424   for (compno = 0; compno < tile->numcomps; compno++) {
1425     tcd_tilecomp_t *tilec = &tile->comps[compno];
1426     if (tcd_tcp->tccps[compno].qmfbid == 1) {
1427       dwt_encode(tilec->data, tilec->x1 - tilec->x0,
1428                  tilec->y1 - tilec->y0, tilec, tilec->numresolutions - 1);
1429     } else if (tcd_tcp->tccps[compno].qmfbid == 0) {
1430       dwt_encode_real(tilec->data, tilec->x1 - tilec->x0,
1431                       tilec->y1 - tilec->y0, tilec,
1432                       tilec->numresolutions - 1);
1433     }
1434   }
1435
1436 /*------------------TIER1-----------------*/
1437
1438   t1_init_luts();
1439   t1_encode_cblks(tile, tcd_tcp);
1440
1441 /*-----------RATE-ALLOCATE------------------*/
1442
1443   info_IM->index_write = 0;     /* INDEX */
1444
1445   if (tcd_cp->disto_alloc || tcd_cp->fixed_quality)     // mod fixed_quality
1446
1447     /* Normal Rate/distortion allocation */
1448
1449     tcd_rateallocate(dest, len, info_IM);
1450
1451   else
1452     /* Fixed layer allocation */
1453
1454     tcd_rateallocate_fixed();
1455
1456 /*--------------TIER2------------------*/
1457   info_IM->index_write = 1;     /* INDEX */
1458
1459   l = t2_encode_packets(tcd_img, tcd_cp, tileno, tile,
1460                         tcd_tcp->numlayers, dest, len, info_IM);
1461
1462  /*---------------CLEAN-------------------*/
1463   time = clock() - time;
1464   fprintf(stdout,"total:     %ld.%.3ld s\n", time / CLOCKS_PER_SEC,
1465          (time % CLOCKS_PER_SEC) * 1000 / CLOCKS_PER_SEC);
1466
1467   for (compno = 0; compno < tile->numcomps; compno++) {
1468     tilec = &tile->comps[compno];
1469     free(tilec->data);
1470   }
1471
1472   return l;
1473 }
1474
1475
1476 int tcd_decode_tile(unsigned char *src, int len, int tileno)
1477 {
1478   int l;
1479   int compno;
1480   int eof = 0;
1481   clock_t time;
1482   tcd_tile_t *tile;
1483
1484   tcd_tileno = tileno;
1485   tcd_tile = &tcd_image.tiles[tileno];
1486   tcd_tcp = &tcd_cp->tcps[tileno];
1487   tile = tcd_tile;
1488
1489   time = clock();
1490
1491   fprintf(stdout, "tile decoding time %d/%d: ", tileno + 1,
1492           tcd_cp->tw * tcd_cp->th);
1493
1494         /*--------------TIER2------------------*/
1495
1496   l = t2_decode_packets(src, len, tcd_img, tcd_cp, tileno, tile);
1497
1498   if (l == -999) {
1499     eof = 1;
1500     fprintf(stderr, "tcd_decode: incomplete bistream\n");
1501   }
1502
1503         /*------------------TIER1-----------------*/
1504   t1_init_luts();
1505   t1_decode_cblks(tile, tcd_tcp);
1506
1507         /*----------------DWT---------------------*/
1508
1509   for (compno = 0; compno < tile->numcomps; compno++) {
1510     tcd_tilecomp_t *tilec = &tile->comps[compno];
1511     if (tcd_cp->reduce_on == 1) {
1512       tcd_img->comps[compno].resno_decoded =
1513         tile->comps[compno].numresolutions - tcd_cp->reduce_value - 1;
1514     }
1515
1516
1517     if (tcd_tcp->tccps[compno].qmfbid == 1) {
1518       dwt_decode(tilec->data, tilec->x1 - tilec->x0,
1519                  tilec->y1 - tilec->y0, tilec,
1520                  tilec->numresolutions - 1,
1521                  tilec->numresolutions - 1 -
1522                  tcd_img->comps[compno].resno_decoded);
1523     } else {
1524       dwt_decode_real(tilec->data, tilec->x1 - tilec->x0,
1525                       tilec->y1 - tilec->y0, tilec,
1526                       tilec->numresolutions - 1,
1527                       tilec->numresolutions - 1 -
1528                       tcd_img->comps[compno].resno_decoded);
1529     }
1530
1531     if (tile->comps[compno].numresolutions > 0)
1532       tcd_img->comps[compno].factor =
1533         tile->comps[compno].numresolutions -
1534         (tcd_img->comps[compno].resno_decoded + 1);
1535   }
1536
1537         /*----------------MCT-------------------*/
1538
1539   if (tcd_tcp->mct) {
1540     if (tcd_tcp->tccps[0].qmfbid == 1) {
1541       mct_decode(tile->comps[0].data, tile->comps[1].data,
1542                  tile->comps[2].data,
1543                  (tile->comps[0].x1 -
1544                   tile->comps[0].x0) * (tile->comps[0].y1 -
1545                                         tile->comps[0].y0));
1546     } else {
1547       mct_decode_real(tile->comps[0].data, tile->comps[1].data,
1548                       tile->comps[2].data,
1549                       (tile->comps[0].x1 -
1550                        tile->comps[0].x0) * (tile->comps[0].y1 -
1551                                              tile->comps[0].y0));
1552     }
1553   }
1554
1555         /*---------------TILE-------------------*/
1556
1557   for (compno = 0; compno < tile->numcomps; compno++) {
1558     tcd_tilecomp_t *tilec = &tile->comps[compno];
1559     tcd_resolution_t *res =
1560       &tilec->resolutions[tcd_img->comps[compno].resno_decoded];
1561     int adjust =
1562       tcd_img->comps[compno].sgnd ? 0 : 1 << (tcd_img->comps[compno].
1563                                               prec - 1);
1564     int min =
1565       tcd_img->comps[compno].
1566       sgnd ? -(1 << (tcd_img->comps[compno].prec - 1)) : 0;
1567     int max =
1568       tcd_img->comps[compno].
1569       sgnd ? (1 << (tcd_img->comps[compno].prec - 1)) -
1570       1 : (1 << tcd_img->comps[compno].prec) - 1;
1571
1572     int tw = tilec->x1 - tilec->x0;
1573     int w = tcd_img->comps[compno].w;
1574
1575     int i, j;
1576     int offset_x = int_ceildivpow2(tcd_img->comps[compno].x0,
1577                                    tcd_img->comps[compno].factor);
1578     int offset_y = int_ceildivpow2(tcd_img->comps[compno].y0,
1579                                    tcd_img->comps[compno].factor);
1580
1581     for (j = res->y0; j < res->y1; j++) {
1582       for (i = res->x0; i < res->x1; i++) {
1583
1584         int v;
1585
1586         double tmp =
1587           (double) tilec->data[i - res->x0 + (j - res->y0) * tw];
1588         if (tcd_tcp->tccps[compno].qmfbid == 1) {
1589           v = (int) tmp;
1590         } else {
1591
1592           //v = (int) tmp >> 13;
1593
1594           //Mod antonin : multbug1
1595           v =
1596             (int) ((fabs(tmp / 8192.0) >=
1597                     floor(fabs(tmp / 8192.0)) +
1598                     0.5) ? fabs(tmp / 8192.0) + 1.0 : fabs(tmp / 8192.0));
1599
1600           v = (tmp < 0) ? -v : v;
1601
1602           //doM
1603         }
1604         v += adjust;
1605
1606         tcd_img->comps[compno].data[(i - offset_x) +
1607                                     (j - offset_y) * w] =
1608           int_clamp(v, min, max);
1609       }
1610     }
1611   }
1612
1613   time = clock() - time;
1614   fprintf(stdout, "total:     %ld.%.3ld s\n", time / CLOCKS_PER_SEC,
1615           (time % CLOCKS_PER_SEC) * 1000 / CLOCKS_PER_SEC);
1616
1617
1618
1619   for (compno = 0; compno < tile->numcomps; compno++) {
1620     free(tcd_image.tiles[tileno].comps[compno].data);
1621   }
1622
1623   if (eof) {
1624     longjmp(j2k_error, 1);
1625   }
1626
1627   return l;
1628 }
1629
1630
1631
1632 void tcd_dec_release()
1633
1634 {
1635
1636   int tileno,compno,resno,bandno,precno;
1637
1638   for (tileno=0;tileno<tcd_image.tw*tcd_image.th;tileno++) {
1639
1640     tcd_tile_t tile=tcd_image.tiles[tileno];
1641
1642     for (compno=0;compno<tile.numcomps;compno++) {
1643
1644       tcd_tilecomp_t tilec=tile.comps[compno];
1645
1646       for (resno=0;resno<tilec.numresolutions;resno++) {
1647
1648         tcd_resolution_t res=tilec.resolutions[resno];
1649
1650         for (bandno=0;bandno<res.numbands;bandno++) {
1651
1652           tcd_band_t band=res.bands[bandno];
1653
1654           for (precno=0;precno<res.ph*res.pw;precno++) {
1655
1656             tcd_precinct_t prec=band.precincts[precno];
1657
1658             if (prec.cblks!=NULL) free(prec.cblks);
1659
1660             if (prec.imsbtree!=NULL) free(prec.imsbtree);
1661
1662             if (prec.incltree!=NULL) free(prec.incltree);
1663
1664           }
1665
1666           if (band.precincts!=NULL) free(band.precincts);
1667
1668         }
1669
1670       }
1671
1672       if (tilec.resolutions!=NULL) free(tilec.resolutions);
1673
1674     }
1675
1676     if (tile.comps!=NULL) free(tile.comps);
1677
1678   }
1679
1680   if (tcd_image.tiles!=NULL) free(tcd_image.tiles);
1681
1682 }