pnt.java

Voronoi图生成算法 ,Voronoi(沃洛诺依)多边形网络常常被用来做为计算机仿真模型.由于按照定义形成Voronoi多边形网络较难用程序实现,人们多采用Delauney(狄洛尼)三角形的方法完

            throw new IllegalArgumentException("Matrix is not square");
        boolean[] columns = new boolean[matrix.length];
        for (int i = 0; i < matrix.length; i++) columns[i] = true;
        try {return determinant(matrix, 0, columns);}
        catch (ArrayIndexOutOfBoundsException e) {
            throw new IllegalArgumentException("Matrix is wrong shape");
        }
    }
    
    /**
     * Compute the determinant of a submatrix specified by starting row
     * and by "active" columns.
     * @param matrix the matrix as an array of Pnts
     * @param row the starting row
     * @param columns a boolean array indicating the "active" columns
     * @return the determinant of the specified submatrix
     * @throws ArrayIndexOutOfBoundsException if dimensions are wrong
     */
    private static double determinant(Pnt[] matrix, int row, boolean[] columns) {
        if (row == matrix.length) return 1;
        double sum = 0;
        int sign = 1;
        for (int col = 0; col < columns.length; col++) {
            if (!columns[col]) continue;
            columns[col] = false;
            sum += sign * matrix[row].coordinates[col] *
                   determinant(matrix, row+1, columns);
            columns[col] = true;
            sign = -sign;
        }
        return sum;
    }
    
    /**
     * Compute generalized cross-product of the rows of a matrix.
     * The result is a Pnt perpendicular (as a vector) to each row of
     * the matrix.  This is not an efficient implementation, but should 
     * be adequate for low dimension.
     * @param matrix the matrix of Pnts (one less row than the Pnt dimension)
     * @return a Pnt perpendicular to each row Pnt
     * @throws IllegalArgumentException if matrix is wrong shape
     */
    public static Pnt cross (Pnt[] matrix) {
        int len = matrix.length + 1;
        if (len != matrix[0].dimension())
            throw new IllegalArgumentException("Dimension mismatch");
        boolean[] columns = new boolean[len];
        for (int i = 0; i < len; i++) columns[i] = true;
        double[] result = new double[len];
        int sign = 1;
        try {
            for (int i = 0; i < len; i++) {
                columns[i] = false;
                result[i] = sign * determinant(matrix, 0, columns);
                columns[i] = true;
                sign = -sign;
            }
        } catch (ArrayIndexOutOfBoundsException e) {
            throw new IllegalArgumentException("Matrix is wrong shape");
        }
        return new Pnt(result);
    }
    
    /* Pnts as simplices */
    
    /**
     * Determine the signed content (i.e., area or volume, etc.) of a simplex.
     * @param simplex the simplex (as an array of Pnts)
     * @return the signed content of the simplex
     */
    public static double content (Pnt[] simplex) {
        Pnt[] matrix = new Pnt[simplex.length];
        for (int i = 0; i < matrix.length; i++)
            matrix[i] = simplex[i].extend(new double[] {1});
        int fact = 1;
        for (int i = 1; i < matrix.length; i++) fact = fact*i;
        return determinant(matrix) / fact;
    }
    
    /**
     * Relation between this Pnt and a simplex (represented as an array of Pnts).
     * Result is an array of signs, one for each vertex of the simplex, indicating
     * the relation between the vertex, the vertex's opposite facet, and this
     * Pnt. <pre>
     *   -1 means Pnt is on same side of facet
     *    0 means Pnt is on the facet
     *   +1 means Pnt is on opposite side of facet</pre>
     * @param simplex an array of Pnts representing a simplex
     * @return an array of signs showing relation between this Pnt and the simplex
     * @throws IllegalArgumentExcpetion if the simplex is degenerate
     */
    public int[] relation (Pnt[] simplex) {
        /* In 2D, we compute the cross of this matrix:
         *    1   1   1   1
         *    p0  a0  b0  c0
         *    p1  a1  b1  c1
         * where (a, b, c) is the simplex and p is this Pnt.  The result
         * is a vector in which the first coordinate is the signed area
         * (all signed areas are off by the same constant factor) of
         * the simplex and the remaining coordinates are the *negated*
         * signed areas for the simplices in which p is substituted for
         * each of the vertices. Analogous results occur in higher dimensions.
         */
        int dim = simplex.length - 1;
        if (this.dimension() != dim)
            throw new IllegalArgumentException("Dimension mismatch");
        
        /* Create and load the matrix */
        Pnt[] matrix = new Pnt[dim+1];
        /* First row */
        double[] coords = new double[dim+2];
        for (int j = 0; j < coords.length; j++) coords[j] = 1;
        matrix[0] = new Pnt(coords);
        /* Other rows */
        for (int i = 0; i < dim; i++) {
            coords[0] = this.coordinates[i];
            for (int j = 0; j < simplex.length; j++) 
                coords[j+1] = simplex[j].coordinates[i];
            matrix[i+1] = new Pnt(coords);
        }
        
        /* Compute and analyze the vector of areas/volumes/contents */
        Pnt vector = cross(matrix);
        double content = vector.coordinates[0];
        int[] result = new int[dim+1];
        for (int i = 0; i < result.length; i++) {
            double value = vector.coordinates[i+1];
            if (Math.abs(value) <= 1.0e-6 * Math.abs(content)) result[i] = 0;
            else if (value < 0) result[i] = -1;
            else result[i] = 1;
        }
        if (content < 0) {
            for (int i = 0; i < result.length; i++) result[i] = -result[i];
        }
        if (content == 0) {
            for (int i = 0; i < result.length; i++) result[i] = Math.abs(result[i]);
        }
        return result;
    }
    
    /**
     * Test if this Pnt is outside of simplex.
     * @param simplex the simplex (an array of Pnts)
     * @return the simplex Pnt that "witnesses" outsideness (or null if not outside)
     */
    public Pnt isOutside (Pnt[] simplex) {
        int[] result = this.relation(simplex);
        for (int i = 0; i < result.length; i++) {
            if (result[i] > 0) return simplex[i];
        }
        return null;
    }
    
    /**
     * Test if this Pnt is on a simplex.
     * @param simplex the simplex (an array of Pnts)
     * @return the simplex Pnt that "witnesses" on-ness (or null if not on)
     */
    public Pnt isOn (Pnt[] simplex) {
        int[] result = this.relation(simplex);
        Pnt witness = null;
        for (int i = 0; i < result.length; i++) {
            if (result[i] == 0) witness = simplex[i];
            else if (result[i] > 0) return null;
        }
        return witness;
    }
    
    /**
     * Test if this Pnt is inside a simplex.
     * @param simplex the simplex (an arary of Pnts)
     * @return true iff this Pnt is inside simplex.
     */
    public boolean isInside (Pnt[] simplex) {
        int[] result = this.relation(simplex);
        for (int i = 0; i < result.length; i++) if (result[i] >= 0) return false;
        return true;
    }
    
    /**
     * Test relation between this Pnt and circumcircle of a simplex.
     * @param simplex the simplex (as an array of Pnts)
     * @return -1, 0, or +1 for inside, on, or outside of circumcircle
     */
    public int vsCircumcircle (Pnt[] simplex) {
        Pnt[] matrix = new Pnt[simplex.length + 1];
        for (int i = 0; i < simplex.length; i++)
            matrix[i] = simplex[i].extend(new double[] {1, simplex[i].dot(simplex[i])});
        matrix[simplex.length] = this.extend(new double[] {1, this.dot(this)});
        double d = determinant(matrix);
        int result = (d < 0)? -1 : ((d > 0)? +1 : 0);
        if (content(simplex) < 0) result = - result;
        return result;
    }
    
    /**
     * Circumcenter of a simplex.
     * @param simplex the simplex (as an array of Pnts)
     * @return the circumcenter (a Pnt) of simplex
     */
    public static Pnt circumcenter (Pnt[] simplex) {
        int dim = simplex[0].dimension();
        if (simplex.length - 1 != dim)
            throw new IllegalArgumentException("Dimension mismatch");
        Pnt[] matrix = new Pnt[dim];
        for (int i = 0; i < dim; i++) 
            matrix[i] = simplex[i].bisector(simplex[i+1]);
        Pnt hCenter = cross(matrix);      // Center in homogeneous coordinates
        double last = hCenter.coordinates[dim];
        double[] result = new double[dim];
        for (int i = 0; i < dim; i++) result[i] = hCenter.coordinates[i] / last;
        return new Pnt(result);
    }
     
    /**
     * Main program (used for testing).
     */
    public static void main (String[] args) {
        Pnt p = new Pnt(1, 2, 3);
        System.out.println("Pnt created: " + p);
        Pnt[] matrix1 = {new Pnt(1,2), new Pnt(3,4)};
        Pnt[] matrix2 = {new Pnt(7,0,5), new Pnt(2,4,6), new Pnt(3,8,1)};
        System.out.print("Results should be -2 and -288: ");
        System.out.println(determinant(matrix1) + " " + determinant(matrix2));
        Pnt p1 = new Pnt(1,1); Pnt p2 = new Pnt(-1,1);
        System.out.println("Angle between " + p1 + " and " + p2 + ": " + p1.angle(p2));
        System.out.println(p1 + " subtract " + p2 + ": " + p1.subtract(p2));
        Pnt v0 = new Pnt(0,0), v1 = new Pnt(1,1), v2 = new Pnt(2,2);
        Pnt[] vs = {v0, new Pnt(0,1), new Pnt(1,0)};
        Pnt vp = new Pnt(.1, .1);
        System.out.println(vp + " isInside " + toString(vs) + ": " + vp.isInside(vs));
        System.out.println(v1 + " isInside " + toString(vs) + ": " + v1.isInside(vs));
        System.out.println(vp + " vsCircumcircle " + toString(vs) + ": " +
                           vp.vsCircumcircle(vs));
        System.out.println(v1 + " vsCircumcircle " + toString(vs) + ": " +
                           v1.vsCircumcircle(vs));
        System.out.println(v2 + " vsCircumcircle " + toString(vs) + ": " +
                           v2.vsCircumcircle(vs));
        System.out.println("Circumcenter of " + toString(vs) + " is " + circumcenter(vs));
    }
}