Sunday, March 8, 2015

Creating order 2 and 3 splines in Java

Wikipedia

Calculating the coefficients:

static public double[] splineCoef(double[] x, double[] y, int degree) {

    if (degree != 2 && degree != 3) {
        throw new IllegalArgumentException();
    }
    int n = (x.length - 1) * (degree + 1);
    double[][] a = new double[n][n];
    double[] b = new double[n];
    double[] coef;
//s
    int k = 0;
    for (int i = 0; i < (x.length - 1) * 2; i += 2, k++) {
        for (int j = 0; j < degree + 1; j++) {
            for (int l = 0; l < 2; l++) {
                a[i + l][j + ((i / 2) * (degree + 1))] = Math.pow(x[k + l], degree - j);
            }
        }
        for (int l = 0; l < 2; l++) {
            b[i + l] = y[k + l];
        }
    }
//s'
    int line = (x.length - 1) * 2;
    int column = 0;
    for (int i = 1; i < x.length - 1; i++) {
        for (int j = 0; j < degree; j++) {
            a[line + i - 1][column + j] = (degree - j) * Math.pow(x[i], degree - 1 - j);
            a[line + i - 1][(column + j) + degree + 1] = -a[line + i - 1][column + j];
        }
        column += degree + 1;
    }
//s"
    if (degree > 2) {
        line += x.length - 2;
        column = 0;
        for (int i = 1; i < x.length - 1; i++) {
            for (int j = 0; j < degree - 1; j++) {
                a[line + i - 1][column + j] = ((degree - j) * (degree - 1 - j)) * Math.pow(x[i], degree - 2 - j);
                a[line + i - 1][(column + j) + degree + 1] = -a[line + i - 1][column + j];
            }
            column += degree + 1;
        }
        if (degree == 3) {
            a[n - 2][0] = 6 * x[0];
            a[n - 2][1] = 2;
            a[n - 1][n - 4] = 6 * x[x.length - 1];
            a[n - 1][n - 3] = 2;
        }

    } else {
        a[n - 1][0] = 1;
    }
    Jama.Matrix A = new Jama.Matrix(a);
    Jama.Matrix c = new Jama.Matrix(b, n);
    Jama.Matrix sol = A.solve(c);
    coef = sol.getColumnPackedCopy();
    return coef;
}

Approximating the value based on the coefficients:

static public double splineVal(double[] coef, double[] x, double val) {

    if (val < x[0]) {
        return Double.NaN;
    }
    int degree = coef.length / x.length;
    for (int i = 1; i < x.length; i++) {
        if (val <= x[i]) {
            int indx = (i - 1) * (degree + 1);
            double s = 0;
            for (int j = 0; j < degree + 1; j++) {
                s += coef[indx + j] * Math.pow(val, degree - j);
            }
            return s;
        }
    }
    return Double.NaN;
}

This is part of the Function Approximations project.

No comments: