`package be.tarsos.dsp.wavelet.lift;`

/**

* <p>

* Line (with slope) wavelet

* </p>

*

* <p>

* The wavelet Lifting Scheme "LineWavelet" wavelet approximates the data set

* using a LineWavelet with with slope (in contrast to the HaarWavelet wavelet

* where a LineWavelet has zero slope is used to approximate the data).

* </p>

*

* <p>

* The predict stage of the LineWavelet wavelet "predicts" that an odd point

* will lie midway between its two neighboring even points. That is, that the

* odd point will lie on a LineWavelet between the two adjacent even points. The

* difference between this "prediction" and the actual odd value replaces the

* odd element.

* </p>

*

* <p>

* The update stage calculates the average of the odd and even element pairs,

* although the method is indirect, since the predict phase has over written the

* odd value.

* </p>

*

* <h4>

* Copyright and Use</h4>

*

* <p>

* You may use this source code without limitation and without fee as long as

* you include:

* </p>

* <blockquote> This software was written and is copyrighted by Ian Kaplan, Bear

* Products International, www.bearcave.com, 2001. </blockquote>

* <p>

* This software is provided "as is", without any warrenty or claim as to its

* usefulness. Anyone who uses this source code uses it at their own risk. Nor

* is any support provided by Ian Kaplan and Bear Products International.

* <p>

* Please send any bug fixes or suggested source changes to:

*

* <pre>

*

iank@bearcave.com

* </pre>

*

* @author Ian Kaplan

*/

public class LineWavelet extends LiftingSchemeBaseWavelet {

/**

* <p>

* Calculate an extra "even" value for the LineWavelet wavelet algorithm at

* the end of the data series. Here we pretend that the last two values in

* the data series are at the x-axis coordinates 0 and 1, respectively. We

* then need to calculate the y-axis value at the x-axis coordinate 2. This

* point lies on a LineWavelet running through the points at 0 and 1.

* </p>

* <p>

* Given two points, x<sub>1</sub>, y<sub>1</sub> and x<sub>2</sub>,

* y<sub>2</sub>, where

* </p>

*

* <pre>

*

x<sub>1</sub> = 0

*

x<sub>2</sub> = 1

* </pre>

* <p>

* calculate the point on the LineWavelet at x<sub>3</sub>, y<sub>3</sub>,

* where

* </p>

*

* <pre>

*

x<sub>3</sub> = 2

* </pre>

* <p>

* The "two-point equation" for a LineWavelet given x<sub>1</sub>,

* y<sub>1</sub> and x<sub>2</sub>, y<sub>2</sub> is

* </p>

*

* <pre>

*

.

y<sub>2</sub> - y<sub>1</sub>

*

(y - y<sub>1</sub>) = -------- (x - x<sub>1</sub>)

*

.

x<sub>2</sub> - x<sub>1</sub>

* </pre>

* <p>

* Solving for y

* </p>

*

* <pre>

*

.

y<sub>2</sub> - y<sub>1</sub>

*

y = -------- (x - x<sub>1</sub>) + y<sub>1</sub>

*

.

x<sub>2</sub> - x<sub>1</sub>

* </pre>

* <p>

* Since x<sub>1</sub> = 0 and x<sub>2</sub> = 1

* </p>

*

* <pre>

*

.

y<sub>2</sub> - y<sub>1</sub>

*

y = -------- (x - 0) + y<sub>1</sub>

*

.

1 - 0

* </pre>

* <p>

* or

* </p>

*

* <pre>

*

y = (y<sub>2</sub> - y<sub>1</sub>)*x + y<sub>1</sub>

* </pre>

* <p>

* We're calculating the value at x<sub>3</sub> = 2, so

* </p>

*

* <pre>

*

y = 2*y<sub>2</sub> - 2*y<sub>1</sub> + y<sub>1</sub>

* </pre>

* <p>

* or

* </p>

*

* <pre>

*

y = 2*y<sub>2</sub> - y<sub>1</sub>

* </pre>

*/

private float new_y(float y1, float y2) {

float y = 2 * y2 - y1;

return y;

}

/**

* <p>

* Predict phase of LineWavelet Lifting Scheme wavelet

* </p>

*

* <p>

* The predict step attempts to "predict" the value of an odd element from

* the even elements. The difference between the prediction and the actual

* element is stored as a wavelet coefficient.

* </p>

* <p>

* The "predict" step takes place after the split step. The split step will

* move the odd elements (b<sub>j</sub>) to the second half of the array,

* leaving the even elements (a<sub>i</sub>) in the first half

* </p>

*

* <pre>

*

a<sub>0</sub>, a<sub>1</sub>, a<sub>1</sub>, a<sub>3</sub>, b<sub>0</sub>, b<sub>1</sub>, b<sub>2</sub>, b<sub>2</sub>,

* </pre>

* <p>

* The predict step of the LineWavelet wavelet "predicts" that the odd

* element will be on a LineWavelet between two even elements.

* </p>

*

* <pre>

*

b<sub>j+1,i</sub> = b<sub>j,i</sub> - (a<sub>j,i</sub> + a<sub>j,i+1</sub>)/2

* </pre>

* <p>

* Note that when we get to the end of the data series the odd element is

* the last element in the data series (remember, wavelet algorithms work on

* data series with 2<sup>n</sup> elements). Here we "predict" that the odd

* element will be on a LineWavelet that runs through the last two even

* elements. This can be calculated by assuming that the last two even

* elements are located at x-axis coordinates 0 and 1, respectively. The odd

* element will be at 2. The <i>new_y()</i> function is called to do this

* simple calculation.

* </p>

*/

protected void predict(float[] vec, int N, int direction) {

int half = N >> 1;

float predictVal;

for (int i = 0; i < half; i++) {

int j = i + half;

if (i < half - 1) {

predictVal = (vec[i] + vec[i + 1]) / 2;

} else if (N == 2) {

predictVal = vec[0];

} else {

// calculate the last "odd" prediction

predictVal = new_y(vec[i - 1], vec[i]);

}

if (direction == forward) {

vec[j] = vec[j] - predictVal;

} else if (direction == inverse) {

vec[j] = vec[j] + predictVal;

} else {

System.out.println("predictline::predict: bad direction value");

}

}

} // predict

/**

* <p>

* The predict phase works on the odd elements in the second half of the

* array. The update phase works on the even elements in the first half of

* the array. The update phase attempts to preserve the average. After the

* update phase is completed the average of the even elements should be

* approximately the same as the average of the input data set from the

* previous iteration. The result of the update phase becomes the input for

* the next iteration.

* </p>

* <p>

* In a HaarWavelet wavelet the average that replaces the even element is

* calculated as the average of the even element and its associated odd

* element (e.g., its odd neighbor before the split). This is not possible

* in the LineWavelet wavelet since the odd element has been replaced by the

* difference between the odd element and the mid-point of its two even

* neighbors. As a result, the odd element cannot be recovered.

* </p>

* <p>

* The value that is added to the even element to preserve the average is

* calculated by the equation shown below. This equation is given in Wim

* Sweldens' journal articles and his tutorial (<i>Building Your Own

* Wavelets at Home</i>) and in <i>Ripples in Mathematics</i>. A somewhat

* more complete derivation of this equation is provided in <i>Ripples in

* Mathematics</i> by A. Jensen and A. la Cour-Harbo, Springer, 2001.

* </p>

* <p>

* The equation used to calculate the average is shown below for a given

* iteratin <i>i</i>. Note that the predict phase has already completed, so

* the odd values belong to iteration <i>i+1</i>.

* </p>

*

* <pre>

*

even<sub>i+1,j</sub> = even<sub>i,j</sub> op (odd<sub>i+1,k-1</sub> + odd<sub>i+1,k</sub>)/4

* </pre>

* <p>

* There is an edge problem here, when i = 0 and k = N/2 (e.g., there is no

* k-1 element). We assume that the odd<sub>i+1,k-1</sub> is the same as

* odd<sub>k</sub>. So for the first element this becomes

*

* <pre>

*

(2 * odd<sub>k</sub>)/4

* </pre>

* <p>

* or

* </p>

*

* <pre>

*

odd<sub>k</sub>/2

* </pre>

*/

protected void update(float[] vec, int N, int direction) {

int half = N >> 1;

for (int i = 0; i < half; i++) {

int j = i + half;

float val;

if (i == 0) {

val = vec[j] / 2.0f;

} else {

val = (vec[j - 1] + vec[j]) / 4.0f;

}

if (direction == forward) {

vec[i] = vec[i] + val;

} else if (direction == inverse) {

vec[i] = vec[i] - val;

} else {

System.out.println("update: bad direction value");

}

} // for

}

} // LineWavelet