Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: libgeotiff
4 : * Purpose: Code to abstract translation between pixel/line and PCS
5 : * coordinates.
6 : * Author: Frank Warmerdam, warmerda@home.com
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 1999, Frank Warmerdam
10 : *
11 : * SPDX-License-Identifier: MIT
12 : *****************************************************************************/
13 :
14 : #include <math.h>
15 : #include <stddef.h>
16 :
17 : #include "geotiff.h"
18 : #include "geo_tiffp.h" /* external TIFF interface */
19 : #include "geo_keyp.h" /* private interface */
20 : #include "geokeys.h"
21 :
22 : /************************************************************************/
23 : /* inv_geotransform() */
24 : /* */
25 : /* Invert a 6 term geotransform style matrix. */
26 : /************************************************************************/
27 :
28 0 : static int inv_geotransform( double *gt_in, double *gt_out )
29 :
30 : {
31 : /* we assume a 3rd row that is [0 0 1] */
32 :
33 : /* Compute determinate */
34 :
35 0 : const double det = gt_in[0] * gt_in[4] - gt_in[1] * gt_in[3];
36 :
37 0 : if( fabs(det) < 0.000000000000001 )
38 0 : return 0;
39 :
40 0 : const double inv_det = 1.0 / det;
41 :
42 : /* compute adjoint, and divide by determinate */
43 :
44 0 : gt_out[0] = gt_in[4] * inv_det;
45 0 : gt_out[3] = -gt_in[3] * inv_det;
46 :
47 0 : gt_out[1] = -gt_in[1] * inv_det;
48 0 : gt_out[4] = gt_in[0] * inv_det;
49 :
50 0 : gt_out[2] = ( gt_in[1] * gt_in[5] - gt_in[2] * gt_in[4]) * inv_det;
51 0 : gt_out[5] = (-gt_in[0] * gt_in[5] + gt_in[2] * gt_in[3]) * inv_det;
52 :
53 0 : return 1;
54 : }
55 :
56 : /************************************************************************/
57 : /* GTIFTiepointTranslate() */
58 : /************************************************************************/
59 :
60 : static
61 0 : int GTIFTiepointTranslate( int gcp_count, double * gcps_in, double * gcps_out,
62 : double x_in, double y_in,
63 : double *x_out, double *y_out )
64 :
65 : {
66 : (void) gcp_count;
67 : (void) gcps_in;
68 : (void) gcps_out;
69 : (void) x_in;
70 : (void) y_in;
71 : (void) x_out;
72 : (void) y_out;
73 :
74 : /* I would appreciate a _brief_ block of code for doing second order
75 : polynomial regression here! */
76 0 : return FALSE;
77 : }
78 :
79 :
80 : /************************************************************************/
81 : /* GTIFImageToPCS() */
82 : /************************************************************************/
83 :
84 : /**
85 : * Translate a pixel/line coordinate to projection coordinates.
86 : *
87 : * At this time this function does not support image to PCS translations for
88 : * tiepoints-only definitions, only pixelscale and transformation matrix
89 : * formulations.
90 : *
91 : * @param gtif The handle from GTIFNew() indicating the target file.
92 : * @param x A pointer to the double containing the pixel offset on input,
93 : * and into which the easting/longitude will be put on completion.
94 : * @param y A pointer to the double containing the line offset on input,
95 : * and into which the northing/latitude will be put on completion.
96 : *
97 : * @return TRUE if the transformation succeeds, or FALSE if it fails. It may
98 : * fail if the file doesn't have properly setup transformation information,
99 : * or it is in a form unsupported by this function.
100 : */
101 :
102 0 : int GTIFImageToPCS( GTIF *gtif, double *x, double *y )
103 :
104 : {
105 0 : tiff_t *tif=gtif->gt_tif;
106 :
107 : int tiepoint_count;
108 0 : double *tiepoints = 0;
109 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_TIEPOINTS,
110 : &tiepoint_count, &tiepoints ))
111 0 : tiepoint_count = 0;
112 :
113 : int count;
114 0 : double *pixel_scale = 0;
115 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_PIXELSCALE, &count, &pixel_scale ))
116 0 : count = 0;
117 :
118 : int transform_count;
119 0 : double *transform = 0;
120 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_TRANSMATRIX,
121 : &transform_count, &transform ))
122 0 : transform_count = 0;
123 :
124 : /* -------------------------------------------------------------------- */
125 : /* If the pixelscale count is zero, but we have tiepoints use */
126 : /* the tiepoint based approach. */
127 : /* -------------------------------------------------------------------- */
128 0 : int res = FALSE;
129 0 : if( tiepoint_count > 6 && count == 0 )
130 : {
131 0 : res = GTIFTiepointTranslate( tiepoint_count / 6,
132 : tiepoints, tiepoints + 3,
133 : *x, *y, x, y );
134 : }
135 :
136 : /* -------------------------------------------------------------------- */
137 : /* If we have a transformation matrix, use it. */
138 : /* -------------------------------------------------------------------- */
139 0 : else if( transform_count == 16 )
140 : {
141 0 : const double x_in = *x;
142 0 : const double y_in = *y;
143 :
144 0 : *x = x_in * transform[0] + y_in * transform[1] + transform[3];
145 0 : *y = x_in * transform[4] + y_in * transform[5] + transform[7];
146 :
147 0 : res = TRUE;
148 : }
149 :
150 : /* -------------------------------------------------------------------- */
151 : /* For now we require one tie point, and a valid pixel scale. */
152 : /* -------------------------------------------------------------------- */
153 0 : else if( count < 3 || tiepoint_count < 6 )
154 : {
155 0 : res = FALSE;
156 : }
157 :
158 : else
159 : {
160 0 : *x = (*x - tiepoints[0]) * pixel_scale[0] + tiepoints[3];
161 0 : *y = (*y - tiepoints[1]) * (-1 * pixel_scale[1]) + tiepoints[4];
162 :
163 0 : res = TRUE;
164 : }
165 :
166 : /* -------------------------------------------------------------------- */
167 : /* Cleanup */
168 : /* -------------------------------------------------------------------- */
169 0 : if(tiepoints)
170 0 : _GTIFFree(tiepoints);
171 0 : if(pixel_scale)
172 0 : _GTIFFree(pixel_scale);
173 0 : if(transform)
174 0 : _GTIFFree(transform);
175 :
176 0 : return res;
177 : }
178 :
179 : /************************************************************************/
180 : /* GTIFPCSToImage() */
181 : /************************************************************************/
182 :
183 : /**
184 : * Translate a projection coordinate to pixel/line coordinates.
185 : *
186 : * At this time this function does not support PCS to image translations for
187 : * tiepoints-only based definitions, only matrix and pixelscale/tiepoints
188 : * formulations are supposed.
189 : *
190 : * @param gtif The handle from GTIFNew() indicating the target file.
191 : * @param x A pointer to the double containing the pixel offset on input,
192 : * and into which the easting/longitude will be put on completion.
193 : * @param y A pointer to the double containing the line offset on input,
194 : * and into which the northing/latitude will be put on completion.
195 : *
196 : * @return TRUE if the transformation succeeds, or FALSE if it fails. It may
197 : * fail if the file doesn't have properly setup transformation information,
198 : * or it is in a form unsupported by this function.
199 : */
200 :
201 0 : int GTIFPCSToImage( GTIF *gtif, double *x, double *y )
202 :
203 : {
204 0 : tiff_t *tif=gtif->gt_tif;
205 0 : int result = FALSE;
206 :
207 : /* -------------------------------------------------------------------- */
208 : /* Fetch tiepoints and pixel scale. */
209 : /* -------------------------------------------------------------------- */
210 0 : double *tiepoints = NULL;
211 : int tiepoint_count;
212 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_TIEPOINTS,
213 : &tiepoint_count, &tiepoints ))
214 0 : tiepoint_count = 0;
215 :
216 : int count;
217 0 : double *pixel_scale = NULL;
218 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_PIXELSCALE, &count, &pixel_scale ))
219 0 : count = 0;
220 :
221 0 : int transform_count = 0;
222 0 : double *transform = NULL;
223 0 : if (!(gtif->gt_methods.get)(tif, GTIFF_TRANSMATRIX,
224 : &transform_count, &transform ))
225 0 : transform_count = 0;
226 :
227 : /* -------------------------------------------------------------------- */
228 : /* If the pixelscale count is zero, but we have tiepoints use */
229 : /* the tiepoint based approach. */
230 : /* -------------------------------------------------------------------- */
231 0 : if( tiepoint_count > 6 && count == 0 )
232 : {
233 0 : result = GTIFTiepointTranslate( tiepoint_count / 6,
234 : tiepoints + 3, tiepoints,
235 : *x, *y, x, y );
236 : }
237 :
238 : /* -------------------------------------------------------------------- */
239 : /* Handle matrix - convert to "geotransform" format, invert and */
240 : /* apply. */
241 : /* -------------------------------------------------------------------- */
242 0 : else if( transform_count == 16 )
243 : {
244 0 : const double x_in = *x;
245 0 : const double y_in = *y;
246 : double gt_in[6];
247 0 : gt_in[0] = transform[0];
248 0 : gt_in[1] = transform[1];
249 0 : gt_in[2] = transform[3];
250 0 : gt_in[3] = transform[4];
251 0 : gt_in[4] = transform[5];
252 0 : gt_in[5] = transform[7];
253 :
254 : double gt_out[6];
255 0 : if( !inv_geotransform( gt_in, gt_out ) )
256 0 : result = FALSE;
257 : else
258 : {
259 0 : *x = x_in * gt_out[0] + y_in * gt_out[1] + gt_out[2];
260 0 : *y = x_in * gt_out[3] + y_in * gt_out[4] + gt_out[5];
261 :
262 0 : result = TRUE;
263 : }
264 : }
265 :
266 : /* -------------------------------------------------------------------- */
267 : /* For now we require one tie point, and a valid pixel scale. */
268 : /* -------------------------------------------------------------------- */
269 0 : else if( count >= 3 && tiepoint_count >= 6 )
270 : {
271 0 : *x = (*x - tiepoints[3]) / pixel_scale[0] + tiepoints[0];
272 0 : *y = (*y - tiepoints[4]) / (-1 * pixel_scale[1]) + tiepoints[1];
273 :
274 0 : result = TRUE;
275 : }
276 :
277 : /* -------------------------------------------------------------------- */
278 : /* Cleanup. */
279 : /* -------------------------------------------------------------------- */
280 0 : if(tiepoints)
281 0 : _GTIFFree(tiepoints);
282 0 : if(pixel_scale)
283 0 : _GTIFFree(pixel_scale);
284 0 : if(transform)
285 0 : _GTIFFree(transform);
286 :
287 0 : return result;
288 : }
|
