diff --git a/wavelet2.c b/wavelet2.c index b0d3b52..b740154 100644 --- a/wavelet2.c +++ b/wavelet2.c @@ -1,741 +1,760 @@ /* ========================================================================= 2-dimensional wavelet transform on (2^l) x (2^l) sized arrays of 3-tuples, where l > 1. Compile with -DWAVELET_TEST to build standalone unit tests. Roland Schregle (roland.schregle@{hslu.ch, gmail.com}) (c) Lucerne University of Applied Sciences and Arts, supported by the Swiss National Science Foundation (SNSF #179067, "Light Fields for Spatio-Temporal Glare Assessment") ========================================================================= $Id$ */ #include "wavelet2.h" #include #include #include #include #define min(a, b) ((a) < (b) ? (a) : (b)) #define coeffAvg(a) (((a) [0] + (a) [1] + (a) [2]) / 3) #define coeffThresh(c, t) ( \ fabs((c) [0]) < (t) || fabs((c) [1]) < (t) || fabs((c) [2]) < (t) \ ) /* The following defs are const, but strict compilers pretend to be dumber than they are, and refuse to init from a func (in this case sqrt(2) and sqrt(3)), or other consts */ #define SQRT2 1.41421356 #define SQRT3 1.73205081 #define H4NORM (0.25 / SQRT2) /* Haar wavelet coeffs */ static const WAVELET_COEFF h2 = 1 / SQRT2; /* Daubechies D4 wavelet coeffs */ static const WAVELET_COEFF h4 [4] = { H4NORM * (1 + SQRT3), H4NORM * (3 + SQRT3), H4NORM * (3 - SQRT3), H4NORM * (1 - SQRT3) }; static const WAVELET_COEFF g4 [4] = { H4NORM * (1 - SQRT3), -H4NORM * (3 - SQRT3), H4NORM * (3 + SQRT3), -H4NORM * (1 + SQRT3) }; /* g4 [4] = { h4 [3], -h4 [2], h4 [1], -h4 [0] }; */ WaveletMatrix allocWaveletMatrix (unsigned l) /* Allocate and return a 2D coefficient array of size (2^l) x (2^l), where l >= 1. Returns NULL if allocation failed. */ { const unsigned len = 1 << l; unsigned i; WaveletMatrix y = NULL; if (l >= 1) { if (!(y = calloc(len, sizeof(WaveletCoeff3*)))) return NULL; for (i = 0; i < len; i++) if (!(y [i] = calloc(len, WAVELET_COEFFSIZE))) return NULL; } return y; } void freeWaveletMatrix (WaveletMatrix y, unsigned l) /* Free previously allocated 2D coefficient array y of size (2^l) x (2^l) */ { unsigned i, j; const unsigned len = 1 << l; if (y) { for (i = 0; i < len; i++) free(y [i]); free(y); } } #ifdef WAVELET_DBG static void zeroCoeffs (WaveletMatrix y, unsigned l) /* Zero output array to facilitate debugging */ { const unsigned len = 1 << l; unsigned i; for (i = 0; i < len; i++) memset(y [i], 0, len * WAVELET_COEFFSIZE); } static void dumpCoeffs (const WaveletMatrix y, const WaveletMatrix yt, unsigned l, float thresh ) /* Dump arrays y and yt side-by-side to stdout (skip yt if NULL) */ { const unsigned len = 1 << l; unsigned i, j, k; for (i = 0; i < len; i++) { for (j = 0; j < len; j++) printf( #ifdef WAVELET_FINALHAAR (i | j) && coeffThresh(y [i][j], thresh) #else (i > 1 || j > 1) && coeffThresh(y [i][j], thresh) #endif - ? "[% 3.2f]\t" : " % 3.2f \t", + ? "[% 5.2f]\t" : " % 5.2f\t", coeffAvg(y [i][j]) ); if (yt) { printf(" -->\t"); for (j = 0; j < len; j++) - printf("% 3.2f \t", coeffAvg(yt [i][j])); + printf("% 5.2f\t", coeffAvg(yt [i][j])); } putchar('\n'); } } static float rmseCoeffs (const WaveletMatrix y1, const WaveletMatrix y2, unsigned l ) /* Calculate RMSE between matrices y and y0 */ { const unsigned len = 1 << l; unsigned i, j; float d, rmse = 0; for (i = 0; i < len; i++) for (j = 0; j < len; j++) { - d = (coeffAvg(y1 [i][j]) - coeffAvg(y2 [i][j])) / + #if 0 + d = (coeffAvg(y1 [i][j]) - coeffAvg(y2 [i][j])) / coeffAvg(y1 [i][j]); + #else + d = coeffAvg(y1 [i][j]) - coeffAvg(y2 [i][j]); + #endif rmse += d * d; } return sqrt(rmse / ((float)len * len)); } #endif static int haarStep (WaveletMatrix y, WaveletMatrix yt, unsigned l) /* Single step of forward 2D Haar wavelet transform on array y of size (2^l) x (2^l) containing 3-tuples, where l >= 1. The transform is performed over y's 2nd axis (i.e. horizontally, assuming row-major addressing as per C convention). Note the triplets per array element are _not_ decorrelated, but transformed independently. The wavelet coefficients are returned in the *TRANSPOSED* output array yt, of identical dimensions to y. The transpose arranges the generated coefficients vertically, and prepares the array for another transform over its 2nd axis in a subsequent call (this time as input y). The result of a subsequent transform restores yt to y's original orientation, in which case both horizontal and vertical axes have been decorrelated. Returns 0 on success, else -1. */ { static unsigned axis = 0; const unsigned len = 1 << l, hlen = 1 << (l - 1); unsigned h, i, j, k; if (len < 2 || !y || !yt) /* Input shorter than wavelet support, or no input/output */ return -1; #ifdef WAVELET_DBG zeroCoeffs(yt, l); #endif /* NOTE: yt is transposed on the fly such that the next function call * transforms over the alternate axis. This is done by simply swapping * the indices during assignment */ for (i = 0; i < len; i++) { for (j = 0; j < len; j += 2) { h = j >> 1; for (k = 0; k < 3; k++) { /* Smooth/approx/avg/lowpass */ yt [h ] [i] [k] = h2 * ( y [i] [j ] [k] + y [i] [j + 1] [k] ); /* Detail/diff/highpass */ yt [hlen + h] [i] [k] = h2 * ( y [i] [j ] [k] - y [i] [j + 1] [k] ); } } } #ifdef WAVELET_DBG printf("%s FWD HAAR (%d x %d)\n", axis ? "VERT" : "HORIZ", len, len); dumpCoeffs(y, yt, l, 0); putchar('\n'); axis ^= 1; #endif return 0; } static int haarInvStep (WaveletMatrix y, WaveletMatrix yt, unsigned l) /* Single step of inverse 2D Haar wavelet transform on coefficient array y of size (2^l) x (2^l) containing 3-tuples, where l >= 1. This reverses the forward transform above. The transform is inverted over y's 2nd axis (i.e. horizontally, assuming row-major addressing as per C convention). The inverted coefficients are returned in the *TRANSPOSED* output array yt, of identical dimensions to y. The transpose arranges the inverted coefficients vertically, and prepares the array for another inverse transform over its 2nd axis in a subsequent call (this time as input y). The result of a subsequent inverse transform restores yt to y's original orientation, in which case both horizontal and vertical axes have been inversely transformed. Returns 0 on success, else -1. */ { static unsigned axis = 1; const unsigned len = 1 << l, hlen = 1 << (l - 1); unsigned h, i, j, k; if (len < 2 || !y || !yt) /* Too few coeffs for reconstruction, or no input/output */ return -1; #ifdef WAVELET_DBG zeroCoeffs(yt, l); #endif /* NOTE: i, j are swapped relative to the forward transform, as axis * order is now reversed. */ /* NOTE: yt is transposed on the fly such that the next function call * inverts over the alternate axis. This is done by simply swapping * the indices during assignment */ for (i = 0; i < len; i++) { for (j = 0; j < len; j += 2) { h = j >> 1; for (k = 0; k < 3; k++) { yt [i] [j ] [k] = h2 * ( y [h ] [i] [k] + /* Avg */ y [hlen + h] [i] [k] /* Diff */ ); yt [i] [j + 1] [k] = h2 * ( y [h ] [i] [k] - /* Avg */ y [hlen + h] [i] [k] /* Diff */ ); } } } #ifdef WAVELET_DBG printf("%s INV HAAR (%d x %d)\n", axis ? "VERT" : "HORIZ", len, len); dumpCoeffs(y, yt, l, 0); putchar('\n'); axis ^= 1; #endif return 0; } static int d4Step (WaveletMatrix y, WaveletMatrix yt, unsigned l) /* Single step of forward 2D Daubechies D4 wavelet transform on array y of size (2^l) x (2^l) containing 3-tuples, where l >= 2. The transform is performed over y's 2nd axis (i.e. horizontally, assuming row-major addressing as per C convention). Note the triplets per array element are _not_ decorrelated, but transformed independently. The wavelet coefficients are returned in the *TRANSPOSED* output array yt, of identical dimensions to y. The transpose arranges the generated coefficients vertically, and prepares the array for another transform over its 2nd axis in a subsequent call (this time as input y). The result of a subsequent transform restores yt to y's original orientation, in which case both horizontal and vertical axes have been decorrelated. Returns 0 on success, else -1. */ { static unsigned axis = 0; const unsigned len = 1 << l, hlen = 1 << (l - 1); unsigned h, i, j, k; if (len < 4 || !y || !yt) /* Input shorter than wavelet support, or no input/output */ return -1; #ifdef WAVELET_DBG zeroCoeffs(yt, l); #endif /* NOTE: yt is transposed on the fly such that the next function call * transforms over the alternate axis. This is done by simply swapping * the indices during assignment */ for (i = 0; i < len; i++) { /* Transform until upper boundary */ for (j = 0; j < len - 2; j += 2) { h = j >> 1; for (k = 0; k < 3; k++) { /* Smooth/approx/avg/lowpass */ yt [h ] [i] [k] = h4 [0] * y [i] [j ] [k] + h4 [1] * y [i] [j + 1] [k] + h4 [2] * y [i] [j + 2] [k] + h4 [3] * y [i] [j + 3] [k]; /* Detail/diff/highpass */ yt [hlen + h] [i] [k] = g4 [0] * y [i] [j ] [k] + g4 [1] * y [i] [j + 1] [k] + g4 [2] * y [i] [j + 2] [k] + g4 [3] * y [i] [j + 3] [k]; } } /* Transform at upper boundary with wraparound. Note j is set to last index from previous loop */ h = j >> 1; for (k = 0; k < 3; k++) { /* Smooth/approx/avg/lowpass */ yt [h ] [i] [k] = h4 [0] * y [i] [j ] [k] + h4 [1] * y [i] [j + 1] [k] + h4 [2] * y [i] [0 ] [k] + h4 [3] * y [i] [1 ] [k]; /* Detail/diff/highpass */ yt [hlen + h] [i] [k] = g4 [0] * y [i] [j ] [k] + g4 [1] * y [i] [j + 1] [k] + g4 [2] * y [i] [0 ] [k] + g4 [3] * y [i] [1 ] [k]; } } #ifdef WAVELET_DBG printf("%s FWD D4 (%d x %d)\n", axis ? "VERT" : "HORIZ", len, len); dumpCoeffs(y, yt, l, 0); putchar('\n'); axis ^= 1; #endif return 0; } static int d4InvStep (WaveletMatrix y, WaveletMatrix yt, unsigned l) /* Single step of inverse 2D Daubechies D4 wavelet transform on coefficient array y of size (2^l) x (2^l) containing 3-tuples, where l >= 2. This reverses the forward transform above. The transform is inverted over y's 2nd axis (i.e. horizontally, assuming row-major addressing as per C convention). The inverted coefficients are returned in the *TRANSPOSED* output array yt, of identical dimensions to y. The transpose arranges the inverted coefficients vertically, and prepares the array for another inverse transform over its 2nd axis in a subsequent call (this time as input y). The result of a subsequent inverse transform restores yt to y's original orientation, in which case both horizontal and vertical axes have been inversely transformed. Returns 0 on success, else -1. */ { static unsigned axis = 1; const unsigned len = 1 << l, hlen = 1 << (l - 1); unsigned h, i, j, k; if (len < 4 || !y || !yt) /* Too few coeffs for reconstruction, or no input/output */ return -1; #ifdef WAVELET_DBG zeroCoeffs(yt, l); #endif /* NOTE: i, j are swapped relative to the forward transform, as axis * order is now reversed. */ /* NOTE: yt is transposed on the fly such that the next function call * inverts over the alternate axis. This is done by simply swapping * the indices during assignment */ for (i = 0; i < len; i++) { /* Invert at lower boundary with wraparound */ for (k = 0; k < 3; k++) { yt [i] [0] [k] = h4 [2] * y [hlen - 1] [i] [k] + /* Last avg */ g4 [2] * y [len - 1] [i] [k] + /* Last diff */ h4 [0] * y [0 ] [i] [k] + /* First avg */ g4 [0] * y [hlen ] [i] [k]; /* First diff */ yt [i] [1] [k] = h4 [3] * y [hlen - 1] [i] [k] + g4 [3] * y [len - 1] [i] [k] + h4 [1] * y [0 ] [i] [k] + g4 [1] * y [hlen ] [i] [k]; } /* Invert until upper boundary */ for (j = 2; j < len; j += 2) { h = (j >> 1) - 1; for (k = 0; k < 3; k++) { yt [i] [j ] [k] = h4 [2] * y [h ] [i] [k] + /* Avg */ g4 [2] * y [hlen + h ] [i] [k] + /* Diff */ h4 [0] * y [h + 1 ] [i] [k] + /* Next avg */ g4 [0] * y [hlen + h + 1] [i] [k]; /* Next diff */ yt [i] [j + 1] [k] = h4 [3] * y [h ] [i] [k] + g4 [3] * y [hlen + h ] [i] [k] + h4 [1] * y [h + 1 ] [i] [k] + g4 [1] * y [hlen + h + 1] [i] [k]; } } } #ifdef WAVELET_DBG printf("%s INV D4 (%d x %d)\n", axis ? "VERT" : "HORIZ", len, len); dumpCoeffs(y, yt, l, 0); putchar('\n'); axis ^= 1; #endif return 0; } int waveletXform2 (WaveletMatrix y, WaveletMatrix yt, unsigned l) /* Perform full 2D multiresolution forward wavelet transform on array y of size (2^l) x (2^l) containing original signal as 3-tuples, where l >= 1. Note no intra-tuple transform occurs. The wavelet coefficients are returned in array y, containing the coarsest approximation in y [0..1][0..1] (or just y [0][0] with WAVELET_FINALHAAR defined) followed by horizontal/vertical details in order of increasing resolution/frequency. A preallocated array yt of identical dimensions to y can be supplied as buffer for intermediate results. If yt == NULL, a buffer is automatically allocated and freed on demand, but this is inefficient for frequent calls. It is recommended to preallocate yt to the maximum expected size. The dimensions of yt are not checked; this is the caller's responsibility. Returns 0 on success, else -1. */ { const unsigned len = 1 << l; unsigned li; WaveletMatrix ytloc = NULL; /* Skip transform if input too short or missing */ if (l < 1 || !y) return -1; if (!yt) /* No buffer supplied; allocate one on demand */ if (!(yt = ytloc = allocWaveletMatrix(l))) { fprintf(stderr, "ERROR - Failed allocating %dx%d buffer array" " in WaveletXform2()", len, len); return -1; } for (li = l; li > 1; li--) { /* Apply horizontal & vertical Daubechies D4 transform, swapping input and transposed output array */ if (d4Step(y, yt, li) || d4Step(yt, y, li)) return -1; } #ifdef WAVELET_FINALHAAR /* Apply horizontal & vertical Haar transform at coarsest resolution (li==1) to obtain single approximation coefficient at y [0][0]; all other coeffs are details. */ if (haarStep(y, yt, li) || haarStep(yt, y, li)) return -1; #endif /* NOTE: All coefficients now in y */ if (ytloc) /* Free yt if allocated on demand */ freeWaveletMatrix(ytloc, l); return 0; } int waveletInvXform2 (WaveletMatrix y, WaveletMatrix yt, unsigned l) /* Perform full 2D multiresolution inverse wavelet transform on array y of size (2^l) x (2^l) containing wavelet coefficients as 3-tuples, where l >= 1. Note no intra-tuple transform occurs. A preallocated array yt of identical dimensions to y can be supplied as buffer for intermediate results. If yt == NULL, a buffer is automatically allocated and freed on demand, but this is inefficient for frequent calls. It is recommended to preallocate yt to the maximum expected size. The dimensions of yt are not checked; this is the caller's responsibility. The reconstructed signal is returned in array y. Returns 0 on success, else -1. */ { const unsigned len = 1 << l; unsigned li; WaveletMatrix ytloc = NULL; /* Skip inverse transform if input too short or missing */ if (l < 1 || !y) return -1; if (!yt) /* No buffer supplied; allocate one on demand */ if (!(yt = ytloc = allocWaveletMatrix(l))) { fprintf(stderr, "ERROR - Failed allocating %dx%d buffer array" " in WaveletInvXform2()", len, len); return -1; } #ifdef WAVELET_FINALHAAR /* Invert horizontal & vertical Haar transform at coarsest level (li==1), swapping input and transposed output array */ if (haarInvStep(y, yt, 1) || haarInvStep(yt, y, 1)) return -1; #endif for (li = 2; li <= l; li++) { /* Invert horizontal & vertical Daubechies D4 transform, swapping input and transposed output arrays */ if (d4InvStep(y, yt, li) || d4InvStep(yt, y, li)) return -1; } /* NOTE: Reconstructed signal now in y */ if (ytloc) /* Free yt if allocated on demand */ freeWaveletMatrix(ytloc, l); return 0; } #ifdef WAVELET_TEST #include int main (int argc, char *argv []) { int i, j, k, l; unsigned len, numThresh = 0; WaveletMatrix y0 = NULL, y = NULL; FILE *dataFile = NULL; float inData, thresh = 0; if (argc < 2) { fprintf(stderr, "%s [threshold] [dataFile]\n", argv [0]); fputs("Missing array resolution l > 1, " "compression threshold >= 0\n", stderr ); return -1; } if (!(l = atoi(argv [1])) || l < 1) { fputs("Invalid array resolution l\n", stderr); return -1; } else len = 1 << l; if (argc > 2 && (thresh = atof(argv [2])) < 0) { fprintf(stderr, "Invalid threshold %.3f\n", thresh); return -1; } /* Allocate arrays for original and reconstruction */ if (!(y0 = allocWaveletMatrix(l)) || !(y = allocWaveletMatrix(l))) { fprintf(stderr, "Failed allocating %dx%d array\n", len, len); return -1; } if (argc > 3) { /* Load data from file; length must not exceed allocated */ if (!(dataFile = fopen(argv [3], "r"))) { fprintf(stderr, "Failed opening data file %s\n", argv [3]); return -1; } for (i = 0; i < len; i++) { for (j = 0; j < len; j++) { if (feof(dataFile)) { fprintf(stderr, "Premature end of file reading data from %s\n", argv [2] ); fclose(dataFile); return -1; } /* Read next float, skipping any leading whitespace */ if (fscanf(dataFile, " %f", &inData)) { y0 [i][j][0] = y0 [i][j][1] = y0 [i][j][2] = y [i][j][0] = y [i][j][1] = y [i][j][2] = (WAVELET_COEFF)inData; } else { fprintf(stderr, "Error reading from data file %s\n", argv [2] ); fclose(dataFile); return -1; } } } fclose(dataFile); } else { - /* Init with random data */ + /* Init input */ srand48(111); for (i = 0; i < len; i++) for (j = 0; j < len; j++) #if 0 + /* Random data, channel-independent */ for (k = 0; k < 3; k++) y0 [i][j][k] = y [i][j][k] = drand48(); - #else + #endif + #if 0 + /* Random data, indentical for all channels */ y0 [i][j][0] = y0 [i][j][1] = y0 [i][j][2] = y [i][j][0] = y [i][j][1] = y [i][j][2] = drand48(); + #endif + #if 0 /* Monotonically increasing along axis 0 */ + y0 [i][j][0] = y0 [i][j][1] = y0 [i][j][2] = + y [i][j][0] = y [i][j][1] = y [i][j][2] = i; + #endif + #if 0 /* Monotonically increasing along axis 1 */ + y0 [i][j][0] = y0 [i][j][1] = y0 [i][j][2] = + y [i][j][0] = y [i][j][1] = y [i][j][2] = j; + #endif + #if 1 /* Monotonically increasing along both axes */ + y0 [i][j][0] = y0 [i][j][1] = y0 [i][j][2] = + y [i][j][0] = y [i][j][1] = y [i][j][2] = i * j; #endif } /* Forward xform */ if (waveletXform2(y, NULL, l)) { fputs("Forward xform failed\n", stderr); return -1; } /* Threshold coefficients; we use hard thresholding as it's easier * to implement than soft thresholding, which requires sorting the * coefficients. * NOTE: y [0..1][0..1] (or just y [0][0] if WAVELET_FINALHAAR is * defned) are omitted as they are the coarsest approximation * coefficients. */ for (i = 0; i < len; i++) for (j = 0; j < len; j++) { #ifdef WAVELET_FINALHAAR if ((i | j) && coeffThresh(y [i][j], thresh)) { #else if ((i > 1 || j > 1) && coeffThresh(y [i][j], thresh)) { #endif y [i][j][0] = y [i][j][1] = y [i][j][2] = 0; numThresh++; #if 0 /* Replace thresholded values with random noise in range [-threshold, threshold] */ y [i][j][0] = y [i][j][1] = y [i][j][2] = thresh * ( 2 * drand48() - 1 ); #endif } } #ifdef WAVELET_DBG /* Dump coefficients */ puts("-----------------------------------------------------------\n"); if (numThresh) printf("%d/%d coefficients thresholded\n", numThresh, len * len ); dumpCoeffs(y, NULL, l, thresh); puts("\n-----------------------------------------------------------\n"); #endif /* Inverse xform */ if (waveletInvXform2(y, NULL, l)) { fputs("Inverse xform failed\n", stderr); return -1; } #ifdef WAVELET_DBG puts("-----------------------------------------------------------\n"); puts("ORIG vs. INV XFORM"); dumpCoeffs(y0, y, l, 0); #endif printf("\nAvg RMSE = %.2f\n", rmseCoeffs(y0, y, l)); freeWaveletMatrix(y0, l); freeWaveletMatrix(y, l); return 0; } #endif diff --git a/wavelet2.h b/wavelet2.h index d189206..aa2b274 100644 --- a/wavelet2.h +++ b/wavelet2.h @@ -1,97 +1,97 @@ /* ========================================================================= 2-dimensional wavelet transform on (2^l) x (2^l) sized arrays of 3-tuples, where l >= 1. Roland Schregle (roland.schregle@{hslu.ch, gmail.com}) (c) Lucerne University of Applied Sciences and Arts, supported by the Swiss National Science Foundation (SNSF #179067, "Light Fields for Spatio-Temporal Glare Assessment") ========================================================================= $Id$ */ #ifndef _WAVELET2_H #define _WAVELET2_H - /* Uncomment to run a Haar wavelet transform on the DB4 coefficients at - * coarsest resolution to obtain a single approximate coefficient. This - * may introduce compression artefacts, however. Use with care! */ - /* #define WAVELET_FINALHAAR */ + /* Run a final Haar wavelet transform on the DB4 coefficients at the + * coarsest resolution to obtain a single approximation coefficient (as + * opposed to 2 x 2 with just DB4. */ + #define WAVELET_FINALHAAR /* Wavelet coefficient type; defaults to double if not already defined. (may be overridden in compiler command line, for example). NOTE: single precision float may improve cache performance during 2D transform with large arrays. */ #ifndef WAVELET_COEFF #define WAVELET_COEFF double #endif /* Wavelet matrix defs; these are stored as arrays of pointers (a.k.a. * Iliffe vectors), as this proved to be more performant than a flattened * representation. Encapsulating the coefficient's triplets in a struct * prevents the compiler from introducing another layer of indirection. * */ typedef WAVELET_COEFF WaveletCoeff3 [3]; #define WAVELET_COEFFSIZE (sizeof(WaveletCoeff3)) typedef WaveletCoeff3 **WaveletMatrix; WaveletMatrix allocWaveletMatrix (unsigned l); /* Allocate and return a 2D coefficient array of size (2^l) x (2^l), where l >= 1. Returns NULL if allocation failed. */ void freeWaveletMatrix (WaveletMatrix y, unsigned l); /* Free previously allocated 2D coefficient array y of size (2^l) x (2^l) */ int waveletXform2 (WaveletMatrix y, WaveletMatrix yt, unsigned l); /* Perform full 2D multiresolution forward wavelet transform on array y of size (2^l) x (2^l) containing original signal as 3-tuples, where l >= 1. Note no intra-tuple transform occurs. The wavelet coefficients are returned in array y, containing the coarsest approximation in y [0][0] followed by horizontal/vertical details in order of increasing resolution/frequency. A preallocated array yt of identical dimensions to y can be supplied as buffer for intermediate results. If yt == NULL, a buffer is automatically allocated and freed on demand, but this is inefficient for frequent calls. Ideally, yt should be preallocated to the maximum expected size. The dimensions of a preallocated yt are not checked; this is the caller's responsibility. Returns 0 on success, else -1. */ int waveletInvXform2 (WaveletMatrix y, WaveletMatrix yt, unsigned l); /* Perform full 2D multiresolution inverse wavelet transform on array y of size (2^l) x (2^l) containing wavelet coefficients as 3-tuples, where l >= 1. Note no intra-tuple transform occurs. A preallocated array yt of identical dimensions to y can be supplied as buffer for intermediate results. If yt == NULL, a buffer is automatically allocated and freed on demand, but this is inefficient for frequent calls. Ideally, yt should be preallocated to the maximum expected size. The dimensions of a preallocated yt are not checked; this is the caller's responsibility. The reconstructed signal is returned in array y. Returns 0 on success, else -1. */ #endif