improve handling of odd sized arrays in DWT/IDWT

[grlee77]

Currently if an input to dwt (or dwt2 or dwtn) is odd in size, the resulting coefficient arrays do not reflect this and calling idwt on the coefficients produced will return an even sized image.

If the original odd image extent is extracted, the reconstruction is perfect as expected.

I am not certain of the proper solution, but I think maybe dwt_buffer_length should return different sizes for the approximation and detail coefficients when the input is of odd size. e.g. for MODE_PERIODIZATION, currently this function returns the size ceil(input_len / 2.0).

At least from what I read in the JPEG 2000 standard, I think the approximation coefficients should use ceil(input_len / 2.0), but the detail coefficients should use floor(input_len / 2.0). Reference: Eqns 7a-b of the JPEG2000 compression standard: http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.136.3993

As a concrete example, dwtn of a 2D array of size 511 x 255 would give the following coefficient array sizes: aa : 256 x 128 da : 255 x 128 ad : 256 x 127 dd : 255 x 127

This would require idwtn to allow unequal shapes of coefficients (see #54 )

[kwohlfahrt]

This proposal looks good to me. It conflicts with current behaviour so it would be a breaking change, but I think this is useful enough to justify it.

[grlee77]

I am not certain if it is worth the trouble, but thought it was worth considering. It is not high priority for me because the signals I work with tend to be even length anyways.

Also, FWIW, I just checked and the behavior we have now seems to match Matlab's dwt2, idwt2 behavior. i.e. for odd size input, all the coefficient arrays are even size.

[kwohlfahrt]

Hm. I've had a quick look, dealing with the edge length functions is no problem. However, the last component(s?) of the approximation coefficients must be adjusted for the fact that the corresponding detail is missing. See below for an example, db3 with MODE_PERIODIZATION, the last F/2 reconstructed values are incorrect. For other modes, only the last reconstructed value is incorrect.

I think I might be missing something in Fig 6 of the paper, I'm familiar with the filter bank looking like this, but don't quite follow their description.

$ clang wt.c convolution.c common.c wavelets.c test_odd.c; and ./a.out 
input:
 0.00  1.00  2.00  3.00  4.00  5.00  6.00  7.00  8.00 
filter (dec_a, dec_d, rec_a, rec_d):
 0.04 -0.09 -0.14  0.46  0.81  0.33 
-0.33  0.81 -0.46 -0.14  0.09  0.04 
 0.33  0.81  0.46 -0.14 -0.09  0.04 
 0.04  0.09 -0.14 -0.46  0.81 -0.33 
approximate:
 8.92  1.16  3.98  6.78 10.28 
detail:
 1.12  0.00  0.00  0.33 
reconstruction:
 3.72 -0.54  2.00  3.00  4.00  5.00  6.19  7.46  7.11

[grlee77]

I don't know that I will be able to take a look this week. I don't have a lot of familiarity with that standard, but came across the example there when looking for info on dealing with odd-sized arrays.

[kwohlfahrt]

Right, the idea proposed here is that the signal must be extrapolated so the detail coefficient of the last element would be 0. This should be possible for all of the wavelets, but he uses the CDF9/7, which is the same used by the JPEG2000 standard, so might rely on some special properties of that wavelet.

[kwohlfahrt]

OK, seems to be doable. However, the question is whether we want to do this, since it adds an arbitrary edge extension to odd-length coefficients that doesn't really fit in with any of the other schemas. I don't really think there should be any objections to doing it with MODE_PERIODIZATION, since that already has a constant-edge padding thrown in, but I'm not so sure about the others since they have well-defined edge extension modes.

(edit) C code for calculating the padding value with MODE_PERIODIZATION:

float padding = 0;
for (size_t i = 0; i < w->dec_len / 2; ++i){
    padding += w->dec_hi_float[w->dec_len / 2 + i] * input[N-1-i];
}
for (size_t i = 0; i < w->dec_len / 2 - 1; ++i){
    padding += w->dec_hi_float[w->dec_len/2-2-i] * input[i];
}
padding /= -w->dec_hi_float[w->dec_len/2 - 1];

[grlee77]

I am closing this because I don't think it is worth the effort that would be required to implement it.