1package com.android.gallery3d.filtershow.filters;
2
3
4public class SplineMath {
5    double[][] mPoints = new double[6][2];
6    double[] mDerivatives;
7    SplineMath(int n) {
8        mPoints = new double[n][2];
9    }
10
11    public void setPoint(int index, double x, double y) {
12        mPoints[index][0] = x;
13        mPoints[index][1] = y;
14        mDerivatives = null;
15    }
16
17    public float[][] calculatetCurve(int n) {
18        float[][] curve = new float[n][2];
19        double[][] points = new double[mPoints.length][2];
20        for (int i = 0; i < mPoints.length; i++) {
21
22            points[i][0] = mPoints[i][0];
23            points[i][1] = mPoints[i][1];
24
25        }
26        double[] derivatives = solveSystem(points);
27        float start = (float) points[0][0];
28        float end = (float) (points[points.length - 1][0]);
29
30        curve[0][0] = (float) (points[0][0]);
31        curve[0][1] = (float) (points[0][1]);
32        int last = curve.length - 1;
33        curve[last][0] = (float) (points[points.length - 1][0]);
34        curve[last][1] = (float) (points[points.length - 1][1]);
35
36        for (int i = 0; i < curve.length; i++) {
37
38            double[] cur = null;
39            double[] next = null;
40            double x = start + i * (end - start) / (curve.length - 1);
41            int pivot = 0;
42            for (int j = 0; j < points.length - 1; j++) {
43                if (x >= points[j][0] && x <= points[j + 1][0]) {
44                    pivot = j;
45                }
46            }
47            cur = points[pivot];
48            next = points[pivot + 1];
49            if (x <= next[0]) {
50                double x1 = cur[0];
51                double x2 = next[0];
52                double y1 = cur[1];
53                double y2 = next[1];
54
55                // Use the second derivatives to apply the cubic spline
56                // equation:
57                double delta = (x2 - x1);
58                double delta2 = delta * delta;
59                double b = (x - x1) / delta;
60                double a = 1 - b;
61                double ta = a * y1;
62                double tb = b * y2;
63                double tc = (a * a * a - a) * derivatives[pivot];
64                double td = (b * b * b - b) * derivatives[pivot + 1];
65                double y = ta + tb + (delta2 / 6) * (tc + td);
66
67                curve[i][0] = (float) (x);
68                curve[i][1] = (float) (y);
69            } else {
70                curve[i][0] = (float) (next[0]);
71                curve[i][1] = (float) (next[1]);
72            }
73        }
74        return curve;
75    }
76
77    public double getValue(double x) {
78        double[] cur = null;
79        double[] next = null;
80        if (mDerivatives == null)
81            mDerivatives = solveSystem(mPoints);
82        int pivot = 0;
83        for (int j = 0; j < mPoints.length - 1; j++) {
84            pivot = j;
85            if (x <= mPoints[j][0]) {
86                break;
87            }
88        }
89        cur = mPoints[pivot];
90        next = mPoints[pivot + 1];
91        double x1 = cur[0];
92        double x2 = next[0];
93        double y1 = cur[1];
94        double y2 = next[1];
95
96        // Use the second derivatives to apply the cubic spline
97        // equation:
98        double delta = (x2 - x1);
99        double delta2 = delta * delta;
100        double b = (x - x1) / delta;
101        double a = 1 - b;
102        double ta = a * y1;
103        double tb = b * y2;
104        double tc = (a * a * a - a) * mDerivatives[pivot];
105        double td = (b * b * b - b) * mDerivatives[pivot + 1];
106        double y = ta + tb + (delta2 / 6) * (tc + td);
107
108        return y;
109
110    }
111
112    double[] solveSystem(double[][] points) {
113        int n = points.length;
114        double[][] system = new double[n][3];
115        double[] result = new double[n]; // d
116        double[] solution = new double[n]; // returned coefficients
117        system[0][1] = 1;
118        system[n - 1][1] = 1;
119        double d6 = 1.0 / 6.0;
120        double d3 = 1.0 / 3.0;
121
122        // let's create a tridiagonal matrix representing the
123        // system, and apply the TDMA algorithm to solve it
124        // (see http://en.wikipedia.org/wiki/Tridiagonal_matrix_algorithm)
125        for (int i = 1; i < n - 1; i++) {
126            double deltaPrevX = points[i][0] - points[i - 1][0];
127            double deltaX = points[i + 1][0] - points[i - 1][0];
128            double deltaNextX = points[i + 1][0] - points[i][0];
129            double deltaNextY = points[i + 1][1] - points[i][1];
130            double deltaPrevY = points[i][1] - points[i - 1][1];
131            system[i][0] = d6 * deltaPrevX; // a_i
132            system[i][1] = d3 * deltaX; // b_i
133            system[i][2] = d6 * deltaNextX; // c_i
134            result[i] = (deltaNextY / deltaNextX) - (deltaPrevY / deltaPrevX); // d_i
135        }
136
137        // Forward sweep
138        for (int i = 1; i < n; i++) {
139            // m = a_i/b_i-1
140            double m = system[i][0] / system[i - 1][1];
141            // b_i = b_i - m(c_i-1)
142            system[i][1] = system[i][1] - m * system[i - 1][2];
143            // d_i = d_i - m(d_i-1)
144            result[i] = result[i] - m * result[i - 1];
145        }
146
147        // Back substitution
148        solution[n - 1] = result[n - 1] / system[n - 1][1];
149        for (int i = n - 2; i >= 0; --i) {
150            solution[i] = (result[i] - system[i][2] * solution[i + 1]) / system[i][1];
151        }
152        return solution;
153    }
154
155    public static void main(String[] args) {
156        SplineMath s = new SplineMath(10);
157        for (int i = 0; i < 10; i++) {
158            s.setPoint(i, i, i);
159        }
160        float[][] curve = s.calculatetCurve(40);
161
162        for (int j = 0; j < curve.length; j++) {
163            System.out.println(curve[j][0] + "," + curve[j][1]);
164        }
165    }
166}
167