/* * GeoTools - The Open Source Java GIS Toolkit * http://geotools.org * * (C) 2002-2008, Open Source Geospatial Foundation (OSGeo) * * This library is free software; you can redistribute it and/or * modify it under the terms of the GNU Lesser General Public * License as published by the Free Software Foundation; * version 2.1 of the License. * * This library is distributed in the hope that it will be useful, * but WITHOUT ANY WARRANTY; without even the implied warranty of * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU * Lesser General Public License for more details. */ package org.geotools.referencing.operation.builder; import org.geotools.referencing.operation.matrix.GeneralMatrix; import org.geotools.referencing.operation.transform.ProjectiveTransform; import org.opengis.referencing.FactoryException; import org.opengis.referencing.cs.CartesianCS; import org.opengis.referencing.operation.MathTransform; import org.opengis.geometry.DirectPosition; import org.opengis.geometry.MismatchedDimensionException; import org.opengis.geometry.MismatchedReferenceSystemException; import java.util.List; import javax.vecmath.MismatchedSizeException; /** * Builds {@linkplain MathTransform * MathTransform} setup as Projective transformation from a list of * {@linkplain org.geotools.referencing.operation.builder.MappedPosition * MappedPosition}. The calculation uses least square method. The Projective * transform equation: (2D). The calculation uses least square method. * Projective transform equation:
[ x'] [ m00 m01 m02 ] [ x ] * [ y'] = [ m10 m11 m12 ] [ y ] * [ 1 ] [ m20 m21 1 ] [ 1 ] x' = m * x *In the case that we have more identical points we can write it * like this (in Matrix):
* [ x'1 ] [ x1 y1 1 0 0 0 -x'x -x'y] [ m00 ] * [ x'2 ] [ x2 y2 1 0 0 0 -x'x -x'y] [ m01 ] * [ . ] [ . ] [ m02 ] * [ . ] [ . ] * [ m10 ] * [ x'n ] = [ xn yn 1 0 0 0 -x'x -x'y] [ m11 ] * [ y'1 ] [ 0 0 0 x1 y1 1 -y'x -y'y] [ m12 ] * [ y'2 ] [ 0 0 0 x2 y2 1 -y'x -y'y] [ m20 ] * [ . ] [ . ] [ m21 ] * [ . ] [ . ] * [ y'n ] [ 0 0 0 xn yn 1 -y'x -y'y] * x' = A*mUsing the least square method we get this result: *
* * @author Jan Jezek * * * @source $URL$ * @version $Id$ * @since 2.4 */ public class ProjectiveTransformBuilder extends MathTransformBuilder { /** Matrix of derivations */ protected GeneralMatrix A; /** Matrix of wights */ protected GeneralMatrix P = null; /** Matrix of target values */ protected GeneralMatrix X; protected ProjectiveTransformBuilder() { } /** * Creates ProjectiveTransformBuilder for the set of properties. * * * @param vectors list of {@linkplain MappedPosition * MappedPosition} * @throws MismatchedSizeException * if the number of properties is not set properly. * @throws MismatchedDimensionException * if the dimension of properties is not set properly. * @throws MismatchedReferenceSystemException * -if there is mismatch in coordinate system in * {@linkplain MappedPosition MappedPosition} */ public ProjectiveTransformBuilder(List* m = (ATPA)-1 ATPx'
* * @return m matrix. */ protected double[] calculateLSM() { fillAMatrix(); // fillPMatrix(); fillXMatrix(); if (P == null) { try { includeWeights(false); } catch (FactoryException e) { // should never reach here - weights are not included } } GeneralMatrix AT = (GeneralMatrix) A.clone(); AT.transpose(); GeneralMatrix ATP = new GeneralMatrix(AT.getNumRow(), P.getNumCol()); GeneralMatrix ATPA = new GeneralMatrix(AT.getNumRow(), A.getNumCol()); GeneralMatrix ATPX = new GeneralMatrix(AT.getNumRow(), 1); GeneralMatrix x = new GeneralMatrix(A.getNumCol(), 1); ATP.mul(AT, P); // ATP ATPA.mul(ATP, A); // ATPA ATPX.mul(ATP, X); // ATPX ATPA.invert(); x.mul(ATPA, ATPX); ATPA.invert(); x.transpose(); return x.getElements()[0]; } /** * Returns the matrix of parameters for Projective transformation. * This method should by override for the special cases like affine or * similar transformation. The M matrix looks like this:* m = (ATA)-1 ATx' *
* * [ m00 m01 m02 ] * [ m10 m11 m12 ] * [ m20 m21 1 ] ** * @return Matrix M */ protected GeneralMatrix getProjectiveMatrix() { GeneralMatrix M = new GeneralMatrix(3, 3); // double[] param = generateMMatrix(); double[] param = calculateLSM(); double[] m0 = { param[0], param[1], param[2] }; double[] m1 = { param[3], param[4], param[5] }; double[] m2 = { param[6], param[7], 1 }; M.setRow(0, m0); M.setRow(1, m1); M.setRow(2, m2); return M; } protected MathTransform computeMathTransform() { return ProjectiveTransform.create(getProjectiveMatrix()); } }