src/gdcmopenjpeg/libopenjpeg/tcd.c

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