001 //$HeadURL: svn+ssh://rbezema@svn.wald.intevation.org/deegree/base/tags/2.1/src/org/deegree/model/csct/ct/MatrixTransform.java $
002 /*---------------- FILE HEADER ------------------------------------------
003
004 This file is part of deegree.
005 Copyright (C) 2001 by:
006 EXSE, Department of Geography, University of Bonn
007 http://www.giub.uni-bonn.de/exse/
008 lat/lon GmbH
009 http://www.lat-lon.de
010
011 It has been implemented within SEAGIS - An OpenSource implementation of OpenGIS specification
012 (C) 2001, Institut de Recherche pour le D�veloppement (http://sourceforge.net/projects/seagis/)
013 SEAGIS Contacts: Surveillance de l'Environnement Assist�e par Satellite
014 Institut de Recherche pour le D�veloppement / US-Espace
015 mailto:seasnet@teledetection.fr
016
017
018 This library is free software; you can redistribute it and/or
019 modify it under the terms of the GNU Lesser General Public
020 License as published by the Free Software Foundation; either
021 version 2.1 of the License, or (at your option) any later version.
022
023 This library is distributed in the hope that it will be useful,
024 but WITHOUT ANY WARRANTY; without even the implied warranty of
025 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
026 Lesser General Public License for more details.
027
028 You should have received a copy of the GNU Lesser General Public
029 License along with this library; if not, write to the Free Software
030 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
031
032 Contact:
033
034 Andreas Poth
035 lat/lon GmbH
036 Aennchenstr. 19
037 53115 Bonn
038 Germany
039 E-Mail: poth@lat-lon.de
040
041 Klaus Greve
042 Department of Geography
043 University of Bonn
044 Meckenheimer Allee 166
045 53115 Bonn
046 Germany
047 E-Mail: klaus.greve@uni-bonn.de
048
049
050 ---------------------------------------------------------------------------*/
051 package org.deegree.model.csct.ct;
052
053 // Geometry
054 import java.awt.geom.Point2D;
055 import java.io.Serializable;
056 import java.util.Arrays;
057
058 import javax.media.jai.ParameterList;
059 import javax.media.jai.PerspectiveTransform;
060 import javax.media.jai.util.Range;
061 import javax.vecmath.GMatrix;
062 import javax.vecmath.SingularMatrixException;
063
064 import org.deegree.model.csct.pt.CoordinatePoint;
065 import org.deegree.model.csct.pt.Matrix;
066 import org.deegree.model.csct.resources.css.ResourceKeys;
067 import org.deegree.model.csct.resources.css.Resources;
068
069 /**
070 * Transforms multi-dimensional coordinate points using a {@link Matrix}.
071 *
072 * @version 1.00
073 * @author OpenGIS (www.opengis.org)
074 * @author Martin Desruisseaux
075 */
076 final class MatrixTransform extends AbstractMathTransform implements Serializable {
077 /**
078 * Serial number for interoperability with different versions.
079 */
080 private static final long serialVersionUID = -2104496465933824935L;
081
082 /**
083 * the number of rows.
084 */
085 private final int numRow;
086
087 /**
088 * the number of columns.
089 */
090 private final int numCol;
091
092 /**
093 * Elements of the matrix. Column indice vary fastest.
094 */
095 private final double[] elt;
096
097 /**
098 * Construct a transform.
099 */
100 protected MatrixTransform( final GMatrix matrix ) {
101 numRow = matrix.getNumRow();
102 numCol = matrix.getNumCol();
103 elt = new double[numRow * numCol];
104 int index = 0;
105 for ( int j = 0; j < numRow; j++ )
106 for ( int i = 0; i < numCol; i++ )
107 elt[index++] = matrix.getElement( j, i );
108 }
109
110 /**
111 * Transforms an array of floating point coordinates by this matrix. Point coordinates must have
112 * a dimension equals to <code>{@link Matrix#getNumCol()}-1</code>. For example, for square
113 * matrix of size 4×4, coordinate points are three-dimensional and stored in the arrays
114 * starting at the specified offset (<code>srcOff</code>) in the order
115 * <code>[x<sub>0</sub>, y<sub>0</sub>, z<sub>0</sub>,
116 * x<sub>1</sub>, y<sub>1</sub>, z<sub>1</sub>...,
117 * x<sub>n</sub>, y<sub>n</sub>, z<sub>n</sub>]</code>.
118 *
119 * The transformed points <code>(x',y',z')</code> are computed as below (note that this
120 * computation is similar to {@link PerspectiveTransform}):
121 *
122 * <blockquote>
123 *
124 * <pre>
125 * [ u ] [ m<sub>00</sub> m<sub>01</sub> m<sub>02</sub> m<sub>03</sub> ] [ x ]
126 * [ v ] = [ m<sub>10</sub> m<sub>11</sub> m<sub>12</sub> m<sub>13</sub> ] [ y ]
127 * [ w ] [ m<sub>20</sub> m<sub>21</sub> m<sub>22</sub> m<sub>23</sub> ] [ z ]
128 * [ t ] [ m<sub>30</sub> m<sub>31</sub> m<sub>32</sub> m<sub>33</sub> ] [ 1 ]
129 *
130 * x' = u/t
131 * y' = v/t
132 * y' = w/t
133 * </pre>
134 *
135 * </blockquote>
136 *
137 * @param srcPts
138 * The array containing the source point coordinates.
139 * @param srcOff
140 * The offset to the first point to be transformed in the source array.
141 * @param dstPts
142 * The array into which the transformed point coordinates are returned.
143 * @param dstOff
144 * The offset to the location of the first transformed point that is stored in the
145 * destination array. The source and destination array sections can be overlaps.
146 * @param numPts
147 * The number of points to be transformed
148 */
149 public void transform( float[] srcPts, int srcOff, final float[] dstPts, int dstOff, int numPts ) {
150 final int inputDimension = numCol - 1; // The last ordinate will be assumed equals to 1.
151 final int outputDimension = numRow - 1;
152 final double[] buffer = new double[numRow];
153 if ( srcPts == dstPts ) {
154 // We are going to write in the source array. Checks if
155 // source and destination sections are going to clash.
156 final int upperSrc = srcOff + numPts * inputDimension;
157 if ( upperSrc > dstOff ) {
158 if ( inputDimension >= outputDimension ? dstOff > srcOff
159 : dstOff + numPts * outputDimension > upperSrc ) {
160 // If source overlaps destination, then the easiest workaround is
161 // to copy source data. This is not the most efficient however...
162 srcPts = new float[numPts * inputDimension];
163 System.arraycopy( dstPts, srcOff, srcPts, 0, srcPts.length );
164 srcOff = 0;
165 }
166 }
167 }
168 while ( --numPts >= 0 ) {
169 int mix = 0;
170 for ( int j = 0; j < numRow; j++ ) {
171 double sum = elt[mix + inputDimension];
172 for ( int i = 0; i < inputDimension; i++ ) {
173 sum += srcPts[srcOff + i] * elt[mix++];
174 }
175 buffer[j] = sum;
176 mix++;
177 }
178 final double w = buffer[outputDimension];
179 for ( int j = 0; j < outputDimension; j++ ) {
180 // 'w' is equals to 1 if the transform is affine.
181 dstPts[dstOff++] = (float) ( buffer[j] / w );
182 }
183 srcOff += inputDimension;
184 }
185 }
186
187 /**
188 * Transforms an array of floating point coordinates by this matrix. Point coordinates must have
189 * a dimension equals to <code>{@link Matrix#getNumCol()}-1</code>. For example, for square
190 * matrix of size 4×4, coordinate points are three-dimensional and stored in the arrays
191 * starting at the specified offset (<code>srcOff</code>) in the order
192 * <code>[x<sub>0</sub>, y<sub>0</sub>, z<sub>0</sub>,
193 * x<sub>1</sub>, y<sub>1</sub>, z<sub>1</sub>...,
194 * x<sub>n</sub>, y<sub>n</sub>, z<sub>n</sub>]</code>.
195 *
196 * The transformed points <code>(x',y',z')</code> are computed as below (note that this
197 * computation is similar to {@link PerspectiveTransform}):
198 *
199 * <blockquote>
200 *
201 * <pre>
202 * [ u ] [ m<sub>00</sub> m<sub>01</sub> m<sub>02</sub> m<sub>03</sub> ] [ x ]
203 * [ v ] = [ m<sub>10</sub> m<sub>11</sub> m<sub>12</sub> m<sub>13</sub> ] [ y ]
204 * [ w ] [ m<sub>20</sub> m<sub>21</sub> m<sub>22</sub> m<sub>23</sub> ] [ z ]
205 * [ t ] [ m<sub>30</sub> m<sub>31</sub> m<sub>32</sub> m<sub>33</sub> ] [ 1 ]
206 *
207 * x' = u/t
208 * y' = v/t
209 * y' = w/t
210 * </pre>
211 *
212 * </blockquote>
213 *
214 * @param srcPts
215 * The array containing the source point coordinates.
216 * @param srcOff
217 * The offset to the first point to be transformed in the source array.
218 * @param dstPts
219 * The array into which the transformed point coordinates are returned.
220 * @param dstOff
221 * The offset to the location of the first transformed point that is stored in the
222 * destination array. The source and destination array sections can be overlaps.
223 * @param numPts
224 * The number of points to be transformed
225 */
226 public void transform( double[] srcPts, int srcOff, final double[] dstPts, int dstOff,
227 int numPts ) {
228 final int inputDimension = numCol - 1; // The last ordinate will be assumed equals to 1.
229 final int outputDimension = numRow - 1;
230 final double[] buffer = new double[numRow];
231 if ( srcPts == dstPts ) {
232 // We are going to write in the source array. Checks if
233 // source and destination sections are going to clash.
234 final int upperSrc = srcOff + numPts * inputDimension;
235 if ( upperSrc > dstOff ) {
236 if ( inputDimension >= outputDimension ? dstOff > srcOff
237 : dstOff + numPts * outputDimension > upperSrc ) {
238 // If source overlaps destination, then the easiest workaround is
239 // to copy source data. This is not the most efficient however...
240 srcPts = new double[numPts * inputDimension];
241 System.arraycopy( dstPts, srcOff, srcPts, 0, srcPts.length );
242 srcOff = 0;
243 }
244 }
245 }
246 while ( --numPts >= 0 ) {
247 int mix = 0;
248 for ( int j = 0; j < numRow; j++ ) {
249 double sum = elt[mix + inputDimension];
250 for ( int i = 0; i < inputDimension; i++ ) {
251 sum += srcPts[srcOff + i] * elt[mix++];
252 }
253 buffer[j] = sum;
254 mix++;
255 }
256 final double w = buffer[outputDimension];
257 for ( int j = 0; j < outputDimension; j++ ) {
258 // 'w' is equals to 1 if the transform is affine.
259 dstPts[dstOff++] = buffer[j] / w;
260 }
261 srcOff += inputDimension;
262 }
263 }
264
265 /**
266 * Gets the derivative of this transform at a point. For a matrix transform, the derivative is
267 * the same everywhere.
268 */
269 public Matrix derivative( final Point2D point ) {
270 return derivative( (CoordinatePoint) null );
271 }
272
273 /**
274 * Gets the derivative of this transform at a point. For a matrix transform, the derivative is
275 * the same everywhere.
276 */
277 public Matrix derivative( final CoordinatePoint point ) {
278 final Matrix matrix = getMatrix();
279 matrix.setSize( numRow - 1, numCol - 1 );
280 return matrix;
281 }
282
283 /**
284 * Returns a copy of the matrix.
285 */
286 public Matrix getMatrix() {
287 return new Matrix( numRow, numCol, elt );
288 }
289
290 /**
291 * Gets the dimension of input points.
292 */
293 public int getDimSource() {
294 return numCol - 1;
295 }
296
297 /**
298 * Gets the dimension of output points.
299 */
300 public int getDimTarget() {
301 return numRow - 1;
302 }
303
304 /**
305 * Tests whether this transform does not move any points.
306 */
307 public boolean isIdentity() {
308 if ( numRow != numCol )
309 return false;
310
311 int index = 0;
312 for ( int j = 0; j < numRow; j++ )
313 for ( int i = 0; i < numCol; i++ )
314 if ( elt[index++] != ( i == j ? 1 : 0 ) )
315 return false;
316 return true;
317 }
318
319 /**
320 * Creates the inverse transform of this object.
321 */
322 public MathTransform inverse()
323 throws NoninvertibleTransformException {
324 if ( isIdentity() )
325 return this;
326 final Matrix matrix = getMatrix();
327 try {
328 matrix.invert();
329 } catch ( SingularMatrixException exception ) {
330 NoninvertibleTransformException e = new NoninvertibleTransformException(
331 Resources.format( ResourceKeys.ERROR_NONINVERTIBLE_TRANSFORM ) );
332 throw e;
333 }
334 return new MatrixTransform( matrix );
335 }
336
337 /**
338 * Returns a hash value for this transform. This value need not remain consistent between
339 * different implementations of the same class.
340 */
341 public int hashCode() {
342 long code = 2563217;
343 for ( int i = elt.length; --i >= 0; ) {
344 code = code * 37 + Double.doubleToLongBits( elt[i] );
345 }
346 return (int) ( code >>> 32 ) ^ (int) code;
347 }
348
349 /**
350 * Compares the specified object with this math transform for equality.
351 */
352 public boolean equals( final Object object ) {
353 if ( object == this )
354 return true; // Slight optimization
355 if ( super.equals( object ) ) {
356 final MatrixTransform that = (MatrixTransform) object;
357 return this.numRow == that.numRow && this.numCol == that.numCol
358 && Arrays.equals( this.elt, that.elt );
359 }
360 return false;
361 }
362
363 /**
364 * Returns the WKT for this math transform.
365 */
366 public String toString() {
367 return toString( getMatrix() );
368 }
369
370 /**
371 * Returns the WKT for an affine transform using the specified matrix.
372 */
373 static String toString( final Matrix matrix ) {
374 final int numRow = matrix.getNumRow();
375 final int numCol = matrix.getNumCol();
376 final StringBuffer buffer = paramMT( "Affine" );
377 final StringBuffer eltBuf = new StringBuffer( "elt_" );
378 addParameter( buffer, "Num_row", numRow );
379 addParameter( buffer, "Num_col", numCol );
380 for ( int j = 0; j < numRow; j++ ) {
381 for ( int i = 0; i < numCol; i++ ) {
382 final double value = matrix.getElement( j, i );
383 if ( value != ( i == j ? 1 : 0 ) ) {
384 eltBuf.setLength( 4 );
385 eltBuf.append( j );
386 eltBuf.append( '_' );
387 eltBuf.append( i );
388 addParameter( buffer, eltBuf.toString(), value );
389 }
390 }
391 }
392 buffer.append( ']' );
393 return buffer.toString();
394 }
395
396 /**
397 * The provider for {@link MatrixTransform}.
398 *
399 * @version 1.0
400 * @author Martin Desruisseaux
401 */
402 static final class Provider extends MathTransformProvider {
403 /**
404 * Range of positives values. Range goes from 1 to the maximum value.
405 */
406 private static final Range POSITIVE_RANGE = new Range( Integer.class, new Integer( 1 ),
407 new Integer( Integer.MAX_VALUE ) );
408
409 /**
410 * Create a provider for affine transforms of the specified dimension. Created affine
411 * transforms will have a size of <code>numRow × numCol</code>.
412 *
413 * @param numRow
414 * The number of matrix's rows.
415 * @param numCol
416 * The number of matrix's columns.
417 */
418 public Provider( final int numRow, final int numCol ) {
419 super( "Affine", ResourceKeys.AFFINE_TRANSFORM, null );
420 putInt( "Num_row", numRow, POSITIVE_RANGE ); // Add integer (not double) parameter
421 putInt( "Num_col", numCol, POSITIVE_RANGE ); // Add integer (not double) parameter
422 final StringBuffer buffer = new StringBuffer( "elt_" );
423 for ( int j = 0; j <= numRow; j++ ) {
424 for ( int i = 0; i <= numCol; i++ ) {
425 buffer.setLength( 4 );
426 buffer.append( j );
427 buffer.append( '_' );
428 buffer.append( i );
429 put( buffer.toString(), ( i == j ) ? 1.0 : 0.0, null );
430 }
431 }
432 }
433
434 /**
435 * Returns a transform for the specified parameters.
436 *
437 * @param parameters
438 * The parameter values in standard units.
439 * @return A {@link MathTransform} object of this classification.
440 */
441 public MathTransform create( final ParameterList parameters ) {
442 return staticCreate( parameters );
443 }
444
445 /**
446 * Static version of {@link #create}, for use by
447 * {@link MathTransformFactory#createParameterizedTransform(String, ParameterList)}.
448 */
449 public static MathTransform staticCreate( final ParameterList parameters ) {
450 final int numRow = parameters.getIntParameter( "Num_row" );
451 final int numCol = parameters.getIntParameter( "Num_col" );
452 final Matrix matrix = new Matrix( numRow, numCol );
453 final String[] names = parameters.getParameterListDescriptor().getParamNames();
454 if ( names != null ) {
455 for ( int i = 0; i < names.length; i++ ) {
456 final String name = names[i];
457 if ( name.regionMatches( true, 0, "elt_", 0, 4 ) ) {
458 final int separator = name.lastIndexOf( '_' );
459 final int row = Integer.parseInt( name.substring( 4, separator ) );
460 final int col = Integer.parseInt( name.substring( separator + 1 ) );
461 matrix.setElement( row, col, parameters.getDoubleParameter( name ) );
462 }
463 }
464 }
465 if ( numRow == 3 && matrix.isAffine() ) {
466 return new AffineTransform2D( matrix.toAffineTransform2D() );
467 }
468 return new MatrixTransform( matrix );
469 }
470 }
471 }