Subtile decoding: only do 9x7 IDWT computations on relevant areas of tile-component...

diff --git a/src/lib/openjp2/dwt.c b/src/lib/openjp2/dwt.c
index e29c1dd6083ef007e8cb95b7b7fcf52527d055d6..96d54df301def8faf74157b5eb06e5f24c057dca 100644 (file)
--- a/src/lib/openjp2/dwt.c
+++ b/src/lib/openjp2/dwt.c
@@ -82,9 +82,9 @@
 
 typedef struct dwt_local {
     OPJ_INT32* mem;
-    OPJ_INT32 dn;
-    OPJ_INT32 sn;
-    OPJ_INT32 cas;
+    OPJ_INT32 dn;   /* number of elements in high pass band */
+    OPJ_INT32 sn;   /* number of elements in low pass band */
+    OPJ_INT32 cas;  /* 0 = start on even coord, 1 = start on odd coord */
 } opj_dwt_t;
 
 typedef union {
@@ -93,9 +93,13 @@ typedef union {
 
 typedef struct v4dwt_local {
     opj_v4_t*   wavelet ;
-    OPJ_INT32       dn ;
-    OPJ_INT32       sn ;
-    OPJ_INT32       cas ;
+    OPJ_INT32       dn ;  /* number of elements in high pass band */
+    OPJ_INT32       sn ;  /* number of elements in low pass band */
+    OPJ_INT32       cas ; /* 0 = start on even coord, 1 = start on odd coord */
+    OPJ_UINT32      win_l_x0; /* start coord in low pass band */
+    OPJ_UINT32      win_l_x1; /* end coord in low pass band */
+    OPJ_UINT32      win_h_x0; /* start coord in high pass band */
+    OPJ_UINT32      win_h_x1; /* end coord in high pass band */
 } opj_v4dwt_t ;
 
 static const OPJ_FLOAT32 opj_dwt_alpha =  1.586134342f; /*  12994 */
@@ -162,25 +166,38 @@ static OPJ_UINT32 opj_dwt_max_resolution(opj_tcd_resolution_t* OPJ_RESTRICT r,
 /* </summary>                            */
 static void opj_v4dwt_decode(opj_v4dwt_t* OPJ_RESTRICT dwt);
 
-static void opj_v4dwt_interleave_h(opj_v4dwt_t* OPJ_RESTRICT w,
-                                   OPJ_FLOAT32* OPJ_RESTRICT a, OPJ_INT32 x, OPJ_INT32 size);
+static void opj_v4dwt_interleave_h(opj_v4dwt_t* OPJ_RESTRICT dwt,
+                                   OPJ_FLOAT32* OPJ_RESTRICT a,
+                                   OPJ_UINT32 width,
+                                   OPJ_UINT32 remaining_height);
 
-static void opj_v4dwt_interleave_v(opj_v4dwt_t* OPJ_RESTRICT v,
-                                   OPJ_FLOAT32* OPJ_RESTRICT a, OPJ_INT32 x, OPJ_INT32 nb_elts_read);
+static void opj_v4dwt_interleave_v(opj_v4dwt_t* OPJ_RESTRICT dwt,
+                                   OPJ_FLOAT32* OPJ_RESTRICT a,
+                                   OPJ_UINT32 width,
+                                   OPJ_UINT32 nb_elts_read);
 
 #ifdef __SSE__
-static void opj_v4dwt_decode_step1_sse(opj_v4_t* w, OPJ_INT32 count,
+static void opj_v4dwt_decode_step1_sse(opj_v4_t* w,
+                                       OPJ_UINT32 start,
+                                       OPJ_UINT32 end,
                                        const __m128 c);
 
-static void opj_v4dwt_decode_step2_sse(opj_v4_t* l, opj_v4_t* w, OPJ_INT32 k,
-                                       OPJ_INT32 m, __m128 c);
+static void opj_v4dwt_decode_step2_sse(opj_v4_t* l, opj_v4_t* w,
+                                       OPJ_UINT32 start,
+                                       OPJ_UINT32 end,
+                                       OPJ_UINT32 m, __m128 c);
 
 #else
-static void opj_v4dwt_decode_step1(opj_v4_t* w, OPJ_INT32 count,
+static void opj_v4dwt_decode_step1(opj_v4_t* w,
+                                   OPJ_UINT32 start,
+                                   OPJ_UINT32 end,
                                    const OPJ_FLOAT32 c);
 
-static void opj_v4dwt_decode_step2(opj_v4_t* l, opj_v4_t* w, OPJ_INT32 k,
-                                   OPJ_INT32 m, OPJ_FLOAT32 c);
+static void opj_v4dwt_decode_step2(opj_v4_t* l, opj_v4_t* w,
+                                   OPJ_UINT32 start,
+                                   OPJ_UINT32 end,
+                                   OPJ_UINT32 m,
+                                   OPJ_FLOAT32 c);
 
 #endif
 
@@ -1839,136 +1856,143 @@ static OPJ_BOOL opj_dwt_decode_partial_tile(opj_tcd_t *tcd,
     return OPJ_TRUE;
 }
 
-static void opj_v4dwt_interleave_h(opj_v4dwt_t* OPJ_RESTRICT w,
-                                   OPJ_FLOAT32* OPJ_RESTRICT a, OPJ_INT32 x, OPJ_INT32 size)
+static void opj_v4dwt_interleave_h(opj_v4dwt_t* OPJ_RESTRICT dwt,
+                                   OPJ_FLOAT32* OPJ_RESTRICT a,
+                                   OPJ_UINT32 width,
+                                   OPJ_UINT32 remaining_height)
 {
-    OPJ_FLOAT32* OPJ_RESTRICT bi = (OPJ_FLOAT32*)(w->wavelet + w->cas);
-    OPJ_INT32 count = w->sn;
-    OPJ_INT32 i, k;
+    OPJ_FLOAT32* OPJ_RESTRICT bi = (OPJ_FLOAT32*)(dwt->wavelet + dwt->cas);
+    OPJ_UINT32 i, k;
+    OPJ_UINT32 x0 = dwt->win_l_x0;
+    OPJ_UINT32 x1 = dwt->win_l_x1;
 
     for (k = 0; k < 2; ++k) {
-        if (count + 3 * x < size && ((size_t) a & 0x0f) == 0 &&
-                ((size_t) bi & 0x0f) == 0 && (x & 0x0f) == 0) {
+        if (remaining_height >= 4 && ((size_t) a & 0x0f) == 0 &&
+                ((size_t) bi & 0x0f) == 0 && (width & 0x0f) == 0) {
             /* Fast code path */
-            for (i = 0; i < count; ++i) {
-                OPJ_INT32 j = i;
+            for (i = x0; i < x1; ++i) {
+                OPJ_UINT32 j = i;
                 bi[i * 8    ] = a[j];
-                j += x;
+                j += width;
                 bi[i * 8 + 1] = a[j];
-                j += x;
+                j += width;
                 bi[i * 8 + 2] = a[j];
-                j += x;
+                j += width;
                 bi[i * 8 + 3] = a[j];
             }
         } else {
             /* Slow code path */
-            for (i = 0; i < count; ++i) {
-                OPJ_INT32 j = i;
+            for (i = x0; i < x1; ++i) {
+                OPJ_UINT32 j = i;
                 bi[i * 8    ] = a[j];
-                j += x;
-                if (j >= size) {
+                j += width;
+                if (remaining_height == 1) {
                     continue;
                 }
                 bi[i * 8 + 1] = a[j];
-                j += x;
-                if (j >= size) {
+                j += width;
+                if (remaining_height == 2) {
                     continue;
                 }
                 bi[i * 8 + 2] = a[j];
-                j += x;
-                if (j >= size) {
+                j += width;
+                if (remaining_height == 3) {
                     continue;
                 }
                 bi[i * 8 + 3] = a[j]; /* This one*/
             }
         }
 
-        bi = (OPJ_FLOAT32*)(w->wavelet + 1 - w->cas);
-        a += w->sn;
-        size -= w->sn;
-        count = w->dn;
+        bi = (OPJ_FLOAT32*)(dwt->wavelet + 1 - dwt->cas);
+        a += dwt->sn;
+        x0 = dwt->win_h_x0;
+        x1 = dwt->win_h_x1;
     }
 }
 
-static void opj_v4dwt_interleave_v(opj_v4dwt_t* OPJ_RESTRICT v,
-                                   OPJ_FLOAT32* OPJ_RESTRICT a, OPJ_INT32 x, OPJ_INT32 nb_elts_read)
+static void opj_v4dwt_interleave_v(opj_v4dwt_t* OPJ_RESTRICT dwt,
+                                   OPJ_FLOAT32* OPJ_RESTRICT a,
+                                   OPJ_UINT32 width,
+                                   OPJ_UINT32 nb_elts_read)
 {
-    opj_v4_t* OPJ_RESTRICT bi = v->wavelet + v->cas;
-    OPJ_INT32 i;
+    opj_v4_t* OPJ_RESTRICT bi = dwt->wavelet + dwt->cas;
+    OPJ_UINT32 i;
 
-    for (i = 0; i < v->sn; ++i) {
-        memcpy(&bi[i * 2], &a[i * x], (size_t)nb_elts_read * sizeof(OPJ_FLOAT32));
+    for (i = dwt->win_l_x0; i < dwt->win_l_x1; ++i) {
+        memcpy(&bi[i * 2], &a[i * width], (size_t)nb_elts_read * sizeof(OPJ_FLOAT32));
     }
 
-    a += v->sn * x;
-    bi = v->wavelet + 1 - v->cas;
+    a += (OPJ_UINT32)dwt->sn * width;
+    bi = dwt->wavelet + 1 - dwt->cas;
 
-    for (i = 0; i < v->dn; ++i) {
-        memcpy(&bi[i * 2], &a[i * x], (size_t)nb_elts_read * sizeof(OPJ_FLOAT32));
+    for (i = dwt->win_h_x0; i < dwt->win_h_x1; ++i) {
+        memcpy(&bi[i * 2], &a[i * width], (size_t)nb_elts_read * sizeof(OPJ_FLOAT32));
     }
 }
 
 #ifdef __SSE__
 
-static void opj_v4dwt_decode_step1_sse(opj_v4_t* w, OPJ_INT32 count,
+static void opj_v4dwt_decode_step1_sse(opj_v4_t* w,
+                                       OPJ_UINT32 start,
+                                       OPJ_UINT32 end,
                                        const __m128 c)
 {
     __m128* OPJ_RESTRICT vw = (__m128*) w;
-    OPJ_INT32 i;
+    OPJ_UINT32 i;
     /* 4x unrolled loop */
-    for (i = 0; i < count >> 2; ++i) {
-        *vw = _mm_mul_ps(*vw, c);
-        vw += 2;
-        *vw = _mm_mul_ps(*vw, c);
-        vw += 2;
-        *vw = _mm_mul_ps(*vw, c);
-        vw += 2;
-        *vw = _mm_mul_ps(*vw, c);
-        vw += 2;
+    for (i = start; i + 3 < end; i += 4) {
+        vw[2 * i] = _mm_mul_ps(vw[2 * i], c);
+        vw[2 * i + 2] = _mm_mul_ps(vw[2 * i + 2], c);
+        vw[2 * i + 4] = _mm_mul_ps(vw[2 * i + 4], c);
+        vw[2 * i + 6] = _mm_mul_ps(vw[2 * i + 6], c);
     }
-    count &= 3;
-    for (i = 0; i < count; ++i) {
-        *vw = _mm_mul_ps(*vw, c);
-        vw += 2;
+    for (; i < end; ++i) {
+        vw[2 * i] = _mm_mul_ps(vw[2 * i], c);
     }
 }
 
-void opj_v4dwt_decode_step2_sse(opj_v4_t* l, opj_v4_t* w, OPJ_INT32 k,
-                                OPJ_INT32 m, __m128 c)
+static void opj_v4dwt_decode_step2_sse(opj_v4_t* l, opj_v4_t* w,
+                                       OPJ_UINT32 start,
+                                       OPJ_UINT32 end,
+                                       OPJ_UINT32 m,
+                                       __m128 c)
 {
     __m128* OPJ_RESTRICT vl = (__m128*) l;
     __m128* OPJ_RESTRICT vw = (__m128*) w;
-    OPJ_INT32 i;
+    OPJ_UINT32 i;
+    OPJ_UINT32 imax = opj_uint_min(end, m);
     __m128 tmp1, tmp2, tmp3;
-    tmp1 = vl[0];
-    for (i = 0; i < m; ++i) {
+    if (start == 0) {
+        tmp1 = vl[0];
+    } else {
+        vw += start * 2;
+        tmp1 = vw[-3];
+    }
+    for (i = start; i < imax; ++i) {
         tmp2 = vw[-1];
         tmp3 = vw[ 0];
         vw[-1] = _mm_add_ps(tmp2, _mm_mul_ps(_mm_add_ps(tmp1, tmp3), c));
         tmp1 = tmp3;
         vw += 2;
     }
-    vl = vw - 2;
-    if (m >= k) {
-        return;
-    }
-    c = _mm_add_ps(c, c);
-    c = _mm_mul_ps(c, vl[0]);
-    for (; m < k; ++m) {
-        __m128 tmp = vw[-1];
-        vw[-1] = _mm_add_ps(tmp, c);
-        vw += 2;
+    if (m < end) {
+        assert(m + 1 == end);
+        c = _mm_add_ps(c, c);
+        c = _mm_mul_ps(c, vw[-2]);
+        vw[-1] = _mm_add_ps(vw[-1], c);
     }
 }
 
 #else
 
-static void opj_v4dwt_decode_step1(opj_v4_t* w, OPJ_INT32 count,
+static void opj_v4dwt_decode_step1(opj_v4_t* w,
+                                   OPJ_UINT32 start,
+                                   OPJ_UINT32 end,
                                    const OPJ_FLOAT32 c)
 {
     OPJ_FLOAT32* OPJ_RESTRICT fw = (OPJ_FLOAT32*) w;
-    OPJ_INT32 i;
-    for (i = 0; i < count; ++i) {
+    OPJ_UINT32 i;
+    for (i = start; i < end; ++i) {
         OPJ_FLOAT32 tmp1 = fw[i * 8    ];
         OPJ_FLOAT32 tmp2 = fw[i * 8 + 1];
         OPJ_FLOAT32 tmp3 = fw[i * 8 + 2];
@@ -1980,13 +2004,21 @@ static void opj_v4dwt_decode_step1(opj_v4_t* w, OPJ_INT32 count,
     }
 }
 
-static void opj_v4dwt_decode_step2(opj_v4_t* l, opj_v4_t* w, OPJ_INT32 k,
-                                   OPJ_INT32 m, OPJ_FLOAT32 c)
+static void opj_v4dwt_decode_step2(opj_v4_t* l, opj_v4_t* w,
+                                   OPJ_UINT32 start,
+                                   OPJ_UINT32 end,
+                                   OPJ_UINT32 m,
+                                   OPJ_FLOAT32 c)
 {
     OPJ_FLOAT32* fl = (OPJ_FLOAT32*) l;
     OPJ_FLOAT32* fw = (OPJ_FLOAT32*) w;
-    OPJ_INT32 i;
-    for (i = 0; i < m; ++i) {
+    OPJ_UINT32 i;
+    OPJ_UINT32 imax = opj_uint_min(end, m);
+    if (start > 0) {
+        fw += 8 * start;
+        fl = fw - 8;
+    }
+    for (i = start; i < imax; ++i) {
         OPJ_FLOAT32 tmp1_1 = fl[0];
         OPJ_FLOAT32 tmp1_2 = fl[1];
         OPJ_FLOAT32 tmp1_3 = fl[2];
@@ -2006,27 +2038,13 @@ static void opj_v4dwt_decode_step2(opj_v4_t* l, opj_v4_t* w, OPJ_INT32 k,
         fl = fw;
         fw += 8;
     }
-    if (m < k) {
-        OPJ_FLOAT32 c1;
-        OPJ_FLOAT32 c2;
-        OPJ_FLOAT32 c3;
-        OPJ_FLOAT32 c4;
+    if (m < end) {
+        assert(m + 1 == end);
         c += c;
-        c1 = fl[0] * c;
-        c2 = fl[1] * c;
-        c3 = fl[2] * c;
-        c4 = fl[3] * c;
-        for (; m < k; ++m) {
-            OPJ_FLOAT32 tmp1 = fw[-4];
-            OPJ_FLOAT32 tmp2 = fw[-3];
-            OPJ_FLOAT32 tmp3 = fw[-2];
-            OPJ_FLOAT32 tmp4 = fw[-1];
-            fw[-4] = tmp1 + c1;
-            fw[-3] = tmp2 + c2;
-            fw[-2] = tmp3 + c3;
-            fw[-1] = tmp4 + c4;
-            fw += 8;
-        }
+        fw[-4] = fw[-4] + fl[0] * c;
+        fw[-3] = fw[-3] + fl[1] * c;
+        fw[-2] = fw[-2] + fl[2] * c;
+        fw[-1] = fw[-1] + fl[3] * c;
     }
 }
 
@@ -2052,27 +2070,47 @@ static void opj_v4dwt_decode(opj_v4dwt_t* OPJ_RESTRICT dwt)
         b = 0;
     }
 #ifdef __SSE__
-    opj_v4dwt_decode_step1_sse(dwt->wavelet + a, dwt->sn, _mm_set1_ps(opj_K));
-    opj_v4dwt_decode_step1_sse(dwt->wavelet + b, dwt->dn, _mm_set1_ps(opj_c13318));
-    opj_v4dwt_decode_step2_sse(dwt->wavelet + b, dwt->wavelet + a + 1, dwt->sn,
-                               opj_int_min(dwt->sn, dwt->dn - a), _mm_set1_ps(opj_dwt_delta));
-    opj_v4dwt_decode_step2_sse(dwt->wavelet + a, dwt->wavelet + b + 1, dwt->dn,
-                               opj_int_min(dwt->dn, dwt->sn - b), _mm_set1_ps(opj_dwt_gamma));
-    opj_v4dwt_decode_step2_sse(dwt->wavelet + b, dwt->wavelet + a + 1, dwt->sn,
-                               opj_int_min(dwt->sn, dwt->dn - a), _mm_set1_ps(opj_dwt_beta));
-    opj_v4dwt_decode_step2_sse(dwt->wavelet + a, dwt->wavelet + b + 1, dwt->dn,
-                               opj_int_min(dwt->dn, dwt->sn - b), _mm_set1_ps(opj_dwt_alpha));
+    opj_v4dwt_decode_step1_sse(dwt->wavelet + a, dwt->win_l_x0, dwt->win_l_x1,
+                               _mm_set1_ps(opj_K));
+    opj_v4dwt_decode_step1_sse(dwt->wavelet + b, dwt->win_h_x0, dwt->win_h_x1,
+                               _mm_set1_ps(opj_c13318));
+    opj_v4dwt_decode_step2_sse(dwt->wavelet + b, dwt->wavelet + a + 1,
+                               dwt->win_l_x0, dwt->win_l_x1,
+                               (OPJ_UINT32)opj_int_min(dwt->sn, dwt->dn - a),
+                               _mm_set1_ps(opj_dwt_delta));
+    opj_v4dwt_decode_step2_sse(dwt->wavelet + a, dwt->wavelet + b + 1,
+                               dwt->win_h_x0, dwt->win_h_x1,
+                               (OPJ_UINT32)opj_int_min(dwt->dn, dwt->sn - b),
+                               _mm_set1_ps(opj_dwt_gamma));
+    opj_v4dwt_decode_step2_sse(dwt->wavelet + b, dwt->wavelet + a + 1,
+                               dwt->win_l_x0, dwt->win_l_x1,
+                               (OPJ_UINT32)opj_int_min(dwt->sn, dwt->dn - a),
+                               _mm_set1_ps(opj_dwt_beta));
+    opj_v4dwt_decode_step2_sse(dwt->wavelet + a, dwt->wavelet + b + 1,
+                               dwt->win_h_x0, dwt->win_h_x1,
+                               (OPJ_UINT32)opj_int_min(dwt->dn, dwt->sn - b),
+                               _mm_set1_ps(opj_dwt_alpha));
 #else
-    opj_v4dwt_decode_step1(dwt->wavelet + a, dwt->sn, opj_K);
-    opj_v4dwt_decode_step1(dwt->wavelet + b, dwt->dn, opj_c13318);
-    opj_v4dwt_decode_step2(dwt->wavelet + b, dwt->wavelet + a + 1, dwt->sn,
-                           opj_int_min(dwt->sn, dwt->dn - a), opj_dwt_delta);
-    opj_v4dwt_decode_step2(dwt->wavelet + a, dwt->wavelet + b + 1, dwt->dn,
-                           opj_int_min(dwt->dn, dwt->sn - b), opj_dwt_gamma);
-    opj_v4dwt_decode_step2(dwt->wavelet + b, dwt->wavelet + a + 1, dwt->sn,
-                           opj_int_min(dwt->sn, dwt->dn - a), opj_dwt_beta);
-    opj_v4dwt_decode_step2(dwt->wavelet + a, dwt->wavelet + b + 1, dwt->dn,
-                           opj_int_min(dwt->dn, dwt->sn - b), opj_dwt_alpha);
+    opj_v4dwt_decode_step1(dwt->wavelet + a, dwt->win_l_x0, dwt->win_l_x1,
+                           opj_K);
+    opj_v4dwt_decode_step1(dwt->wavelet + b, dwt->win_h_x0, dwt->win_h_x1,
+                           opj_c13318);
+    opj_v4dwt_decode_step2(dwt->wavelet + b, dwt->wavelet + a + 1,
+                           dwt->win_l_x0, dwt->win_l_x1,
+                           (OPJ_UINT32)opj_int_min(dwt->sn, dwt->dn - a),
+                           opj_dwt_delta);
+    opj_v4dwt_decode_step2(dwt->wavelet + a, dwt->wavelet + b + 1,
+                           dwt->win_h_x0, dwt->win_h_x1,
+                           (OPJ_UINT32)opj_int_min(dwt->dn, dwt->sn - b),
+                           opj_dwt_gamma);
+    opj_v4dwt_decode_step2(dwt->wavelet + b, dwt->wavelet + a + 1,
+                           dwt->win_l_x0, dwt->win_l_x1,
+                           (OPJ_UINT32)opj_int_min(dwt->sn, dwt->dn - a),
+                           opj_dwt_beta);
+    opj_v4dwt_decode_step2(dwt->wavelet + a, dwt->wavelet + b + 1,
+                           dwt->win_h_x0, dwt->win_h_x1,
+                           (OPJ_UINT32)opj_int_min(dwt->dn, dwt->sn - b),
+                           opj_dwt_alpha);
 #endif
 }
 
@@ -2080,8 +2118,9 @@ static void opj_v4dwt_decode(opj_v4dwt_t* OPJ_RESTRICT dwt)
 /* <summary>                             */
 /* Inverse 9-7 wavelet transform in 2-D. */
 /* </summary>                            */
-OPJ_BOOL opj_dwt_decode_real(opj_tcd_tilecomp_t* OPJ_RESTRICT tilec,
-                             OPJ_UINT32 numres)
+static
+OPJ_BOOL opj_dwt_decode_tile_97(opj_tcd_tilecomp_t* OPJ_RESTRICT tilec,
+                                OPJ_UINT32 numres)
 {
     opj_v4dwt_t h;
     opj_v4dwt_t v;
@@ -2118,9 +2157,7 @@ OPJ_BOOL opj_dwt_decode_real(opj_tcd_tilecomp_t* OPJ_RESTRICT tilec,
 
     while (--numres) {
         OPJ_FLOAT32 * OPJ_RESTRICT aj = (OPJ_FLOAT32*) tilec->data;
-        OPJ_UINT32 bufsize = (OPJ_UINT32)((tilec->x1 - tilec->x0) *
-                                          (tilec->y1 - tilec->y0));
-        OPJ_INT32 j;
+        OPJ_UINT32 j;
 
         h.sn = (OPJ_INT32)rw;
         v.sn = (OPJ_INT32)rh;
@@ -2135,49 +2172,55 @@ OPJ_BOOL opj_dwt_decode_real(opj_tcd_tilecomp_t* OPJ_RESTRICT tilec,
         h.dn = (OPJ_INT32)(rw - (OPJ_UINT32)h.sn);
         h.cas = res->x0 % 2;
 
-        for (j = (OPJ_INT32)rh; j > 3; j -= 4) {
-            OPJ_INT32 k;
-            opj_v4dwt_interleave_h(&h, aj, (OPJ_INT32)w, (OPJ_INT32)bufsize);
+        h.win_l_x0 = 0;
+        h.win_l_x1 = (OPJ_UINT32)h.sn;
+        h.win_h_x0 = 0;
+        h.win_h_x1 = (OPJ_UINT32)h.dn;
+        for (j = 0; j + 3 < rh; j += 4) {
+            OPJ_UINT32 k;
+            opj_v4dwt_interleave_h(&h, aj, w, rh - j);
             opj_v4dwt_decode(&h);
 
-            for (k = (OPJ_INT32)rw; --k >= 0;) {
-                aj[k               ] = h.wavelet[k].f[0];
-                aj[k + (OPJ_INT32)w  ] = h.wavelet[k].f[1];
-                aj[k + (OPJ_INT32)w * 2] = h.wavelet[k].f[2];
-                aj[k + (OPJ_INT32)w * 3] = h.wavelet[k].f[3];
+            for (k = 0; k < rw; k++) {
+                aj[k      ] = h.wavelet[k].f[0];
+                aj[k + w  ] = h.wavelet[k].f[1];
+                aj[k + w * 2] = h.wavelet[k].f[2];
+                aj[k + w * 3] = h.wavelet[k].f[3];
             }
 
             aj += w * 4;
-            bufsize -= w * 4;
         }
 
-        if (rh & 0x03) {
-            OPJ_INT32 k;
-            j = rh & 0x03;
-            opj_v4dwt_interleave_h(&h, aj, (OPJ_INT32)w, (OPJ_INT32)bufsize);
+        if (j < rh) {
+            OPJ_UINT32 k;
+            opj_v4dwt_interleave_h(&h, aj, w, rh - j);
             opj_v4dwt_decode(&h);
-            for (k = (OPJ_INT32)rw; --k >= 0;) {
-                switch (j) {
+            for (k = 0; k < rw; k++) {
+                switch (rh - j) {
                 case 3:
-                    aj[k + (OPJ_INT32)w * 2] = h.wavelet[k].f[2];
+                    aj[k + w * 2] = h.wavelet[k].f[2];
                 /* FALLTHRU */
                 case 2:
-                    aj[k + (OPJ_INT32)w  ] = h.wavelet[k].f[1];
+                    aj[k + w  ] = h.wavelet[k].f[1];
                 /* FALLTHRU */
                 case 1:
-                    aj[k               ] = h.wavelet[k].f[0];
+                    aj[k] = h.wavelet[k].f[0];
                 }
             }
         }
 
         v.dn = (OPJ_INT32)(rh - (OPJ_UINT32)v.sn);
         v.cas = res->y0 % 2;
+        v.win_l_x0 = 0;
+        v.win_l_x1 = (OPJ_UINT32)v.sn;
+        v.win_h_x0 = 0;
+        v.win_h_x1 = (OPJ_UINT32)v.dn;
 
         aj = (OPJ_FLOAT32*) tilec->data;
-        for (j = (OPJ_INT32)rw; j > 3; j -= 4) {
+        for (j = rw; j > 3; j -= 4) {
             OPJ_UINT32 k;
 
-            opj_v4dwt_interleave_v(&v, aj, (OPJ_INT32)w, 4);
+            opj_v4dwt_interleave_v(&v, aj, w, 4);
             opj_v4dwt_decode(&v);
 
             for (k = 0; k < rh; ++k) {
@@ -2191,7 +2234,7 @@ OPJ_BOOL opj_dwt_decode_real(opj_tcd_tilecomp_t* OPJ_RESTRICT tilec,
 
             j = rw & 0x03;
 
-            opj_v4dwt_interleave_v(&v, aj, (OPJ_INT32)w, j);
+            opj_v4dwt_interleave_v(&v, aj, w, j);
             opj_v4dwt_decode(&v);
 
             for (k = 0; k < rh; ++k) {
@@ -2203,3 +2246,220 @@ OPJ_BOOL opj_dwt_decode_real(opj_tcd_tilecomp_t* OPJ_RESTRICT tilec,
     opj_aligned_free(h.wavelet);
     return OPJ_TRUE;
 }
+
+static
+OPJ_BOOL opj_dwt_decode_partial_97(opj_tcd_t *tcd,
+                                   opj_tcd_tilecomp_t* OPJ_RESTRICT tilec,
+                                   OPJ_UINT32 numres)
+{
+    opj_v4dwt_t h;
+    opj_v4dwt_t v;
+    OPJ_UINT32 resno;
+    const OPJ_UINT32 filter_width = 4U;
+
+    opj_tcd_resolution_t* tr = tilec->resolutions;
+
+    OPJ_UINT32 rw = (OPJ_UINT32)(tr->x1 -
+                                 tr->x0);    /* width of the resolution level computed */
+    OPJ_UINT32 rh = (OPJ_UINT32)(tr->y1 -
+                                 tr->y0);    /* height of the resolution level computed */
+
+    OPJ_UINT32 w = (OPJ_UINT32)(tilec->x1 - tilec->x0);
+
+    size_t l_data_size;
+
+    opj_image_comp_t* image_comp = &(tcd->image->comps[tilec->compno]);
+    /* Compute the intersection of the area of interest, expressed in tile coordinates */
+    /* with the tile coordinates */
+    OPJ_UINT32 win_tcx0 = opj_uint_max(
+                              (OPJ_UINT32)tilec->x0,
+                              opj_uint_ceildiv(tcd->decoded_x0, image_comp->dx));
+    OPJ_UINT32 win_tcy0 = opj_uint_max(
+                              (OPJ_UINT32)tilec->y0,
+                              opj_uint_ceildiv(tcd->decoded_y0, image_comp->dy));
+    OPJ_UINT32 win_tcx1 = opj_uint_min(
+                              (OPJ_UINT32)tilec->x1,
+                              opj_uint_ceildiv(tcd->decoded_x1, image_comp->dx));
+    OPJ_UINT32 win_tcy1 = opj_uint_min(
+                              (OPJ_UINT32)tilec->y1,
+                              opj_uint_ceildiv(tcd->decoded_y1, image_comp->dy));
+
+    l_data_size = opj_dwt_max_resolution(tr, numres);
+    /* overflow check */
+    if (l_data_size > (SIZE_MAX - 5U)) {
+        /* FIXME event manager error callback */
+        return OPJ_FALSE;
+    }
+    l_data_size += 5U;
+    /* overflow check */
+    if (l_data_size > (SIZE_MAX / sizeof(opj_v4_t))) {
+        /* FIXME event manager error callback */
+        return OPJ_FALSE;
+    }
+    h.wavelet = (opj_v4_t*) opj_aligned_malloc(l_data_size * sizeof(opj_v4_t));
+    if (!h.wavelet) {
+        /* FIXME event manager error callback */
+        return OPJ_FALSE;
+    }
+    v.wavelet = h.wavelet;
+
+    for (resno = 1; --numres; resno++) {
+        OPJ_FLOAT32 * OPJ_RESTRICT aj = (OPJ_FLOAT32*) tilec->data;
+        OPJ_UINT32 j;
+        /* Window of interest subband-based coordinates */
+        OPJ_UINT32 win_ll_x0, win_ll_y0, win_ll_x1, win_ll_y1;
+        OPJ_UINT32 win_hl_x0, win_hl_x1;
+        OPJ_UINT32 win_lh_y0, win_lh_y1;
+        /* Window of interest tile-resolution-based coordinates */
+        OPJ_UINT32 win_tr_x0, win_tr_x1, win_tr_y0, win_tr_y1;
+        /* Tile-resolution subband-based coordinates */
+        OPJ_UINT32 tr_ll_x0, tr_ll_y0, tr_hl_x0, tr_lh_y0;
+
+        ++tr;
+
+        h.sn = (OPJ_INT32)rw;
+        v.sn = (OPJ_INT32)rh;
+
+        rw = (OPJ_UINT32)(tr->x1 - tr->x0);
+        rh = (OPJ_UINT32)(tr->y1 - tr->y0);
+
+        h.dn = (OPJ_INT32)(rw - (OPJ_UINT32)h.sn);
+        h.cas = tr->x0 % 2;
+
+        v.dn = (OPJ_INT32)(rh - (OPJ_UINT32)v.sn);
+        v.cas = tr->y0 % 2;
+
+        /* Get the subband coordinates for the window of interest */
+        /* LL band */
+        opj_dwt_get_band_coordinates(tilec, resno, 0,
+                                     win_tcx0, win_tcy0, win_tcx1, win_tcy1,
+                                     &win_ll_x0, &win_ll_y0,
+                                     &win_ll_x1, &win_ll_y1);
+
+        /* HL band */
+        opj_dwt_get_band_coordinates(tilec, resno, 1,
+                                     win_tcx0, win_tcy0, win_tcx1, win_tcy1,
+                                     &win_hl_x0, NULL, &win_hl_x1, NULL);
+
+        /* LH band */
+        opj_dwt_get_band_coordinates(tilec, resno, 2,
+                                     win_tcx0, win_tcy0, win_tcx1, win_tcy1,
+                                     NULL, &win_lh_y0, NULL, &win_lh_y1);
+
+        /* Beware: band index for non-LL0 resolution are 0=HL, 1=LH and 2=HH */
+        tr_ll_x0 = (OPJ_UINT32)tr->bands[1].x0;
+        tr_ll_y0 = (OPJ_UINT32)tr->bands[0].y0;
+        tr_hl_x0 = (OPJ_UINT32)tr->bands[0].x0;
+        tr_lh_y0 = (OPJ_UINT32)tr->bands[1].y0;
+
+        /* Substract the origin of the bands for this tile, to the subwindow */
+        /* of interest band coordinates, so as to get them relative to the */
+        /* tile */
+        win_ll_x0 = opj_uint_subs(win_ll_x0, tr_ll_x0);
+        win_ll_y0 = opj_uint_subs(win_ll_y0, tr_ll_y0);
+        win_ll_x1 = opj_uint_subs(win_ll_x1, tr_ll_x0);
+        win_ll_y1 = opj_uint_subs(win_ll_y1, tr_ll_y0);
+        win_hl_x0 = opj_uint_subs(win_hl_x0, tr_hl_x0);
+        win_hl_x1 = opj_uint_subs(win_hl_x1, tr_hl_x0);
+        win_lh_y0 = opj_uint_subs(win_lh_y0, tr_lh_y0);
+        win_lh_y1 = opj_uint_subs(win_lh_y1, tr_lh_y0);
+
+        opj_dwt_segment_grow(filter_width, (OPJ_UINT32)h.sn, &win_ll_x0, &win_ll_x1);
+        opj_dwt_segment_grow(filter_width, (OPJ_UINT32)h.dn, &win_hl_x0, &win_hl_x1);
+
+        opj_dwt_segment_grow(filter_width, (OPJ_UINT32)v.sn, &win_ll_y0, &win_ll_y1);
+        opj_dwt_segment_grow(filter_width, (OPJ_UINT32)v.dn, &win_lh_y0, &win_lh_y1);
+
+        /* Compute the tile-resolution-based coordinates for the window of interest */
+        if (h.cas == 0) {
+            win_tr_x0 = opj_uint_min(2 * win_ll_x0, 2 * win_hl_x0 + 1);
+            win_tr_x1 = opj_uint_min(opj_uint_max(2 * win_ll_x1, 2 * win_hl_x1 + 1), rw);
+        } else {
+            win_tr_x0 = opj_uint_min(2 * win_hl_x0, 2 * win_ll_x0 + 1);
+            win_tr_x1 = opj_uint_min(opj_uint_max(2 * win_hl_x1, 2 * win_ll_x1 + 1), rw);
+        }
+
+        if (v.cas == 0) {
+            win_tr_y0 = opj_uint_min(2 * win_ll_y0, 2 * win_lh_y0 + 1);
+            win_tr_y1 = opj_uint_min(opj_uint_max(2 * win_ll_y1, 2 * win_lh_y1 + 1), rh);
+        } else {
+            win_tr_y0 = opj_uint_min(2 * win_lh_y0, 2 * win_ll_y0 + 1);
+            win_tr_y1 = opj_uint_min(opj_uint_max(2 * win_lh_y1, 2 * win_ll_y1 + 1), rh);
+        }
+
+        h.win_l_x0 = win_ll_x0;
+        h.win_l_x1 = win_ll_x1;
+        h.win_h_x0 = win_hl_x0;
+        h.win_h_x1 = win_hl_x1;
+        for (j = 0; j + 3 < rh; j += 4, aj += w * 4) {
+            if ((j + 3 >= win_ll_y0 && j < win_ll_y1) ||
+                    (j + 3 >= win_lh_y0 + (OPJ_UINT32)v.sn &&
+                     j < win_lh_y1 + (OPJ_UINT32)v.sn)) {
+                OPJ_UINT32 k;
+                opj_v4dwt_interleave_h(&h, aj, w, rh - j);
+                opj_v4dwt_decode(&h);
+
+                for (k = 0; k < rw; k++) {
+                    aj[k        ] = h.wavelet[k].f[0];
+                    aj[k + w    ] = h.wavelet[k].f[1];
+                    aj[k + w * 2] = h.wavelet[k].f[2];
+                    aj[k + w * 3] = h.wavelet[k].f[3];
+                }
+            }
+        }
+
+        if (j < rh &&
+                ((j + 3 >= win_ll_y0 && j < win_ll_y1) ||
+                 (j + 3 >= win_lh_y0 + (OPJ_UINT32)v.sn &&
+                  j < win_lh_y1 + (OPJ_UINT32)v.sn))) {
+            OPJ_UINT32 k;
+            opj_v4dwt_interleave_h(&h, aj, w, rh - j);
+            opj_v4dwt_decode(&h);
+            for (k = 0; k < rw; k++) {
+                switch (rh - j) {
+                case 3:
+                    aj[k + w * 2] = h.wavelet[k].f[2];
+                /* FALLTHRU */
+                case 2:
+                    aj[k + w    ] = h.wavelet[k].f[1];
+                /* FALLTHRU */
+                case 1:
+                    aj[k        ] = h.wavelet[k].f[0];
+                }
+            }
+        }
+
+        v.win_l_x0 = win_ll_y0;
+        v.win_l_x1 = win_ll_y1;
+        v.win_h_x0 = win_lh_y0;
+        v.win_h_x1 = win_lh_y1;
+        aj = (OPJ_FLOAT32*) tilec->data;
+        aj += win_tr_x0;
+        for (j = win_tr_x0; j < win_tr_x1; j += 4, aj += 4) {
+            OPJ_UINT32 nb_elts = opj_uint_min(4U, win_tr_x1 - j);
+            OPJ_UINT32 k;
+
+            opj_v4dwt_interleave_v(&v, aj, w, nb_elts);
+            opj_v4dwt_decode(&v);
+
+            for (k = win_tr_y0; k < win_tr_y1; ++k) {
+                memcpy(&aj[k * w], &v.wavelet[k], nb_elts * sizeof(OPJ_FLOAT32));
+            }
+        }
+    }
+
+    opj_aligned_free(h.wavelet);
+    return OPJ_TRUE;
+}
+
+
+OPJ_BOOL opj_dwt_decode_real(opj_tcd_t *p_tcd,
+                             opj_tcd_tilecomp_t* OPJ_RESTRICT tilec,
+                             OPJ_UINT32 numres)
+{
+    if (opj_dwt_is_whole_tile_decoding(p_tcd, tilec)) {
+        return opj_dwt_decode_tile_97(tilec, numres);
+    } else {
+        return opj_dwt_decode_partial_97(p_tcd, tilec, numres);
+    }
+}

diff --git a/src/lib/openjp2/dwt.h b/src/lib/openjp2/dwt.h
index ac1e157d0ea2662e5ba88bef9af34297f0437d51..a66ac71e0b588498630ae45b554509d2963b2dd2 100644 (file)
--- a/src/lib/openjp2/dwt.h
+++ b/src/lib/openjp2/dwt.h
@@ -93,10 +93,12 @@ OPJ_BOOL opj_dwt_encode_real(opj_tcd_tilecomp_t * tilec);
 /**
 Inverse 9-7 wavelet transform in 2-D.
 Apply an irreversible inverse DWT transform to a component of an image.
+@param tcd TCD handle
 @param tilec Tile component information (current tile)
 @param numres Number of resolution levels to decode
 */
-OPJ_BOOL opj_dwt_decode_real(opj_tcd_tilecomp_t* OPJ_RESTRICT tilec,
+OPJ_BOOL opj_dwt_decode_real(opj_tcd_t *p_tcd,
+                             opj_tcd_tilecomp_t* OPJ_RESTRICT tilec,
                              OPJ_UINT32 numres);
 
 /**

diff --git a/src/lib/openjp2/tcd.c b/src/lib/openjp2/tcd.c
index f61c65d408e6ecb83b6d97c733a163565cd39762..5e4681c4543598b8c92e386ddfb66e7bcd1ac60d 100644 (file)
--- a/src/lib/openjp2/tcd.c
+++ b/src/lib/openjp2/tcd.c
@@ -1783,7 +1783,8 @@ static OPJ_BOOL opj_tcd_dwt_decode(opj_tcd_t *p_tcd)
                 return OPJ_FALSE;
             }
         } else {
-            if (! opj_dwt_decode_real(l_tile_comp, l_img_comp->resno_decoded + 1)) {
+            if (! opj_dwt_decode_real(p_tcd, l_tile_comp,
+                                      l_img_comp->resno_decoded + 1)) {
                 return OPJ_FALSE;
             }
         }

diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt
index c5e9b22dde07342a5cf5e2b124661b4d70015b6d..579e066b80846260fc6cb52986593fa6c7c571c4 100644 (file)
--- a/tests/CMakeLists.txt
+++ b/tests/CMakeLists.txt
@@ -99,6 +99,10 @@ add_test(NAME tda_prep_irreversible_no_precinct COMMAND test_tile_encoder 1 256
 add_test(NAME tda_irreversible_no_precinct COMMAND test_decode_area -q irreversible_no_precinct.j2k)
 set_property(TEST tda_irreversible_no_precinct APPEND PROPERTY DEPENDS tda_prep_irreversible_no_precinct)
 
+add_test(NAME tda_prep_irreversible_203_201_17_19_no_precinct COMMAND test_tile_encoder 1 203 201 17 19 8 1 irreversible_203_201_17_19_no_precinct.j2k 4 4 3 0 0 1)
+add_test(NAME tda_irreversible_203_201_17_19_no_precinct COMMAND test_decode_area -q irreversible_203_201_17_19_no_precinct.j2k)
+set_property(TEST tda_irreversible_203_201_17_19_no_precinct APPEND PROPERTY DEPENDS tda_prep_irreversible_203_201_17_19_no_precinct)
+
 add_executable(include_openjpeg include_openjpeg.c)
 
 # No image send to the dashboard if lib PNG is not available.