src/lib/openjp3d/pi.c

   1 /*
   2  * The copyright in this software is being made available under the 2-clauses
   3  * BSD License, included below. This software may be subject to other third
   4  * party and contributor rights, including patent rights, and no such rights
   5  * are granted under this license.
   6  *
   7  * Copyright (c) 2001-2003, David Janssens
   8  * Copyright (c) 2002-2003, Yannick Verschueren
   9  * Copyright (c) 2003-2005, Francois Devaux and Antonin Descampe
  10  * Copyright (c) 2005, Herve Drolon, FreeImage Team
  11  * Copyright (c) 2002-2005, Communications and remote sensing Laboratory, Universite catholique de Louvain, Belgium
  12  * Copyright (c) 2006, M�nica D�ez, LPI-UVA, Spain
  13  * All rights reserved.
  14  *
  15  * Redistribution and use in source and binary forms, with or without
  16  * modification, are permitted provided that the following conditions
  17  * are met:
  18  * 1. Redistributions of source code must retain the above copyright
  19  *    notice, this list of conditions and the following disclaimer.
  20  * 2. Redistributions in binary form must reproduce the above copyright
  21  *    notice, this list of conditions and the following disclaimer in the
  22  *    documentation and/or other materials provided with the distribution.
  23  *
  24  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS `AS IS'
  25  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  26  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
  27  * ARE DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
  28  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
  29  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
  30  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
  31  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
  32  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
  33  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
  34  * POSSIBILITY OF SUCH DAMAGE.
  35  */
  36
  37 #include "opj_includes.h"
  38
  39 /** @defgroup PI PI - Implementation of a packet iterator */
  40 /*@{*/
  41
  42 /** @name Funciones locales */
  43 /*@{*/
  44
  45 /**
  46 Get next packet in layer-resolution-component-precinct order.
  47 @param pi packet iterator to modify
  48 @return returns false if pi pointed to the last packet or else returns true
  49 */
  50 static bool pi_next_lrcp(opj_pi_iterator_t * pi);
  51 /**
  52 Get next packet in resolution-layer-component-precinct order.
  53 @param pi packet iterator to modify
  54 @return returns false if pi pointed to the last packet or else returns true
  55 */
  56 static bool pi_next_rlcp(opj_pi_iterator_t * pi);
  57 /**
  58 Get next packet in resolution-precinct-component-layer order.
  59 @param pi packet iterator to modify
  60 @return returns false if pi pointed to the last packet or else returns true
  61 */
  62 static bool pi_next_rpcl(opj_pi_iterator_t * pi);
  63 /**
  64 Get next packet in precinct-component-resolution-layer order.
  65 @param pi packet iterator to modify
  66 @return returns false if pi pointed to the last packet or else returns true
  67 */
  68 static bool pi_next_pcrl(opj_pi_iterator_t * pi);
  69 /**
  70 Get next packet in component-precinct-resolution-layer order.
  71 @param pi packet iterator to modify
  72 @return returns false if pi pointed to the last packet or else returns true
  73 */
  74 static bool pi_next_cprl(opj_pi_iterator_t * pi);
  75
  76 /*@}*/
  77
  78 /*@}*/
  79
  80 /*
  81 ==========================================================
  82    local functions
  83 ==========================================================
  84 */
  85
  86 static bool pi_next_lrcp(opj_pi_iterator_t * pi) {
  87         opj_pi_comp_t *comp = NULL;
  88         opj_pi_resolution_t *res = NULL;
  89         long index = 0;
  90
  91         if (!pi->first) {
  92                 comp = &pi->comps[pi->compno];
  93                 res = &comp->resolutions[pi->resno];
  94                 goto LABEL_SKIP;
  95         } else {
  96                 pi->first = 0;
  97         }
  98
  99         for (pi->layno = 0; pi->layno < pi->poc.layno1; pi->layno++) {
 100                 for (pi->resno = pi->poc.resno0; pi->resno < pi->poc.resno1; pi->resno++) {
 101                         for (pi->compno = pi->poc.compno0; pi->compno < pi->poc.compno1; pi->compno++) {
 102                                 comp = &pi->comps[pi->compno];
 103                                 if (pi->resno >= comp->numresolution[0]) {
 104                                         continue;
 105                                 }
 106                                 res = &comp->resolutions[pi->resno];
 107                                 /*for (pi->precno = 0; pi->precno < (res->prctno[0] * res->prctno[1]); pi->precno++) {*/
 108                                 for (pi->precno = 0; pi->precno < (res->prctno[0] * res->prctno[1] * res->prctno[2]); pi->precno++) {
 109                                         index = pi->layno * pi->step_l
 110                                                 + pi->resno * pi->step_r
 111                                                 + pi->compno * pi->step_c
 112                                                 + pi->precno * pi->step_p;
 113                                         if (!pi->include[index]) {
 114                                                 pi->include[index] = 1;
 115                                                 return true;
 116                                         }
 117 LABEL_SKIP:;
 118
 119                                 }
 120                         }
 121                 }
 122         }
 123
 124         return false;
 125 }
 126
 127 static bool pi_next_rlcp(opj_pi_iterator_t * pi) {
 128         opj_pi_comp_t *comp = NULL;
 129         opj_pi_resolution_t *res = NULL;
 130         long index = 0;
 131
 132         if (!pi->first) {
 133                 comp = &pi->comps[pi->compno];
 134                 res = &comp->resolutions[pi->resno];
 135                 goto LABEL_SKIP;
 136         } else {
 137                 pi->first = 0;
 138         }
 139
 140         for (pi->resno = pi->poc.resno0; pi->resno < pi->poc.resno1; pi->resno++) {
 141                 for (pi->layno = 0; pi->layno < pi->poc.layno1; pi->layno++) {
 142                         for (pi->compno = pi->poc.compno0; pi->compno < pi->poc.compno1; pi->compno++) {
 143                                 comp = &pi->comps[pi->compno];
 144                                 if (pi->resno >= comp->numresolution[0]) {
 145                                         continue;
 146                                 }
 147                                 res = &comp->resolutions[pi->resno];
 148                                 /*for (pi->precno = 0; pi->precno < (res->prctno[0] * res->prctno[1]); pi->precno++) {*/
 149                                 for (pi->precno = 0; pi->precno < (res->prctno[0] * res->prctno[1] * res->prctno[2]); pi->precno++) {
 150                                         index = pi->layno * pi->step_l + pi->resno * pi->step_r + pi->compno * pi->step_c + pi->precno * pi->step_p;
 151                                         if (!pi->include[index]) {
 152                                                 pi->include[index] = 1;
 153                                                 return true;
 154                                         }
 155 LABEL_SKIP:;
 156                                 }
 157                         }
 158                 }
 159         }
 160
 161         return false;
 162 }
 163
 164 static bool pi_next_rpcl(opj_pi_iterator_t * pi) {
 165         opj_pi_comp_t *comp = NULL;
 166         opj_pi_resolution_t *res = NULL;
 167         long index = 0;
 168
 169         if (!pi->first) {
 170                 goto LABEL_SKIP;
 171         } else {
 172                 int compno, resno;
 173                 pi->first = 0;
 174                 pi->dx = 0;
 175                 pi->dy = 0;
 176                 for (compno = 0; compno < pi->numcomps; compno++) {
 177                         comp = &pi->comps[compno];
 178                         for (resno = 0; resno < comp->numresolution[0]; resno++) {
 179                                 int dx, dy,dz;
 180                                 res = &comp->resolutions[resno];
 181                                 dx = comp->dx * (1 << (res->pdx + comp->numresolution[0] - 1 - resno));
 182                                 dy = comp->dy * (1 << (res->pdy + comp->numresolution[1] - 1 - resno));
 183                                 dz = comp->dz * (1 << (res->pdz + comp->numresolution[2] - 1 - resno));
 184                                 pi->dx = !pi->dx ? dx : int_min(pi->dx, dx);
 185                                 pi->dy = !pi->dy ? dy : int_min(pi->dy, dy);
 186                                 pi->dz = !pi->dz ? dz : int_min(pi->dz, dz);
 187                         }
 188                 }
 189         }
 190
 191         for (pi->resno = pi->poc.resno0; pi->resno < pi->poc.resno1; pi->resno++) {
 192                 for (pi->z = pi->tz0; pi->z < pi->tz1; pi->z += pi->dz - (pi->z % pi->dz)) {
 193                         for (pi->y = pi->ty0; pi->y < pi->ty1; pi->y += pi->dy - (pi->y % pi->dy)) {
 194                                 for (pi->x = pi->tx0; pi->x < pi->tx1; pi->x += pi->dx - (pi->x % pi->dx)) {
 195                                         for (pi->compno = pi->poc.compno0; pi->compno < pi->poc.compno1; pi->compno++) {
 196                                                 int levelnox, levelnoy, levelnoz;
 197                                                 int trx0, try0, trz0;
 198                                                 int trx1, try1, trz1;
 199                                                 int rpx, rpy, rpz;
 200                                                 int prci, prcj, prck;
 201                                                 comp = &pi->comps[pi->compno];
 202                                                 if (pi->resno >= comp->numresolution[0]) {
 203                                                         continue;
 204                                                 }
 205                                                 res = &comp->resolutions[pi->resno];
 206                                                 levelnox = comp->numresolution[0] - 1 - pi->resno;
 207                                                 levelnoy = comp->numresolution[1] - 1 - pi->resno;
 208                                                 levelnoz = comp->numresolution[2] - 1 - pi->resno;
 209                                                 trx0 = int_ceildiv(pi->tx0, comp->dx << levelnox);
 210                                                 try0 = int_ceildiv(pi->ty0, comp->dy << levelnoy);
 211                                                 trz0 = int_ceildiv(pi->tz0, comp->dz << levelnoz);
 212                                                 trx1 = int_ceildiv(pi->tx1, comp->dx << levelnox);
 213                                                 try1 = int_ceildiv(pi->ty1, comp->dy << levelnoy);
 214                                                 trz1 = int_ceildiv(pi->tz1, comp->dz << levelnoz);
 215                                                 rpx = res->pdx + levelnox;
 216                                                 rpy = res->pdy + levelnoy;
 217                                                 rpz = res->pdz + levelnoz;
 218                                                 if ((!(pi->x % (comp->dx << rpx) == 0) || (pi->x == pi->tx0 && (trx0 << levelnox) % (1 << rpx)))) {
 219                                                         continue;
 220                                                 }
 221                                                 if ((!(pi->y % (comp->dy << rpy) == 0) || (pi->y == pi->ty0 && (try0 << levelnoy) % (1 << rpx)))) {
 222                                                         continue;
 223                                                 }
 224                                                 if ((!(pi->z % (comp->dz << rpz) == 0) || (pi->z == pi->tz0 && (trz0 << levelnoz) % (1 << rpx)))) {
 225                                                         continue;
 226                                                 }
 227                                                 if ((res->prctno[0]==0)||(res->prctno[1]==0)||(res->prctno[2]==0)) continue;
 228
 229                                                 if ((trx0==trx1)||(try0==try1)||(trz0==trz1)) continue;
 230
 231                                                 prci = int_floordivpow2(int_ceildiv(pi->x, comp->dx << levelnox), res->pdx)
 232                                                         - int_floordivpow2(trx0, res->pdx);
 233                                                 prcj = int_floordivpow2(int_ceildiv(pi->y, comp->dy << levelnoy), res->pdy)
 234                                                         - int_floordivpow2(try0, res->pdy);
 235                                                 prck = int_floordivpow2(int_ceildiv(pi->z, comp->dz << levelnoz), res->pdz)
 236                                                         - int_floordivpow2(trz0, res->pdz);
 237                                                 pi->precno = prci + prcj * res->prctno[0] + prck * res->prctno[0] * res->prctno[1];
 238                                                 for (pi->layno = 0; pi->layno < pi->poc.layno1; pi->layno++) {
 239                                                         index = pi->layno * pi->step_l + pi->resno * pi->step_r + pi->compno * pi->step_c + pi->precno * pi->step_p;
 240                                                         if (!pi->include[index]) {
 241                                                                 pi->include[index] = 1;
 242                                                                 return true;
 243                                                         }
 244         LABEL_SKIP:;
 245                                                 }
 246                                         }
 247                                 }
 248                         }
 249                 }
 250         }
 251
 252         return false;
 253 }
 254
 255 static bool pi_next_pcrl(opj_pi_iterator_t * pi) {
 256         opj_pi_comp_t *comp = NULL;
 257         opj_pi_resolution_t *res = NULL;
 258         long index = 0;
 259
 260         if (!pi->first) {
 261                 comp = &pi->comps[pi->compno];
 262                 goto LABEL_SKIP;
 263         } else {
 264                 int compno, resno;
 265                 pi->first = 0;
 266                 pi->dx = 0;
 267                 pi->dy = 0;
 268                 pi->dz = 0;
 269                 for (compno = 0; compno < pi->numcomps; compno++) {
 270                         comp = &pi->comps[compno];
 271                         for (resno = 0; resno < comp->numresolution[0]; resno++) {
 272                                 int dx, dy, dz;
 273                                 res = &comp->resolutions[resno];
 274                                 dx = comp->dx * (1 << (res->pdx + comp->numresolution[0] - 1 - resno));
 275                                 dy = comp->dy * (1 << (res->pdy + comp->numresolution[1] - 1 - resno));
 276                                 dz = comp->dz * (1 << (res->pdy + comp->numresolution[2] - 1 - resno));
 277                                 pi->dx = !pi->dx ? dx : int_min(pi->dx, dx);
 278                                 pi->dy = !pi->dy ? dy : int_min(pi->dy, dy);
 279                                 pi->dz = !pi->dz ? dz : int_min(pi->dz, dz);
 280                         }
 281                 }
 282         }
 283
 284 for (pi->z = pi->tz0; pi->z < pi->tz1; pi->z += pi->dz - (pi->z % pi->dz)) {
 285         for (pi->y = pi->ty0; pi->y < pi->ty1; pi->y += pi->dy - (pi->y % pi->dy)) {
 286                 for (pi->x = pi->tx0; pi->x < pi->tx1; pi->x += pi->dx - (pi->x % pi->dx)) {
 287                         for (pi->compno = pi->poc.compno0; pi->compno < pi->poc.compno1; pi->compno++) {
 288                                 comp = &pi->comps[pi->compno];
 289                                 for (pi->resno = pi->poc.resno0; pi->resno < int_min(pi->poc.resno1, comp->numresolution[0]); pi->resno++) {
 290                                                 int levelnox, levelnoy, levelnoz;
 291                                                 int trx0, try0, trz0;
 292                                                 int trx1, try1, trz1;
 293                                                 int rpx, rpy, rpz;
 294                                                 int prci, prcj, prck;
 295                                                 comp = &pi->comps[pi->compno];
 296                                                 if (pi->resno >= comp->numresolution[0]) {
 297                                                         continue;
 298                                                 }
 299                                                 res = &comp->resolutions[pi->resno];
 300                                                 levelnox = comp->numresolution[0] - 1 - pi->resno;
 301                                                 levelnoy = comp->numresolution[1] - 1 - pi->resno;
 302                                                 levelnoz = comp->numresolution[2] - 1 - pi->resno;
 303                                                 trx0 = int_ceildiv(pi->tx0, comp->dx << levelnox);
 304                                                 try0 = int_ceildiv(pi->ty0, comp->dy << levelnoy);
 305                                                 trz0 = int_ceildiv(pi->tz0, comp->dz << levelnoz);
 306                                                 trx1 = int_ceildiv(pi->tx1, comp->dx << levelnox);
 307                                                 try1 = int_ceildiv(pi->ty1, comp->dy << levelnoy);
 308                                                 trz1 = int_ceildiv(pi->tz1, comp->dz << levelnoz);
 309                                                 rpx = res->pdx + levelnox;
 310                                                 rpy = res->pdy + levelnoy;
 311                                                 rpz = res->pdz + levelnoz;
 312                                                 if ((!(pi->x % (comp->dx << rpx) == 0) || (pi->x == pi->tx0 && (trx0 << levelnox) % (1 << rpx)))) {
 313                                                         continue;
 314                                                 }
 315                                                 if ((!(pi->y % (comp->dy << rpy) == 0) || (pi->y == pi->ty0 && (try0 << levelnoy) % (1 << rpx)))) {
 316                                                         continue;
 317                                                 }
 318                                                 if ((!(pi->z % (comp->dz << rpz) == 0) || (pi->z == pi->tz0 && (trz0 << levelnoz) % (1 << rpx)))) {
 319                                                         continue;
 320                                                 }
 321                                                 if ((res->prctno[0]==0)||(res->prctno[1]==0)||(res->prctno[2]==0)) continue;
 322
 323                                                 if ((trx0==trx1)||(try0==try1)||(trz0==trz1)) continue;
 324
 325                                                 prci = int_floordivpow2(int_ceildiv(pi->x, comp->dx << levelnox), res->pdx)
 326                                                         - int_floordivpow2(trx0, res->pdx);
 327                                                 prcj = int_floordivpow2(int_ceildiv(pi->y, comp->dy << levelnoy), res->pdy)
 328                                                         - int_floordivpow2(try0, res->pdy);
 329                                                 prck = int_floordivpow2(int_ceildiv(pi->z, comp->dz << levelnoz), res->pdz)
 330                                                         - int_floordivpow2(trz0, res->pdz);
 331                                                 pi->precno = prci + prcj * res->prctno[0] + prck * res->prctno[0] * res->prctno[1];
 332                                                 for (pi->layno = 0; pi->layno < pi->poc.layno1; pi->layno++) {
 333                                                         index = pi->layno * pi->step_l + pi->resno * pi->step_r + pi->compno * pi->step_c + pi->precno * pi->step_p;
 334                                                         if (!pi->include[index]) {
 335                                                                 pi->include[index] = 1;
 336                                                                 return true;
 337                                                         }
 338 LABEL_SKIP:;
 339                                         }
 340                                 }
 341                         }
 342                 }
 343         }
 344 }
 345
 346         return false;
 347 }
 348
 349 static bool pi_next_cprl(opj_pi_iterator_t * pi) {
 350         opj_pi_comp_t *comp = NULL;
 351         opj_pi_resolution_t *res = NULL;
 352         long index = 0;
 353
 354         if (!pi->first) {
 355                 comp = &pi->comps[pi->compno];
 356                 goto LABEL_SKIP;
 357         } else {
 358                 pi->first = 0;
 359         }
 360
 361         for (pi->compno = pi->poc.compno0; pi->compno < pi->poc.compno1; pi->compno++) {
 362                 int resno;
 363                 comp = &pi->comps[pi->compno];
 364                 pi->dx = 0;
 365                 pi->dy = 0;
 366                 pi->dz = 0;
 367                 for (resno = 0; resno < comp->numresolution[0]; resno++) {
 368                         int dx, dy, dz;
 369                         res = &comp->resolutions[resno];
 370                         dx = comp->dx * (1 << (res->pdx + comp->numresolution[0] - 1 - resno));
 371                         dy = comp->dy * (1 << (res->pdy + comp->numresolution[1] - 1 - resno));
 372                         dz = comp->dz * (1 << (res->pdz + comp->numresolution[2] - 1 - resno));
 373                         pi->dx = !pi->dx ? dx : int_min(pi->dx, dx);
 374                         pi->dy = !pi->dy ? dy : int_min(pi->dy, dy);
 375                         pi->dz = !pi->dz ? dz : int_min(pi->dz, dz);
 376                 }
 377         for (pi->z = pi->tz0; pi->z < pi->tz1; pi->z += pi->dz - (pi->z % pi->dz)) {
 378                 for (pi->y = pi->ty0; pi->y < pi->ty1; pi->y += pi->dy - (pi->y % pi->dy)) {
 379                         for (pi->x = pi->tx0; pi->x < pi->tx1; pi->x += pi->dx - (pi->x % pi->dx)) {
 380                                 for (pi->resno = pi->poc.resno0; pi->resno < int_min(pi->poc.resno1, comp->numresolution[0]); pi->resno++) {
 381                                                 int levelnox, levelnoy, levelnoz;
 382                                                 int trx0, try0, trz0;
 383                                                 int trx1, try1, trz1;
 384                                                 int rpx, rpy, rpz;
 385                                                 int prci, prcj, prck;
 386                                                 comp = &pi->comps[pi->compno];
 387                                                 if (pi->resno >= comp->numresolution[0]) {
 388                                                         continue;
 389                                                 }
 390                                                 res = &comp->resolutions[pi->resno];
 391                                                 levelnox = comp->numresolution[0] - 1 - pi->resno;
 392                                                 levelnoy = comp->numresolution[1] - 1 - pi->resno;
 393                                                 levelnoz = comp->numresolution[2] - 1 - pi->resno;
 394                                                 trx0 = int_ceildiv(pi->tx0, comp->dx << levelnox);
 395                                                 try0 = int_ceildiv(pi->ty0, comp->dy << levelnoy);
 396                                                 trz0 = int_ceildiv(pi->tz0, comp->dz << levelnoz);
 397                                                 trx1 = int_ceildiv(pi->tx1, comp->dx << levelnox);
 398                                                 try1 = int_ceildiv(pi->ty1, comp->dy << levelnoy);
 399                                                 trz1 = int_ceildiv(pi->tz1, comp->dz << levelnoz);
 400                                                 rpx = res->pdx + levelnox;
 401                                                 rpy = res->pdy + levelnoy;
 402                                                 rpz = res->pdz + levelnoz;
 403                                                 if ((!(pi->x % (comp->dx << rpx) == 0) || (pi->x == pi->tx0 && (trx0 << levelnox) % (1 << rpx)))) {
 404                                                         continue;
 405                                                 }
 406                                                 if ((!(pi->y % (comp->dy << rpy) == 0) || (pi->y == pi->ty0 && (try0 << levelnoy) % (1 << rpx)))) {
 407                                                         continue;
 408                                                 }
 409                                                 if ((!(pi->z % (comp->dz << rpz) == 0) || (pi->z == pi->tz0 && (trz0 << levelnoz) % (1 << rpx)))) {
 410                                                         continue;
 411                                                 }
 412                                                 if ((res->prctno[0]==0)||(res->prctno[1]==0)||(res->prctno[2]==0)) continue;
 413
 414                                                 if ((trx0==trx1)||(try0==try1)||(trz0==trz1)) continue;
 415
 416                                                 prci = int_floordivpow2(int_ceildiv(pi->x, comp->dx << levelnox), res->pdx)
 417                                                         - int_floordivpow2(trx0, res->pdx);
 418                                                 prcj = int_floordivpow2(int_ceildiv(pi->y, comp->dy << levelnoy), res->pdy)
 419                                                         - int_floordivpow2(try0, res->pdy);
 420                                                 prck = int_floordivpow2(int_ceildiv(pi->z, comp->dz << levelnoz), res->pdz)
 421                                                         - int_floordivpow2(trz0, res->pdz);
 422                                                 pi->precno = prci + prcj * res->prctno[0] + prck * res->prctno[0] * res->prctno[1];
 423                                                 for (pi->layno = 0; pi->layno < pi->poc.layno1; pi->layno++) {
 424                                                         index = pi->layno * pi->step_l + pi->resno * pi->step_r + pi->compno * pi->step_c + pi->precno * pi->step_p;
 425                                                         if (!pi->include[index]) {
 426                                                                 pi->include[index] = 1;
 427                                                                 return true;
 428                                                         }
 429         LABEL_SKIP:;
 430                                                 }
 431                                         }
 432                                 }
 433                         }
 434                 }
 435         }
 436
 437         return false;
 438 }
 439
 440 /*
 441 ==========================================================
 442    Packet iterator interface
 443 ==========================================================
 444 */
 445
 446 opj_pi_iterator_t *pi_create(opj_volume_t *volume, opj_cp_t *cp, int tileno) {
 447         int p, q, r;
 448         int compno, resno, pino;
 449         opj_pi_iterator_t *pi = NULL;
 450         opj_tcp_t *tcp = NULL;
 451         opj_tccp_t *tccp = NULL;
 452         size_t array_size;
 453
 454         tcp = &cp->tcps[tileno];
 455
 456         array_size = (tcp->numpocs + 1) * sizeof(opj_pi_iterator_t);
 457         pi = (opj_pi_iterator_t *) opj_malloc(array_size);
 458         if(!pi) {
 459                 fprintf(stdout,"[ERROR] Malloc of opj_pi_iterator failed \n");
 460                 return NULL;
 461         }
 462
 463         for (pino = 0; pino < tcp->numpocs + 1; pino++) {       /* change */
 464                 int maxres = 0;
 465                 int maxprec = 0;
 466                 p = tileno % cp->tw;
 467                 q = tileno / cp->tw;
 468                 r = tileno / (cp->tw * cp->th);
 469
 470                 pi[pino].tx0 = int_max(cp->tx0 + p * cp->tdx, volume->x0);
 471                 pi[pino].ty0 = int_max(cp->ty0 + q * cp->tdy, volume->y0);
 472                 pi[pino].tz0 = int_max(cp->tz0 + r * cp->tdz, volume->z0);
 473                 pi[pino].tx1 = int_min(cp->tx0 + (p + 1) * cp->tdx, volume->x1);
 474                 pi[pino].ty1 = int_min(cp->ty0 + (q + 1) * cp->tdy, volume->y1);
 475                 pi[pino].tz1 = int_min(cp->tz0 + (r + 1) * cp->tdz, volume->z1);
 476                 pi[pino].numcomps = volume->numcomps;
 477
 478                 array_size = volume->numcomps * sizeof(opj_pi_comp_t);
 479                 pi[pino].comps = (opj_pi_comp_t *) opj_malloc(array_size);
 480                 if(!pi[pino].comps) {
 481                         fprintf(stdout,"[ERROR] Malloc of opj_pi_comp failed \n");
 482                         pi_destroy(pi, cp, tileno);
 483                         return NULL;
 484                 }
 485                 memset(pi[pino].comps, 0, array_size);
 486
 487                 for (compno = 0; compno < pi->numcomps; compno++) {
 488                         int tcx0, tcx1, tcy0, tcy1, tcz0, tcz1;
 489                         int i;
 490                         opj_pi_comp_t *comp = &pi[pino].comps[compno];
 491                         tccp = &tcp->tccps[compno];
 492
 493                         comp->dx = volume->comps[compno].dx;
 494                         comp->dy = volume->comps[compno].dy;
 495                         comp->dz = volume->comps[compno].dz;
 496                         for (i = 0; i < 3; i++) {
 497                                 comp->numresolution[i] = tccp->numresolution[i];
 498                                 if (comp->numresolution[i] > maxres) {
 499                                         maxres = comp->numresolution[i];
 500                                 }
 501                         }
 502                         array_size = comp->numresolution[0] * sizeof(opj_pi_resolution_t);
 503                         comp->resolutions =     (opj_pi_resolution_t *) opj_malloc(array_size);
 504                         if(!comp->resolutions) {
 505                                 fprintf(stdout,"[ERROR] Malloc of opj_pi_resolution failed \n");
 506                                 pi_destroy(pi, cp, tileno);
 507                                 return NULL;
 508                         }
 509
 510                         tcx0 = int_ceildiv(pi->tx0, comp->dx);
 511                         tcy0 = int_ceildiv(pi->ty0, comp->dy);
 512                         tcz0 = int_ceildiv(pi->tz0, comp->dz);
 513                         tcx1 = int_ceildiv(pi->tx1, comp->dx);
 514                         tcy1 = int_ceildiv(pi->ty1, comp->dy);
 515                         tcz1 = int_ceildiv(pi->tz1, comp->dz);
 516
 517                         for (resno = 0; resno < comp->numresolution[0]; resno++) {
 518                                 int levelnox, levelnoy, levelnoz, diff;
 519                                 int rx0, ry0, rz0, rx1, ry1, rz1;
 520                                 int px0, py0, pz0, px1, py1, pz1;
 521                                 opj_pi_resolution_t *res = &comp->resolutions[resno];
 522                                 if (tccp->csty & J3D_CCP_CSTY_PRT) {
 523                                         res->pdx = tccp->prctsiz[0][resno];
 524                                         res->pdy = tccp->prctsiz[1][resno];
 525                                         res->pdz = tccp->prctsiz[2][resno];
 526                                 } else {
 527                                         res->pdx = 15;
 528                                         res->pdy = 15;
 529                                         res->pdz = 15;
 530                                 }
 531                                 levelnox = comp->numresolution[0] - 1 - resno;
 532                                 levelnoy = comp->numresolution[1] - 1 - resno;
 533         levelnoz = comp->numresolution[2] - 1 - resno;
 534                                 if (levelnoz < 0) levelnoz = 0;
 535                                 diff = comp->numresolution[0] - comp->numresolution[2];
 536
 537                                 rx0 = int_ceildivpow2(tcx0, levelnox);
 538                                 ry0 = int_ceildivpow2(tcy0, levelnoy);
 539                                 rz0 = int_ceildivpow2(tcz0, levelnoz);
 540                                 rx1 = int_ceildivpow2(tcx1, levelnox);
 541                                 ry1 = int_ceildivpow2(tcy1, levelnoy);
 542                                 rz1 = int_ceildivpow2(tcz1, levelnoz);
 543                                 px0 = int_floordivpow2(rx0, res->pdx) << res->pdx;
 544                                 py0 = int_floordivpow2(ry0, res->pdy) << res->pdy;
 545                                 pz0 = int_floordivpow2(rz0, res->pdz) << res->pdz;
 546                                 px1 = int_ceildivpow2(rx1, res->pdx) << res->pdx;
 547                                 py1 = int_ceildivpow2(ry1, res->pdy) << res->pdy;
 548                                 pz1 = int_ceildivpow2(rz1, res->pdz) << res->pdz;
 549                                 res->prctno[0] = (rx0==rx1)? 0 : ((px1 - px0) >> res->pdx);
 550                                 res->prctno[1] = (ry0==ry1)? 0 : ((py1 - py0) >> res->pdy);
 551                                 res->prctno[2] = (rz0==rz1)? 0 : ((pz1 - pz0) >> res->pdz);
 552
 553                                 if (res->prctno[0]*res->prctno[1]*res->prctno[2] > maxprec) {
 554                                         maxprec = res->prctno[0]*res->prctno[1]*res->prctno[2];
 555                                 }
 556                         }
 557                 }
 558
 559                 tccp = &tcp->tccps[0];
 560                 pi[pino].step_p = 1;
 561                 pi[pino].step_c = maxprec * pi[pino].step_p;
 562                 pi[pino].step_r = volume->numcomps * pi[pino].step_c;
 563                 pi[pino].step_l = maxres * pi[pino].step_r;
 564
 565                 if (pino == 0) {
 566                         array_size = volume->numcomps * maxres * tcp->numlayers * maxprec * sizeof(short int);
 567                         pi[pino].include = (short int *) opj_malloc(array_size);
 568                         if(!pi[pino].include) {
 569                                 fprintf(stdout,"[ERROR] Malloc of pi[pino].include failed \n");
 570                                 pi_destroy(pi, cp, tileno);
 571                                 return NULL;
 572                         }
 573                 }
 574                 else {
 575                         pi[pino].include = pi[pino - 1].include;
 576                 }
 577
 578                 if (tcp->POC == 0) {
 579                         pi[pino].first = 1;
 580                         pi[pino].poc.resno0 = 0;
 581                         pi[pino].poc.compno0 = 0;
 582                         pi[pino].poc.layno1 = tcp->numlayers;
 583                         pi[pino].poc.resno1 = maxres;
 584                         pi[pino].poc.compno1 = volume->numcomps;
 585                         pi[pino].poc.prg = tcp->prg;
 586                 } else {
 587                         pi[pino].first = 1;
 588                         pi[pino].poc.resno0 = tcp->pocs[pino].resno0;
 589                         pi[pino].poc.compno0 = tcp->pocs[pino].compno0;
 590                         pi[pino].poc.layno1 = tcp->pocs[pino].layno1;
 591                         pi[pino].poc.resno1 = tcp->pocs[pino].resno1;
 592                         pi[pino].poc.compno1 = tcp->pocs[pino].compno1;
 593                         pi[pino].poc.prg = tcp->pocs[pino].prg;
 594                 }
 595         }
 596
 597         return pi;
 598 }
 599
 600 void pi_destroy(opj_pi_iterator_t *pi, opj_cp_t *cp, int tileno) {
 601         int compno, pino;
 602         opj_tcp_t *tcp = &cp->tcps[tileno];
 603         if(pi) {
 604                 for (pino = 0; pino < tcp->numpocs + 1; pino++) {
 605                         if(pi[pino].comps) {
 606                                 for (compno = 0; compno < pi->numcomps; compno++) {
 607                                         opj_pi_comp_t *comp = &pi[pino].comps[compno];
 608                                         if(comp->resolutions) {
 609                                                 opj_free(comp->resolutions);
 610                                         }
 611                                 }
 612                                 opj_free(pi[pino].comps);
 613                         }
 614                 }
 615                 if(pi->include) {
 616                         opj_free(pi->include);
 617                 }
 618                 opj_free(pi);
 619         }
 620 }
 621
 622 bool pi_next(opj_pi_iterator_t * pi) {
 623         switch (pi->poc.prg) {
 624                 case LRCP:
 625                         return pi_next_lrcp(pi);
 626                 case RLCP:
 627                         return pi_next_rlcp(pi);
 628                 case RPCL:
 629                         return pi_next_rpcl(pi);
 630                 case PCRL:
 631                         return pi_next_pcrl(pi);
 632                 case CPRL:
 633                         return pi_next_cprl(pi);
 634         }
 635
 636         return false;
 637 }
 638