1 : /******************************************************************************
2 : * $Id: geo_trans.c 1568 2009-04-22 21:10:55Z warmerdam $
3 : *
4 : * Project: libgeotiff
5 : * Purpose: Code to abstract translation between pixel/line and PCS
6 : * coordinates.
7 : * Author: Frank Warmerdam, warmerda@home.com
8 : *
9 : ******************************************************************************
10 : * Copyright (c) 1999, Frank Warmerdam
11 : *
12 : * Permission is hereby granted, free of charge, to any person obtaining a
13 : * copy of this software and associated documentation files (the "Software"),
14 : * to deal in the Software without restriction, including without limitation
15 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
16 : * and/or sell copies of the Software, and to permit persons to whom the
17 : * Software is furnished to do so, subject to the following conditions:
18 : *
19 : * The above copyright notice and this permission notice shall be included
20 : * in all copies or substantial portions of the Software.
21 : *
22 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
23 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
24 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
25 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
26 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
27 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
28 : * DEALINGS IN THE SOFTWARE.
29 : *****************************************************************************/
30 :
31 : #include "geotiff.h"
32 : #include "geo_tiffp.h" /* external TIFF interface */
33 : #include "geo_keyp.h" /* private interface */
34 : #include "geokeys.h"
35 :
36 : /************************************************************************/
37 : /* inv_geotransform() */
38 : /* */
39 : /* Invert a 6 term geotransform style matrix. */
40 : /************************************************************************/
41 :
42 0 : static int inv_geotransform( double *gt_in, double *gt_out )
43 :
44 : {
45 : double det, inv_det;
46 :
47 : /* we assume a 3rd row that is [0 0 1] */
48 :
49 : /* Compute determinate */
50 :
51 0 : det = gt_in[0] * gt_in[4] - gt_in[1] * gt_in[3];
52 :
53 0 : if( fabs(det) < 0.000000000000001 )
54 0 : return 0;
55 :
56 0 : inv_det = 1.0 / det;
57 :
58 : /* compute adjoint, and devide by determinate */
59 :
60 0 : gt_out[0] = gt_in[4] * inv_det;
61 0 : gt_out[3] = -gt_in[3] * inv_det;
62 :
63 0 : gt_out[1] = -gt_in[1] * inv_det;
64 0 : gt_out[4] = gt_in[0] * inv_det;
65 :
66 0 : gt_out[2] = ( gt_in[1] * gt_in[5] - gt_in[2] * gt_in[4]) * inv_det;
67 0 : gt_out[5] = (-gt_in[0] * gt_in[5] + gt_in[2] * gt_in[3]) * inv_det;
68 :
69 0 : return 1;
70 : }
71 :
72 : /************************************************************************/
73 : /* GTIFTiepointTranslate() */
74 : /************************************************************************/
75 :
76 0 : int GTIFTiepointTranslate( int gcp_count, double * gcps_in, double * gcps_out,
77 : double x_in, double y_in,
78 : double *x_out, double *y_out )
79 :
80 : {
81 : (void) gcp_count;
82 : (void) gcps_in;
83 : (void) gcps_out;
84 : (void) x_in;
85 : (void) y_in;
86 : (void) x_out;
87 : (void) y_out;
88 :
89 : /* I would appreciate a _brief_ block of code for doing second order
90 : polynomial regression here! */
91 0 : return FALSE;
92 : }
93 :
94 :
95 : /************************************************************************/
96 : /* GTIFImageToPCS() */
97 : /************************************************************************/
98 :
99 : /**
100 : * Translate a pixel/line coordinate to projection coordinates.
101 : *
102 : * At this time this function does not support image to PCS translations for
103 : * tiepoints-only definitions, only pixelscale and transformation matrix
104 : * formulations.
105 : *
106 : * @param gtif The handle from GTIFNew() indicating the target file.
107 : * @param x A pointer to the double containing the pixel offset on input,
108 : * and into which the easting/longitude will be put on completion.
109 : * @param y A pointer to the double containing the line offset on input,
110 : * and into which the northing/latitude will be put on completion.
111 : *
112 : * @return TRUE if the transformation succeeds, or FALSE if it fails. It may
113 : * fail if the file doesn't have properly setup transformation information,
114 : * or it is in a form unsupported by this function.
115 : */
116 :
117 0 : int GTIFImageToPCS( GTIF *gtif, double *x, double *y )
118 :
119 : {
120 0 : int res = FALSE;
121 : int tiepoint_count, count, transform_count;
122 0 : tiff_t *tif=gtif->gt_tif;
123 0 : double *tiepoints = 0;
124 0 : double *pixel_scale = 0;
125 0 : double *transform = 0;
126 :
127 :
128 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_TIEPOINTS,
129 : &tiepoint_count, &tiepoints ))
130 0 : tiepoint_count = 0;
131 :
132 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_PIXELSCALE, &count, &pixel_scale ))
133 0 : count = 0;
134 :
135 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_TRANSMATRIX,
136 : &transform_count, &transform ))
137 0 : transform_count = 0;
138 :
139 : /* -------------------------------------------------------------------- */
140 : /* If the pixelscale count is zero, but we have tiepoints use */
141 : /* the tiepoint based approach. */
142 : /* -------------------------------------------------------------------- */
143 0 : if( tiepoint_count > 6 && count == 0 )
144 : {
145 0 : res = GTIFTiepointTranslate( tiepoint_count / 6,
146 : tiepoints, tiepoints + 3,
147 : *x, *y, x, y );
148 : }
149 :
150 : /* -------------------------------------------------------------------- */
151 : /* If we have a transformation matrix, use it. */
152 : /* -------------------------------------------------------------------- */
153 0 : else if( transform_count == 16 )
154 : {
155 0 : double x_in = *x, y_in = *y;
156 :
157 0 : *x = x_in * transform[0] + y_in * transform[1] + transform[3];
158 0 : *y = x_in * transform[4] + y_in * transform[5] + transform[7];
159 :
160 0 : res = TRUE;
161 : }
162 :
163 : /* -------------------------------------------------------------------- */
164 : /* For now we require one tie point, and a valid pixel scale. */
165 : /* -------------------------------------------------------------------- */
166 0 : else if( count < 3 || tiepoint_count < 6 )
167 : {
168 0 : res = FALSE;
169 : }
170 :
171 : else
172 : {
173 0 : *x = (*x - tiepoints[0]) * pixel_scale[0] + tiepoints[3];
174 0 : *y = (*y - tiepoints[1]) * (-1 * pixel_scale[1]) + tiepoints[4];
175 :
176 0 : res = TRUE;
177 : }
178 :
179 : /* -------------------------------------------------------------------- */
180 : /* Cleanup */
181 : /* -------------------------------------------------------------------- */
182 0 : if(tiepoints)
183 0 : _GTIFFree(tiepoints);
184 0 : if(pixel_scale)
185 0 : _GTIFFree(pixel_scale);
186 0 : if(transform)
187 0 : _GTIFFree(transform);
188 :
189 0 : return res;
190 : }
191 :
192 : /************************************************************************/
193 : /* GTIFPCSToImage() */
194 : /************************************************************************/
195 :
196 : /**
197 : * Translate a projection coordinate to pixel/line coordinates.
198 : *
199 : * At this time this function does not support PCS to image translations for
200 : * tiepoints-only based definitions, only matrix and pixelscale/tiepoints
201 : * formulations are supposed.
202 : *
203 : * @param gtif The handle from GTIFNew() indicating the target file.
204 : * @param x A pointer to the double containing the pixel offset on input,
205 : * and into which the easting/longitude will be put on completion.
206 : * @param y A pointer to the double containing the line offset on input,
207 : * and into which the northing/latitude will be put on completion.
208 : *
209 : * @return TRUE if the transformation succeeds, or FALSE if it fails. It may
210 : * fail if the file doesn't have properly setup transformation information,
211 : * or it is in a form unsupported by this function.
212 : */
213 :
214 0 : int GTIFPCSToImage( GTIF *gtif, double *x, double *y )
215 :
216 : {
217 0 : double *tiepoints = NULL;
218 0 : int tiepoint_count, count, transform_count = 0;
219 0 : double *pixel_scale = NULL;
220 0 : double *transform = NULL;
221 0 : tiff_t *tif=gtif->gt_tif;
222 0 : int result = FALSE;
223 :
224 : /* -------------------------------------------------------------------- */
225 : /* Fetch tiepoints and pixel scale. */
226 : /* -------------------------------------------------------------------- */
227 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_TIEPOINTS,
228 : &tiepoint_count, &tiepoints ))
229 0 : tiepoint_count = 0;
230 :
231 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_PIXELSCALE, &count, &pixel_scale ))
232 0 : count = 0;
233 :
234 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_TRANSMATRIX,
235 : &transform_count, &transform ))
236 0 : transform_count = 0;
237 :
238 : /* -------------------------------------------------------------------- */
239 : /* If the pixelscale count is zero, but we have tiepoints use */
240 : /* the tiepoint based approach. */
241 : /* -------------------------------------------------------------------- */
242 0 : if( tiepoint_count > 6 && count == 0 )
243 : {
244 0 : result = GTIFTiepointTranslate( tiepoint_count / 6,
245 : tiepoints + 3, tiepoints,
246 : *x, *y, x, y );
247 : }
248 :
249 : /* -------------------------------------------------------------------- */
250 : /* Handle matrix - convert to "geotransform" format, invert and */
251 : /* apply. */
252 : /* -------------------------------------------------------------------- */
253 0 : else if( transform_count == 16 )
254 : {
255 0 : double x_in = *x, y_in = *y;
256 : double gt_in[6], gt_out[6];
257 :
258 0 : gt_in[0] = transform[0];
259 0 : gt_in[1] = transform[1];
260 0 : gt_in[2] = transform[3];
261 0 : gt_in[3] = transform[4];
262 0 : gt_in[4] = transform[5];
263 0 : gt_in[5] = transform[7];
264 :
265 0 : if( !inv_geotransform( gt_in, gt_out ) )
266 0 : result = FALSE;
267 : else
268 : {
269 0 : *x = x_in * gt_out[0] + y_in * gt_out[1] + gt_out[2];
270 0 : *y = x_in * gt_out[3] + y_in * gt_out[4] + gt_out[5];
271 :
272 0 : result = TRUE;
273 : }
274 : }
275 :
276 : /* -------------------------------------------------------------------- */
277 : /* For now we require one tie point, and a valid pixel scale. */
278 : /* -------------------------------------------------------------------- */
279 0 : else if( count >= 3 && tiepoint_count >= 6 )
280 : {
281 0 : *x = (*x - tiepoints[3]) / pixel_scale[0] + tiepoints[0];
282 0 : *y = (*y - tiepoints[4]) / (-1 * pixel_scale[1]) + tiepoints[1];
283 :
284 0 : result = TRUE;
285 : }
286 :
287 : /* -------------------------------------------------------------------- */
288 : /* Cleanup. */
289 : /* -------------------------------------------------------------------- */
290 0 : if(tiepoints)
291 0 : _GTIFFree(tiepoints);
292 0 : if(pixel_scale)
293 0 : _GTIFFree(pixel_scale);
294 0 : if(transform)
295 0 : _GTIFFree(transform);
296 :
297 0 : return result;
298 : }
299 :
|
