18921c28c7333ad2b4d34f013904ad4737044f366John Hofordpackage com.android.gallery3d.filtershow.filters; 28921c28c7333ad2b4d34f013904ad4737044f366John Hoford 38921c28c7333ad2b4d34f013904ad4737044f366John Hoford 48921c28c7333ad2b4d34f013904ad4737044f366John Hofordpublic class SplineMath { 58921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[][] mPoints = new double[6][2]; 68921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[] mDerivatives; 78921c28c7333ad2b4d34f013904ad4737044f366John Hoford SplineMath(int n) { 88921c28c7333ad2b4d34f013904ad4737044f366John Hoford mPoints = new double[n][2]; 98921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 108921c28c7333ad2b4d34f013904ad4737044f366John Hoford 118921c28c7333ad2b4d34f013904ad4737044f366John Hoford public void setPoint(int index, double x, double y) { 128921c28c7333ad2b4d34f013904ad4737044f366John Hoford mPoints[index][0] = x; 138921c28c7333ad2b4d34f013904ad4737044f366John Hoford mPoints[index][1] = y; 148921c28c7333ad2b4d34f013904ad4737044f366John Hoford mDerivatives = null; 158921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 168921c28c7333ad2b4d34f013904ad4737044f366John Hoford 178921c28c7333ad2b4d34f013904ad4737044f366John Hoford public float[][] calculatetCurve(int n) { 188921c28c7333ad2b4d34f013904ad4737044f366John Hoford float[][] curve = new float[n][2]; 198921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[][] points = new double[mPoints.length][2]; 208921c28c7333ad2b4d34f013904ad4737044f366John Hoford for (int i = 0; i < mPoints.length; i++) { 218921c28c7333ad2b4d34f013904ad4737044f366John Hoford 228921c28c7333ad2b4d34f013904ad4737044f366John Hoford points[i][0] = mPoints[i][0]; 238921c28c7333ad2b4d34f013904ad4737044f366John Hoford points[i][1] = mPoints[i][1]; 248921c28c7333ad2b4d34f013904ad4737044f366John Hoford 258921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 268921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[] derivatives = solveSystem(points); 278921c28c7333ad2b4d34f013904ad4737044f366John Hoford float start = (float) points[0][0]; 288921c28c7333ad2b4d34f013904ad4737044f366John Hoford float end = (float) (points[points.length - 1][0]); 298921c28c7333ad2b4d34f013904ad4737044f366John Hoford 308921c28c7333ad2b4d34f013904ad4737044f366John Hoford curve[0][0] = (float) (points[0][0]); 318921c28c7333ad2b4d34f013904ad4737044f366John Hoford curve[0][1] = (float) (points[0][1]); 328921c28c7333ad2b4d34f013904ad4737044f366John Hoford int last = curve.length - 1; 338921c28c7333ad2b4d34f013904ad4737044f366John Hoford curve[last][0] = (float) (points[points.length - 1][0]); 348921c28c7333ad2b4d34f013904ad4737044f366John Hoford curve[last][1] = (float) (points[points.length - 1][1]); 358921c28c7333ad2b4d34f013904ad4737044f366John Hoford 368921c28c7333ad2b4d34f013904ad4737044f366John Hoford for (int i = 0; i < curve.length; i++) { 378921c28c7333ad2b4d34f013904ad4737044f366John Hoford 388921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[] cur = null; 398921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[] next = null; 408921c28c7333ad2b4d34f013904ad4737044f366John Hoford double x = start + i * (end - start) / (curve.length - 1); 418921c28c7333ad2b4d34f013904ad4737044f366John Hoford int pivot = 0; 428921c28c7333ad2b4d34f013904ad4737044f366John Hoford for (int j = 0; j < points.length - 1; j++) { 438921c28c7333ad2b4d34f013904ad4737044f366John Hoford if (x >= points[j][0] && x <= points[j + 1][0]) { 448921c28c7333ad2b4d34f013904ad4737044f366John Hoford pivot = j; 458921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 468921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 478921c28c7333ad2b4d34f013904ad4737044f366John Hoford cur = points[pivot]; 488921c28c7333ad2b4d34f013904ad4737044f366John Hoford next = points[pivot + 1]; 498921c28c7333ad2b4d34f013904ad4737044f366John Hoford if (x <= next[0]) { 508921c28c7333ad2b4d34f013904ad4737044f366John Hoford double x1 = cur[0]; 518921c28c7333ad2b4d34f013904ad4737044f366John Hoford double x2 = next[0]; 528921c28c7333ad2b4d34f013904ad4737044f366John Hoford double y1 = cur[1]; 538921c28c7333ad2b4d34f013904ad4737044f366John Hoford double y2 = next[1]; 548921c28c7333ad2b4d34f013904ad4737044f366John Hoford 558921c28c7333ad2b4d34f013904ad4737044f366John Hoford // Use the second derivatives to apply the cubic spline 568921c28c7333ad2b4d34f013904ad4737044f366John Hoford // equation: 578921c28c7333ad2b4d34f013904ad4737044f366John Hoford double delta = (x2 - x1); 588921c28c7333ad2b4d34f013904ad4737044f366John Hoford double delta2 = delta * delta; 598921c28c7333ad2b4d34f013904ad4737044f366John Hoford double b = (x - x1) / delta; 608921c28c7333ad2b4d34f013904ad4737044f366John Hoford double a = 1 - b; 618921c28c7333ad2b4d34f013904ad4737044f366John Hoford double ta = a * y1; 628921c28c7333ad2b4d34f013904ad4737044f366John Hoford double tb = b * y2; 638921c28c7333ad2b4d34f013904ad4737044f366John Hoford double tc = (a * a * a - a) * derivatives[pivot]; 648921c28c7333ad2b4d34f013904ad4737044f366John Hoford double td = (b * b * b - b) * derivatives[pivot + 1]; 658921c28c7333ad2b4d34f013904ad4737044f366John Hoford double y = ta + tb + (delta2 / 6) * (tc + td); 668921c28c7333ad2b4d34f013904ad4737044f366John Hoford 678921c28c7333ad2b4d34f013904ad4737044f366John Hoford curve[i][0] = (float) (x); 688921c28c7333ad2b4d34f013904ad4737044f366John Hoford curve[i][1] = (float) (y); 698921c28c7333ad2b4d34f013904ad4737044f366John Hoford } else { 708921c28c7333ad2b4d34f013904ad4737044f366John Hoford curve[i][0] = (float) (next[0]); 718921c28c7333ad2b4d34f013904ad4737044f366John Hoford curve[i][1] = (float) (next[1]); 728921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 738921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 748921c28c7333ad2b4d34f013904ad4737044f366John Hoford return curve; 758921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 768921c28c7333ad2b4d34f013904ad4737044f366John Hoford 778921c28c7333ad2b4d34f013904ad4737044f366John Hoford public double getValue(double x) { 788921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[] cur = null; 798921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[] next = null; 808921c28c7333ad2b4d34f013904ad4737044f366John Hoford if (mDerivatives == null) 818921c28c7333ad2b4d34f013904ad4737044f366John Hoford mDerivatives = solveSystem(mPoints); 828921c28c7333ad2b4d34f013904ad4737044f366John Hoford int pivot = 0; 838921c28c7333ad2b4d34f013904ad4737044f366John Hoford for (int j = 0; j < mPoints.length - 1; j++) { 848921c28c7333ad2b4d34f013904ad4737044f366John Hoford pivot = j; 858921c28c7333ad2b4d34f013904ad4737044f366John Hoford if (x <= mPoints[j][0]) { 868921c28c7333ad2b4d34f013904ad4737044f366John Hoford break; 878921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 888921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 898921c28c7333ad2b4d34f013904ad4737044f366John Hoford cur = mPoints[pivot]; 908921c28c7333ad2b4d34f013904ad4737044f366John Hoford next = mPoints[pivot + 1]; 918921c28c7333ad2b4d34f013904ad4737044f366John Hoford double x1 = cur[0]; 928921c28c7333ad2b4d34f013904ad4737044f366John Hoford double x2 = next[0]; 938921c28c7333ad2b4d34f013904ad4737044f366John Hoford double y1 = cur[1]; 948921c28c7333ad2b4d34f013904ad4737044f366John Hoford double y2 = next[1]; 958921c28c7333ad2b4d34f013904ad4737044f366John Hoford 968921c28c7333ad2b4d34f013904ad4737044f366John Hoford // Use the second derivatives to apply the cubic spline 978921c28c7333ad2b4d34f013904ad4737044f366John Hoford // equation: 988921c28c7333ad2b4d34f013904ad4737044f366John Hoford double delta = (x2 - x1); 998921c28c7333ad2b4d34f013904ad4737044f366John Hoford double delta2 = delta * delta; 1008921c28c7333ad2b4d34f013904ad4737044f366John Hoford double b = (x - x1) / delta; 1018921c28c7333ad2b4d34f013904ad4737044f366John Hoford double a = 1 - b; 1028921c28c7333ad2b4d34f013904ad4737044f366John Hoford double ta = a * y1; 1038921c28c7333ad2b4d34f013904ad4737044f366John Hoford double tb = b * y2; 1048921c28c7333ad2b4d34f013904ad4737044f366John Hoford double tc = (a * a * a - a) * mDerivatives[pivot]; 1058921c28c7333ad2b4d34f013904ad4737044f366John Hoford double td = (b * b * b - b) * mDerivatives[pivot + 1]; 1068921c28c7333ad2b4d34f013904ad4737044f366John Hoford double y = ta + tb + (delta2 / 6) * (tc + td); 1078921c28c7333ad2b4d34f013904ad4737044f366John Hoford 1088921c28c7333ad2b4d34f013904ad4737044f366John Hoford return y; 1098921c28c7333ad2b4d34f013904ad4737044f366John Hoford 1108921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 1118921c28c7333ad2b4d34f013904ad4737044f366John Hoford 1128921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[] solveSystem(double[][] points) { 1138921c28c7333ad2b4d34f013904ad4737044f366John Hoford int n = points.length; 1148921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[][] system = new double[n][3]; 1158921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[] result = new double[n]; // d 1168921c28c7333ad2b4d34f013904ad4737044f366John Hoford double[] solution = new double[n]; // returned coefficients 1178921c28c7333ad2b4d34f013904ad4737044f366John Hoford system[0][1] = 1; 1188921c28c7333ad2b4d34f013904ad4737044f366John Hoford system[n - 1][1] = 1; 1198921c28c7333ad2b4d34f013904ad4737044f366John Hoford double d6 = 1.0 / 6.0; 1208921c28c7333ad2b4d34f013904ad4737044f366John Hoford double d3 = 1.0 / 3.0; 1218921c28c7333ad2b4d34f013904ad4737044f366John Hoford 1228921c28c7333ad2b4d34f013904ad4737044f366John Hoford // let's create a tridiagonal matrix representing the 1238921c28c7333ad2b4d34f013904ad4737044f366John Hoford // system, and apply the TDMA algorithm to solve it 1248921c28c7333ad2b4d34f013904ad4737044f366John Hoford // (see http://en.wikipedia.org/wiki/Tridiagonal_matrix_algorithm) 1258921c28c7333ad2b4d34f013904ad4737044f366John Hoford for (int i = 1; i < n - 1; i++) { 1268921c28c7333ad2b4d34f013904ad4737044f366John Hoford double deltaPrevX = points[i][0] - points[i - 1][0]; 1278921c28c7333ad2b4d34f013904ad4737044f366John Hoford double deltaX = points[i + 1][0] - points[i - 1][0]; 1288921c28c7333ad2b4d34f013904ad4737044f366John Hoford double deltaNextX = points[i + 1][0] - points[i][0]; 1298921c28c7333ad2b4d34f013904ad4737044f366John Hoford double deltaNextY = points[i + 1][1] - points[i][1]; 1308921c28c7333ad2b4d34f013904ad4737044f366John Hoford double deltaPrevY = points[i][1] - points[i - 1][1]; 1318921c28c7333ad2b4d34f013904ad4737044f366John Hoford system[i][0] = d6 * deltaPrevX; // a_i 1328921c28c7333ad2b4d34f013904ad4737044f366John Hoford system[i][1] = d3 * deltaX; // b_i 1338921c28c7333ad2b4d34f013904ad4737044f366John Hoford system[i][2] = d6 * deltaNextX; // c_i 1348921c28c7333ad2b4d34f013904ad4737044f366John Hoford result[i] = (deltaNextY / deltaNextX) - (deltaPrevY / deltaPrevX); // d_i 1358921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 1368921c28c7333ad2b4d34f013904ad4737044f366John Hoford 1378921c28c7333ad2b4d34f013904ad4737044f366John Hoford // Forward sweep 1388921c28c7333ad2b4d34f013904ad4737044f366John Hoford for (int i = 1; i < n; i++) { 1398921c28c7333ad2b4d34f013904ad4737044f366John Hoford // m = a_i/b_i-1 1408921c28c7333ad2b4d34f013904ad4737044f366John Hoford double m = system[i][0] / system[i - 1][1]; 1418921c28c7333ad2b4d34f013904ad4737044f366John Hoford // b_i = b_i - m(c_i-1) 1428921c28c7333ad2b4d34f013904ad4737044f366John Hoford system[i][1] = system[i][1] - m * system[i - 1][2]; 1438921c28c7333ad2b4d34f013904ad4737044f366John Hoford // d_i = d_i - m(d_i-1) 1448921c28c7333ad2b4d34f013904ad4737044f366John Hoford result[i] = result[i] - m * result[i - 1]; 1458921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 1468921c28c7333ad2b4d34f013904ad4737044f366John Hoford 1478921c28c7333ad2b4d34f013904ad4737044f366John Hoford // Back substitution 1488921c28c7333ad2b4d34f013904ad4737044f366John Hoford solution[n - 1] = result[n - 1] / system[n - 1][1]; 1498921c28c7333ad2b4d34f013904ad4737044f366John Hoford for (int i = n - 2; i >= 0; --i) { 1508921c28c7333ad2b4d34f013904ad4737044f366John Hoford solution[i] = (result[i] - system[i][2] * solution[i + 1]) / system[i][1]; 1518921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 1528921c28c7333ad2b4d34f013904ad4737044f366John Hoford return solution; 1538921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 1548921c28c7333ad2b4d34f013904ad4737044f366John Hoford 1558921c28c7333ad2b4d34f013904ad4737044f366John Hoford public static void main(String[] args) { 1568921c28c7333ad2b4d34f013904ad4737044f366John Hoford SplineMath s = new SplineMath(10); 1578921c28c7333ad2b4d34f013904ad4737044f366John Hoford for (int i = 0; i < 10; i++) { 1588921c28c7333ad2b4d34f013904ad4737044f366John Hoford s.setPoint(i, i, i); 1598921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 1608921c28c7333ad2b4d34f013904ad4737044f366John Hoford float[][] curve = s.calculatetCurve(40); 1618921c28c7333ad2b4d34f013904ad4737044f366John Hoford 1628921c28c7333ad2b4d34f013904ad4737044f366John Hoford for (int j = 0; j < curve.length; j++) { 1638921c28c7333ad2b4d34f013904ad4737044f366John Hoford System.out.println(curve[j][0] + "," + curve[j][1]); 1648921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 1658921c28c7333ad2b4d34f013904ad4737044f366John Hoford } 1668921c28c7333ad2b4d34f013904ad4737044f366John Hoford} 167