natcubicspline.java

openmap java写的开源数字地图程序. 用applet实现,可以像google map 那样放大缩小地图.

//**********************************************************************
//
//<copyright>
//
//BBN Technologies
//10 Moulton Street
//Cambridge, MA 02138
//[phone]
//
//Copyright (C) BBNT Solutions LLC. All rights reserved.
//
//</copyright>
//**********************************************************************
//
//$Source: /cvs/distapps/openmap/src/openmap/com/bbn/openmap/omGraphics/NatCubicSpline.java,v $
//$RCSfile: NatCubicSpline.java,v $
//$Revision: 1.3.2.1 $
//$Date: 2004/10/14 18:27:25 $
//$Author: dietrick $
//
//**********************************************************************

package com.bbn.openmap.omGraphics;

import java.awt.Polygon;

/**
 * A natural cubic spline calculation.
 * 
 * @author Eric LEPICIER
 * @see <a href="http://www.cse.unsw.edu.au/~lambert/splines/">Splines
 *      </a>
 * @version 21 juil. 2002
 */
public class NatCubicSpline {

    /**
     * Calculates the natural cubic spline that interpolates y[0],
     * y[1], ... y[n]. The first segment is returned as C[0].a +
     * C[0].b*u + C[0].c*u^2 + C[0].d*u^3 0 <=u <1 the other segments
     * are in C[1], C[2], ... C[n-1]
     * 
     * @param n
     * @param x
     * @return Cubic[]
     */
    Cubic[] calcNaturalCubic(int n, int[] x) {
        float[] gamma = new float[n + 1];
        float[] delta = new float[n + 1];
        float[] D = new float[n + 1];
        int i;
        /*
         * We solve the equation [2 1 ] [D[0]] [3(x[1] - x[0]) ] |1 4
         * 1 | |D[1]| |3(x[2] - x[0]) | | 1 4 1 | | . | = | . | |
         * ..... | | . | | . | | 1 4 1| | . | |3(x[n] - x[n-2])| [ 1
         * 2] [D[n]] [3(x[n] - x[n-1])]
         * 
         * by using row operations to convert the matrix to upper
         * triangular and then back sustitution. The D[i] are the
         * derivatives at the knots.
         */

        gamma[0] = 1.0f / 2.0f;
        for (i = 1; i < n; i++) {
            gamma[i] = 1 / (4 - gamma[i - 1]);
        }
        gamma[n] = 1 / (2 - gamma[n - 1]);

        delta[0] = 3 * (x[1] - x[0]) * gamma[0];
        for (i = 1; i < n; i++) {
            delta[i] = (3 * (x[i + 1] - x[i - 1]) - delta[i - 1]) * gamma[i];
        }
        delta[n] = (3 * (x[n] - x[n - 1]) - delta[n - 1]) * gamma[n];

        D[n] = delta[n];
        for (i = n - 1; i >= 0; i--) {
            D[i] = delta[i] - gamma[i] * D[i + 1];
        }

        /* now compute the coefficients of the cubics */
        Cubic[] C = new Cubic[n];
        for (i = 0; i < n; i++) {
            C[i] = new Cubic((float) x[i], D[i], 3 * (x[i + 1] - x[i]) - 2
                    * D[i] - D[i + 1], 2 * (x[i] - x[i + 1]) + D[i] + D[i + 1]);
        }
        return C;
    }

    /**
     * Calculates a cubic spline polyline
     * 
     * @param xpoints
     * @param ypoints
     * @return int[][]
     */
    public int[][] calc(int[] xpoints, int[] ypoints) {
        int[][] res = new int[2][0];
        if (xpoints.length > 2) {
            Cubic[] X = calcNaturalCubic(xpoints.length - 1, xpoints);
            Cubic[] Y = calcNaturalCubic(ypoints.length - 1, ypoints);

            /*
             * very crude technique just break each segment up into
             * steps lines
             */
            Polygon p = new Polygon();
            p.addPoint((int) Math.round(X[0].eval(0)),
                    (int) Math.round(Y[0].eval(0)));
            for (int i = 0; i < X.length; i++) {
                for (int j = 1; j <= steps; j++) {
                    float u = j / (float) steps;
                    p.addPoint(Math.round(X[i].eval(u)),
                            Math.round(Y[i].eval(u)));
                }
            }

            // copy polygon points to the return array
            res[0] = new int[p.npoints];
            System.arraycopy(p.xpoints, 0, res[0], 0, p.npoints);
            res[1] = new int[p.npoints];
            System.arraycopy(p.ypoints, 0, res[1], 0, p.npoints);

            p = null;
        } else {
            res[0] = xpoints;
            res[1] = ypoints;
        }
        return res;
    }

    /**
     * Calculates a float lat/lon cubic spline
     * 
     * @param llpoints
     * @param precision for dividing floating coordinates to become
     *        int, e.g 0.01 means spline to be calculated with
     *        coordinates * 100
     * @return float[]
     */
    public float[] calc(float[] llpoints, float precision) {
        float[] res;
        if (llpoints.length > 4) { // 2 points

            int[] xpoints = new int[(int) (llpoints.length / 2)];
            int[] ypoints = new int[xpoints.length];
            for (int i = 0, j = 0; i < llpoints.length; i += 2, j++) {
                xpoints[j] = (int) (llpoints[i] / precision);
                ypoints[j] = (int) (llpoints[i + 1] / precision);
            }

            Cubic[] X = calcNaturalCubic(xpoints.length - 1, xpoints);
            Cubic[] Y = calcNaturalCubic(ypoints.length - 1, ypoints);

            /*
             * very crude technique just break each segment up into
             * steps lines
             */
            Polygon p = new Polygon();
            p.addPoint((int) Math.round(X[0].eval(0)),
                    (int) Math.round(Y[0].eval(0)));
            for (int i = 0; i < X.length; i++) {
                for (int j = 1; j <= steps; j++) {
                    float u = j / (float) steps;
                    p.addPoint(Math.round(X[i].eval(u)),
                            Math.round(Y[i].eval(u)));
                }
            }

            res = new float[p.npoints * 2];
            for (int i = 0, j = 0; i < p.npoints; i++, j += 2) {
                res[j] = (float) p.xpoints[i] * precision;
                res[j + 1] = (float) p.ypoints[i] * precision;
            }

            p = null;
        } else {
            res = llpoints;
        }
        return res;
    }

    /**
     * Returns the steps.
     * 
     * @return int
     */
    public int getSteps() {
        return steps;
    }

    /**
     * Sets the number of points (steps) interpolated on the curve
     * between the original points to draw it as a polyline.
     * 
     * @param steps The steps to set
     */
    public void setSteps(int steps) {
        this.steps = steps > 0 ? steps : 12;
    }

    private int steps = 12;
}