1 /******************************************************************************
2 * $Id: gdal_misc.cpp 15744 2008-11-16 18:19:00Z rouault $
3 *
4 * Project: GDAL Core
5 * Purpose: Free standing functions for GDAL.
6 * Author: Frank Warmerdam, warmerdam@pobox.com
7 *
8 ******************************************************************************
9 * Copyright (c) 1999, Frank Warmerdam
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 "../port/cpl_string.h"
32 #include "../port/cpl_minixml.h"
33 #include <ctype.h>
34 #include <string>
35
36 CPL_CVSID("$Id: gdal_misc.cpp 15744 2008-11-16 18:19:00Z rouault $");
37
38 #include "../ogr/ogr_spatialref.h"
39
40 /************************************************************************/
41 /* __pure_virtual() */
42 /* */
43 /* The following is a gross hack to remove the last remaining */
44 /* dependency on the GNU C++ standard library. */
45 /************************************************************************/
46
47 #ifdef __GNUC__
48
49 extern "C"
__pure_virtual()50 void __pure_virtual()
51
52 {
53 }
54
55 #endif
56
57 /************************************************************************/
58 /* GDALDataTypeUnion() */
59 /************************************************************************/
60
61 /**
62 * Return the smallest data type that can fully express both input data
63 * types.
64 *
65 * @param eType1
66 * @param eType2
67 *
68 * @return a data type able to express eType1 and eType2.
69 */
70
71 GDALDataType CPL_STDCALL
GDALDataTypeUnion(GDALDataType eType1,GDALDataType eType2)72 GDALDataTypeUnion( GDALDataType eType1, GDALDataType eType2 )
73
74 {
75 int bFloating, bComplex, nBits, bSigned;
76
77 bComplex = GDALDataTypeIsComplex(eType1) | GDALDataTypeIsComplex(eType2);
78
79 switch( eType1 )
80 {
81 case GDT_Byte:
82 nBits = 8;
83 bSigned = FALSE;
84 bFloating = FALSE;
85 break;
86
87 case GDT_Int16:
88 case GDT_CInt16:
89 nBits = 16;
90 bSigned = TRUE;
91 bFloating = FALSE;
92 break;
93
94 case GDT_UInt16:
95 nBits = 16;
96 bSigned = FALSE;
97 bFloating = FALSE;
98 break;
99
100 case GDT_Int32:
101 case GDT_CInt32:
102 nBits = 32;
103 bSigned = TRUE;
104 bFloating = FALSE;
105 break;
106
107 case GDT_UInt32:
108 nBits = 32;
109 bSigned = FALSE;
110 bFloating = FALSE;
111 break;
112
113 case GDT_Float32:
114 case GDT_CFloat32:
115 nBits = 32;
116 bSigned = TRUE;
117 bFloating = TRUE;
118 break;
119
120 case GDT_Float64:
121 case GDT_CFloat64:
122 nBits = 64;
123 bSigned = TRUE;
124 bFloating = TRUE;
125 break;
126
127 default:
128 CPLAssert( FALSE );
129 return GDT_Unknown;
130 }
131
132 switch( eType2 )
133 {
134 case GDT_Byte:
135 break;
136
137 case GDT_Int16:
138 case GDT_CInt16:
139 nBits = MAX(nBits,16);
140 bSigned = TRUE;
141 break;
142
143 case GDT_UInt16:
144 nBits = MAX(nBits,16);
145 break;
146
147 case GDT_Int32:
148 case GDT_CInt32:
149 nBits = MAX(nBits,32);
150 bSigned = TRUE;
151 break;
152
153 case GDT_UInt32:
154 nBits = MAX(nBits,32);
155 break;
156
157 case GDT_Float32:
158 case GDT_CFloat32:
159 nBits = MAX(nBits,32);
160 bSigned = TRUE;
161 bFloating = TRUE;
162 break;
163
164 case GDT_Float64:
165 case GDT_CFloat64:
166 nBits = MAX(nBits,64);
167 bSigned = TRUE;
168 bFloating = TRUE;
169 break;
170
171 default:
172 CPLAssert( FALSE );
173 return GDT_Unknown;
174 }
175
176 if( nBits == 8 )
177 return GDT_Byte;
178 else if( nBits == 16 && bComplex )
179 return GDT_CInt16;
180 else if( nBits == 16 && bSigned )
181 return GDT_Int16;
182 else if( nBits == 16 && !bSigned )
183 return GDT_UInt16;
184 else if( nBits == 32 && bFloating && bComplex )
185 return GDT_CFloat32;
186 else if( nBits == 32 && bFloating )
187 return GDT_Float32;
188 else if( nBits == 32 && bComplex )
189 return GDT_CInt32;
190 else if( nBits == 32 && bSigned )
191 return GDT_Int32;
192 else if( nBits == 32 && !bSigned )
193 return GDT_UInt32;
194 else if( nBits == 64 && bComplex )
195 return GDT_CFloat64;
196 else
197 return GDT_Float64;
198 }
199
200
201 /************************************************************************/
202 /* GDALGetDataTypeSize() */
203 /************************************************************************/
204
205 /**
206 * Get data type size in bits.
207 *
208 * Returns the size of a a GDT_* type in bits, <b>not bytes</b>!
209 *
210 * @param data type, such as GDT_Byte.
211 * @return the number of bits or zero if it is not recognised.
212 */
213
GDALGetDataTypeSize(GDALDataType eDataType)214 int CPL_STDCALL GDALGetDataTypeSize( GDALDataType eDataType )
215
216 {
217 switch( eDataType )
218 {
219 case GDT_Byte:
220 return 8;
221
222 case GDT_UInt16:
223 case GDT_Int16:
224 return 16;
225
226 case GDT_UInt32:
227 case GDT_Int32:
228 case GDT_Float32:
229 case GDT_CInt16:
230 return 32;
231
232 case GDT_Float64:
233 case GDT_CInt32:
234 case GDT_CFloat32:
235 return 64;
236
237 case GDT_CFloat64:
238 return 128;
239
240 default:
241 CPLAssert( FALSE );
242 return 0;
243 }
244 }
245
246 /************************************************************************/
247 /* GDALDataTypeIsComplex() */
248 /************************************************************************/
249
250 /**
251 * Is data type complex?
252 *
253 * @return TRUE if the passed type is complex (one of GDT_CInt16, GDT_CInt32,
254 * GDT_CFloat32 or GDT_CFloat64), that is it consists of a real and imaginary
255 * component.
256 */
257
GDALDataTypeIsComplex(GDALDataType eDataType)258 int CPL_STDCALL GDALDataTypeIsComplex( GDALDataType eDataType )
259
260 {
261 switch( eDataType )
262 {
263 case GDT_CInt16:
264 case GDT_CInt32:
265 case GDT_CFloat32:
266 case GDT_CFloat64:
267 return TRUE;
268
269 default:
270 return FALSE;
271 }
272 }
273
274 /************************************************************************/
275 /* GDALGetDataTypeName() */
276 /************************************************************************/
277
278 /**
279 * Get name of data type.
280 *
281 * Returns a symbolic name for the data type. This is essentially the
282 * the enumerated item name with the GDT_ prefix removed. So GDT_Byte returns
283 * "Byte". The returned strings are static strings and should not be modified
284 * or freed by the application. These strings are useful for reporting
285 * datatypes in debug statements, errors and other user output.
286 *
287 * @param eDataType type to get name of.
288 * @return string corresponding to type.
289 */
290
GDALGetDataTypeName(GDALDataType eDataType)291 const char * CPL_STDCALL GDALGetDataTypeName( GDALDataType eDataType )
292
293 {
294 switch( eDataType )
295 {
296 case GDT_Unknown:
297 return "Unknown";
298
299 case GDT_Byte:
300 return "Byte";
301
302 case GDT_UInt16:
303 return "UInt16";
304
305 case GDT_Int16:
306 return "Int16";
307
308 case GDT_UInt32:
309 return "UInt32";
310
311 case GDT_Int32:
312 return "Int32";
313
314 case GDT_Float32:
315 return "Float32";
316
317 case GDT_Float64:
318 return "Float64";
319
320 case GDT_CInt16:
321 return "CInt16";
322
323 case GDT_CInt32:
324 return "CInt32";
325
326 case GDT_CFloat32:
327 return "CFloat32";
328
329 case GDT_CFloat64:
330 return "CFloat64";
331
332 default:
333 return NULL;
334 }
335 }
336
337 /************************************************************************/
338 /* GDALGetDataTypeByName() */
339 /************************************************************************/
340
341 /**
342 * Get data type by symbolic name.
343 *
344 * Returns a data type corresponding to the given symbolic name. This
345 * function is opposite to the GDALGetDataTypeName().
346 *
347 * @param pszName string containing the symbolic name of the type.
348 *
349 * @return GDAL data type.
350 */
351
GDALGetDataTypeByName(const char * pszName)352 GDALDataType CPL_STDCALL GDALGetDataTypeByName( const char *pszName )
353
354 {
355 VALIDATE_POINTER1( pszName, "GDALGetDataTypeByName", GDT_Unknown );
356
357 int iType;
358
359 for( iType = 1; iType < GDT_TypeCount; iType++ )
360 {
361 if( GDALGetDataTypeName((GDALDataType)iType) != NULL
362 && EQUAL(GDALGetDataTypeName((GDALDataType)iType), pszName) )
363 {
364 return (GDALDataType)iType;
365 }
366 }
367
368 return GDT_Unknown;
369 }
370
371 /************************************************************************/
372 /* GDALGetPaletteInterpretationName() */
373 /************************************************************************/
374
GDALGetPaletteInterpretationName(GDALPaletteInterp eInterp)375 const char *GDALGetPaletteInterpretationName( GDALPaletteInterp eInterp )
376
377 {
378 switch( eInterp )
379 {
380 case GPI_Gray:
381 return "Gray";
382
383 case GPI_RGB:
384 return "RGB";
385
386 case GPI_CMYK:
387 return "CMYK";
388
389 case GPI_HLS:
390 return "HLS";
391
392 default:
393 return "Unknown";
394 }
395 }
396
397 /************************************************************************/
398 /* GDALGetColorInterpretationName() */
399 /************************************************************************/
400
GDALGetColorInterpretationName(GDALColorInterp eInterp)401 const char *GDALGetColorInterpretationName( GDALColorInterp eInterp )
402
403 {
404 switch( eInterp )
405 {
406 case GCI_Undefined:
407 return "Undefined";
408
409 case GCI_GrayIndex:
410 return "Gray";
411
412 case GCI_PaletteIndex:
413 return "Palette";
414
415 case GCI_RedBand:
416 return "Red";
417
418 case GCI_GreenBand:
419 return "Green";
420
421 case GCI_BlueBand:
422 return "Blue";
423
424 case GCI_AlphaBand:
425 return "Alpha";
426
427 case GCI_HueBand:
428 return "Hue";
429
430 case GCI_SaturationBand:
431 return "Saturation";
432
433 case GCI_LightnessBand:
434 return "Lightness";
435
436 case GCI_CyanBand:
437 return "Cyan";
438
439 case GCI_MagentaBand:
440 return "Magenta";
441
442 case GCI_YellowBand:
443 return "Yellow";
444
445 case GCI_BlackBand:
446 return "Black";
447
448 case GCI_YCbCr_YBand:
449 return "YCbCr_Y";
450
451 case GCI_YCbCr_CbBand:
452 return "YCbCr_Cb";
453
454 case GCI_YCbCr_CrBand:
455 return "YCbCr_Cr";
456
457 default:
458 return "Unknown";
459 }
460 }
461
462 /************************************************************************/
463 /* GDALDummyProgress() */
464 /************************************************************************/
465
466 /**
467 * Stub progress function.
468 *
469 * This is a stub (does nothing) implementation of the GDALProgressFunc()
470 * semantics. It is primarily useful for passing to functions that take
471 * a GDALProgressFunc() argument but for which the application does not want
472 * to use one of the other progress functions that actually do something.
473 */
474
GDALDummyProgress(double,const char *,void *)475 int CPL_STDCALL GDALDummyProgress( double, const char *, void * )
476
477 {
478 return TRUE;
479 }
480
481 /************************************************************************/
482 /* GDALScaledProgress() */
483 /************************************************************************/
484 typedef struct {
485 GDALProgressFunc pfnProgress;
486 void *pData;
487 double dfMin;
488 double dfMax;
489 } GDALScaledProgressInfo;
490
491 /**
492 * Scaled progress transformer.
493 *
494 * This is the progress function that should be passed along with the
495 * callback data returned by GDALCreateScaledProgress().
496 */
497
GDALScaledProgress(double dfComplete,const char * pszMessage,void * pData)498 int CPL_STDCALL GDALScaledProgress( double dfComplete, const char *pszMessage,
499 void *pData )
500
501 {
502 GDALScaledProgressInfo *psInfo = (GDALScaledProgressInfo *) pData;
503
504 return psInfo->pfnProgress( dfComplete * (psInfo->dfMax - psInfo->dfMin)
505 + psInfo->dfMin,
506 pszMessage, psInfo->pData );
507 }
508
509 /************************************************************************/
510 /* GDALCreateScaledProgress() */
511 /************************************************************************/
512
513 /**
514 * Create scaled progress transformer.
515 *
516 * Sometimes when an operations wants to report progress it actually
517 * invokes several subprocesses which also take GDALProgressFunc()s,
518 * and it is desirable to map the progress of each sub operation into
519 * a portion of 0.0 to 1.0 progress of the overall process. The scaled
520 * progress function can be used for this.
521 *
522 * For each subsection a scaled progress function is created and
523 * instead of passing the overall progress func down to the sub functions,
524 * the GDALScaledProgress() function is passed instead.
525 *
526 * @param dfMin the value to which 0.0 in the sub operation is mapped.
527 * @param dfMax the value to which 1.0 is the sub operation is mapped.
528 * @param pfnProgress the overall progress function.
529 * @param pData the overall progress function callback data.
530 *
531 * @return pointer to pass as pProgressArg to sub functions. Should be freed
532 * with GDALDestroyScaledProgress().
533 *
534 * Example:
535 *
536 * \code
537 * int MyOperation( ..., GDALProgressFunc pfnProgress, void *pProgressData );
538 *
539 * {
540 * void *pScaledProgress;
541 *
542 * pScaledProgress = GDALCreateScaledProgress( 0.0, 0.5, pfnProgress,
543 * pProgressData );
544 * GDALDoLongSlowOperation( ..., GDALScaledProgress, pScaledProgress );
545 * GDALDestroyScaledProgress( pScaledProgress );
546 *
547 * pScaledProgress = GDALCreateScaledProgress( 0.5, 1.0, pfnProgress,
548 * pProgressData );
549 * GDALDoAnotherOperation( ..., GDALScaledProgress, pScaledProgress );
550 * GDALDestroyScaledProgress( pScaledProgress );
551 *
552 * return ...;
553 * }
554 * \endcode
555 */
556
GDALCreateScaledProgress(double dfMin,double dfMax,GDALProgressFunc pfnProgress,void * pData)557 void * CPL_STDCALL GDALCreateScaledProgress( double dfMin, double dfMax,
558 GDALProgressFunc pfnProgress,
559 void * pData )
560
561 {
562 GDALScaledProgressInfo *psInfo;
563
564 psInfo = (GDALScaledProgressInfo *)
565 CPLCalloc(sizeof(GDALScaledProgressInfo),1);
566
567 if( ABS(dfMin-dfMax) < 0.0000001 )
568 dfMax = dfMin + 0.01;
569
570 psInfo->pData = pData;
571 psInfo->pfnProgress = pfnProgress;
572 psInfo->dfMin = dfMin;
573 psInfo->dfMax = dfMax;
574
575 return (void *) psInfo;
576 }
577
578 /************************************************************************/
579 /* GDALDestroyScaledProgress() */
580 /************************************************************************/
581
582 /**
583 * Cleanup scaled progress handle.
584 *
585 * This function cleans up the data associated with a scaled progress function
586 * as returned by GADLCreateScaledProgress().
587 *
588 * @param pData scaled progress handle returned by GDALCreateScaledProgress().
589 */
590
GDALDestroyScaledProgress(void * pData)591 void CPL_STDCALL GDALDestroyScaledProgress( void * pData )
592
593 {
594 CPLFree( pData );
595 }
596
597 /************************************************************************/
598 /* GDALTermProgress() */
599 /************************************************************************/
600
601 /**
602 * Simple progress report to terminal.
603 *
604 * This progress reporter prints simple progress report to the
605 * terminal window. The progress report generally looks something like
606 * this:
607
608 \verbatim
609 0...10...20...30...40...50...60...70...80...90...100 - done.
610 \endverbatim
611
612 * Every 2.5% of progress another number or period is emitted. Note that
613 * GDALTermProgress() uses internal static data to keep track of the last
614 * percentage reported and will get confused if two terminal based progress
615 * reportings are active at the same time.
616 *
617 * The GDALTermProgress() function maintains an internal memory of the
618 * last percentage complete reported in a static variable, and this makes
619 * it unsuitable to have multiple GDALTermProgress()'s active eithin a
620 * single thread or across multiple threads.
621 *
622 * @param dfComplete completion ratio from 0.0 to 1.0.
623 * @param pszMessage optional message.
624 * @param pProgressArg ignored callback data argument.
625 *
626 * @return Always returns TRUE indicating the process should continue.
627 */
628
GDALTermProgress(double dfComplete,const char * pszMessage,void * pProgressArg)629 int CPL_STDCALL GDALTermProgress( double dfComplete, const char *pszMessage,
630 void * pProgressArg )
631
632 {
633 static int nLastTick = -1;
634 int nThisTick = (int) (dfComplete * 40.0);
635
636 (void) pProgressArg;
637
638 nThisTick = MIN(40,MAX(0,nThisTick));
639
640 // Have we started a new progress run?
641 if( nThisTick < nLastTick && nLastTick >= 39 )
642 nLastTick = -1;
643
644 if( nThisTick <= nLastTick )
645 return TRUE;
646
647 while( nThisTick > nLastTick )
648 {
649 nLastTick++;
650 if( nLastTick % 4 == 0 )
651 fprintf( stdout, "%d", (nLastTick / 4) * 10 );
652 else
653 fprintf( stdout, "." );
654 }
655
656 if( nThisTick == 40 )
657 fprintf( stdout, " - done.\n" );
658 else
659 fflush( stdout );
660
661 return TRUE;
662 }
663
664 /************************************************************************/
665 /* GDALGetRandomRasterSample() */
666 /************************************************************************/
667
668 int CPL_STDCALL
GDALGetRandomRasterSample(GDALRasterBandH hBand,int nSamples,float * pafSampleBuf)669 GDALGetRandomRasterSample( GDALRasterBandH hBand, int nSamples,
670 float *pafSampleBuf )
671
672 {
673 VALIDATE_POINTER1( hBand, "GDALGetRandomRasterSample", 0 );
674
675 GDALRasterBand *poBand;
676
677 poBand = (GDALRasterBand *) GDALGetRasterSampleOverview( hBand, nSamples );
678 CPLAssert( NULL != poBand );
679
680 /* -------------------------------------------------------------------- */
681 /* Figure out the ratio of blocks we will read to get an */
682 /* approximate value. */
683 /* -------------------------------------------------------------------- */
684 int nBlockXSize, nBlockYSize;
685 int nBlocksPerRow, nBlocksPerColumn;
686 int nSampleRate;
687 int bGotNoDataValue;
688 double dfNoDataValue;
689 int nActualSamples = 0;
690 int nBlockSampleRate;
691 int nBlockPixels, nBlockCount;
692
693 dfNoDataValue = poBand->GetNoDataValue( &bGotNoDataValue );
694
695 poBand->GetBlockSize( &nBlockXSize, &nBlockYSize );
696
697 nBlocksPerRow = (poBand->GetXSize() + nBlockXSize - 1) / nBlockXSize;
698 nBlocksPerColumn = (poBand->GetYSize() + nBlockYSize - 1) / nBlockYSize;
699
700 nBlockPixels = nBlockXSize * nBlockYSize;
701 nBlockCount = nBlocksPerRow * nBlocksPerColumn;
702
703 if( nBlocksPerRow == 0 || nBlocksPerColumn == 0 || nBlockPixels == 0
704 || nBlockCount == 0 )
705 {
706 CPLError( CE_Failure, CPLE_AppDefined,
707 "GDALGetRandomRasterSample(): returning because band"
708 " appears degenerate." );
709
710 return FALSE;
711 }
712
713 nSampleRate = (int) MAX(1,sqrt((double) nBlockCount)-2.0);
714
715 if( nSampleRate == nBlocksPerRow && nSampleRate > 1 )
716 nSampleRate--;
717
718 while( nSampleRate > 1
719 && ((nBlockCount-1) / nSampleRate + 1) * nBlockPixels < nSamples )
720 nSampleRate--;
721
722 if ((nSamples / ((nBlockCount-1) / nSampleRate + 1)) == 0)
723 nBlockSampleRate = 1;
724 else
725 nBlockSampleRate =
726 MAX(1,nBlockPixels / (nSamples / ((nBlockCount-1) / nSampleRate + 1)));
727
728 for( int iSampleBlock = 0;
729 iSampleBlock < nBlockCount;
730 iSampleBlock += nSampleRate )
731 {
732 double dfValue = 0.0, dfReal, dfImag;
733 int iXBlock, iYBlock, iX, iY, iXValid, iYValid, iRemainder = 0;
734 GDALRasterBlock *poBlock;
735
736 iYBlock = iSampleBlock / nBlocksPerRow;
737 iXBlock = iSampleBlock - nBlocksPerRow * iYBlock;
738
739 poBlock = poBand->GetLockedBlockRef( iXBlock, iYBlock );
740 if( poBlock == NULL )
741 continue;
742
743 if( (iXBlock + 1) * nBlockXSize > poBand->GetXSize() )
744 iXValid = poBand->GetXSize() - iXBlock * nBlockXSize;
745 else
746 iXValid = nBlockXSize;
747
748 if( (iYBlock + 1) * nBlockYSize > poBand->GetYSize() )
749 iYValid = poBand->GetYSize() - iYBlock * nBlockYSize;
750 else
751 iYValid = nBlockYSize;
752
753 for( iY = 0; iY < iYValid; iY++ )
754 {
755 for( iX = iRemainder; iX < iXValid; iX += nBlockSampleRate )
756 {
757 int iOffset;
758
759 iOffset = iX + iY * nBlockXSize;
760 switch( poBlock->GetDataType() )
761 {
762 case GDT_Byte:
763 dfValue = ((GByte *) poBlock->GetDataRef())[iOffset];
764 break;
765 case GDT_UInt16:
766 dfValue = ((GUInt16 *) poBlock->GetDataRef())[iOffset];
767 break;
768 case GDT_Int16:
769 dfValue = ((GInt16 *) poBlock->GetDataRef())[iOffset];
770 break;
771 case GDT_UInt32:
772 dfValue = ((GUInt32 *) poBlock->GetDataRef())[iOffset];
773 break;
774 case GDT_Int32:
775 dfValue = ((GInt32 *) poBlock->GetDataRef())[iOffset];
776 break;
777 case GDT_Float32:
778 dfValue = ((float *) poBlock->GetDataRef())[iOffset];
779 break;
780 case GDT_Float64:
781 dfValue = ((double *) poBlock->GetDataRef())[iOffset];
782 break;
783 case GDT_CInt16:
784 dfReal = ((GInt16 *) poBlock->GetDataRef())[iOffset*2];
785 dfImag = ((GInt16 *) poBlock->GetDataRef())[iOffset*2+1];
786 dfValue = sqrt(dfReal*dfReal + dfImag*dfImag);
787 break;
788 case GDT_CInt32:
789 dfReal = ((GInt32 *) poBlock->GetDataRef())[iOffset*2];
790 dfImag = ((GInt32 *) poBlock->GetDataRef())[iOffset*2+1];
791 dfValue = sqrt(dfReal*dfReal + dfImag*dfImag);
792 break;
793 case GDT_CFloat32:
794 dfReal = ((float *) poBlock->GetDataRef())[iOffset*2];
795 dfImag = ((float *) poBlock->GetDataRef())[iOffset*2+1];
796 dfValue = sqrt(dfReal*dfReal + dfImag*dfImag);
797 break;
798 case GDT_CFloat64:
799 dfReal = ((double *) poBlock->GetDataRef())[iOffset*2];
800 dfImag = ((double *) poBlock->GetDataRef())[iOffset*2+1];
801 dfValue = sqrt(dfReal*dfReal + dfImag*dfImag);
802 break;
803 default:
804 CPLAssert( FALSE );
805 }
806
807 if( bGotNoDataValue && dfValue == dfNoDataValue )
808 continue;
809
810 if( nActualSamples < nSamples )
811 pafSampleBuf[nActualSamples++] = (float) dfValue;
812 }
813
814 iRemainder = iX - iXValid;
815 }
816
817 poBlock->DropLock();
818 }
819
820 return nActualSamples;
821 }
822
823 /************************************************************************/
824 /* GDALInitGCPs() */
825 /************************************************************************/
826
GDALInitGCPs(int nCount,GDAL_GCP * psGCP)827 void CPL_STDCALL GDALInitGCPs( int nCount, GDAL_GCP * psGCP )
828
829 {
830 if( nCount > 0 )
831 {
832 VALIDATE_POINTER0( psGCP, "GDALInitGCPs" );
833 }
834
835 for( int iGCP = 0; iGCP < nCount; iGCP++ )
836 {
837 memset( psGCP, 0, sizeof(GDAL_GCP) );
838 psGCP->pszId = CPLStrdup("");
839 psGCP->pszInfo = CPLStrdup("");
840 psGCP++;
841 }
842 }
843
844 /************************************************************************/
845 /* GDALDeinitGCPs() */
846 /************************************************************************/
847
GDALDeinitGCPs(int nCount,GDAL_GCP * psGCP)848 void CPL_STDCALL GDALDeinitGCPs( int nCount, GDAL_GCP * psGCP )
849
850 {
851 if ( nCount > 0 )
852 {
853 VALIDATE_POINTER0( psGCP, "GDALDeinitGCPs" );
854 }
855
856 for( int iGCP = 0; iGCP < nCount; iGCP++ )
857 {
858 CPLFree( psGCP->pszId );
859 CPLFree( psGCP->pszInfo );
860 psGCP++;
861 }
862 }
863
864 /************************************************************************/
865 /* GDALDuplicateGCPs() */
866 /************************************************************************/
867
868 GDAL_GCP * CPL_STDCALL
GDALDuplicateGCPs(int nCount,const GDAL_GCP * pasGCPList)869 GDALDuplicateGCPs( int nCount, const GDAL_GCP *pasGCPList )
870
871 {
872 GDAL_GCP *pasReturn;
873
874 pasReturn = (GDAL_GCP *) CPLMalloc(sizeof(GDAL_GCP) * nCount);
875 GDALInitGCPs( nCount, pasReturn );
876
877 for( int iGCP = 0; iGCP < nCount; iGCP++ )
878 {
879 CPLFree( pasReturn[iGCP].pszId );
880 pasReturn[iGCP].pszId = CPLStrdup( pasGCPList[iGCP].pszId );
881
882 CPLFree( pasReturn[iGCP].pszInfo );
883 pasReturn[iGCP].pszInfo = CPLStrdup( pasGCPList[iGCP].pszInfo );
884
885 pasReturn[iGCP].dfGCPPixel = pasGCPList[iGCP].dfGCPPixel;
886 pasReturn[iGCP].dfGCPLine = pasGCPList[iGCP].dfGCPLine;
887 pasReturn[iGCP].dfGCPX = pasGCPList[iGCP].dfGCPX;
888 pasReturn[iGCP].dfGCPY = pasGCPList[iGCP].dfGCPY;
889 pasReturn[iGCP].dfGCPZ = pasGCPList[iGCP].dfGCPZ;
890 }
891
892 return pasReturn;
893 }
894
895 /************************************************************************/
896 /* GDALLoadTabFile() */
897 /* */
898 /* Helper function for translator implementators wanting */
899 /* support for MapInfo .tab-files. */
900 /************************************************************************/
901
902 #define MAX_GCP 256
903
GDALLoadTabFile(const char * pszFilename,double * padfGeoTransform,char ** ppszWKT,int * pnGCPCount,GDAL_GCP ** ppasGCPs)904 int CPL_STDCALL GDALLoadTabFile( const char *pszFilename,
905 double *padfGeoTransform, char **ppszWKT,
906 int *pnGCPCount, GDAL_GCP **ppasGCPs )
907
908
909 {
910 char **papszLines;
911 char **papszTok=NULL;
912 int bTypeRasterFound = FALSE;
913 int bInsideTableDef = FALSE;
914 int iLine, numLines=0;
915 int nCoordinateCount = 0;
916 GDAL_GCP asGCPs[MAX_GCP];
917
918 papszLines = CSLLoad( pszFilename );
919
920 if ( !papszLines )
921 return FALSE;
922
923 numLines = CSLCount(papszLines);
924
925 // Iterate all lines in the TAB-file
926 for(iLine=0; iLine<numLines; iLine++)
927 {
928 CSLDestroy(papszTok);
929 papszTok = CSLTokenizeStringComplex(papszLines[iLine], " \t(),;",
930 TRUE, FALSE);
931
932 if (CSLCount(papszTok) < 2)
933 continue;
934
935 // Did we find table definition
936 if (EQUAL(papszTok[0], "Definition") && EQUAL(papszTok[1], "Table") )
937 {
938 bInsideTableDef = TRUE;
939 }
940 else if (bInsideTableDef && (EQUAL(papszTok[0], "Type")) )
941 {
942 // Only RASTER-type will be handled
943 if (EQUAL(papszTok[1], "RASTER"))
944 {
945 bTypeRasterFound = TRUE;
946 }
947 else
948 {
949 CSLDestroy(papszTok);
950 CSLDestroy(papszLines);
951 return FALSE;
952 }
953 }
954 else if (bTypeRasterFound && bInsideTableDef
955 && CSLCount(papszTok) > 4
956 && EQUAL(papszTok[4], "Label")
957 && nCoordinateCount < MAX_GCP )
958 {
959 GDALInitGCPs( 1, asGCPs + nCoordinateCount );
960
961 asGCPs[nCoordinateCount].dfGCPPixel = CPLAtofM(papszTok[2]);
962 asGCPs[nCoordinateCount].dfGCPLine = CPLAtofM(papszTok[3]);
963 asGCPs[nCoordinateCount].dfGCPX = CPLAtofM(papszTok[0]);
964 asGCPs[nCoordinateCount].dfGCPY = CPLAtofM(papszTok[1]);
965 if( papszTok[5] != NULL )
966 {
967 CPLFree( asGCPs[nCoordinateCount].pszId );
968 asGCPs[nCoordinateCount].pszId = CPLStrdup(papszTok[5]);
969 }
970
971 nCoordinateCount++;
972 }
973 else if( bTypeRasterFound && bInsideTableDef
974 && EQUAL(papszTok[0],"CoordSys")
975 && ppszWKT != NULL )
976 {
977 OGRSpatialReference oSRS;
978
979 if( oSRS.importFromMICoordSys( papszLines[iLine] ) == OGRERR_NONE )
980 oSRS.exportToWkt( ppszWKT );
981 }
982 else if( EQUAL(papszTok[0],"Units")
983 && CSLCount(papszTok) > 1
984 && EQUAL(papszTok[1],"degree") )
985 {
986 /*
987 ** If we have units of "degree", but a projected coordinate
988 ** system we need to convert it to geographic. See to01_02.TAB.
989 */
990 if( ppszWKT != NULL && *ppszWKT != NULL
991 && EQUALN(*ppszWKT,"PROJCS",6) )
992 {
993 OGRSpatialReference oSRS, oSRSGeogCS;
994 char *pszSrcWKT = *ppszWKT;
995
996 oSRS.importFromWkt( &pszSrcWKT );
997 oSRSGeogCS.CopyGeogCSFrom( &oSRS );
998 CPLFree( *ppszWKT );
999 oSRSGeogCS.exportToWkt( ppszWKT );
1000 }
1001 }
1002
1003 }
1004
1005 CSLDestroy(papszTok);
1006 CSLDestroy(papszLines);
1007
1008 if( nCoordinateCount == 0 )
1009 {
1010 CPLDebug( "GDAL", "GDALLoadTabFile(%s) did not get any GCPs.",
1011 pszFilename );
1012 return FALSE;
1013 }
1014
1015 /* -------------------------------------------------------------------- */
1016 /* Try to convert the GCPs into a geotransform definition, if */
1017 /* possible. Otherwise we will need to use them as GCPs. */
1018 /* -------------------------------------------------------------------- */
1019 if( !GDALGCPsToGeoTransform( nCoordinateCount, asGCPs, padfGeoTransform,
1020 FALSE ) )
1021 {
1022 if (pnGCPCount && ppasGCPs)
1023 {
1024 CPLDebug( "GDAL",
1025 "GDALLoadTabFile(%s) found file, wasn't able to derive a\n"
1026 "first order geotransform. Using points as GCPs.",
1027 pszFilename );
1028
1029 *ppasGCPs = (GDAL_GCP *)
1030 CPLCalloc( sizeof(GDAL_GCP),nCoordinateCount );
1031 memcpy( *ppasGCPs, asGCPs, sizeof(GDAL_GCP) * nCoordinateCount );
1032 *pnGCPCount = nCoordinateCount;
1033 }
1034 }
1035 else
1036 {
1037 GDALDeinitGCPs( nCoordinateCount, asGCPs );
1038 }
1039
1040 return TRUE;
1041 }
1042
1043 /************************************************************************/
1044 /* GDALReadTabFile() */
1045 /* */
1046 /* Helper function for translator implementators wanting */
1047 /* support for MapInfo .tab-files. */
1048 /************************************************************************/
1049
GDALReadTabFile(const char * pszBaseFilename,double * padfGeoTransform,char ** ppszWKT,int * pnGCPCount,GDAL_GCP ** ppasGCPs)1050 int CPL_STDCALL GDALReadTabFile( const char * pszBaseFilename,
1051 double *padfGeoTransform, char **ppszWKT,
1052 int *pnGCPCount, GDAL_GCP **ppasGCPs )
1053
1054
1055 {
1056 const char *pszTAB;
1057 FILE *fpTAB;
1058
1059 /* -------------------------------------------------------------------- */
1060 /* Try lower case, then upper case. */
1061 /* -------------------------------------------------------------------- */
1062 pszTAB = CPLResetExtension( pszBaseFilename, "tab" );
1063
1064 fpTAB = VSIFOpen( pszTAB, "rt" );
1065
1066 #ifndef WIN32
1067 if( fpTAB == NULL )
1068 {
1069 pszTAB = CPLResetExtension( pszBaseFilename, "TAB" );
1070 fpTAB = VSIFOpen( pszTAB, "rt" );
1071 }
1072 #endif
1073
1074 if( fpTAB == NULL )
1075 return FALSE;
1076
1077 VSIFClose( fpTAB );
1078
1079 /* -------------------------------------------------------------------- */
1080 /* We found the file, now load and parse it. */
1081 /* -------------------------------------------------------------------- */
1082 return GDALLoadTabFile( pszTAB, padfGeoTransform, ppszWKT,
1083 pnGCPCount, ppasGCPs );
1084 }
1085
1086 /************************************************************************/
1087 /* GDALLoadWorldFile() */
1088 /************************************************************************/
1089
1090 /**
1091 * Read ESRI world file.
1092 *
1093 * This function reads an ESRI style world file, and formats a geotransform
1094 * from its contents.
1095 *
1096 * The world file contains an affine transformation with the parameters
1097 * in a different order than in a geotransform array.
1098 *
1099 * <ul>
1100 * <li> geotransform[1] : width of pixel
1101 * <li> geotransform[4] : rotational coefficient, zero for north up images.
1102 * <li> geotransform[2] : rotational coefficient, zero for north up images.
1103 * <li> geotransform[5] : height of pixel (but negative)
1104 * <li> geotransform[0] + 0.5 * geotransform[1] + 0.5 * geotransform[2] : x offset to center of top left pixel.
1105 * <li> geotransform[3] + 0.5 * geotransform[4] + 0.5 * geotransform[5] : y offset to center of top left pixel.
1106 * </ul>
1107 *
1108 * @param pszFilename the world file name.
1109 * @param padfGeoTransform the six double array into which the
1110 * geotransformation should be placed.
1111 *
1112 * @return TRUE on success or FALSE on failure.
1113 */
1114
1115 int CPL_STDCALL
GDALLoadWorldFile(const char * pszFilename,double * padfGeoTransform)1116 GDALLoadWorldFile( const char *pszFilename, double *padfGeoTransform )
1117
1118 {
1119 int nLinesCount = 0;
1120 char **papszLines = NULL;
1121 bool bValid = false;
1122
1123 VALIDATE_POINTER1( pszFilename, "GDALLoadWorldFile", FALSE );
1124 VALIDATE_POINTER1( padfGeoTransform, "GDALLoadWorldFile", FALSE );
1125
1126 papszLines = CSLLoad( pszFilename );
1127
1128 if ( !papszLines )
1129 return FALSE;
1130
1131 /* Handling of blank lines is not supported. */
1132 nLinesCount = CSLCount(papszLines);
1133 if( nLinesCount >= 6 )
1134 {
1135 int i;
1136 bValid = true;
1137
1138 /* First six lines of a world file can not be empty. */
1139 for ( i = 0; i < 6; i++ )
1140 {
1141 CPLString line(papszLines[i]);
1142 if( line.Trim().empty() )
1143 {
1144 bValid = false;
1145 break;
1146 }
1147 }
1148 }
1149
1150 if( bValid
1151 && (CPLAtofM(papszLines[0]) != 0.0 || CPLAtofM(papszLines[2]) != 0.0)
1152 && (CPLAtofM(papszLines[3]) != 0.0 || CPLAtofM(papszLines[1]) != 0.0) )
1153 {
1154 padfGeoTransform[0] = CPLAtofM(papszLines[4]);
1155 padfGeoTransform[1] = CPLAtofM(papszLines[0]);
1156 padfGeoTransform[2] = CPLAtofM(papszLines[2]);
1157 padfGeoTransform[3] = CPLAtofM(papszLines[5]);
1158 padfGeoTransform[4] = CPLAtofM(papszLines[1]);
1159 padfGeoTransform[5] = CPLAtofM(papszLines[3]);
1160
1161 // correct for center of pixel vs. top left of pixel
1162 padfGeoTransform[0] -= 0.5 * padfGeoTransform[1];
1163 padfGeoTransform[0] -= 0.5 * padfGeoTransform[2];
1164 padfGeoTransform[3] -= 0.5 * padfGeoTransform[4];
1165 padfGeoTransform[3] -= 0.5 * padfGeoTransform[5];
1166
1167 CSLDestroy(papszLines);
1168
1169 return TRUE;
1170 }
1171 else
1172 {
1173 CPLDebug( "GDAL",
1174 "GDALLoadWorldFile(%s) found file, but it was corrupt.",
1175 pszFilename );
1176 CSLDestroy(papszLines);
1177 return FALSE;
1178 }
1179 }
1180
1181 /************************************************************************/
1182 /* GDALReadWorldFile() */
1183 /************************************************************************/
1184
1185 /**
1186 * Read ESRI world file.
1187 *
1188 * This function reads an ESRI style world file, and formats a geotransform
1189 * from it's contents. It does the same as GDALLoadWorldFile() function, but
1190 * it will form the filename for the worldfile from the filename of the raster
1191 * file referred and the suggested extension. If no extension is provided,
1192 * the code will internally try the unix style and windows style world file
1193 * extensions (eg. for .tif these would be .tfw and .tifw).
1194 *
1195 * The world file contains an affine transformation with the parameters
1196 * in a different order than in a geotransform array.
1197 *
1198 * <ul>
1199 * <li> geotransform[1] : width of pixel
1200 * <li> geotransform[4] : rotational coefficient, zero for north up images.
1201 * <li> geotransform[2] : rotational coefficient, zero for north up images.
1202 * <li> geotransform[5] : height of pixel (but negative)
1203 * <li> geotransform[0] + 0.5 * geotransform[1] + 0.5 * geotransform[2] : x offset to center of top left pixel.
1204 * <li> geotransform[3] + 0.5 * geotransform[4] + 0.5 * geotransform[5] : y offset to center of top left pixel.
1205 * </ul>
1206 *
1207 * @param pszBaseFilename the target raster file.
1208 * @param pszExtension the extension to use (ie. ".wld") or NULL to derive it
1209 * from the pszBaseFilename
1210 * @param padfGeoTransform the six double array into which the
1211 * geotransformation should be placed.
1212 *
1213 * @return TRUE on success or FALSE on failure.
1214 */
1215
1216 int CPL_STDCALL
GDALReadWorldFile(const char * pszBaseFilename,const char * pszExtension,double * padfGeoTransform)1217 GDALReadWorldFile( const char *pszBaseFilename, const char *pszExtension,
1218 double *padfGeoTransform )
1219
1220 {
1221 const char *pszTFW;
1222 char szExtUpper[32], szExtLower[32];
1223 int i;
1224
1225 VALIDATE_POINTER1( pszBaseFilename, "GDALReadWorldFile", FALSE );
1226 VALIDATE_POINTER1( padfGeoTransform, "GDALReadWorldFile", FALSE );
1227
1228 /* -------------------------------------------------------------------- */
1229 /* If we aren't given an extension, try both the unix and */
1230 /* windows style extensions. */
1231 /* -------------------------------------------------------------------- */
1232 if( pszExtension == NULL )
1233 {
1234 char szDerivedExtension[100];
1235 std::string oBaseExt = CPLGetExtension( pszBaseFilename );
1236
1237 if( oBaseExt.length() < 2 )
1238 return FALSE;
1239
1240 // windows version - first + last + 'w'
1241 szDerivedExtension[0] = oBaseExt[0];
1242 szDerivedExtension[1] = oBaseExt[oBaseExt.length()-1];
1243 szDerivedExtension[2] = 'w';
1244 szDerivedExtension[3] = '\0';
1245
1246 if( GDALReadWorldFile( pszBaseFilename, szDerivedExtension,
1247 padfGeoTransform ) )
1248 return TRUE;
1249
1250 // unix version - extension + 'w'
1251 if( oBaseExt.length() > sizeof(szDerivedExtension)-2 )
1252 return FALSE;
1253
1254 strcpy( szDerivedExtension, oBaseExt.c_str() );
1255 strcat( szDerivedExtension, "w" );
1256 return GDALReadWorldFile( pszBaseFilename, szDerivedExtension,
1257 padfGeoTransform );
1258 }
1259
1260 /* -------------------------------------------------------------------- */
1261 /* Skip the leading period in the extension if there is one. */
1262 /* -------------------------------------------------------------------- */
1263 if( *pszExtension == '.' )
1264 pszExtension++;
1265
1266 /* -------------------------------------------------------------------- */
1267 /* Generate upper and lower case versions of the extension. */
1268 /* -------------------------------------------------------------------- */
1269 strncpy( szExtUpper, pszExtension, 32 );
1270 strncpy( szExtLower, pszExtension, 32 );
1271
1272 for( i = 0; szExtUpper[i] != '\0' && i < 32; i++ )
1273 {
1274 szExtUpper[i] = (char) toupper(szExtUpper[i]);
1275 szExtLower[i] = (char) tolower(szExtLower[i]);
1276 }
1277
1278 /* -------------------------------------------------------------------- */
1279 /* Try lower case, then upper case. */
1280 /* -------------------------------------------------------------------- */
1281 VSIStatBufL sStatBuf;
1282 int bGotTFW;
1283
1284 pszTFW = CPLResetExtension( pszBaseFilename, szExtLower );
1285
1286 bGotTFW = VSIStatL( pszTFW, &sStatBuf ) == 0;
1287
1288 #ifndef WIN32
1289 if( !bGotTFW )
1290 {
1291 pszTFW = CPLResetExtension( pszBaseFilename, szExtUpper );
1292 bGotTFW = VSIStatL( pszTFW, &sStatBuf ) == 0;
1293 }
1294 #endif
1295
1296 if( !bGotTFW )
1297 return FALSE;
1298
1299 /* -------------------------------------------------------------------- */
1300 /* We found the file, now load and parse it. */
1301 /* -------------------------------------------------------------------- */
1302 return GDALLoadWorldFile( pszTFW, padfGeoTransform );
1303 }
1304
1305 /************************************************************************/
1306 /* GDALWriteWorldFile() */
1307 /* */
1308 /* Helper function for translator implementators wanting */
1309 /* support for ESRI world files. */
1310 /************************************************************************/
1311
1312 /**
1313 * Write ESRI world file.
1314 *
1315 * This function writes an ESRI style world file from the passed geotransform.
1316 *
1317 * The world file contains an affine transformation with the parameters
1318 * in a different order than in a geotransform array.
1319 *
1320 * <ul>
1321 * <li> geotransform[1] : width of pixel
1322 * <li> geotransform[4] : rotational coefficient, zero for north up images.
1323 * <li> geotransform[2] : rotational coefficient, zero for north up images.
1324 * <li> geotransform[5] : height of pixel (but negative)
1325 * <li> geotransform[0] + 0.5 * geotransform[1] + 0.5 * geotransform[2] : x offset to center of top left pixel.
1326 * <li> geotransform[3] + 0.5 * geotransform[4] + 0.5 * geotransform[5] : y offset to center of top left pixel.
1327 * </ul>
1328 *
1329 * @param pszBaseFilename the target raster file.
1330 * @param pszExtension the extension to use (ie. ".wld"). Must not be NULL
1331 * @param padfGeoTransform the six double array from which the
1332 * geotransformation should be read.
1333 *
1334 * @return TRUE on success or FALSE on failure.
1335 */
1336
1337 int CPL_STDCALL
GDALWriteWorldFile(const char * pszBaseFilename,const char * pszExtension,double * padfGeoTransform)1338 GDALWriteWorldFile( const char * pszBaseFilename, const char *pszExtension,
1339 double *padfGeoTransform )
1340
1341 {
1342 VALIDATE_POINTER1( pszBaseFilename, "GDALWriteWorldFile", FALSE );
1343 VALIDATE_POINTER1( pszExtension, "GDALWriteWorldFile", FALSE );
1344 VALIDATE_POINTER1( padfGeoTransform, "GDALWriteWorldFile", FALSE );
1345
1346 /* -------------------------------------------------------------------- */
1347 /* Prepare the text to write to the file. */
1348 /* -------------------------------------------------------------------- */
1349 CPLString osTFWText;
1350
1351 osTFWText.Printf( "%.10f\n%.10f\n%.10f\n%.10f\n%.10f\n%.10f\n",
1352 padfGeoTransform[1],
1353 padfGeoTransform[4],
1354 padfGeoTransform[2],
1355 padfGeoTransform[5],
1356 padfGeoTransform[0]
1357 + 0.5 * padfGeoTransform[1]
1358 + 0.5 * padfGeoTransform[2],
1359 padfGeoTransform[3]
1360 + 0.5 * padfGeoTransform[4]
1361 + 0.5 * padfGeoTransform[5] );
1362
1363 /* -------------------------------------------------------------------- */
1364 /* Update extention, and write to disk. */
1365 /* -------------------------------------------------------------------- */
1366 const char *pszTFW;
1367 FILE *fpTFW;
1368
1369 pszTFW = CPLResetExtension( pszBaseFilename, pszExtension );
1370 fpTFW = VSIFOpenL( pszTFW, "wt" );
1371 if( fpTFW == NULL )
1372 return FALSE;
1373
1374 VSIFWriteL( (void *) osTFWText.c_str(), 1, osTFWText.size(), fpTFW );
1375 VSIFCloseL( fpTFW );
1376
1377 return TRUE;
1378 }
1379
1380 /************************************************************************/
1381 /* GDALVersionInfo() */
1382 /************************************************************************/
1383
1384 /**
1385 * Get runtime version information.
1386 *
1387 * Available pszRequest values:
1388 * <ul>
1389 * <li> "VERSION_NUM": Returns GDAL_VERSION_NUM formatted as a string. ie. "1170"
1390 * <li> "RELEASE_DATE": Returns GDAL_RELEASE_DATE formatted as a string.
1391 * ie. "20020416".
1392 * <li> "RELEASE_NAME": Returns the GDAL_RELEASE_NAME. ie. "1.1.7"
1393 * <li> "--version": Returns one line version message suitable for use in
1394 * response to --version requests. ie. "GDAL 1.1.7, released 2002/04/16"
1395 * <li> "LICENCE": Returns the content of the LICENSE.TXT file from the GDAL_DATA directory.
1396 * Currently, the returned string is leaking memory but applications are discouraged from
1397 * deallocating the returned memory in this case since in the
1398 * future we may resolve the leak issue internally.
1399 * </ul>
1400 *
1401 * @param pszRequest the type of version info desired, as listed above.
1402 *
1403 * @return an internal string containing the requested information.
1404 */
1405
GDALVersionInfo(const char * pszRequest)1406 const char * CPL_STDCALL GDALVersionInfo( const char *pszRequest )
1407
1408 {
1409 static char szResult[128];
1410
1411 /* -------------------------------------------------------------------- */
1412 /* LICENSE is a special case. We try to find and read the */
1413 /* LICENSE.TXT file from the GDAL_DATA directory and return */
1414 /* it. We leak memory. Applications are discouraged from */
1415 /* deallocating the returned memory in this case since in the */
1416 /* future we may resolve the leak issue internally. */
1417 /* -------------------------------------------------------------------- */
1418 if( EQUAL(pszRequest,"LICENSE") )
1419 {
1420 const char *pszFilename = CPLFindFile( "etc", "LICENSE.TXT" );
1421 FILE *fp = NULL;
1422 int nLength;
1423 char *pszLICENSE;
1424
1425 if( pszFilename != NULL )
1426 fp = VSIFOpenL( pszFilename, "r" );
1427
1428 if( fp == NULL )
1429 {
1430 sprintf( szResult,
1431 "GDAL/OGR is released under the MIT/X license.\n"
1432 "The LICENSE.TXT distributed with GDAL/OGR should\n"
1433 "contain additional details.\n" );
1434 return szResult;
1435 }
1436
1437 VSIFSeekL( fp, 0, SEEK_END );
1438 nLength = VSIFTellL( fp ) + 1;
1439 VSIFSeekL( fp, SEEK_SET, 0 );
1440
1441 pszLICENSE = (char *) CPLCalloc(1,nLength);
1442 VSIFReadL( pszLICENSE, 1, nLength-1, fp );
1443
1444 VSIFCloseL( fp );
1445
1446 return pszLICENSE;
1447 }
1448
1449 // NOTE: There is a slight risk of a multithreaded race condition if
1450 // one thread is in the process of sprintf()ing into this buffer while
1451 // another is using it but that seems pretty low risk. All threads
1452 // want the same value in the buffer.
1453
1454 if( pszRequest == NULL || EQUAL(pszRequest,"VERSION_NUM") )
1455 sprintf( szResult, "%d", GDAL_VERSION_NUM );
1456 else if( EQUAL(pszRequest,"RELEASE_DATE") )
1457 sprintf( szResult, "%d", GDAL_RELEASE_DATE );
1458 else if( EQUAL(pszRequest,"RELEASE_NAME") )
1459 sprintf( szResult, "%s", GDAL_RELEASE_NAME );
1460 else // --version
1461 sprintf( szResult, "GDAL %s, released %d/%02d/%02d",
1462 GDAL_RELEASE_NAME,
1463 GDAL_RELEASE_DATE / 10000,
1464 (GDAL_RELEASE_DATE % 10000) / 100,
1465 GDAL_RELEASE_DATE % 100 );
1466
1467 return szResult;
1468 }
1469
1470 /************************************************************************/
1471 /* GDALCheckVersion() */
1472 /************************************************************************/
1473
1474 /** Return TRUE if GDAL library version at runtime matches nVersionMajor.nVersionMinor.
1475
1476 The purpose of this method is to ensure that calling code will run with the GDAL
1477 version it is compiled for. It is primarly intented for external plugins.
1478
1479 @param nVersionMajor Major version to be tested against
1480 @param nVersionMinor Minor version to be tested against
1481 @param pszCallingComponentName If not NULL, in case of version mismatch, the method
1482 will issue a failure mentionning the name of
1483 the calling component.
1484 */
GDALCheckVersion(int nVersionMajor,int nVersionMinor,const char * pszCallingComponentName)1485 int CPL_STDCALL GDALCheckVersion( int nVersionMajor, int nVersionMinor,
1486 const char* pszCallingComponentName)
1487 {
1488 if (nVersionMajor == GDAL_VERSION_MAJOR &&
1489 nVersionMinor == GDAL_VERSION_MINOR)
1490 return TRUE;
1491
1492 if (pszCallingComponentName)
1493 {
1494 CPLError( CE_Failure, CPLE_AppDefined,
1495 "%s was compiled against GDAL %d.%d but current library version is %d.%d\n",
1496 pszCallingComponentName, nVersionMajor, nVersionMinor,
1497 GDAL_VERSION_MAJOR, GDAL_VERSION_MINOR);
1498 }
1499 return FALSE;
1500 }
1501
1502 /************************************************************************/
1503 /* GDALDecToDMS() */
1504 /* */
1505 /* Translate a decimal degrees value to a DMS string with */
1506 /* hemisphere. */
1507 /************************************************************************/
1508
GDALDecToDMS(double dfAngle,const char * pszAxis,int nPrecision)1509 const char * CPL_STDCALL GDALDecToDMS( double dfAngle, const char * pszAxis,
1510 int nPrecision )
1511
1512 {
1513 return CPLDecToDMS( dfAngle, pszAxis, nPrecision );
1514 }
1515
1516 /************************************************************************/
1517 /* GDALPackedDMSToDec() */
1518 /************************************************************************/
1519
1520 /**
1521 * Convert a packed DMS value (DDDMMMSSS.SS) into decimal degrees.
1522 *
1523 * See CPLPackedDMSToDec().
1524 */
1525
GDALPackedDMSToDec(double dfPacked)1526 double CPL_STDCALL GDALPackedDMSToDec( double dfPacked )
1527
1528 {
1529 return CPLPackedDMSToDec( dfPacked );
1530 }
1531
1532 /************************************************************************/
1533 /* GDALDecToPackedDMS() */
1534 /************************************************************************/
1535
1536 /**
1537 * Convert decimal degrees into packed DMS value (DDDMMMSSS.SS).
1538 *
1539 * See CPLDecToPackedDMS().
1540 */
1541
GDALDecToPackedDMS(double dfDec)1542 double CPL_STDCALL GDALDecToPackedDMS( double dfDec )
1543
1544 {
1545 return CPLDecToPackedDMS( dfDec );
1546 }
1547
1548 /************************************************************************/
1549 /* GDALGCPsToGeoTransform() */
1550 /************************************************************************/
1551
1552 /**
1553 * Generate Geotransform from GCPs.
1554 *
1555 * Given a set of GCPs perform first order fit as a geotransform.
1556 *
1557 * Due to imprecision in the calculations the fit algorithm will often
1558 * return non-zero rotational coefficients even if given perfectly non-rotated
1559 * inputs. A special case has been implemented for corner corner coordinates
1560 * given in TL, TR, BR, BL order. So when using this to get a geotransform
1561 * from 4 corner coordinates, pass them in this order.
1562 *
1563 * @param nGCPCount the number of GCPs being passed in.
1564 * @param pasGCPs the list of GCP structures.
1565 * @param padfGeoTransform the six double array in which the affine
1566 * geotransformation will be returned.
1567 * @param bApproxOK If FALSE the function will fail if the geotransform is not
1568 * essentially an exact fit (within 0.25 pixel) for all GCPs.
1569 *
1570 * @return TRUE on success or FALSE if there aren't enough points to prepare a
1571 * geotransform, the pointers are ill-determined or if bApproxOK is FALSE
1572 * and the fit is poor.
1573 */
1574
1575 int CPL_STDCALL
GDALGCPsToGeoTransform(int nGCPCount,const GDAL_GCP * pasGCPs,double * padfGeoTransform,int bApproxOK)1576 GDALGCPsToGeoTransform( int nGCPCount, const GDAL_GCP *pasGCPs,
1577 double *padfGeoTransform, int bApproxOK )
1578
1579 {
1580 int i;
1581
1582 /* -------------------------------------------------------------------- */
1583 /* Recognise a few special cases. */
1584 /* -------------------------------------------------------------------- */
1585 if( nGCPCount < 2 )
1586 return FALSE;
1587
1588 if( nGCPCount == 2 )
1589 {
1590 if( pasGCPs[1].dfGCPPixel == pasGCPs[0].dfGCPPixel
1591 || pasGCPs[1].dfGCPLine == pasGCPs[0].dfGCPLine )
1592 return FALSE;
1593
1594 padfGeoTransform[1] = (pasGCPs[1].dfGCPX - pasGCPs[0].dfGCPX)
1595 / (pasGCPs[1].dfGCPPixel - pasGCPs[0].dfGCPPixel);
1596 padfGeoTransform[2] = 0.0;
1597
1598 padfGeoTransform[4] = 0.0;
1599 padfGeoTransform[5] = (pasGCPs[1].dfGCPY - pasGCPs[0].dfGCPY)
1600 / (pasGCPs[1].dfGCPLine - pasGCPs[0].dfGCPLine);
1601
1602 padfGeoTransform[0] = pasGCPs[0].dfGCPX
1603 - pasGCPs[0].dfGCPPixel * padfGeoTransform[1]
1604 - pasGCPs[0].dfGCPLine * padfGeoTransform[2];
1605
1606 padfGeoTransform[3] = pasGCPs[0].dfGCPY
1607 - pasGCPs[0].dfGCPPixel * padfGeoTransform[4]
1608 - pasGCPs[0].dfGCPLine * padfGeoTransform[5];
1609
1610 return TRUE;
1611 }
1612
1613 /* -------------------------------------------------------------------- */
1614 /* Special case of 4 corner coordinates of a non-rotated */
1615 /* image. The points must be in TL-TR-BR-BL order for now. */
1616 /* This case helps avoid some imprecision in the general */
1617 /* calcuations. */
1618 /* -------------------------------------------------------------------- */
1619 if( nGCPCount == 4
1620 && pasGCPs[0].dfGCPLine == pasGCPs[1].dfGCPLine
1621 && pasGCPs[2].dfGCPLine == pasGCPs[3].dfGCPLine
1622 && pasGCPs[0].dfGCPPixel == pasGCPs[3].dfGCPPixel
1623 && pasGCPs[1].dfGCPPixel == pasGCPs[2].dfGCPPixel
1624 && pasGCPs[0].dfGCPLine != pasGCPs[2].dfGCPLine
1625 && pasGCPs[0].dfGCPPixel != pasGCPs[1].dfGCPPixel
1626 && pasGCPs[0].dfGCPY == pasGCPs[1].dfGCPY
1627 && pasGCPs[2].dfGCPY == pasGCPs[3].dfGCPY
1628 && pasGCPs[0].dfGCPX == pasGCPs[3].dfGCPX
1629 && pasGCPs[1].dfGCPX == pasGCPs[2].dfGCPX
1630 && pasGCPs[0].dfGCPY != pasGCPs[2].dfGCPY
1631 && pasGCPs[0].dfGCPX != pasGCPs[1].dfGCPX )
1632 {
1633 padfGeoTransform[1] = (pasGCPs[1].dfGCPX - pasGCPs[0].dfGCPX)
1634 / (pasGCPs[1].dfGCPPixel - pasGCPs[0].dfGCPPixel);
1635 padfGeoTransform[2] = 0.0;
1636 padfGeoTransform[4] = 0.0;
1637 padfGeoTransform[5] = (pasGCPs[2].dfGCPY - pasGCPs[1].dfGCPY)
1638 / (pasGCPs[2].dfGCPLine - pasGCPs[1].dfGCPLine);
1639
1640 padfGeoTransform[0] =
1641 pasGCPs[0].dfGCPX - pasGCPs[0].dfGCPPixel * padfGeoTransform[1];
1642 padfGeoTransform[3] =
1643 pasGCPs[0].dfGCPY - pasGCPs[0].dfGCPLine * padfGeoTransform[5];
1644 return TRUE;
1645 }
1646
1647 /* -------------------------------------------------------------------- */
1648 /* In the general case, do a least squares error approximation by */
1649 /* solving the equation Sum[(A - B*x + C*y - Lon)^2] = minimum */
1650 /* -------------------------------------------------------------------- */
1651
1652 double sum_x = 0.0, sum_y = 0.0, sum_xy = 0.0, sum_xx = 0.0, sum_yy = 0.0;
1653 double sum_Lon = 0.0, sum_Lonx = 0.0, sum_Lony = 0.0;
1654 double sum_Lat = 0.0, sum_Latx = 0.0, sum_Laty = 0.0;
1655 double divisor;
1656
1657 for (i = 0; i < nGCPCount; ++i) {
1658 sum_x += pasGCPs[i].dfGCPPixel;
1659 sum_y += pasGCPs[i].dfGCPLine;
1660 sum_xy += pasGCPs[i].dfGCPPixel * pasGCPs[i].dfGCPLine;
1661 sum_xx += pasGCPs[i].dfGCPPixel * pasGCPs[i].dfGCPPixel;
1662 sum_yy += pasGCPs[i].dfGCPLine * pasGCPs[i].dfGCPLine;
1663 sum_Lon += pasGCPs[i].dfGCPX;
1664 sum_Lonx += pasGCPs[i].dfGCPX * pasGCPs[i].dfGCPPixel;
1665 sum_Lony += pasGCPs[i].dfGCPX * pasGCPs[i].dfGCPLine;
1666 sum_Lat += pasGCPs[i].dfGCPY;
1667 sum_Latx += pasGCPs[i].dfGCPY * pasGCPs[i].dfGCPPixel;
1668 sum_Laty += pasGCPs[i].dfGCPY * pasGCPs[i].dfGCPLine;
1669 }
1670
1671 divisor = nGCPCount * (sum_xx * sum_yy - sum_xy * sum_xy)
1672 + 2 * sum_x * sum_y * sum_xy - sum_y * sum_y * sum_xx
1673 - sum_x * sum_x * sum_yy;
1674
1675 /* -------------------------------------------------------------------- */
1676 /* If the divisor is zero, there is no valid solution. */
1677 /* -------------------------------------------------------------------- */
1678 if (divisor == 0.0)
1679 return FALSE;
1680
1681 /* -------------------------------------------------------------------- */
1682 /* Compute top/left origin. */
1683 /* -------------------------------------------------------------------- */
1684
1685 padfGeoTransform[0] = (sum_Lon * (sum_xx * sum_yy - sum_xy * sum_xy)
1686 + sum_Lonx * (sum_y * sum_xy - sum_x * sum_yy)
1687 + sum_Lony * (sum_x * sum_xy - sum_y * sum_xx))
1688 / divisor;
1689
1690 padfGeoTransform[3] = (sum_Lat * (sum_xx * sum_yy - sum_xy * sum_xy)
1691 + sum_Latx * (sum_y * sum_xy - sum_x * sum_yy)
1692 + sum_Laty * (sum_x * sum_xy - sum_y * sum_xx))
1693 / divisor;
1694
1695 /* -------------------------------------------------------------------- */
1696 /* Compute X related coefficients. */
1697 /* -------------------------------------------------------------------- */
1698 padfGeoTransform[1] = (sum_Lon * (sum_y * sum_xy - sum_x * sum_yy)
1699 + sum_Lonx * (nGCPCount * sum_yy - sum_y * sum_y)
1700 + sum_Lony * (sum_x * sum_y - sum_xy * nGCPCount))
1701 / divisor;
1702
1703 padfGeoTransform[2] = (sum_Lon * (sum_x * sum_xy - sum_y * sum_xx)
1704 + sum_Lonx * (sum_x * sum_y - nGCPCount * sum_xy)
1705 + sum_Lony * (nGCPCount * sum_xx - sum_x * sum_x))
1706 / divisor;
1707
1708 /* -------------------------------------------------------------------- */
1709 /* Compute Y related coefficients. */
1710 /* -------------------------------------------------------------------- */
1711 padfGeoTransform[4] = (sum_Lat * (sum_y * sum_xy - sum_x * sum_yy)
1712 + sum_Latx * (nGCPCount * sum_yy - sum_y * sum_y)
1713 + sum_Laty * (sum_x * sum_y - sum_xy * nGCPCount))
1714 / divisor;
1715
1716 padfGeoTransform[5] = (sum_Lat * (sum_x * sum_xy - sum_y * sum_xx)
1717 + sum_Latx * (sum_x * sum_y - nGCPCount * sum_xy)
1718 + sum_Laty * (nGCPCount * sum_xx - sum_x * sum_x))
1719 / divisor;
1720
1721 /* -------------------------------------------------------------------- */
1722 /* Now check if any of the input points fit this poorly. */
1723 /* -------------------------------------------------------------------- */
1724 if( !bApproxOK )
1725 {
1726 double dfPixelSize = ABS(padfGeoTransform[1])
1727 + ABS(padfGeoTransform[2])
1728 + ABS(padfGeoTransform[4])
1729 + ABS(padfGeoTransform[5]);
1730
1731 for( i = 0; i < nGCPCount; i++ )
1732 {
1733 double dfErrorX, dfErrorY;
1734
1735 dfErrorX =
1736 (pasGCPs[i].dfGCPPixel * padfGeoTransform[1]
1737 + pasGCPs[i].dfGCPLine * padfGeoTransform[2]
1738 + padfGeoTransform[0])
1739 - pasGCPs[i].dfGCPX;
1740 dfErrorY =
1741 (pasGCPs[i].dfGCPPixel * padfGeoTransform[4]
1742 + pasGCPs[i].dfGCPLine * padfGeoTransform[5]
1743 + padfGeoTransform[3])
1744 - pasGCPs[i].dfGCPY;
1745
1746 if( ABS(dfErrorX) > 0.25 * dfPixelSize
1747 || ABS(dfErrorY) > 0.25 * dfPixelSize )
1748 return FALSE;
1749 }
1750 }
1751
1752 return TRUE;
1753 }
1754
1755 /************************************************************************/
1756 /* GDALGeneralCmdLineProcessor() */
1757 /************************************************************************/
1758
1759 /**
1760 * General utility option processing.
1761 *
1762 * This function is intended to provide a variety of generic commandline
1763 * options for all GDAL commandline utilities. It takes care of the following
1764 * commandline options:
1765 *
1766 * --version: report version of GDAL in use.
1767 * --license: report GDAL license info.
1768 * --formats: report all format drivers configured.
1769 * --format [format]: report details of one format driver.
1770 * --optfile filename: expand an option file into the argument list.
1771 * --config key value: set system configuration option.
1772 * --debug [on/off/value]: set debug level.
1773 * --mempreload dir: preload directory contents into /vsimem
1774 * --help-general: report detailed help on general options.
1775 *
1776 * The argument array is replaced "in place" and should be freed with
1777 * CSLDestroy() when no longer needed. The typical usage looks something
1778 * like the following. Note that the formats should be registered so that
1779 * the --formats and --format options will work properly.
1780 *
1781 * int main( int argc, char ** argv )
1782 * {
1783 * GDALAllRegister();
1784 *
1785 * argc = GDALGeneralCmdLineProcessor( argc, &argv, 0 );
1786 * if( argc < 1 )
1787 * exit( -argc );
1788 *
1789 * @param nArgc number of values in the argument list.
1790 * @param Pointer to the argument list array (will be updated in place).
1791 *
1792 * @return updated nArgc argument count. Return of 0 requests terminate
1793 * without error, return of -1 requests exit with error code.
1794 */
1795
1796 int CPL_STDCALL
GDALGeneralCmdLineProcessor(int nArgc,char *** ppapszArgv,int nOptions)1797 GDALGeneralCmdLineProcessor( int nArgc, char ***ppapszArgv, int nOptions )
1798
1799 {
1800 char **papszReturn = NULL;
1801 int iArg;
1802 char **papszArgv = *ppapszArgv;
1803
1804 (void) nOptions;
1805
1806 /* -------------------------------------------------------------------- */
1807 /* Preserve the program name. */
1808 /* -------------------------------------------------------------------- */
1809 papszReturn = CSLAddString( papszReturn, papszArgv[0] );
1810
1811 /* ==================================================================== */
1812 /* Loop over all arguments. */
1813 /* ==================================================================== */
1814 for( iArg = 1; iArg < nArgc; iArg++ )
1815 {
1816 /* -------------------------------------------------------------------- */
1817 /* --version */
1818 /* -------------------------------------------------------------------- */
1819 if( EQUAL(papszArgv[iArg],"--version") )
1820 {
1821 printf( "%s\n", GDALVersionInfo( "--version" ) );
1822 CSLDestroy( papszReturn );
1823 return 0;
1824 }
1825
1826 /* -------------------------------------------------------------------- */
1827 /* --license */
1828 /* -------------------------------------------------------------------- */
1829 else if( EQUAL(papszArgv[iArg],"--license") )
1830 {
1831 printf( "%s\n", GDALVersionInfo( "LICENSE" ) );
1832 CSLDestroy( papszReturn );
1833 return 0;
1834 }
1835
1836 /* -------------------------------------------------------------------- */
1837 /* --config */
1838 /* -------------------------------------------------------------------- */
1839 else if( EQUAL(papszArgv[iArg],"--config") )
1840 {
1841 if( iArg + 2 >= nArgc )
1842 {
1843 CPLError( CE_Failure, CPLE_AppDefined,
1844 "--config option given without a key and value argument." );
1845 CSLDestroy( papszReturn );
1846 return -1;
1847 }
1848
1849 CPLSetConfigOption( papszArgv[iArg+1], papszArgv[iArg+2] );
1850
1851 iArg += 2;
1852 }
1853
1854 /* -------------------------------------------------------------------- */
1855 /* --mempreload */
1856 /* -------------------------------------------------------------------- */
1857 else if( EQUAL(papszArgv[iArg],"--mempreload") )
1858 {
1859 int i;
1860
1861 if( iArg + 1 >= nArgc )
1862 {
1863 CPLError( CE_Failure, CPLE_AppDefined,
1864 "--mempreload option given without directory path.");
1865 CSLDestroy( papszReturn );
1866 return -1;
1867 }
1868
1869 char **papszFiles = CPLReadDir( papszArgv[iArg+1] );
1870 if( CSLCount(papszFiles) == 0 )
1871 {
1872 CPLError( CE_Failure, CPLE_AppDefined,
1873 "--mempreload given invalid or empty directory.");
1874 CSLDestroy( papszReturn );
1875 return -1;
1876 }
1877
1878 for( i = 0; papszFiles[i] != NULL; i++ )
1879 {
1880 CPLString osOldPath, osNewPath;
1881
1882 if( EQUAL(papszFiles[i],".") || EQUAL(papszFiles[i],"..") )
1883 continue;
1884
1885 osOldPath = CPLFormFilename( papszArgv[iArg+1],
1886 papszFiles[i], NULL );
1887 osNewPath.Printf( "/vsimem/%s", papszFiles[i] );
1888
1889 CPLDebug( "VSI", "Preloading %s to %s.",
1890 osOldPath.c_str(), osNewPath.c_str() );
1891
1892 if( CPLCopyFile( osNewPath, osOldPath ) != 0 )
1893 return -1;
1894 }
1895
1896 CSLDestroy( papszFiles );
1897 iArg += 1;
1898 }
1899
1900 /* -------------------------------------------------------------------- */
1901 /* --debug */
1902 /* -------------------------------------------------------------------- */
1903 else if( EQUAL(papszArgv[iArg],"--debug") )
1904 {
1905 if( iArg + 1 >= nArgc )
1906 {
1907 CPLError( CE_Failure, CPLE_AppDefined,
1908 "--debug option given without debug level." );
1909 CSLDestroy( papszReturn );
1910 return -1;
1911 }
1912
1913 CPLSetConfigOption( "CPL_DEBUG", papszArgv[iArg+1] );
1914 iArg += 1;
1915 }
1916
1917 /* -------------------------------------------------------------------- */
1918 /* --optfile */
1919 /* */
1920 /* Annoyingly the options inserted by --optfile will *not* be */
1921 /* processed properly if they are general options. */
1922 /* -------------------------------------------------------------------- */
1923 else if( EQUAL(papszArgv[iArg],"--optfile") )
1924 {
1925 const char *pszLine;
1926 FILE *fpOptFile;
1927
1928 if( iArg + 1 >= nArgc )
1929 {
1930 CPLError( CE_Failure, CPLE_AppDefined,
1931 "--optfile option given without filename." );
1932 CSLDestroy( papszReturn );
1933 return -1;
1934 }
1935
1936 fpOptFile = VSIFOpen( papszArgv[iArg+1], "rb" );
1937
1938 if( fpOptFile == NULL )
1939 {
1940 CPLError( CE_Failure, CPLE_AppDefined,
1941 "Unable to open optfile '%s'.\n%s",
1942 papszArgv[iArg+1], VSIStrerror( errno ) );
1943 CSLDestroy( papszReturn );
1944 return -1;
1945 }
1946
1947 while( (pszLine = CPLReadLine( fpOptFile )) != NULL )
1948 {
1949 char **papszTokens;
1950 int i;
1951
1952 if( pszLine[0] == '#' || strlen(pszLine) == 0 )
1953 continue;
1954
1955 papszTokens = CSLTokenizeString( pszLine );
1956 for( i = 0; papszTokens != NULL && papszTokens[i] != NULL; i++)
1957 papszReturn = CSLAddString( papszReturn, papszTokens[i] );
1958 CSLDestroy( papszTokens );
1959 }
1960
1961 VSIFClose( fpOptFile );
1962
1963 iArg += 1;
1964 }
1965
1966 /* -------------------------------------------------------------------- */
1967 /* --formats */
1968 /* -------------------------------------------------------------------- */
1969 else if( EQUAL(papszArgv[iArg], "--formats") )
1970 {
1971 int iDr;
1972
1973 printf( "Supported Formats:\n" );
1974 for( iDr = 0; iDr < GDALGetDriverCount(); iDr++ )
1975 {
1976 GDALDriverH hDriver = GDALGetDriver(iDr);
1977 const char *pszRWFlag;
1978
1979 if( GDALGetMetadataItem( hDriver, GDAL_DCAP_CREATE, NULL ) )
1980 pszRWFlag = "rw+";
1981 else if( GDALGetMetadataItem( hDriver, GDAL_DCAP_CREATECOPY,
1982 NULL ) )
1983 pszRWFlag = "rw";
1984 else
1985 pszRWFlag = "ro";
1986
1987 printf( " %s (%s): %s\n",
1988 GDALGetDriverShortName( hDriver ),
1989 pszRWFlag,
1990 GDALGetDriverLongName( hDriver ) );
1991 }
1992
1993 CSLDestroy( papszReturn );
1994 return 0;
1995 }
1996
1997 /* -------------------------------------------------------------------- */
1998 /* --format */
1999 /* -------------------------------------------------------------------- */
2000 else if( EQUAL(papszArgv[iArg], "--format") )
2001 {
2002 GDALDriverH hDriver;
2003 char **papszMD;
2004
2005 CSLDestroy( papszReturn );
2006
2007 if( iArg + 1 >= nArgc )
2008 {
2009 CPLError( CE_Failure, CPLE_AppDefined,
2010 "--format option given without a format code." );
2011 return -1;
2012 }
2013
2014 hDriver = GDALGetDriverByName( papszArgv[iArg+1] );
2015 if( hDriver == NULL )
2016 {
2017 CPLError( CE_Failure, CPLE_AppDefined,
2018 "--format option given with format '%s', but that format not\n"
2019 "recognised. Use the --formats option to get a list of available formats,\n"
2020 "and use the short code (ie. GTiff or HFA) as the format identifier.\n",
2021 papszArgv[iArg+1] );
2022 return -1;
2023 }
2024
2025 printf( "Format Details:\n" );
2026 printf( " Short Name: %s\n", GDALGetDriverShortName( hDriver ) );
2027 printf( " Long Name: %s\n", GDALGetDriverLongName( hDriver ) );
2028
2029 papszMD = GDALGetMetadata( hDriver, NULL );
2030
2031 if( CSLFetchNameValue( papszMD, GDAL_DMD_EXTENSION ) )
2032 printf( " Extension: %s\n",
2033 CSLFetchNameValue( papszMD, GDAL_DMD_EXTENSION ) );
2034 if( CSLFetchNameValue( papszMD, GDAL_DMD_MIMETYPE ) )
2035 printf( " Mime Type: %s\n",
2036 CSLFetchNameValue( papszMD, GDAL_DMD_MIMETYPE ) );
2037 if( CSLFetchNameValue( papszMD, GDAL_DMD_HELPTOPIC ) )
2038 printf( " Help Topic: %s\n",
2039 CSLFetchNameValue( papszMD, GDAL_DMD_HELPTOPIC ) );
2040
2041 if( CSLFetchBoolean( papszMD, GDAL_DCAP_CREATE, FALSE ) )
2042 printf( " Supports: Create() - Create writeable dataset.\n" );
2043 if( CSLFetchBoolean( papszMD, GDAL_DCAP_CREATECOPY, FALSE ) )
2044 printf( " Supports: CreateCopy() - Create dataset by copying another.\n" );
2045 if( CSLFetchBoolean( papszMD, GDAL_DCAP_VIRTUALIO, FALSE ) )
2046 printf( " Supports: Virtual IO - eg. /vsimem/\n" );
2047 if( CSLFetchNameValue( papszMD, GDAL_DMD_CREATIONDATATYPES ) )
2048 printf( " Creation Datatypes: %s\n",
2049 CSLFetchNameValue( papszMD, GDAL_DMD_CREATIONDATATYPES ) );
2050 if( CSLFetchNameValue( papszMD, GDAL_DMD_CREATIONOPTIONLIST ) )
2051 {
2052 CPLXMLNode *psCOL =
2053 CPLParseXMLString(
2054 CSLFetchNameValue( papszMD,
2055 GDAL_DMD_CREATIONOPTIONLIST ) );
2056 char *pszFormattedXML =
2057 CPLSerializeXMLTree( psCOL );
2058
2059 CPLDestroyXMLNode( psCOL );
2060
2061 printf( "\n%s\n", pszFormattedXML );
2062 CPLFree( pszFormattedXML );
2063 }
2064 return 0;
2065 }
2066
2067 /* -------------------------------------------------------------------- */
2068 /* --help-general */
2069 /* -------------------------------------------------------------------- */
2070 else if( EQUAL(papszArgv[iArg],"--help-general") )
2071 {
2072 printf( "Generic GDAL utility command options:\n" );
2073 printf( " --version: report version of GDAL in use.\n" );
2074 printf( " --license: report GDAL license info.\n" );
2075 printf( " --formats: report all configured format drivers.\n" );
2076 printf( " --format [format]: details of one format.\n" );
2077 printf( " --optfile filename: expand an option file into the argument list.\n" );
2078 printf( " --config key value: set system configuration option.\n" );
2079 printf( " --debug [on/off/value]: set debug level.\n" );
2080 printf( " --help-general: report detailed help on general options.\n" );
2081 CSLDestroy( papszReturn );
2082 return 0;
2083 }
2084
2085 /* -------------------------------------------------------------------- */
2086 /* --locale */
2087 /* -------------------------------------------------------------------- */
2088 else if( EQUAL(papszArgv[iArg],"--locale") && iArg < nArgc-1 )
2089 {
2090 setlocale( LC_ALL, papszArgv[++iArg] );
2091 }
2092
2093 /* -------------------------------------------------------------------- */
2094 /* carry through unrecognised options. */
2095 /* -------------------------------------------------------------------- */
2096 else
2097 {
2098 papszReturn = CSLAddString( papszReturn, papszArgv[iArg] );
2099 }
2100 }
2101
2102 *ppapszArgv = papszReturn;
2103
2104 return CSLCount( *ppapszArgv );
2105 }
2106
2107
2108 /************************************************************************/
2109 /* _FetchDblFromMD() */
2110 /************************************************************************/
2111
_FetchDblFromMD(char ** papszMD,const char * pszKey,double * padfTarget,int nCount,double dfDefault)2112 static int _FetchDblFromMD( char **papszMD, const char *pszKey,
2113 double *padfTarget, int nCount, double dfDefault )
2114
2115 {
2116 char szFullKey[200];
2117
2118 sprintf( szFullKey, "%s", pszKey );
2119
2120 const char *pszValue = CSLFetchNameValue( papszMD, szFullKey );
2121 int i;
2122
2123 for( i = 0; i < nCount; i++ )
2124 padfTarget[i] = dfDefault;
2125
2126 if( pszValue == NULL )
2127 return FALSE;
2128
2129 if( nCount == 1 )
2130 {
2131 *padfTarget = CPLAtofM( pszValue );
2132 return TRUE;
2133 }
2134
2135 char **papszTokens = CSLTokenizeStringComplex( pszValue, " ,",
2136 FALSE, FALSE );
2137
2138 if( CSLCount( papszTokens ) != nCount )
2139 {
2140 CSLDestroy( papszTokens );
2141 return FALSE;
2142 }
2143
2144 for( i = 0; i < nCount; i++ )
2145 padfTarget[i] = CPLAtofM(papszTokens[i]);
2146
2147 CSLDestroy( papszTokens );
2148
2149 return TRUE;
2150 }
2151
2152 /************************************************************************/
2153 /* GDALExtractRPCInfo() */
2154 /* */
2155 /* Extract RPC info from metadata, and apply to an RPCInfo */
2156 /* structure. The inverse of this function is RPCInfoToMD() in */
2157 /* alg/gdal_rpc.cpp (should it be needed). */
2158 /************************************************************************/
2159
GDALExtractRPCInfo(char ** papszMD,GDALRPCInfo * psRPC)2160 int CPL_STDCALL GDALExtractRPCInfo( char **papszMD, GDALRPCInfo *psRPC )
2161
2162 {
2163 if( CSLFetchNameValue( papszMD, "LINE_NUM_COEFF" ) == NULL )
2164 return FALSE;
2165
2166 if( CSLFetchNameValue( papszMD, "LINE_NUM_COEFF" ) == NULL
2167 || CSLFetchNameValue( papszMD, "LINE_DEN_COEFF" ) == NULL
2168 || CSLFetchNameValue( papszMD, "SAMP_NUM_COEFF" ) == NULL
2169 || CSLFetchNameValue( papszMD, "SAMP_DEN_COEFF" ) == NULL )
2170 {
2171 CPLError( CE_Failure, CPLE_AppDefined,
2172 "Some required RPC metadata missing in GDALExtractRPCInfo()");
2173 return FALSE;
2174 }
2175
2176 _FetchDblFromMD( papszMD, "LINE_OFF", &(psRPC->dfLINE_OFF), 1, 0.0 );
2177 _FetchDblFromMD( papszMD, "LINE_SCALE", &(psRPC->dfLINE_SCALE), 1, 1.0 );
2178 _FetchDblFromMD( papszMD, "SAMP_OFF", &(psRPC->dfSAMP_OFF), 1, 0.0 );
2179 _FetchDblFromMD( papszMD, "SAMP_SCALE", &(psRPC->dfSAMP_SCALE), 1, 1.0 );
2180 _FetchDblFromMD( papszMD, "HEIGHT_OFF", &(psRPC->dfHEIGHT_OFF), 1, 0.0 );
2181 _FetchDblFromMD( papszMD, "HEIGHT_SCALE", &(psRPC->dfHEIGHT_SCALE),1, 1.0);
2182 _FetchDblFromMD( papszMD, "LAT_OFF", &(psRPC->dfLAT_OFF), 1, 0.0 );
2183 _FetchDblFromMD( papszMD, "LAT_SCALE", &(psRPC->dfLAT_SCALE), 1, 1.0 );
2184 _FetchDblFromMD( papszMD, "LONG_OFF", &(psRPC->dfLONG_OFF), 1, 0.0 );
2185 _FetchDblFromMD( papszMD, "LONG_SCALE", &(psRPC->dfLONG_SCALE), 1, 1.0 );
2186
2187 _FetchDblFromMD( papszMD, "LINE_NUM_COEFF", psRPC->adfLINE_NUM_COEFF,
2188 20, 0.0 );
2189 _FetchDblFromMD( papszMD, "LINE_DEN_COEFF", psRPC->adfLINE_DEN_COEFF,
2190 20, 0.0 );
2191 _FetchDblFromMD( papszMD, "SAMP_NUM_COEFF", psRPC->adfSAMP_NUM_COEFF,
2192 20, 0.0 );
2193 _FetchDblFromMD( papszMD, "SAMP_DEN_COEFF", psRPC->adfSAMP_DEN_COEFF,
2194 20, 0.0 );
2195
2196 _FetchDblFromMD( papszMD, "MIN_LONG", &(psRPC->dfMIN_LONG), 1, -180.0 );
2197 _FetchDblFromMD( papszMD, "MIN_LAT", &(psRPC->dfMIN_LAT), 1, -90.0 );
2198 _FetchDblFromMD( papszMD, "MAX_LONG", &(psRPC->dfMAX_LONG), 1, 180.0 );
2199 _FetchDblFromMD( papszMD, "MAX_LAT", &(psRPC->dfMAX_LAT), 1, 90.0 );
2200
2201 return TRUE;
2202 }
2203
2204 /************************************************************************/
2205 /* GDALFindAssociatedAuxFile() */
2206 /************************************************************************/
2207
GDALFindAssociatedAuxFile(const char * pszBasename,GDALAccess eAccess,GDALDataset * poDependentDS)2208 GDALDataset *GDALFindAssociatedAuxFile( const char *pszBasename,
2209 GDALAccess eAccess,
2210 GDALDataset *poDependentDS )
2211
2212 {
2213 const char *pszAuxSuffixLC = "aux";
2214 const char *pszAuxSuffixUC = "AUX";
2215
2216 if( EQUAL(CPLGetExtension(pszBasename), pszAuxSuffixLC) )
2217 return NULL;
2218
2219 /* -------------------------------------------------------------------- */
2220 /* Don't even try to look for an .aux file if we don't have a */
2221 /* path of any kind. */
2222 /* -------------------------------------------------------------------- */
2223 if( strlen(pszBasename) == 0 )
2224 return NULL;
2225
2226 /* -------------------------------------------------------------------- */
2227 /* We didn't find that, so try and find a corresponding aux */
2228 /* file. Check that we are the dependent file of the aux */
2229 /* file, or if we aren't verify that the dependent file does */
2230 /* not exist, likely mean it is us but some sort of renaming */
2231 /* has occured. */
2232 /* -------------------------------------------------------------------- */
2233 CPLString osJustFile = CPLGetFilename(pszBasename); // without dir
2234 CPLString osAuxFilename = CPLResetExtension(pszBasename, pszAuxSuffixLC);
2235 GDALDataset *poODS = NULL;
2236 GByte abyHeader[32];
2237 FILE *fp;
2238
2239 fp = VSIFOpenL( osAuxFilename, "rb" );
2240
2241
2242 if ( fp == NULL )
2243 {
2244 // Can't found file with lower case suffix. Try the upper case one.
2245 // no point in doing this on Win32 with case insensitive filenames.
2246 #ifndef WIN32
2247 osAuxFilename = CPLResetExtension(pszBasename, pszAuxSuffixUC);
2248 fp = VSIFOpenL( osAuxFilename, "rb" );
2249 #endif
2250 }
2251
2252 if( fp != NULL )
2253 {
2254 VSIFReadL( abyHeader, 1, 32, fp );
2255 if( EQUALN((char *) abyHeader,"EHFA_HEADER_TAG",15) )
2256 poODS = (GDALDataset *) GDALOpenShared( osAuxFilename, eAccess );
2257 VSIFCloseL( fp );
2258 }
2259
2260 /* -------------------------------------------------------------------- */
2261 /* Try replacing extension with .aux */
2262 /* -------------------------------------------------------------------- */
2263 if( poODS != NULL )
2264 {
2265 const char *pszDep
2266 = poODS->GetMetadataItem( "HFA_DEPENDENT_FILE", "HFA" );
2267 if( pszDep == NULL )
2268 {
2269 CPLDebug( "AUX",
2270 "Found %s but it has no dependent file, ignoring.",
2271 osAuxFilename.c_str() );
2272 GDALClose( poODS );
2273 poODS = NULL;
2274 }
2275 else if( !EQUAL(pszDep,osJustFile) )
2276 {
2277 VSIStatBufL sStatBuf;
2278
2279 if( VSIStatL( pszDep, &sStatBuf ) == 0 )
2280 {
2281 CPLDebug( "AUX", "%s is for file %s, not %s, ignoring.",
2282 osAuxFilename.c_str(),
2283 pszDep, osJustFile.c_str() );
2284 GDALClose( poODS );
2285 poODS = NULL;
2286 }
2287 else
2288 {
2289 CPLDebug( "AUX", "%s is for file %s, not %s, but since\n"
2290 "%s does not exist, we will use .aux file as our own.",
2291 osAuxFilename.c_str(),
2292 pszDep, osJustFile.c_str(),
2293 pszDep );
2294 }
2295 }
2296
2297 /* -------------------------------------------------------------------- */
2298 /* Confirm that the aux file matches the configuration of the */
2299 /* dependent dataset. */
2300 /* -------------------------------------------------------------------- */
2301 if( poODS != NULL && poDependentDS != NULL
2302 && (poODS->GetRasterCount() != poDependentDS->GetRasterCount()
2303 || poODS->GetRasterXSize() != poDependentDS->GetRasterXSize()
2304 || poODS->GetRasterYSize() != poDependentDS->GetRasterYSize()) )
2305 {
2306 CPLDebug( "AUX",
2307 "Ignoring aux file %s as it's raster configuration\n"
2308 "(%dP x %dL x %dB) does not match master file (%dP x %dL x %dB)",
2309 osAuxFilename.c_str(),
2310 poODS->GetRasterXSize(),
2311 poODS->GetRasterYSize(),
2312 poODS->GetRasterCount(),
2313 poDependentDS->GetRasterXSize(),
2314 poDependentDS->GetRasterYSize(),
2315 poDependentDS->GetRasterCount() );
2316
2317 GDALClose( poODS );
2318 poODS = NULL;
2319 }
2320 }
2321
2322 /* -------------------------------------------------------------------- */
2323 /* Try appending .aux to the end of the filename. */
2324 /* -------------------------------------------------------------------- */
2325 if( poODS == NULL )
2326 {
2327 osAuxFilename = pszBasename;
2328 osAuxFilename += ".";
2329 osAuxFilename += pszAuxSuffixLC;
2330 fp = VSIFOpenL( osAuxFilename, "rb" );
2331 #ifndef WIN32
2332 if ( fp == NULL )
2333 {
2334 // Can't found file with lower case suffix. Try the upper case one.
2335 osAuxFilename = pszBasename;
2336 osAuxFilename += ".";
2337 osAuxFilename += pszAuxSuffixUC;
2338 fp = VSIFOpenL( osAuxFilename, "rb" );
2339 }
2340 #endif
2341
2342 if( fp != NULL )
2343 {
2344 VSIFReadL( abyHeader, 1, 32, fp );
2345 if( EQUALN((char *) abyHeader,"EHFA_HEADER_TAG",15) )
2346 poODS = (GDALDataset *) GDALOpenShared( osAuxFilename, eAccess );
2347 VSIFCloseL( fp );
2348 }
2349
2350 if( poODS != NULL )
2351 {
2352 const char *pszDep
2353 = poODS->GetMetadataItem( "HFA_DEPENDENT_FILE", "HFA" );
2354 if( pszDep == NULL )
2355 {
2356 CPLDebug( "AUX",
2357 "Found %s but it has no dependent file, ignoring.",
2358 osAuxFilename.c_str() );
2359 GDALClose( poODS );
2360 poODS = NULL;
2361 }
2362 else if( !EQUAL(pszDep,osJustFile) )
2363 {
2364 VSIStatBufL sStatBuf;
2365
2366 if( VSIStatL( pszDep, &sStatBuf ) == 0 )
2367 {
2368 CPLDebug( "AUX", "%s is for file %s, not %s, ignoring.",
2369 osAuxFilename.c_str(),
2370 pszDep, osJustFile.c_str() );
2371 GDALClose( poODS );
2372 poODS = NULL;
2373 }
2374 else
2375 {
2376 CPLDebug( "AUX", "%s is for file %s, not %s, but since\n"
2377 "%s does not exist, we will use .aux file as our own.",
2378 osAuxFilename.c_str(),
2379 pszDep, osJustFile.c_str(),
2380 pszDep );
2381 }
2382 }
2383 }
2384 }
2385
2386 /* -------------------------------------------------------------------- */
2387 /* Confirm that the aux file matches the configuration of the */
2388 /* dependent dataset. */
2389 /* -------------------------------------------------------------------- */
2390 if( poODS != NULL && poDependentDS != NULL
2391 && (poODS->GetRasterCount() != poDependentDS->GetRasterCount()
2392 || poODS->GetRasterXSize() != poDependentDS->GetRasterXSize()
2393 || poODS->GetRasterYSize() != poDependentDS->GetRasterYSize()) )
2394 {
2395 CPLDebug( "AUX",
2396 "Ignoring aux file %s as it's raster configuration\n"
2397 "(%dP x %dL x %dB) does not match master file (%dP x %dL x %dB)",
2398 osAuxFilename.c_str(),
2399 poODS->GetRasterXSize(),
2400 poODS->GetRasterYSize(),
2401 poODS->GetRasterCount(),
2402 poDependentDS->GetRasterXSize(),
2403 poDependentDS->GetRasterYSize(),
2404 poDependentDS->GetRasterCount() );
2405
2406 GDALClose( poODS );
2407 poODS = NULL;
2408 }
2409
2410 return poODS;
2411 }
2412
2413 /************************************************************************/
2414 /* -------------------------------------------------------------------- */
2415 /* The following stubs are present to ensure that older GDAL */
2416 /* bridges don't fail with newer libraries. */
2417 /* -------------------------------------------------------------------- */
2418 /************************************************************************/
2419
2420 CPL_C_START
2421
GDALCreateProjDef(const char *)2422 void * CPL_STDCALL GDALCreateProjDef( const char * )
2423 {
2424 CPLDebug( "GDAL", "GDALCreateProjDef no longer supported." );
2425 return NULL;
2426 }
2427
GDALReprojectToLongLat(void *,double *,double *)2428 CPLErr CPL_STDCALL GDALReprojectToLongLat( void *, double *, double * )
2429 {
2430 CPLDebug( "GDAL", "GDALReprojectToLatLong no longer supported." );
2431 return CE_Failure;
2432 }
2433
GDALReprojectFromLongLat(void *,double *,double *)2434 CPLErr CPL_STDCALL GDALReprojectFromLongLat( void *, double *, double * )
2435 {
2436 CPLDebug( "GDAL", "GDALReprojectFromLatLong no longer supported." );
2437 return CE_Failure;
2438 }
2439
GDALDestroyProjDef(void *)2440 void CPL_STDCALL GDALDestroyProjDef( void * )
2441
2442 {
2443 CPLDebug( "GDAL", "GDALDestroyProjDef no longer supported." );
2444 }
2445
2446 CPL_C_END
2447