1 : /******************************************************************************
2 : * $Id: gdalgeoloc.cpp 18948 2010-02-27 21:27:53Z rouault $
3 : *
4 : * Project: GDAL
5 : * Purpose: Implements Geolocation array based transformer.
6 : * Author: Frank Warmerdam, warmerdam@pobox.com
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 2006, Frank Warmerdam <warmerdam@pobox.com>
10 : *
11 : * Permission is hereby granted, free of charge, to any person obtaining a
12 : * copy of this software and associated documentation files (the "Software"),
13 : * to deal in the Software without restriction, including without limitation
14 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
15 : * and/or sell copies of the Software, and to permit persons to whom the
16 : * Software is furnished to do so, subject to the following conditions:
17 : *
18 : * The above copyright notice and this permission notice shall be included
19 : * in all copies or substantial portions of the Software.
20 : *
21 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
22 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
23 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
24 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
25 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
26 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
27 : * DEALINGS IN THE SOFTWARE.
28 : ****************************************************************************/
29 :
30 : #include "gdal_priv.h"
31 : #include "gdal_alg.h"
32 :
33 : #ifdef SHAPE_DEBUG
34 : #include "/u/pkg/shapelib/shapefil.h"
35 :
36 : SHPHandle hSHP = NULL;
37 : DBFHandle hDBF = NULL;
38 : #endif
39 :
40 : CPL_CVSID("$Id: gdalgeoloc.cpp 18948 2010-02-27 21:27:53Z rouault $");
41 :
42 : CPL_C_START
43 : CPLXMLNode *GDALSerializeGeoLocTransformer( void *pTransformArg );
44 : void *GDALDeserializeGeoLocTransformer( CPLXMLNode *psTree );
45 : CPL_C_END
46 :
47 : /************************************************************************/
48 : /* ==================================================================== */
49 : /* GDALGeoLocTransformer */
50 : /* ==================================================================== */
51 : /************************************************************************/
52 :
53 : typedef struct {
54 :
55 : GDALTransformerInfo sTI;
56 :
57 : int bReversed;
58 :
59 : // Map from target georef coordinates back to geolocation array
60 : // pixel line coordinates. Built only if needed.
61 :
62 : int nBackMapWidth;
63 : int nBackMapHeight;
64 : double adfBackMapGeoTransform[6]; // maps georef to pixel/line.
65 : float *pafBackMapX;
66 : float *pafBackMapY;
67 :
68 : // geolocation bands.
69 :
70 : GDALDatasetH hDS_X;
71 : GDALRasterBandH hBand_X;
72 : GDALDatasetH hDS_Y;
73 : GDALRasterBandH hBand_Y;
74 :
75 : // Located geolocation data.
76 : int nGeoLocXSize;
77 : int nGeoLocYSize;
78 : double *padfGeoLocX;
79 : double *padfGeoLocY;
80 :
81 : double dfNoDataX;
82 : double dfNoDataY;
83 :
84 : // geolocation <-> base image mapping.
85 : double dfPIXEL_OFFSET;
86 : double dfPIXEL_STEP;
87 : double dfLINE_OFFSET;
88 : double dfLINE_STEP;
89 :
90 : char ** papszGeolocationInfo;
91 :
92 : } GDALGeoLocTransformInfo;
93 :
94 : /************************************************************************/
95 : /* GeoLocLoadFullData() */
96 : /************************************************************************/
97 :
98 2 : static int GeoLocLoadFullData( GDALGeoLocTransformInfo *psTransform )
99 :
100 : {
101 : int nXSize, nYSize;
102 :
103 2 : int nXSize_XBand = GDALGetRasterXSize( psTransform->hDS_X );
104 2 : int nYSize_XBand = GDALGetRasterYSize( psTransform->hDS_X );
105 2 : int nXSize_YBand = GDALGetRasterXSize( psTransform->hDS_Y );
106 2 : int nYSize_YBand = GDALGetRasterYSize( psTransform->hDS_Y );
107 2 : if (nYSize_XBand == 1 && nYSize_YBand == 1)
108 : {
109 0 : nXSize = nXSize_XBand;
110 0 : nYSize = nXSize_YBand;
111 : }
112 : else
113 : {
114 2 : nXSize = nXSize_XBand;
115 2 : nYSize = nYSize_XBand;
116 : }
117 :
118 2 : psTransform->nGeoLocXSize = nXSize;
119 2 : psTransform->nGeoLocYSize = nYSize;
120 :
121 : psTransform->padfGeoLocY = (double *)
122 2 : VSIMalloc3(sizeof(double), nXSize, nYSize);
123 : psTransform->padfGeoLocX = (double *)
124 2 : VSIMalloc3(sizeof(double), nXSize, nYSize);
125 :
126 2 : if( psTransform->padfGeoLocX == NULL ||
127 : psTransform->padfGeoLocY == NULL )
128 : {
129 : CPLError(CE_Failure, CPLE_OutOfMemory,
130 0 : "GeoLocLoadFullData : Out of memory");
131 0 : return FALSE;
132 : }
133 :
134 2 : if (nYSize_XBand == 1 && nYSize_YBand == 1)
135 : {
136 : /* Case of regular grid */
137 : /* The XBAND contains the x coordinates for all lines */
138 : /* The YBAND contains the y coordinates for all columns */
139 :
140 0 : double* padfTempX = (double*)VSIMalloc2(nXSize, sizeof(double));
141 0 : double* padfTempY = (double*)VSIMalloc2(nYSize, sizeof(double));
142 0 : if (padfTempX == NULL || padfTempY == NULL)
143 : {
144 0 : CPLFree(padfTempX);
145 0 : CPLFree(padfTempY);
146 : CPLError(CE_Failure, CPLE_OutOfMemory,
147 0 : "GeoLocLoadFullData : Out of memory");
148 0 : return FALSE;
149 : }
150 :
151 0 : CPLErr eErr = CE_None;
152 :
153 : eErr = GDALRasterIO( psTransform->hBand_X, GF_Read,
154 : 0, 0, nXSize, 1,
155 : padfTempX, nXSize, 1,
156 0 : GDT_Float64, 0, 0 );
157 :
158 : int i,j;
159 0 : for(j=0;j<nYSize;j++)
160 : {
161 : memcpy( psTransform->padfGeoLocX + j * nXSize,
162 : padfTempX,
163 0 : nXSize * sizeof(double) );
164 : }
165 :
166 0 : if (eErr == CE_None)
167 : {
168 : eErr = GDALRasterIO( psTransform->hBand_Y, GF_Read,
169 : 0, 0, nYSize, 1,
170 : padfTempY, nYSize, 1,
171 0 : GDT_Float64, 0, 0 );
172 :
173 0 : for(j=0;j<nYSize;j++)
174 : {
175 0 : for(i=0;i<nXSize;i++)
176 : {
177 0 : psTransform->padfGeoLocY[j * nXSize + i] = padfTempY[j];
178 : }
179 : }
180 : }
181 :
182 0 : CPLFree(padfTempX);
183 0 : CPLFree(padfTempY);
184 :
185 0 : if (eErr != CE_None)
186 0 : return FALSE;
187 : }
188 : else
189 : {
190 2 : if( GDALRasterIO( psTransform->hBand_X, GF_Read,
191 : 0, 0, nXSize, nYSize,
192 : psTransform->padfGeoLocX, nXSize, nYSize,
193 : GDT_Float64, 0, 0 ) != CE_None
194 : || GDALRasterIO( psTransform->hBand_Y, GF_Read,
195 : 0, 0, nXSize, nYSize,
196 : psTransform->padfGeoLocY, nXSize, nYSize,
197 : GDT_Float64, 0, 0 ) != CE_None )
198 0 : return FALSE;
199 : }
200 :
201 : psTransform->dfNoDataX = GDALGetRasterNoDataValue( psTransform->hBand_X,
202 2 : NULL );
203 : psTransform->dfNoDataY = GDALGetRasterNoDataValue( psTransform->hBand_Y,
204 2 : NULL );
205 :
206 2 : return TRUE;
207 : }
208 :
209 : /************************************************************************/
210 : /* GeoLocGenerateBackMap() */
211 : /************************************************************************/
212 :
213 2 : static int GeoLocGenerateBackMap( GDALGeoLocTransformInfo *psTransform )
214 :
215 : {
216 2 : int nXSize = psTransform->nGeoLocXSize;
217 2 : int nYSize = psTransform->nGeoLocYSize;
218 2 : int nMaxIter = 3;
219 :
220 : /* -------------------------------------------------------------------- */
221 : /* Scan forward map for lat/long extents. */
222 : /* -------------------------------------------------------------------- */
223 2 : double dfMinX=0, dfMaxX=0, dfMinY=0, dfMaxY=0;
224 2 : int i, bInit = FALSE;
225 :
226 4682 : for( i = nXSize * nYSize - 1; i >= 0; i-- )
227 : {
228 4680 : if( psTransform->padfGeoLocX[i] != psTransform->dfNoDataX )
229 : {
230 4680 : if( bInit )
231 : {
232 4678 : dfMinX = MIN(dfMinX,psTransform->padfGeoLocX[i]);
233 4678 : dfMaxX = MAX(dfMaxX,psTransform->padfGeoLocX[i]);
234 4678 : dfMinY = MIN(dfMinY,psTransform->padfGeoLocY[i]);
235 4678 : dfMaxY = MAX(dfMaxY,psTransform->padfGeoLocY[i]);
236 : }
237 : else
238 : {
239 2 : bInit = TRUE;
240 2 : dfMinX = dfMaxX = psTransform->padfGeoLocX[i];
241 2 : dfMinY = dfMaxY = psTransform->padfGeoLocY[i];
242 : }
243 : }
244 : }
245 :
246 : /* -------------------------------------------------------------------- */
247 : /* Decide on resolution for backmap. We aim for slightly */
248 : /* higher resolution than the source but we can't easily */
249 : /* establish how much dead space there is in the backmap, so it */
250 : /* is approximate. */
251 : /* -------------------------------------------------------------------- */
252 2 : double dfTargetPixels = (nXSize * nYSize * 1.3);
253 : double dfPixelSize = sqrt((dfMaxX - dfMinX) * (dfMaxY - dfMinY)
254 2 : / dfTargetPixels);
255 : int nBMXSize, nBMYSize;
256 :
257 : nBMYSize = psTransform->nBackMapHeight =
258 2 : (int) ((dfMaxY - dfMinY) / dfPixelSize + 1);
259 : nBMXSize= psTransform->nBackMapWidth =
260 2 : (int) ((dfMaxX - dfMinX) / dfPixelSize + 1);
261 :
262 2 : if (nBMXSize > INT_MAX / nBMYSize)
263 : {
264 : CPLError(CE_Failure, CPLE_AppDefined, "Int overflow : %d x %d",
265 0 : nBMXSize, nBMYSize);
266 0 : return FALSE;
267 : }
268 :
269 2 : dfMinX -= dfPixelSize/2.0;
270 2 : dfMaxY += dfPixelSize/2.0;
271 :
272 2 : psTransform->adfBackMapGeoTransform[0] = dfMinX;
273 2 : psTransform->adfBackMapGeoTransform[1] = dfPixelSize;
274 2 : psTransform->adfBackMapGeoTransform[2] = 0.0;
275 2 : psTransform->adfBackMapGeoTransform[3] = dfMaxY;
276 2 : psTransform->adfBackMapGeoTransform[4] = 0.0;
277 2 : psTransform->adfBackMapGeoTransform[5] = -dfPixelSize;
278 :
279 : /* -------------------------------------------------------------------- */
280 : /* Allocate backmap, and initialize to nodata value (-1.0). */
281 : /* -------------------------------------------------------------------- */
282 : GByte *pabyValidFlag;
283 :
284 : pabyValidFlag = (GByte *)
285 2 : VSICalloc(nBMXSize, nBMYSize);
286 :
287 : psTransform->pafBackMapX = (float *)
288 2 : VSIMalloc3(nBMXSize, nBMYSize, sizeof(float));
289 : psTransform->pafBackMapY = (float *)
290 2 : VSIMalloc3(nBMXSize, nBMYSize, sizeof(float));
291 :
292 2 : if( pabyValidFlag == NULL ||
293 : psTransform->pafBackMapX == NULL ||
294 : psTransform->pafBackMapY == NULL )
295 : {
296 : CPLError( CE_Failure, CPLE_OutOfMemory,
297 : "Unable to allocate %dx%d back-map for geolocation array transformer.",
298 0 : nBMXSize, nBMYSize );
299 0 : CPLFree( pabyValidFlag );
300 0 : return FALSE;
301 : }
302 :
303 6258 : for( i = nBMXSize * nBMYSize - 1; i >= 0; i-- )
304 : {
305 6256 : psTransform->pafBackMapX[i] = -1.0;
306 6256 : psTransform->pafBackMapY[i] = -1.0;
307 : }
308 :
309 : /* -------------------------------------------------------------------- */
310 : /* Run through the whole geoloc array forward projecting and */
311 : /* pushing into the backmap. */
312 : /* Initialise to the nMaxIter+1 value so we can spot genuinely */
313 : /* valid pixels in the hole-filling loop. */
314 : /* -------------------------------------------------------------------- */
315 : int iBMX, iBMY;
316 : int iX, iY;
317 :
318 80 : for( iY = 0; iY < nYSize; iY++ )
319 : {
320 4758 : for( iX = 0; iX < nXSize; iX++ )
321 : {
322 4680 : if( psTransform->padfGeoLocX[iX + iY * nXSize]
323 : == psTransform->dfNoDataX )
324 0 : continue;
325 :
326 4680 : i = iX + iY * nXSize;
327 :
328 4680 : iBMX = (int) ((psTransform->padfGeoLocX[i] - dfMinX) / dfPixelSize);
329 4680 : iBMY = (int) ((dfMaxY - psTransform->padfGeoLocY[i]) / dfPixelSize);
330 :
331 4680 : if( iBMX < 0 || iBMY < 0 || iBMX >= nBMXSize || iBMY >= nBMYSize )
332 2 : continue;
333 :
334 : psTransform->pafBackMapX[iBMX + iBMY * nBMXSize] =
335 4678 : (float)(iX * psTransform->dfPIXEL_STEP + psTransform->dfPIXEL_OFFSET);
336 : psTransform->pafBackMapY[iBMX + iBMY * nBMXSize] =
337 4678 : (float)(iY * psTransform->dfLINE_STEP + psTransform->dfLINE_OFFSET);
338 :
339 4678 : pabyValidFlag[iBMX + iBMY * nBMXSize] = nMaxIter+1;
340 :
341 : }
342 : }
343 :
344 : /* -------------------------------------------------------------------- */
345 : /* Now, loop over the backmap trying to fill in holes with */
346 : /* nearby values. */
347 : /* -------------------------------------------------------------------- */
348 : int iIter;
349 : int nNumValid, nExpectedValid;
350 :
351 8 : for( iIter = 0; iIter < nMaxIter; iIter++ )
352 : {
353 6 : nNumValid = 0;
354 6 : nExpectedValid = (nBMYSize - 2) * (nBMXSize - 2);
355 270 : for( iBMY = 1; iBMY < nBMYSize-1; iBMY++ )
356 : {
357 17688 : for( iBMX = 1; iBMX < nBMXSize-1; iBMX++ )
358 : {
359 : // if this point is already set, ignore it.
360 17424 : if( pabyValidFlag[iBMX + iBMY*nBMXSize] )
361 : {
362 12860 : nNumValid++;
363 12860 : continue;
364 : }
365 :
366 4564 : int nCount = 0;
367 4564 : double dfXSum = 0.0, dfYSum = 0.0;
368 4564 : int nMarkedAsGood = nMaxIter - iIter;
369 :
370 : // left?
371 4564 : if( pabyValidFlag[iBMX-1+iBMY*nBMXSize] > nMarkedAsGood )
372 : {
373 788 : dfXSum += psTransform->pafBackMapX[iBMX-1+iBMY*nBMXSize];
374 788 : dfYSum += psTransform->pafBackMapY[iBMX-1+iBMY*nBMXSize];
375 788 : nCount++;
376 : }
377 : // right?
378 4564 : if( pabyValidFlag[iBMX+1+iBMY*nBMXSize] > nMarkedAsGood )
379 : {
380 812 : dfXSum += psTransform->pafBackMapX[iBMX+1+iBMY*nBMXSize];
381 812 : dfYSum += psTransform->pafBackMapY[iBMX+1+iBMY*nBMXSize];
382 812 : nCount++;
383 : }
384 : // top?
385 4564 : if( pabyValidFlag[iBMX+(iBMY-1)*nBMXSize] > nMarkedAsGood )
386 : {
387 756 : dfXSum += psTransform->pafBackMapX[iBMX+(iBMY-1)*nBMXSize];
388 756 : dfYSum += psTransform->pafBackMapY[iBMX+(iBMY-1)*nBMXSize];
389 756 : nCount++;
390 : }
391 : // bottom?
392 4564 : if( pabyValidFlag[iBMX+(iBMY+1)*nBMXSize] > nMarkedAsGood )
393 : {
394 736 : dfXSum += psTransform->pafBackMapX[iBMX+(iBMY+1)*nBMXSize];
395 736 : dfYSum += psTransform->pafBackMapY[iBMX+(iBMY+1)*nBMXSize];
396 736 : nCount++;
397 : }
398 : // top-left?
399 4564 : if( pabyValidFlag[iBMX-1+(iBMY-1)*nBMXSize] > nMarkedAsGood )
400 : {
401 908 : dfXSum += psTransform->pafBackMapX[iBMX-1+(iBMY-1)*nBMXSize];
402 908 : dfYSum += psTransform->pafBackMapY[iBMX-1+(iBMY-1)*nBMXSize];
403 908 : nCount++;
404 : }
405 : // top-right?
406 4564 : if( pabyValidFlag[iBMX+1+(iBMY-1)*nBMXSize] > nMarkedAsGood )
407 : {
408 778 : dfXSum += psTransform->pafBackMapX[iBMX+1+(iBMY-1)*nBMXSize];
409 778 : dfYSum += psTransform->pafBackMapY[iBMX+1+(iBMY-1)*nBMXSize];
410 778 : nCount++;
411 : }
412 : // bottom-left?
413 4564 : if( pabyValidFlag[iBMX-1+(iBMY+1)*nBMXSize] > nMarkedAsGood )
414 : {
415 734 : dfXSum += psTransform->pafBackMapX[iBMX-1+(iBMY+1)*nBMXSize];
416 734 : dfYSum += psTransform->pafBackMapY[iBMX-1+(iBMY+1)*nBMXSize];
417 734 : nCount++;
418 : }
419 : // bottom-right?
420 4564 : if( pabyValidFlag[iBMX+1+(iBMY+1)*nBMXSize] > nMarkedAsGood )
421 : {
422 916 : dfXSum += psTransform->pafBackMapX[iBMX+1+(iBMY+1)*nBMXSize];
423 916 : dfYSum += psTransform->pafBackMapY[iBMX+1+(iBMY+1)*nBMXSize];
424 916 : nCount++;
425 : }
426 :
427 4564 : if( nCount > 0 )
428 : {
429 1702 : psTransform->pafBackMapX[iBMX + iBMY * nBMXSize] = (float)(dfXSum/nCount);
430 1702 : psTransform->pafBackMapY[iBMX + iBMY * nBMXSize] = (float)(dfYSum/nCount);
431 : // genuinely valid points will have value iMaxIter+1
432 : // On each iteration mark newly valid points with a
433 : // descending value so that it will not be used on the
434 : // current iteration only on subsequent ones.
435 1702 : pabyValidFlag[iBMX+iBMY*nBMXSize] = nMaxIter - iIter;
436 : }
437 : }
438 : }
439 6 : if (nNumValid == nExpectedValid)
440 0 : break;
441 : }
442 :
443 : #ifdef notdef
444 : GDALDatasetH hBMDS = GDALCreate( GDALGetDriverByName( "GTiff" ),
445 : "backmap.tif", nBMXSize, nBMYSize, 2,
446 : GDT_Float32, NULL );
447 : GDALSetGeoTransform( hBMDS, psTransform->adfBackMapGeoTransform );
448 : GDALRasterIO( GDALGetRasterBand(hBMDS,1), GF_Write,
449 : 0, 0, nBMXSize, nBMYSize,
450 : psTransform->pafBackMapX, nBMXSize, nBMYSize,
451 : GDT_Float32, 0, 0 );
452 : GDALRasterIO( GDALGetRasterBand(hBMDS,2), GF_Write,
453 : 0, 0, nBMXSize, nBMYSize,
454 : psTransform->pafBackMapY, nBMXSize, nBMYSize,
455 : GDT_Float32, 0, 0 );
456 : GDALClose( hBMDS );
457 : #endif
458 :
459 2 : CPLFree( pabyValidFlag );
460 :
461 2 : return TRUE;
462 : }
463 :
464 : /************************************************************************/
465 : /* FindGeoLocPosition() */
466 : /************************************************************************/
467 :
468 : #ifdef notdef
469 :
470 : /*
471 : This searching approach has been abandoned because it is too sensitive
472 : to discontinuities in the data. Left in case it might be revived in
473 : the future.
474 : */
475 :
476 : static int FindGeoLocPosition( GDALGeoLocTransformInfo *psTransform,
477 : double dfGeoX, double dfGeoY,
478 : int nStartX, int nStartY,
479 : double *pdfFoundX, double *pdfFoundY )
480 :
481 : {
482 : double adfPathX[5000], adfPathY[5000];
483 :
484 : if( psTransform->padfGeoLocX == NULL )
485 : return FALSE;
486 :
487 : int nXSize = psTransform->nGeoLocXSize;
488 : int nYSize = psTransform->nGeoLocYSize;
489 : int nStepCount = 0;
490 :
491 : // Start in center if we don't have any provided info.
492 : if( nStartX < 0 || nStartY < 0
493 : || nStartX >= nXSize || nStartY >= nYSize )
494 : {
495 : nStartX = nXSize / 2;
496 : nStartY = nYSize / 2;
497 : }
498 :
499 : nStartX = MIN(nStartX,nXSize-2);
500 : nStartY = MIN(nStartY,nYSize-2);
501 :
502 : int iX = nStartX, iY = nStartY;
503 : int iLastX = -1, iLastY = -1;
504 : int iSecondLastX = -1, iSecondLastY = -1;
505 :
506 : while( nStepCount < MAX(nXSize,nYSize) )
507 : {
508 : int iXNext = -1, iYNext = -1;
509 : double dfDeltaXRight, dfDeltaYRight, dfDeltaXDown, dfDeltaYDown;
510 :
511 : double *padfThisX = psTransform->padfGeoLocX + iX + iY * nXSize;
512 : double *padfThisY = psTransform->padfGeoLocY + iX + iY * nXSize;
513 :
514 : double dfDeltaX = dfGeoX - *padfThisX;
515 : double dfDeltaY = dfGeoY - *padfThisY;
516 :
517 : if( iX == nXSize-1 )
518 : {
519 : dfDeltaXRight = *(padfThisX) - *(padfThisX-1);
520 : dfDeltaYRight = *(padfThisY) - *(padfThisY-1);
521 : }
522 : else
523 : {
524 : dfDeltaXRight = *(padfThisX+1) - *padfThisX;
525 : dfDeltaYRight = *(padfThisY+1) - *padfThisY;
526 : }
527 :
528 : if( iY == nYSize - 1 )
529 : {
530 : dfDeltaXDown = *(padfThisX) - *(padfThisX-nXSize);
531 : dfDeltaYDown = *(padfThisY) - *(padfThisY-nXSize);
532 : }
533 : else
534 : {
535 : dfDeltaXDown = *(padfThisX+nXSize) - *padfThisX;
536 : dfDeltaYDown = *(padfThisY+nXSize) - *padfThisY;
537 : }
538 :
539 : double dfRightProjection =
540 : (dfDeltaXRight * dfDeltaX + dfDeltaYRight * dfDeltaY)
541 : / (dfDeltaXRight*dfDeltaXRight + dfDeltaYRight*dfDeltaYRight);
542 :
543 : double dfDownProjection =
544 : (dfDeltaXDown * dfDeltaX + dfDeltaYDown * dfDeltaY)
545 : / (dfDeltaXDown*dfDeltaXDown + dfDeltaYDown*dfDeltaYDown);
546 :
547 : // Are we in our target cell?
548 : if( dfRightProjection >= 0.0 && dfRightProjection < 1.0
549 : && dfDownProjection >= 0.0 && dfDownProjection < 1.0 )
550 : {
551 : *pdfFoundX = iX + dfRightProjection;
552 : *pdfFoundY = iY + dfDownProjection;
553 :
554 : return TRUE;
555 : }
556 :
557 : if( ABS(dfRightProjection) > ABS(dfDownProjection) )
558 : {
559 : // Do we want to move right?
560 : if( dfRightProjection > 1.0 && iX < nXSize-1 )
561 : {
562 : iXNext = iX + MAX(1,(int)(dfRightProjection - nStepCount)/2);
563 : iYNext = iY;
564 : }
565 :
566 : // Do we want to move left?
567 : else if( dfRightProjection < 0.0 && iX > 0 )
568 : {
569 : iXNext = iX - MAX(1,(int)(ABS(dfRightProjection) - nStepCount)/2);
570 : iYNext = iY;
571 : }
572 :
573 : // Do we want to move down.
574 : else if( dfDownProjection > 1.0 && iY < nYSize-1 )
575 : {
576 : iXNext = iX;
577 : iYNext = iY + MAX(1,(int)(dfDownProjection - nStepCount)/2);
578 : }
579 :
580 : // Do we want to move up?
581 : else if( dfDownProjection < 0.0 && iY > 0 )
582 : {
583 : iXNext = iX;
584 : iYNext = iY - MAX(1,(int)(ABS(dfDownProjection) - nStepCount)/2);
585 : }
586 :
587 : // We aren't there, and we have no where to go
588 : else
589 : {
590 : return FALSE;
591 : }
592 : }
593 : else
594 : {
595 : // Do we want to move down.
596 : if( dfDownProjection > 1.0 && iY < nYSize-1 )
597 : {
598 : iXNext = iX;
599 : iYNext = iY + MAX(1,(int)(dfDownProjection - nStepCount)/2);
600 : }
601 :
602 : // Do we want to move up?
603 : else if( dfDownProjection < 0.0 && iY > 0 )
604 : {
605 : iXNext = iX;
606 : iYNext = iY - MAX(1,(int)(ABS(dfDownProjection) - nStepCount)/2);
607 : }
608 :
609 : // Do we want to move right?
610 : else if( dfRightProjection > 1.0 && iX < nXSize-1 )
611 : {
612 : iXNext = iX + MAX(1,(int)(dfRightProjection - nStepCount)/2);
613 : iYNext = iY;
614 : }
615 :
616 : // Do we want to move left?
617 : else if( dfRightProjection < 0.0 && iX > 0 )
618 : {
619 : iXNext = iX - MAX(1,(int)(ABS(dfRightProjection) - nStepCount)/2);
620 : iYNext = iY;
621 : }
622 :
623 : // We aren't there, and we have no where to go
624 : else
625 : {
626 : return FALSE;
627 : }
628 : }
629 : adfPathX[nStepCount] = iX;
630 : adfPathY[nStepCount] = iY;
631 :
632 : nStepCount++;
633 : iX = MAX(0,MIN(iXNext,nXSize-1));
634 : iY = MAX(0,MIN(iYNext,nYSize-1));
635 :
636 : if( iX == iSecondLastX && iY == iSecondLastY )
637 : {
638 : // Are we *near* our target cell?
639 : if( dfRightProjection >= -1.0 && dfRightProjection < 2.0
640 : && dfDownProjection >= -1.0 && dfDownProjection < 2.0 )
641 : {
642 : *pdfFoundX = iX + dfRightProjection;
643 : *pdfFoundY = iY + dfDownProjection;
644 :
645 : return TRUE;
646 : }
647 :
648 : #ifdef SHAPE_DEBUG
649 : if( hSHP != NULL )
650 : {
651 : SHPObject *hObj;
652 :
653 : hObj = SHPCreateSimpleObject( SHPT_ARC, nStepCount,
654 : adfPathX, adfPathY, NULL );
655 : SHPWriteObject( hSHP, -1, hObj );
656 : SHPDestroyObject( hObj );
657 :
658 : int iShape = DBFGetRecordCount( hDBF );
659 : DBFWriteDoubleAttribute( hDBF, iShape, 0, dfGeoX );
660 : DBFWriteDoubleAttribute( hDBF, iShape, 1, dfGeoY );
661 : }
662 : #endif
663 : //CPLDebug( "GeoL", "Looping at step (%d) on search for %g,%g.",
664 : // nStepCount, dfGeoX, dfGeoY );
665 : return FALSE;
666 : }
667 :
668 : iSecondLastX = iLastX;
669 : iSecondLastY = iLastY;
670 :
671 : iLastX = iX;
672 : iLastY = iY;
673 :
674 : }
675 :
676 : //CPLDebug( "GeoL", "Exceeded step count max (%d) on search for %g,%g.",
677 : // MAX(nXSize,nYSize),
678 : // dfGeoX, dfGeoY );
679 :
680 : #ifdef SHAPE_DEBUG
681 : if( hSHP != NULL )
682 : {
683 : SHPObject *hObj;
684 :
685 : hObj = SHPCreateSimpleObject( SHPT_ARC, nStepCount,
686 : adfPathX, adfPathY, NULL );
687 : SHPWriteObject( hSHP, -1, hObj );
688 : SHPDestroyObject( hObj );
689 :
690 : int iShape = DBFGetRecordCount( hDBF );
691 : DBFWriteDoubleAttribute( hDBF, iShape, 0, dfGeoX );
692 : DBFWriteDoubleAttribute( hDBF, iShape, 1, dfGeoY );
693 : }
694 : #endif
695 :
696 : return FALSE;
697 : }
698 : #endif /* def notdef */
699 :
700 :
701 :
702 : /************************************************************************/
703 : /* GDALCreateGeoLocTransformer() */
704 : /************************************************************************/
705 :
706 : void *GDALCreateGeoLocTransformer( GDALDatasetH hBaseDS,
707 : char **papszGeolocationInfo,
708 2 : int bReversed )
709 :
710 : {
711 : GDALGeoLocTransformInfo *psTransform;
712 :
713 2 : if( CSLFetchNameValue(papszGeolocationInfo,"PIXEL_OFFSET") == NULL
714 : || CSLFetchNameValue(papszGeolocationInfo,"LINE_OFFSET") == NULL
715 : || CSLFetchNameValue(papszGeolocationInfo,"PIXEL_STEP") == NULL
716 : || CSLFetchNameValue(papszGeolocationInfo,"LINE_STEP") == NULL
717 : || CSLFetchNameValue(papszGeolocationInfo,"X_BAND") == NULL
718 : || CSLFetchNameValue(papszGeolocationInfo,"Y_BAND") == NULL )
719 : {
720 : CPLError( CE_Failure, CPLE_AppDefined,
721 0 : "Missing some geolocation fields in GDALCreateGeoLocTransformer()" );
722 0 : return NULL;
723 : }
724 :
725 : /* -------------------------------------------------------------------- */
726 : /* Initialize core info. */
727 : /* -------------------------------------------------------------------- */
728 : psTransform = (GDALGeoLocTransformInfo *)
729 2 : CPLCalloc(sizeof(GDALGeoLocTransformInfo),1);
730 :
731 2 : psTransform->bReversed = bReversed;
732 :
733 2 : strcpy( psTransform->sTI.szSignature, "GTI" );
734 2 : psTransform->sTI.pszClassName = "GDALGeoLocTransformer";
735 2 : psTransform->sTI.pfnTransform = GDALGeoLocTransform;
736 2 : psTransform->sTI.pfnCleanup = GDALDestroyGeoLocTransformer;
737 2 : psTransform->sTI.pfnSerialize = GDALSerializeGeoLocTransformer;
738 2 : psTransform->papszGeolocationInfo = CSLDuplicate( papszGeolocationInfo );
739 :
740 : /* -------------------------------------------------------------------- */
741 : /* Pull geolocation info from the options/metadata. */
742 : /* -------------------------------------------------------------------- */
743 : psTransform->dfPIXEL_OFFSET = atof(CSLFetchNameValue( papszGeolocationInfo,
744 2 : "PIXEL_OFFSET" ));
745 : psTransform->dfLINE_OFFSET = atof(CSLFetchNameValue( papszGeolocationInfo,
746 2 : "LINE_OFFSET" ));
747 : psTransform->dfPIXEL_STEP = atof(CSLFetchNameValue( papszGeolocationInfo,
748 2 : "PIXEL_STEP" ));
749 : psTransform->dfLINE_STEP = atof(CSLFetchNameValue( papszGeolocationInfo,
750 2 : "LINE_STEP" ));
751 :
752 : /* -------------------------------------------------------------------- */
753 : /* Establish access to geolocation dataset(s). */
754 : /* -------------------------------------------------------------------- */
755 : const char *pszDSName = CSLFetchNameValue( papszGeolocationInfo,
756 2 : "X_DATASET" );
757 2 : if( pszDSName != NULL )
758 : {
759 2 : psTransform->hDS_X = GDALOpenShared( pszDSName, GA_ReadOnly );
760 : }
761 : else
762 : {
763 0 : psTransform->hDS_X = hBaseDS;
764 0 : GDALReferenceDataset( psTransform->hDS_X );
765 : psTransform->papszGeolocationInfo =
766 : CSLSetNameValue( psTransform->papszGeolocationInfo,
767 : "X_DATASET",
768 0 : GDALGetDescription( hBaseDS ) );
769 : }
770 :
771 2 : pszDSName = CSLFetchNameValue( papszGeolocationInfo, "Y_DATASET" );
772 2 : if( pszDSName != NULL )
773 : {
774 2 : psTransform->hDS_Y = GDALOpenShared( pszDSName, GA_ReadOnly );
775 : }
776 : else
777 : {
778 0 : psTransform->hDS_Y = hBaseDS;
779 0 : GDALReferenceDataset( psTransform->hDS_Y );
780 : psTransform->papszGeolocationInfo =
781 : CSLSetNameValue( psTransform->papszGeolocationInfo,
782 : "Y_DATASET",
783 0 : GDALGetDescription( hBaseDS ) );
784 : }
785 :
786 2 : if (psTransform->hDS_X == NULL ||
787 : psTransform->hDS_Y == NULL)
788 : {
789 0 : GDALDestroyGeoLocTransformer( psTransform );
790 0 : return NULL;
791 : }
792 :
793 : /* -------------------------------------------------------------------- */
794 : /* Get the band handles. */
795 : /* -------------------------------------------------------------------- */
796 : int nBand;
797 :
798 2 : nBand = MAX(1,atoi(CSLFetchNameValue( papszGeolocationInfo, "X_BAND" )));
799 2 : psTransform->hBand_X = GDALGetRasterBand( psTransform->hDS_X, nBand );
800 :
801 2 : nBand = MAX(1,atoi(CSLFetchNameValue( papszGeolocationInfo, "Y_BAND" )));
802 2 : psTransform->hBand_Y = GDALGetRasterBand( psTransform->hDS_Y, nBand );
803 :
804 2 : if (psTransform->hBand_X == NULL ||
805 : psTransform->hBand_Y == NULL)
806 : {
807 0 : GDALDestroyGeoLocTransformer( psTransform );
808 0 : return NULL;
809 : }
810 :
811 : /* -------------------------------------------------------------------- */
812 : /* Check that X and Y bands have the same dimensions */
813 : /* -------------------------------------------------------------------- */
814 2 : int nXSize_XBand = GDALGetRasterXSize( psTransform->hDS_X );
815 2 : int nYSize_XBand = GDALGetRasterYSize( psTransform->hDS_X );
816 2 : int nXSize_YBand = GDALGetRasterXSize( psTransform->hDS_Y );
817 2 : int nYSize_YBand = GDALGetRasterYSize( psTransform->hDS_Y );
818 2 : if (nYSize_XBand == 1 || nYSize_YBand == 1)
819 : {
820 0 : if (nYSize_XBand != 1 || nYSize_YBand != 1)
821 : {
822 : CPLError(CE_Failure, CPLE_AppDefined,
823 0 : "X_BAND and Y_BAND should have both nYSize == 1");
824 0 : GDALDestroyGeoLocTransformer( psTransform );
825 0 : return NULL;
826 : }
827 : }
828 2 : else if (nXSize_XBand != nXSize_YBand ||
829 : nYSize_XBand != nYSize_YBand )
830 : {
831 : CPLError(CE_Failure, CPLE_AppDefined,
832 0 : "X_BAND and Y_BAND do not have the same dimensions");
833 0 : GDALDestroyGeoLocTransformer( psTransform );
834 0 : return NULL;
835 : }
836 :
837 2 : if (nXSize_XBand > INT_MAX / nYSize_XBand)
838 : {
839 : CPLError(CE_Failure, CPLE_AppDefined, "Int overflow : %d x %d",
840 0 : nXSize_XBand, nYSize_XBand);
841 0 : GDALDestroyGeoLocTransformer( psTransform );
842 0 : return NULL;
843 : }
844 :
845 : /* -------------------------------------------------------------------- */
846 : /* Load the geolocation array. */
847 : /* -------------------------------------------------------------------- */
848 2 : if( !GeoLocLoadFullData( psTransform )
849 : || !GeoLocGenerateBackMap( psTransform ) )
850 : {
851 0 : GDALDestroyGeoLocTransformer( psTransform );
852 0 : return NULL;
853 : }
854 :
855 2 : return psTransform;
856 : }
857 :
858 : /************************************************************************/
859 : /* GDALDestroyGeoLocTransformer() */
860 : /************************************************************************/
861 :
862 2 : void GDALDestroyGeoLocTransformer( void *pTransformAlg )
863 :
864 : {
865 : GDALGeoLocTransformInfo *psTransform =
866 2 : (GDALGeoLocTransformInfo *) pTransformAlg;
867 :
868 2 : CPLFree( psTransform->pafBackMapX );
869 2 : CPLFree( psTransform->pafBackMapY );
870 2 : CSLDestroy( psTransform->papszGeolocationInfo );
871 2 : CPLFree( psTransform->padfGeoLocX );
872 2 : CPLFree( psTransform->padfGeoLocY );
873 :
874 2 : if( psTransform->hDS_X != NULL
875 : && GDALDereferenceDataset( psTransform->hDS_X ) == 0 )
876 0 : GDALClose( psTransform->hDS_X );
877 :
878 2 : if( psTransform->hDS_Y != NULL
879 : && GDALDereferenceDataset( psTransform->hDS_Y ) == 0 )
880 2 : GDALClose( psTransform->hDS_Y );
881 :
882 2 : CPLFree( pTransformAlg );
883 2 : }
884 :
885 : /************************************************************************/
886 : /* GDALGeoLocTransform() */
887 : /************************************************************************/
888 :
889 : int GDALGeoLocTransform( void *pTransformArg, int bDstToSrc,
890 : int nPointCount,
891 : double *padfX, double *padfY, double *padfZ,
892 260 : int *panSuccess )
893 :
894 : {
895 : GDALGeoLocTransformInfo *psTransform =
896 260 : (GDALGeoLocTransformInfo *) pTransformArg;
897 :
898 260 : if( psTransform->bReversed )
899 0 : bDstToSrc = !bDstToSrc;
900 :
901 : /* -------------------------------------------------------------------- */
902 : /* Do original pixel line to target geox/geoy. */
903 : /* -------------------------------------------------------------------- */
904 260 : if( !bDstToSrc )
905 : {
906 1 : int i, nXSize = psTransform->nGeoLocXSize;
907 :
908 2 : for( i = 0; i < nPointCount; i++ )
909 : {
910 1 : if( !panSuccess[i] )
911 0 : continue;
912 :
913 1 : if( padfX[i] == HUGE_VAL || padfY[i] == HUGE_VAL )
914 : {
915 0 : panSuccess[i] = FALSE;
916 0 : continue;
917 : }
918 :
919 : double dfGeoLocPixel = (padfX[i] - psTransform->dfPIXEL_OFFSET)
920 1 : / psTransform->dfPIXEL_STEP;
921 : double dfGeoLocLine = (padfY[i] - psTransform->dfLINE_OFFSET)
922 1 : / psTransform->dfLINE_STEP;
923 :
924 : int iX, iY;
925 :
926 1 : iX = MAX(0,(int) dfGeoLocPixel);
927 1 : iX = MIN(iX,psTransform->nGeoLocXSize-2);
928 1 : iY = MAX(0,(int) dfGeoLocLine);
929 1 : iY = MIN(iY,psTransform->nGeoLocYSize-2);
930 :
931 1 : double *padfGLX = psTransform->padfGeoLocX + iX + iY * nXSize;
932 1 : double *padfGLY = psTransform->padfGeoLocY + iX + iY * nXSize;
933 :
934 : // This assumes infinite extension beyond borders of available
935 : // data based on closest grid square.
936 :
937 : padfX[i] = padfGLX[0]
938 : + (dfGeoLocPixel-iX) * (padfGLX[1] - padfGLX[0])
939 1 : + (dfGeoLocLine -iY) * (padfGLX[nXSize] - padfGLX[0]);
940 : padfY[i] = padfGLY[0]
941 : + (dfGeoLocPixel-iX) * (padfGLY[1] - padfGLY[0])
942 1 : + (dfGeoLocLine -iY) * (padfGLY[nXSize] - padfGLY[0]);
943 :
944 1 : panSuccess[i] = TRUE;
945 : }
946 : }
947 :
948 : /* -------------------------------------------------------------------- */
949 : /* geox/geoy to pixel/line using backmap. */
950 : /* -------------------------------------------------------------------- */
951 : else
952 : {
953 : int i;
954 :
955 66705 : for( i = 0; i < nPointCount; i++ )
956 : {
957 66446 : if( !panSuccess[i] )
958 0 : continue;
959 :
960 66446 : if( padfX[i] == HUGE_VAL || padfY[i] == HUGE_VAL )
961 : {
962 0 : panSuccess[i] = FALSE;
963 0 : continue;
964 : }
965 :
966 : int iBMX, iBMY;
967 :
968 : iBMX = (int) ((padfX[i] - psTransform->adfBackMapGeoTransform[0])
969 66446 : / psTransform->adfBackMapGeoTransform[1]);
970 : iBMY = (int) ((padfY[i] - psTransform->adfBackMapGeoTransform[3])
971 66446 : / psTransform->adfBackMapGeoTransform[5]);
972 :
973 66446 : int iBM = iBMX + iBMY * psTransform->nBackMapWidth;
974 :
975 66446 : if( iBMX < 0 || iBMY < 0
976 : || iBMX >= psTransform->nBackMapWidth
977 : || iBMY >= psTransform->nBackMapHeight
978 : || psTransform->pafBackMapX[iBM] < 0 )
979 : {
980 63969 : panSuccess[i] = FALSE;
981 63969 : padfX[i] = HUGE_VAL;
982 63969 : padfY[i] = HUGE_VAL;
983 63969 : continue;
984 : }
985 :
986 2477 : padfX[i] = psTransform->pafBackMapX[iBM];
987 2477 : padfY[i] = psTransform->pafBackMapY[iBM];
988 2477 : panSuccess[i] = TRUE;
989 : }
990 : }
991 :
992 : /* -------------------------------------------------------------------- */
993 : /* geox/geoy to pixel/line using search algorithm. */
994 : /* -------------------------------------------------------------------- */
995 : #ifdef notdef
996 : else
997 : {
998 : int i;
999 : int nStartX = -1, nStartY = -1;
1000 :
1001 : #ifdef SHAPE_DEBUG
1002 : hSHP = SHPCreate( "tracks.shp", SHPT_ARC );
1003 : hDBF = DBFCreate( "tracks.dbf" );
1004 : DBFAddField( hDBF, "GEOX", FTDouble, 10, 4 );
1005 : DBFAddField( hDBF, "GEOY", FTDouble, 10, 4 );
1006 : #endif
1007 : for( i = 0; i < nPointCount; i++ )
1008 : {
1009 : double dfGeoLocX, dfGeoLocY;
1010 :
1011 : if( !panSuccess[i] )
1012 : continue;
1013 :
1014 : if( padfX[i] == HUGE_VAL || padfY[i] == HUGE_VAL )
1015 : {
1016 : panSuccess[i] = FALSE;
1017 : continue;
1018 : }
1019 :
1020 : if( !FindGeoLocPosition( psTransform, padfX[i], padfY[i],
1021 : -1, -1, &dfGeoLocX, &dfGeoLocY ) )
1022 : {
1023 : padfX[i] = HUGE_VAL;
1024 : padfY[i] = HUGE_VAL;
1025 :
1026 : panSuccess[i] = FALSE;
1027 : continue;
1028 : }
1029 : nStartX = (int) dfGeoLocX;
1030 : nStartY = (int) dfGeoLocY;
1031 :
1032 : padfX[i] = dfGeoLocX * psTransform->dfPIXEL_STEP
1033 : + psTransform->dfPIXEL_OFFSET;
1034 : padfY[i] = dfGeoLocY * psTransform->dfLINE_STEP
1035 : + psTransform->dfLINE_OFFSET;
1036 :
1037 : panSuccess[i] = TRUE;
1038 : }
1039 :
1040 : #ifdef SHAPE_DEBUG
1041 : if( hSHP != NULL )
1042 : {
1043 : DBFClose( hDBF );
1044 : hDBF = NULL;
1045 :
1046 : SHPClose( hSHP );
1047 : hSHP = NULL;
1048 : }
1049 : #endif
1050 : }
1051 : #endif
1052 :
1053 260 : return TRUE;
1054 : }
1055 :
1056 : /************************************************************************/
1057 : /* GDALSerializeGeoLocTransformer() */
1058 : /************************************************************************/
1059 :
1060 0 : CPLXMLNode *GDALSerializeGeoLocTransformer( void *pTransformArg )
1061 :
1062 : {
1063 0 : VALIDATE_POINTER1( pTransformArg, "GDALSerializeGeoLocTransformer", NULL );
1064 :
1065 : CPLXMLNode *psTree;
1066 : GDALGeoLocTransformInfo *psInfo =
1067 0 : (GDALGeoLocTransformInfo *)(pTransformArg);
1068 :
1069 0 : psTree = CPLCreateXMLNode( NULL, CXT_Element, "GeoLocTransformer" );
1070 :
1071 : /* -------------------------------------------------------------------- */
1072 : /* Serialize bReversed. */
1073 : /* -------------------------------------------------------------------- */
1074 : CPLCreateXMLElementAndValue(
1075 : psTree, "Reversed",
1076 0 : CPLString().Printf( "%d", psInfo->bReversed ) );
1077 :
1078 : /* -------------------------------------------------------------------- */
1079 : /* geoloc metadata. */
1080 : /* -------------------------------------------------------------------- */
1081 0 : char **papszMD = psInfo->papszGeolocationInfo;
1082 : CPLXMLNode *psMD= CPLCreateXMLNode( psTree, CXT_Element,
1083 0 : "Metadata" );
1084 :
1085 0 : for( int i = 0; papszMD != NULL && papszMD[i] != NULL; i++ )
1086 : {
1087 : const char *pszRawValue;
1088 : char *pszKey;
1089 : CPLXMLNode *psMDI;
1090 :
1091 0 : pszRawValue = CPLParseNameValue( papszMD[i], &pszKey );
1092 :
1093 0 : psMDI = CPLCreateXMLNode( psMD, CXT_Element, "MDI" );
1094 0 : CPLSetXMLValue( psMDI, "#key", pszKey );
1095 0 : CPLCreateXMLNode( psMDI, CXT_Text, pszRawValue );
1096 :
1097 0 : CPLFree( pszKey );
1098 : }
1099 :
1100 0 : return psTree;
1101 : }
1102 :
1103 : /************************************************************************/
1104 : /* GDALDeserializeGeoLocTransformer() */
1105 : /************************************************************************/
1106 :
1107 1 : void *GDALDeserializeGeoLocTransformer( CPLXMLNode *psTree )
1108 :
1109 : {
1110 : void *pResult;
1111 : int bReversed;
1112 1 : char **papszMD = NULL;
1113 : CPLXMLNode *psMDI, *psMetadata;
1114 :
1115 : /* -------------------------------------------------------------------- */
1116 : /* Collect metadata. */
1117 : /* -------------------------------------------------------------------- */
1118 1 : psMetadata = CPLGetXMLNode( psTree, "Metadata" );
1119 :
1120 1 : if( psMetadata == NULL ||
1121 : psMetadata->eType != CXT_Element
1122 : || !EQUAL(psMetadata->pszValue,"Metadata") )
1123 0 : return NULL;
1124 :
1125 10 : for( psMDI = psMetadata->psChild; psMDI != NULL;
1126 : psMDI = psMDI->psNext )
1127 : {
1128 9 : if( !EQUAL(psMDI->pszValue,"MDI")
1129 : || psMDI->eType != CXT_Element
1130 : || psMDI->psChild == NULL
1131 : || psMDI->psChild->psNext == NULL
1132 : || psMDI->psChild->eType != CXT_Attribute
1133 : || psMDI->psChild->psChild == NULL )
1134 0 : continue;
1135 :
1136 : papszMD =
1137 : CSLSetNameValue( papszMD,
1138 : psMDI->psChild->psChild->pszValue,
1139 9 : psMDI->psChild->psNext->pszValue );
1140 : }
1141 :
1142 : /* -------------------------------------------------------------------- */
1143 : /* Get other flags. */
1144 : /* -------------------------------------------------------------------- */
1145 1 : bReversed = atoi(CPLGetXMLValue(psTree,"Reversed","0"));
1146 :
1147 : /* -------------------------------------------------------------------- */
1148 : /* Generate transformation. */
1149 : /* -------------------------------------------------------------------- */
1150 1 : pResult = GDALCreateGeoLocTransformer( NULL, papszMD, bReversed );
1151 :
1152 : /* -------------------------------------------------------------------- */
1153 : /* Cleanup GCP copy. */
1154 : /* -------------------------------------------------------------------- */
1155 1 : CSLDestroy( papszMD );
1156 :
1157 1 : return pResult;
1158 : }
|
