plectrum

PolynomialInterpolation.java (4871B)

---

 package be.tarsos.dsp.wavelet.lift;
 
 /**
 * 
 * @author Ian Kaplan
 */
 class PolynomialInterpolation {
 	/** number of polynomial interpolation ponts */
 	private final static int numPts = 4;
 
 	/** Table for 4-point interpolation coefficients */
 	private float fourPointTable[][];
 
 	/** Table for 2-point interpolation coefficients */
 	private float twoPointTable[][];
 
 	/**
 	 * <p>
 	 * The polynomial interpolation algorithm assumes that the known points are
 	 * located at x-coordinates 0, 1,.. N-1. An interpolated point is calculated
 	 * at <b><i>x</i></b>, using N coefficients. The polynomial coefficients for
 	 * the point <b><i>x</i></b> can be calculated staticly, using the Lagrange
 	 * method.
 	 * </p>
 	 * 
 	 * @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 lagrange(float x, int N, float c[]) {
 		float 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
 	} // lagrange
 
 	/**
 	 * <p>
 	 * For a given N-point polynomial interpolation, fill the coefficient table,
 	 * for points 0.5 ... (N-0.5).
 	 * </p>
 	 */
 	private void fillTable(int N, float table[][]) {
 		float x;
 		float n = N;
 		int i = 0;
 
 		for (x = 0.5f; x < n; x = x + 1.0f) {
 			lagrange(x, N, table[i]);
 			i++;
 		}
 	} // fillTable
 
 	/**
 	 * <p>
 	 * PolynomialWavelets constructor
 	 * </p>
 	 * <p>
 	 * Build the 4-point and 2-point polynomial coefficient tables.
 	 * </p>
 	 */
 	public PolynomialInterpolation() {
 
 		// Fill in the 4-point polynomial interplation table
 		// for the points 0.5, 1.5, 2.5, 3.5
 		fourPointTable = new float[numPts][numPts];
 
 		fillTable(numPts, fourPointTable);
 
 		// Fill in the 2-point polynomial interpolation table
 		// for 0.5 and 1.5
 		twoPointTable = new float[2][2];
 
 		fillTable(2, twoPointTable);
 	} // PolynomialWavelets constructor
 
 	/**
 	 * Print an N x N table polynomial coefficient table
 	 */
 	private void printTable(float 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
 
 	/**
 	 * <p>
 	 * For the polynomial interpolation point x-coordinate <b><i>x</i></b>,
 	 * return the associated polynomial interpolation coefficients.
 	 * </p>
 	 * 
 	 * @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(float x, int n, float c[]) {
 		float table[][] = null;
 
 		int j = (int) x;
 		if (j < 0 || j >= n) {
 			System.out.println("PolynomialWavelets::getCoef: n = " + n
 					+ ", bad x value");
 		}
 
 		if (n == numPts) {
 			table = fourPointTable;
 		} else if (n == 2) {
 			table = twoPointTable;
 			c[2] = 0.0f;
 			c[3] = 0.0f;
 		} else {
 			System.out.println("PolynomialWavelets::getCoef: bad value for N");
 		}
 
 		if (table != null) {
 			for (int i = 0; i < n; i++) {
 				c[i] = table[j][i];
 			}
 		}
 	} // getCoef
 
 	/**
 	 * <p>
 	 * Given four points at the x,y coordinates {0,d<sub>0</sub>},
 	 * {1,d<sub>1</sub>}, {2,d<sub>2</sub>}, {3,d<sub>3</sub>} return the
 	 * y-coordinate value for the polynomial interpolated point at
 	 * <b><i>x</i></b>.
 	 * </p>
 	 * 
 	 * @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).
 	 * @return the y-coordinate value for the polynomial interpolated point at
 	 *         <b><i>x</i></b>.
 	 */
 	public float interpPoint(float x, int N, float d[]) {
 		float c[] = new float[numPts];
 		float 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
 
 } // PolynomialInterpolation