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