/*
*
* Two-dimensional best basis and custom basis development.
*
* Copyright (C) 1991--94 Wickerhauser Consulting. All Rights Reserved.
* May be freely copied for personal, noncommercial use by owners of
* ``Adapted Wavelet Analysis from Theory to Software'' ISBN 1-56881-041-5
* by Mladen Victor Wickerhauser [AK Peters, Ltd., Wellesley, Mass., 1994]
*/
#include
#include "real.h"
#include "qf.h"
#include "hedge.h"
#include "infocost.h"
#include "cd.h"
#include "cdp2d.h"
#include "wpt2.h"
/*******************************************************************
* bbwp2()
*
* Compute and write the best basis of 2-dimensional isotropic
* wavelet packets to an output hedge, using the space-saving
* periodic discrete wavelet packet analysis algorithm.
*
* Calling sequence and basic algorithm:
*
* bbwp2( GRAPH, KD, IN, X,Y, S, L, WK, H, G ):
* Let XY = X*Y
* Let COST = infocost( IN, 0, XY-1 )
* If Scontents[0]' is replaced with its coefficients
* in the 2-dimensional wavelet packet best basis defined by the cost
* function `infocost()', and the other pointers in `graph->contents[]'
* are set to their corresponding wavelet packet subspaces.
*
* Assumptions:
* 1. x>0; y>0
* 2. L >= s >= 0
* 3. in != NULL
* 4. wk != NULL
* 5. kd != NULL
* 6. graph != NULL
*
* External functions called:
* scdpe2(), infocost(), assert()
*
*/
extern real
bbwp2(
hedge *graph, /* Preallocated struct for the best basis */
real *kd, /* Array to hold one branch of 2-D DWPAP */
const real *in, /* Input array, will not be changed */
int x, /* Length of one column of `in[]' */
int y, /* Length of one row of `in[]' */
int s, /* Current level in the DWPA tree */
int L, /* Maximum depth of the decomposition */
real *wk, /* Scratch array for transposition */
pqf *H, /* Low-pass periodic quadrature filter */
pqf *G) /* High-pass periodic quadrature filter */
{
int block, i, j, xy;
real *out, *k[4];
real cost, kcost;
assert(x>0);
assert(y>0);
assert(s<=L);
assert(s>=0);
assert(in);
assert(wk);
assert(kd);
assert(graph);
xy = x*y;
cost = infocost( in, 0, xy-1 );
if( sblocks;
for( j=0; j<4; j++ )
k[j] = kd + j*xy/4;
kd = k[3] + xy/4;
scdpe2( k[0], k[1], k[2], k[3], in, x,y, wk, H, G);
kcost = 0;
for( j=0; j<4; j++ )
kcost += bbwp2( graph, kd, k[j], x/2, y/2, s+1, L, wk, H, G);
if( kcostlevels[block] = s;
out = (real *)graph->contents[block];
for( i=0; iblocks = block+1;
graph->contents[graph->blocks] = out + xy;
}
}
else
{
graph->levels[graph->blocks] = L;
out = graph->contents[graph->blocks];
for( i=0; iblocks += 1;
graph->contents[graph->blocks] = out + xy;
}
return(cost);
}
/*******************************************************************
* cbwp2()
*
* Compute and write a custom basis of 2-dimensional isotropic
* wavelet packets to an output hedge, doing all development
* into periodic discrete wavelet packets in place.
*
* Calling sequence and basic algorithm:
*
* cbwp2( GRAPH, LEVEL, IX, IY, WORK, H, G ):
* If LEVEL < GRAPH.LEVELS[GRAPH.BLOCKS] then
* scdpi2( GRAPH.CONTENTS[GRAPH.BLOCKS], IX,IY, WORK, H, G )
* For K = 0 to 3
* cbwp2( GRAPH, LEVEL+1, IX/2, IY/2, WORK, H, G )
* Else
* Let GRAPH.CONTENTS[GRAPH.BLOCKS+1] =
* GRAPH.CONTENTS[GRAPH.BLOCKS] + IX*IY
* GRAPH.BLOCKS += 1
* Return
*
*
* Inputs:
*
* (hedge *)graph This preallocated HEDGE must contain enough
* space for the desired graph basis, up to
* `1<<(2*L)' elements in its
* `contents' and `levels' arrays.
*
* (int)level This small integer is the current level.
*
* (int)ix These positive integers are the numbers of
* (int)iy rows and columns of the input.
*
* (real *)work This is a scratch array for transposition.
*
* (pqf *)H Use these low-pass and high-pass periodized
* (pqf *)G filters for convolution-decimation.
*
* Output:
* The array at `graph->contents[0]' is replaced with its coefficients
* in the 2-dimensional wavelet packet basis defined by the encounter
* list `graph->levels[]', and the other pointers in `graph->contents[]'
* are set to their corresponding wavelet packet subspaces.
*
* Assumptions:
* 1. ix>0; iy>0
* 2. level >= 0
* 3. graph != NULL
* 4. work != NULL
*
* External functions called:
* assert(), scdpi2()
*
*/
extern void
cbwp2(
hedge *graph, /* Preallocated struct for the best basis */
int level, /* Current level in the DWPA tree */
int ix, /* Length of one column of `in[]' */
int iy, /* Length of one row of `in[]' */
real *work, /* Scratch array for transposition */
pqf *H, /* Low-pass periodic quadrature filter */
pqf *G) /* High-pass periodic quadrature filter */
{
real *out;
int k;
assert(ix>0);
assert(iy>0);
assert(level>=0);
assert(work);
assert(graph);
if( levellevels[graph->blocks] )
{
out = (real *)graph->contents[graph->blocks];
scdpi2( out, ix,iy, work, H, G );
for( k=0; k<4; k++ )
cbwp2( graph, level+1, ix/2, iy/2, work, H, G );
}
else
{
graph->contents[graph->blocks+1]
= (real *)graph->contents[graph->blocks] + ix*iy;
graph->blocks += 1;
}
return;
}