predicttest.java100644 1040 1001 4670 7375114437 12470 0ustar everyone import predict.*; class predicttest { static double[] copy( double vec[] ) { double newVec[] = null; if (vec != null) { final int len = vec.length; newVec = new double[ len ]; for (int i = 0; i < len; i++) { newVec[i] = vec[i]; } } return newVec; } // copy static boolean equal( double v1[], double v2[] ) { boolean rslt = true; int len = v1.length; for (int i = 0; i < len; i++) { if (v1[i] != v2[i]) { rslt = false; break; } } // for return rslt; } // equal static void pr( double vec[] ) { int cnt = 0; final int len = vec.length; for (int i = 0; i < len; i++) { System.out.print( vec[i] + " "); cnt++; if (cnt == 6) { cnt = 0; System.out.println(); } } System.out.println(); } // pr public static void main( String[] args ) { double vals1[] = { 32.0, 10.0, 20.0, 38.0, 37.0, 28.0, 38.0, 34.0, 18.0, 24.0, 18.0, 9.0, 23.0, 24.0, 28.0, 34.0 }; double vals2[] = { 77.6875, 78.1875, 82.0625, 85.5625, 86.7500, 82.4375, 82.2500, 82.7500, 81.2500, 79.5625, 80.2813, 79.8750, 77.7500, 74.7500, 78.5000, 79.1875, 78.8125, 80.3125, 80.1250, 79.3125, 83.7500, 89.8125, 87.7500, 91.1250, 94.4375, 92.7500, 98.0000, 97.1875, 99.4375, 101.7500, 108.5000, 109.0000, 105.2500, 105.5000, 110.0000, 107.0000, 107.2500, 103.3125, 102.8750, 102.4375, 102.0000, 101.3125, 97.4375, 100.5000, 107.7500, 110.2500, 114.3125, 111.2500, 114.8125, 112.6875, 109.4375, 108.0625, 104.5625, 103.2500, 110.5625, 110.7500, 116.3125, 123.6250, 120.9375, 121.6250, 127.6875, 126.0625, 126.3750, 124.3750 }; double save[] = copy( vals2 ); // predicthaar haar = new predicthaar(); // haar.forwardTrans( vals2 ); // haar.inverseTrans( vals2 ); // predictline line = new predictline(); // line.forwardTrans( vals2 ); // line.inverseTrans( vals2 ); predictpoly poly = new predictpoly(); poly.forwardTrans( vals2 ); poly.pr_ordered( vals2 ); System.out.println(); poly.inverseTrans( vals2 ); if (equal(save, vals2)) { System.out.println("Arrays are the same"); } else { System.out.println("Arrays are NOT the same"); pr( vals2 ); } } // main } predict/ 40755 1040 1001 0 7375117454 10623 5ustar everyonepredict/predictbase.java100644 1040 1001 14664 7375117102 14071 0ustar everyone package predict; /**

class predictbase: base class for ur-Lifting Scheme wavelets using split/predict or predict/merge steps.

Predict wavelets are a prototype for Lifting Scheme wavlets. Like the Lifting Scheme wavelets, predict wavelets can be thought of as a set of steps that perform the in-place wavelet calculation. The building blocks consist of

The split and merge methods are shared by all Lifting Scheme wavelet algorithms. This base class provides the transform and inverse transform methods (forwardTrans and inverseTrans). The predict method is an abstract method that is defined for a particular predict wavelet sub-class.

References:

Copyright and Use

You may use this source code without limitation and without fee as long as you include:

This software was written and is copyrighted by Ian Kaplan, Bear Products International, www.bearcave.com, 2001.

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.

Please send any bug fixes or suggested source changes to:

     iank@bearcave.com
@author Ian Kaplan */ public abstract class predictbase { /** "enumeration" for forward wavelet transform */ protected final int forward = 1; /** "enumeration" for inverse wavelet transform */ protected final int inverse = 2; /** Move the odd elements to the second half of the N element region in the vec array. */ protected void split( double[] vec, int N ) { int start = 1; int end = N - 1; while (start < end) { for (int i = start; i < end; i = i + 2) { double tmp = vec[i]; vec[i] = vec[i+1]; vec[i+1] = tmp; } start = start + 1; end = end - 1; } } /** Merge the odd elements from the second half of the N element region in the array with the even elements in the first half of the N element region. The result will be the combination of the odd and even elements in a region of length N. */ protected void merge( double[] vec, int N ) { int half = N >> 1; int start = half-1; int end = half; while (start > 0) { for (int i = start; i < end; i = i + 2) { double tmp = vec[i]; vec[i] = vec[i+1]; vec[i+1] = tmp; } start = start - 1; end = end + 1; } } /** Round a double to four decimal places. */ private double round4( double d ) { final double multiplier = 10000.0; double rounded; rounded = Math.round(d * multiplier)/multiplier; return rounded; } // round4 /**

Print the result of an ordered wavelet transform. the result is printed to show the coefficient (difference) bands, where the bands are in increasing frequency.

@param vec[] result of a wavelet transform */ public void pr_ordered( double[] vec ) { if (vec != null) { int len = vec.length; if (len > 0) { System.out.println(round4(vec[0])); int num_in_freq = 1; int cnt = 0; for (int i = 1; i < len; i++) { System.out.print(round4(vec[i]) + " "); cnt++; if (cnt == num_in_freq) { System.out.println(); cnt = 0; num_in_freq = num_in_freq * 2; } } } } } // pr_ordered /** Predict step, to be defined by the subclass @param vec input array @param N size of region to act on (from 0..N-1) @param direction forward or inverse transform */ protected abstract void predict( double[] vec, int N, int direction ); /**

Predict wavelet forward transform

forwardTrans is passed an array of doubles. The array size must be a power of two. The predict wavelet transform is calculated in-place and the result is returned in the argument array.

The result of forwardTrans is a set of differences, ordered by powers of two. Element vec[0] will contain the element that started the data set. Element vec[1] will contain the difference between the prediction for vec[1] and the original value at vec[1]. The next set of differences will be at vec[2..3]... the differences at vec[N/2..N-1] will contain the largest set of differences.

*/ public void forwardTrans( double[] vec ) { final int N = vec.length; for (int n = N; n > 1; n = n >> 1) { split( vec, n ); predict( vec, n, forward ); } } // forwardTrans /**

Default two step Lifting Scheme inverse wavelet transform

inverseTrans is passed the result of an ordered set of predict wavelet coefficients. The inverse transform is calculated in-place and the result is returned in the argument array.

*/ public void inverseTrans( double[] vec ) { final int N = vec.length; for (int n = 2; n <= N; n = n << 1) { predict( vec, n, inverse ); merge( vec, n ); } } // inverseTrans } // liftbase predict/predicthaar.java100644 1040 1001 4074 7375117056 14054 0ustar everyone package predict; /**

Haar (flat line) predict wavelets.

The wavelet Lifting Scheme predict wavelets are a proto-wavelet, which does not include an averaging step (referred to as an update step in the Lifting Scheme).

The predicthaar class has a Haar wavelet predict step. The Haar predict step "predicts" that an odd element will have the same value as it preceeding even element. Stated another way, the odd element is "predicted" to be on a float (zero slope line) shared with the odd point. The difference between this "prediction" and the actual odd value replaces the odd element.

Before the predict step is invoked, the split step has run. The split step moves the odd elements to the second half of the N element region.

Copyright and Use

You may use this source code without limitation and without fee as long as you include:

This software was written and is copyrighted by Ian Kaplan, Bear Products International, www.bearcave.com, 2001.

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.

Please send any bug fixes or suggested source changes to:

     iank@bearcave.com
@author Ian Kaplan */ public class predicthaar extends predictbase { /** Haar predict step */ protected void predict( double[] vec, int N, int direction ) { int half = N >> 1; int cnt = 0; for (int i = 0; i < half; i++) { double predictVal = vec[i]; int j = i + half; if (direction == forward) { vec[j] = vec[j] - predictVal; } else if (direction == inverse) { vec[j] = vec[j] + predictVal; } else { System.out.println("predicthaar::predict: bad direction value"); } } } } // predicthaar predict/predictline.java100644 1040 1001 12545 7375117120 14102 0ustar everyone package predict; /**

Line (with slope) predict wavelets

The wavelet Lifting Scheme predict wavelets are a proto-wavelet, which does not include an averaging step (referred to as an update step in the Lifting Scheme).

The line wavelet "predicts" that an even point will like midway between its two neighboring even points. That is, that the odd point will lie on a line between the two adjacent even points. The difference between this "prediction" and the actual odd value replaces the odd element.

Copyright and Use

You may use this source code without limitation and without fee as long as you include:

This software was written and is copyrighted by Ian Kaplan, Bear Products International, www.bearcave.com, 2001.

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.

Please send any bug fixes or suggested source changes to:

     iank@bearcave.com
@author Ian Kaplan */ public class predictline extends predictbase { /**

Calculate an extra value for the line 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 data series (y-axis value) at the x-axis coordinate 2.

Given two points, x1, y1 and x2, y2, where

        x1 = 0
        x2 = 1
     

calculate the point on the line at x3, y3, where

        x3 = 2
     

The "two-point equation" for a line given x1, y1 and x2, y2 is

     .          y2 - y1
     (y - y1) = -------- (x - x1)
     .          x2 - x1
     

Solving for y

     .    y2 - y1
     y = -------- (x - x1) + y1
     .    x2 - x1
     

Since x1 = 0 and x2 = 1

     .    y2 - y1
     y = -------- (x - 0) + y1
     .    1 - 0
     

or

     y = (y2 - y1)*x + y1
     

We're calculating the value at x3 = 2, so

     y = 2*y2 - 2*y1 + y1
     

or

     y = 2*y2 - y1
     
*/ private double new_y( double y1, double y2) { double y = 2 * y2 - y1; return y; } /**

Predict phase of line Lifting Scheme wavelet

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.

The "predict" step takes place after the split step. The split step will move the odd elements (bj) to the second half of the array, leaving the even elements (ai) in the first half

    a0, a1, a1, a3, b0, b1, b2, b2, 
    

The predict step of the line wavelet "predicts" that the odd element will be on a line between two even elements.

    bj+1,i = bj,i - (aj,i + aj,i+1)/2
    

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 2n elements). Here we "predict" that the odd element will be on a line 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 new_y() function is called to do this simple calculation.

*/ protected void predict( double[] vec, int N, int direction ) { int half = N >> 1; double 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 } predict/predictpoly.java100644 1040 1001 22273 7375116700 14140 0ustar everyone package predict; /**

Polynomial predict wavelets

The wavelet Lifting Scheme predict wavelets are a proto-wavelet, which does not include an averaging step (referred to as an update step in the Lifting Scheme).

The polynomial wavelet uses 4-point polynomial interpolation to predict an even point from the odd points. The difference between the "prediction" and the actual value of the odd point replaces the odd point. For complete documentation see

  
  http://www.bearcave.com/misl/misl_tech/wavelets/lifting/index.html
  

Copyright and Use

You may use this source code without limitation and without fee as long as you include:

This software was written and is copyrighted by Ian Kaplan, Bear Products International, www.bearcave.com, 2001.

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.

Please send any bug fixes or suggested source changes to:

     iank@bearcave.com
@author Ian Kaplan */ public class predictpoly extends predictbase { /** number of polynomial interpolation ponts */ final int numPts = 4; /** Table for 4-point interpolation coefficients */ private double fourPointTable[][]; /** Table for 2-point interpolation coefficients */ private double twoPointTable[][]; /**

The polynomial interpolation algorithm assumes that the known points are located at x-coordinates 0, 1,.. N-1. An interpolated point is calculated at x, using N coefficients. The polynomial coefficients for the point x can be calculated staticly, using the Legrange method.

@param x the x-coordinate of the interpolated point @param N the number of polynomial points. @param c[] an array for returning the coefficients */ private void legrange( double x, int N, double c[] ) { double num, denom; for (int i = 0; i < N; i++) { num = 1; denom = 1; for (int k = 0; k < N; k++) { if (i != k) { num = num * (x - k); denom = denom * (i - k); } } // for k c[i] = num / denom; } // for i } // legrange /**

For a given N-point polynomial interpolation, fill the coefficient table, for points 0.5 ... (N-0.5).

*/ private void fillTable( int N, double table[][] ) { double x; double n = N; int i = 0; for (x = 0.5; x < n; x = x + 1.0) { legrange( x, N, table[i] ); i++; } } // fillTable /**

predictpoly constructor

Build the 4-point and 2-point polynomial coefficient tables.

*/ public predictpoly() { // Fill in the 4-point polynomial interplation table // for the points 0.5, 1.5, 2.5, 3.5 fourPointTable = new double[numPts][numPts]; fillTable( numPts, fourPointTable ); // Fill in the 2-point polynomial interpolation table // for 0.5 and 1.5 twoPointTable = new double[2][2]; fillTable( 2, twoPointTable ); } // predictpoly constructor /** Print an N x N table polynomial coefficient table */ private void printTable( double table[][], int N) { System.out.println(N + "-point interpolation table:"); double x = 0.5; for (int i = 0; i < N; i++) { System.out.print(x + ": "); for (int j = 0; j < N; j++) { System.out.print( table[i][j] ); if (j < N-1) System.out.print(", "); } System.out.println(); x = x + 1.0; } } /** Print the 4-point and 2-point polynomial coefficient tables. */ public void printTables() { printTable( fourPointTable, numPts ); printTable( twoPointTable, 2 ); } // printTables /**

For the polynomial interpolation point x-coordinate x, return the associated polynomial interpolation coefficients.

@param x the x-coordinate for the interpolated pont @param n the number of polynomial interpolation points @param c[] an array to return the polynomial coefficients */ private void getCoef( double x, int n, double c[]) { double table[][] = null; int j = (int)x; if (j < 0 || j >= n) { System.out.println("predictpoly::getCoef: n = " + n + ", bad x value"); } if (n == numPts) { table = fourPointTable; } else if (n == 2) { table = twoPointTable; c[2] = 0.0; c[3] = 0.0; } else { System.out.println("predictpoly::getCoef: bad value for N"); } if (table != null) { for (int i = 0; i < n; i++) { c[i] = table[j][i]; } } } // getCoef /**

Given four points at the x,y coordinates {0,d0}, {1,d1}, {2,d2}, {3,d3} return the y-coordinate value for the polynomial interpolated point at x.

@param x the x-coordinate for the point to be interpolated @param N the number of interpolation points @param d[] an array containing the y-coordinate values for the known points (which are located at x-coordinates 0..N-1). */ private double interpPoint( double x, int N, double d[]) { double c[] = new double[ numPts ]; double point = 0; int n = numPts; if (N < numPts) n = N; getCoef( x, n, c ); if (n == numPts) { point = c[0]*d[0] + c[1]*d[1] + c[2]*d[2] + c[3]*d[3]; } else if (n == 2) { point = c[0]*d[0] + c[1]*d[1]; } return point; } // interpPoint /**

Copy N data points from vec into d These points are the "known" points used in the polynomial interpolation.

@param vec[] the input data set on which the wavelet is calculated @param d[] an array into which N data points, starting at start are copied. @param N the number of polynomial interpolation points @param start the index in vec from which copying starts */ private void fill( double vec[], double d[], int N, int start ) { int n = numPts; if (n > N) n = N; int end = start + n; int j = 0; for (int i = start; i < end; i++) { d[j] = vec[i]; j++; } } // fill /**

Predict an odd point from the even points, using 4-point polynomial interpolation.

The four points used in the polynomial interpolation are the even points. We pretend that these four points are located at the x-coordinates 0,1,2,3. The first odd point interpolated will be located between the first and second even point, at 0.5. The next N-3 points are located at 1.5 (in the middle of the four points). The last two points are located at 2.5 and 3.5. For complete documentation see

  
  http://www.bearcave.com/misl/misl_tech/wavelets/lifting/index.html
    

The difference between the predicted (interpolated) value and the actual odd value replaces the odd value in the forward transform.

As the recursive steps proceed, N will eventually be 4 and then 2. When N = 4, linear interpolation is used. When N = 2, Haar interpolation is used (the prediction for the odd value is that it is equal to the even value).

@param vec the input data on which the forward or inverse transform is calculated. @param N the area of vec over which the transform is calculated @param direction forward or inverse transform */ protected void predict( double[] vec, int N, int direction ) { int half = N >> 1; double d[] = new double[ numPts ]; int k = 42; for (int i = 0; i < half; i++) { double predictVal; if (i == 0) { if (half == 1) { // e.g., N == 2, and we use Haar interpolation predictVal = vec[0]; } else { fill( vec, d, N, 0 ); predictVal = interpPoint( 0.5, half, d ); } } else if (i == 1) { predictVal = interpPoint( 1.5, half, d ); } else if (i == half-2) { predictVal = interpPoint( 2.5, half, d ); } else if (i == half-1) { predictVal = interpPoint( 3.5, half, d ); } else { fill( vec, d, N, i-1); predictVal = interpPoint( 1.5, half, d ); } int j = i + half; if (direction == forward) { vec[j] = vec[j] - predictVal; } else if (direction == inverse) { vec[j] = vec[j] + predictVal; } else { System.out.println("predictpoly::predict: bad direction value"); } } } // predict } // predictpoly