Page MenuHomec4science

wavelet2.h
No OneTemporary

File Metadata

Created
Sat, Apr 27, 16:30

wavelet2.h

/*
=========================================================================
Definitions for 2D wavelet transforms on arrays of 3-tuples.
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
/* 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
/* Perform final Haar wavelet transform on coarsest two detail levels.
* This has the advantage of resulting in a single coarsest level
* approximation coefficient, possibly at the expense of poorer
* compression */
#define WAVELET_FINAL_HAAR
/* Boundary extension modes for padded Daubechies D4 wavelet transform */
#define WAVELET_EXTEND_CONST 0 /* Constant (default if undef) */
#define WAVELET_EXTEND_GRAD1 1 /* 1st order gradient */
#define WAVELET_EXTEND_GRAD2 2 /* 2nd order gradient */
#if 0
/* Not implemented (yet) */
#define WAVELET_EXTEND_PERIOD 3 /* Periodic (wraparound) */
#define WAVELET_EXTEND_SYMM 4 /* Symmetric (reflection) */
#endif
#ifndef WAVELET_EXTEND
#define WAVELET_EXTEND_MODE WAVELET_EXTEND_GRAD2
#endif
/* Convenience macros for handling coefficients */
#define min(a, b) ((a) < (b) ? (a) : (b))
#define max(a, b) ((a) > (b) ? (a) : (b))
#define coeffAvg(a) (((a) [0] + (a) [1] + (a) [2]) / 3)
#define coeffIsEmpty(c) ( \
isnan((c) [0]) || isnan((c) [1]) || isnan((c) [2]) \
)
#define coeffIsZero(c) ((c) [0] == 0 && (c) [1] == 0 && (c) [2] == 0)
/* Determine if coeff tuple c is thresholded, also including empty coeffs
* in debugging mode, since these will be unconditionally dropped. Note
* that the coarsest approximation coefficients are not thresholded,
* since they are essential for reconstruction! */
#ifdef WAVELET_FINAL_HAAR
#define coeffThresh(c, i, j, t) ( \
fabs(t) > 0 && ((i) || (j)) && ( \
coeffIsEmpty((c) [i][j]) || fabs((c) [i][j][0]) < (t) || \
fabs((c) [i][j][1]) < (t) || fabs((c) [i][j][2]) < (t) \
) \
)
#else
#define coeffThresh(c, i, j, t) ( \
fabs(t) > 0 && ((i) > 2 || (j) > 2) && ( \
coeffIsEmpty((c) [i][j]) || fabs((c) [i][j][0]) < (t) || \
fabs((c) [i][j][1]) < (t) || fabs((c) [i][j][2]) < (t) \
) \
)
#endif
/* Haar / Daubechies D4 coefficients */
#define SQRT2 1.41421356
#define SQRT3 1.73205081
#define H4NORM (0.25 / SQRT2)
/* Haar wavelet coeffs */
extern const WAVELET_COEFF h2;
/* Daubechies D4 wavelet coeffs */
extern const WAVELET_COEFF h4 [4];
extern const WAVELET_COEFF g4 [4];
/* 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 **WaveletMatrix2;
/*
----------- 2D COEFFICIENT MATRIX ALLOC/INIT ROUTINES -----------
*/
WaveletMatrix2 allocWaveletMatrix2 (unsigned len);
/* Allocate and return a 2D coefficient array of size (len x len).
Returns NULL if allocation failed. */
void freeWaveletMatrix2 (WaveletMatrix2 y, unsigned len);
/* Free previously allocated 2D coeff array y of size (len x len) */
void zeroCoeffs2 (WaveletMatrix2 y, unsigned len);
/* Set 2D array coefficients to zero */
/*
---------- 2D WAVELET TRANSFORM OPTIMISED FOR SIZES 2^l ----------
------ THESE USE PERIODIC (WRAP-AROUND) BOUNDARY EXTENSION -------
*/
int waveletXform2 (WaveletMatrix2 y, WaveletMatrix2 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. 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.
*/
int waveletInvXform2 (WaveletMatrix2 y, WaveletMatrix2 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.
*/
/*
------------ 2D WAVELET TRANSFORM FOR ARBITRARY SIZES -------------
---- THESE USE CONSTANT BOUNDARY EXTENSION WITH PADDING COEFFS-----
*/
unsigned padWaveletXform2 (WaveletMatrix2 y, WaveletMatrix2 yt,
unsigned yLen0, unsigned *numNonZero
);
/* Perform full 2D multiresolution forward wavelet transform with padded
boundary coefficients on array y containing (yLen0 x yLen0) samples of
3-tuples. Note no intra-tuple transform occurs. An additional array
yt of identical dimension to y is required as buffer for intermediate
results.
The wavelet coefficients are returned in array y, containing the
coarsest approximation coefficients in y [0][0] followed by
horizontal/vertical detail coefficients in order of increasing
resolution/frequency.
NOTE: Both y and yt must be large enough to accommodate the additional
padding coefficients introduced at the array boundaries during the
transform. It is the caller's responsibility to preallocate and
dimension these arrays accordingly with allocWaveletMatrix2(). The
dimension of the output array (which generally exceeds the input
dimension yLen0) must be interrogated beforehand by calling the
function with y or yt set to NULL.
The returned value is the output dimension, or 0 if the transform
failed. In addition, the number of nonzero coefficients is returned
in numNonZero (if not NULL), which is a subset of the number of output
coefficients, i.e. (output dimension)^2.
*/
unsigned padWaveletInvXform2 (WaveletMatrix2 y, WaveletMatrix2 yt,
unsigned yLen0, unsigned yLenOrg
);
/* Perform full 2D multiresolution inverse wavelet transform on array y
containing (yLen0 x yLen0) boundary padded wavelet coefficients as
3-tuples. The original number of samples (yLenOrg x yLenOrg) is
required for the reconstruction. An additional array yt of identical
dimension to y is required as buffer for intermediate results.
The reconstructed samples are returned in y. This is smaller than
(yLen0 x yLen0) since the boundary coefficients from the forward
transform are collapsed.
The returned value is the reconstructed number of samples per
dimension (yLenOrg), or 0 if the inverse transform failed.
*/
#ifdef WAVELET_DBG
/*
----------------- SUPPORT ROUTINES FOR DEBUGGING -----------------
*/
void clearCoeffs2 (WaveletMatrix2 y, unsigned len);
/* Clear 2D array to facilitate debugging */
void dumpCoeffs2 (const WaveletMatrix2 y1, const WaveletMatrix2 y2,
unsigned y1Len, unsigned y2Len, float thresh
);
/*
Dump 2D arrays y1 and y2 of dims y1Len resp. y2Len side-by-side to
stdout (skip y2 if NULL). Coefficients with absolute magnitude
less than thresh are marked with brackets ('[]') as thresholded.
*/
float rmseCoeffs2 (const WaveletMatrix2 y1, const WaveletMatrix2 y2,
unsigned len
);
/* Calculate RMSE between 2D matrices y1 and y2 */
#endif
#endif

Event Timeline