1 : /******************************************************************************
2 : * $Id: hfadataset.cpp 18157 2009-12-02 16:53:53Z warmerdam $
3 : *
4 : * Name: hfadataset.cpp
5 : * Project: Erdas Imagine Driver
6 : * Purpose: Main driver for Erdas Imagine format.
7 : * Author: Frank Warmerdam, warmerdam@pobox.com
8 : *
9 : ******************************************************************************
10 : * Copyright (c) 1999, Frank Warmerdam
11 : *
12 : * Permission is hereby granted, free of charge, to any person obtaining a
13 : * copy of this software and associated documentation files (the "Software"),
14 : * to deal in the Software without restriction, including without limitation
15 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
16 : * and/or sell copies of the Software, and to permit persons to whom the
17 : * Software is furnished to do so, subject to the following conditions:
18 : *
19 : * The above copyright notice and this permission notice shall be included
20 : * in all copies or substantial portions of the Software.
21 : *
22 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
23 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
24 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
25 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
26 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
27 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
28 : * DEALINGS IN THE SOFTWARE.
29 : ****************************************************************************/
30 :
31 : #include "gdal_pam.h"
32 : #include "gdal_rat.h"
33 : #include "hfa_p.h"
34 : #include "ogr_spatialref.h"
35 :
36 : CPL_CVSID("$Id: hfadataset.cpp 18157 2009-12-02 16:53:53Z warmerdam $");
37 :
38 : CPL_C_START
39 : void GDALRegister_HFA(void);
40 : CPL_C_END
41 :
42 : #ifndef PI
43 : # define PI 3.14159265358979323846
44 : #endif
45 :
46 : #ifndef R2D
47 : # define R2D (180/PI)
48 : #endif
49 : #ifndef D2R
50 : # define D2R (PI/180)
51 : #endif
52 :
53 : OGRBoolean WritePeStringIfNeeded(OGRSpatialReference* poSRS, HFAHandle hHFA);
54 : void ClearSR(HFAHandle hHFA);
55 :
56 : static const char *apszDatumMap[] = {
57 : /* Imagine name, WKT name */
58 : "NAD27", "North_American_Datum_1927",
59 : "NAD83", "North_American_Datum_1983",
60 : "WGS 84", "WGS_1984",
61 : "WGS 1972", "WGS_1972",
62 : NULL, NULL
63 : };
64 :
65 : static const char *apszUnitMap[] = {
66 : "meters", "1.0",
67 : "meter", "1.0",
68 : "m", "1.0",
69 : "centimeters", "0.01",
70 : "centimeter", "0.01",
71 : "cm", "0.01",
72 : "millimeters", "0.001",
73 : "millimeter", "0.001",
74 : "mm", "0.001",
75 : "kilometers", "1000.0",
76 : "kilometer", "1000.0",
77 : "km", "1000.0",
78 : "us_survey_feet", "0.3048006096012192",
79 : "us_survey_foot", "0.3048006096012192",
80 : "feet", "0.3048006096012192",
81 : "foot", "0.3048006096012192",
82 : "ft", "0.3048006096012192",
83 : "international_feet", "0.3048",
84 : "international_foot", "0.3048",
85 : "inches", "0.0254000508001",
86 : "inch", "0.0254000508001",
87 : "in", "0.0254000508001",
88 : "yards", "0.9144",
89 : "yard", "0.9144",
90 : "yd", "0.9144",
91 : "miles", "1304.544",
92 : "mile", "1304.544",
93 : "mi", "1304.544",
94 : "modified_american_feet", "0.3048122530",
95 : "modified_american_foot", "0.3048122530",
96 : "clarke_feet", "0.3047972651",
97 : "clarke_foot", "0.3047972651",
98 : "indian_feet", "0.3047995142",
99 : "indian_foot", "0.3047995142",
100 : NULL, NULL
101 : };
102 :
103 :
104 : /* ==================================================================== */
105 : /* Table relating USGS and ESRI state plane zones. */
106 : /* ==================================================================== */
107 : static const int anUsgsEsriZones[] =
108 : {
109 : 101, 3101,
110 : 102, 3126,
111 : 201, 3151,
112 : 202, 3176,
113 : 203, 3201,
114 : 301, 3226,
115 : 302, 3251,
116 : 401, 3276,
117 : 402, 3301,
118 : 403, 3326,
119 : 404, 3351,
120 : 405, 3376,
121 : 406, 3401,
122 : 407, 3426,
123 : 501, 3451,
124 : 502, 3476,
125 : 503, 3501,
126 : 600, 3526,
127 : 700, 3551,
128 : 901, 3601,
129 : 902, 3626,
130 : 903, 3576,
131 : 1001, 3651,
132 : 1002, 3676,
133 : 1101, 3701,
134 : 1102, 3726,
135 : 1103, 3751,
136 : 1201, 3776,
137 : 1202, 3801,
138 : 1301, 3826,
139 : 1302, 3851,
140 : 1401, 3876,
141 : 1402, 3901,
142 : 1501, 3926,
143 : 1502, 3951,
144 : 1601, 3976,
145 : 1602, 4001,
146 : 1701, 4026,
147 : 1702, 4051,
148 : 1703, 6426,
149 : 1801, 4076,
150 : 1802, 4101,
151 : 1900, 4126,
152 : 2001, 4151,
153 : 2002, 4176,
154 : 2101, 4201,
155 : 2102, 4226,
156 : 2103, 4251,
157 : 2111, 6351,
158 : 2112, 6376,
159 : 2113, 6401,
160 : 2201, 4276,
161 : 2202, 4301,
162 : 2203, 4326,
163 : 2301, 4351,
164 : 2302, 4376,
165 : 2401, 4401,
166 : 2402, 4426,
167 : 2403, 4451,
168 : 2500, 0,
169 : 2501, 4476,
170 : 2502, 4501,
171 : 2503, 4526,
172 : 2600, 0,
173 : 2601, 4551,
174 : 2602, 4576,
175 : 2701, 4601,
176 : 2702, 4626,
177 : 2703, 4651,
178 : 2800, 4676,
179 : 2900, 4701,
180 : 3001, 4726,
181 : 3002, 4751,
182 : 3003, 4776,
183 : 3101, 4801,
184 : 3102, 4826,
185 : 3103, 4851,
186 : 3104, 4876,
187 : 3200, 4901,
188 : 3301, 4926,
189 : 3302, 4951,
190 : 3401, 4976,
191 : 3402, 5001,
192 : 3501, 5026,
193 : 3502, 5051,
194 : 3601, 5076,
195 : 3602, 5101,
196 : 3701, 5126,
197 : 3702, 5151,
198 : 3800, 5176,
199 : 3900, 0,
200 : 3901, 5201,
201 : 3902, 5226,
202 : 4001, 5251,
203 : 4002, 5276,
204 : 4100, 5301,
205 : 4201, 5326,
206 : 4202, 5351,
207 : 4203, 5376,
208 : 4204, 5401,
209 : 4205, 5426,
210 : 4301, 5451,
211 : 4302, 5476,
212 : 4303, 5501,
213 : 4400, 5526,
214 : 4501, 5551,
215 : 4502, 5576,
216 : 4601, 5601,
217 : 4602, 5626,
218 : 4701, 5651,
219 : 4702, 5676,
220 : 4801, 5701,
221 : 4802, 5726,
222 : 4803, 5751,
223 : 4901, 5776,
224 : 4902, 5801,
225 : 4903, 5826,
226 : 4904, 5851,
227 : 5001, 6101,
228 : 5002, 6126,
229 : 5003, 6151,
230 : 5004, 6176,
231 : 5005, 6201,
232 : 5006, 6226,
233 : 5007, 6251,
234 : 5008, 6276,
235 : 5009, 6301,
236 : 5010, 6326,
237 : 5101, 5876,
238 : 5102, 5901,
239 : 5103, 5926,
240 : 5104, 5951,
241 : 5105, 5976,
242 : 5201, 6001,
243 : 5200, 6026,
244 : 5200, 6076,
245 : 5201, 6051,
246 : 5202, 6051,
247 : 5300, 0,
248 : 5400, 0
249 : };
250 :
251 :
252 : /************************************************************************/
253 : /* ==================================================================== */
254 : /* HFADataset */
255 : /* ==================================================================== */
256 : /************************************************************************/
257 :
258 : class HFARasterBand;
259 :
260 : class CPL_DLL HFADataset : public GDALPamDataset
261 : {
262 : friend class HFARasterBand;
263 :
264 : HFAHandle hHFA;
265 :
266 : int bMetadataDirty;
267 :
268 : int bGeoDirty;
269 : double adfGeoTransform[6];
270 : char *pszProjection;
271 :
272 : int bIgnoreUTM;
273 :
274 : CPLErr ReadProjection();
275 : CPLErr WriteProjection();
276 : int bForceToPEString;
277 :
278 : int nGCPCount;
279 : GDAL_GCP asGCPList[36];
280 :
281 : void UseXFormStack( int nStepCount,
282 : Efga_Polynomial *pasPolyListForward,
283 : Efga_Polynomial *pasPolyListReverse );
284 :
285 : protected:
286 : virtual CPLErr IRasterIO( GDALRWFlag, int, int, int, int,
287 : void *, int, int, GDALDataType,
288 : int, int *, int, int, int );
289 :
290 : public:
291 : HFADataset();
292 : ~HFADataset();
293 :
294 : static int Identify( GDALOpenInfo * );
295 : static GDALDataset *Open( GDALOpenInfo * );
296 : static GDALDataset *Create( const char * pszFilename,
297 : int nXSize, int nYSize, int nBands,
298 : GDALDataType eType, char ** papszParmList );
299 : static GDALDataset *CreateCopy( const char * pszFilename,
300 : GDALDataset *poSrcDS,
301 : int bStrict, char ** papszOptions,
302 : GDALProgressFunc pfnProgress,
303 : void * pProgressData );
304 : static CPLErr Delete( const char *pszFilename );
305 :
306 : virtual char **GetFileList(void);
307 :
308 : virtual const char *GetProjectionRef(void);
309 : virtual CPLErr SetProjection( const char * );
310 :
311 : virtual CPLErr GetGeoTransform( double * );
312 : virtual CPLErr SetGeoTransform( double * );
313 :
314 : virtual int GetGCPCount();
315 : virtual const char *GetGCPProjection();
316 : virtual const GDAL_GCP *GetGCPs();
317 :
318 : virtual CPLErr SetMetadata( char **, const char * = "" );
319 : virtual CPLErr SetMetadataItem( const char *, const char *, const char * = "" );
320 :
321 : virtual void FlushCache( void );
322 : virtual CPLErr IBuildOverviews( const char *pszResampling,
323 : int nOverviews, int *panOverviewList,
324 : int nListBands, int *panBandList,
325 : GDALProgressFunc pfnProgress,
326 : void * pProgressData );
327 : };
328 :
329 : /************************************************************************/
330 : /* ==================================================================== */
331 : /* HFARasterBand */
332 : /* ==================================================================== */
333 : /************************************************************************/
334 :
335 : class HFARasterBand : public GDALPamRasterBand
336 : {
337 : friend class HFADataset;
338 :
339 : GDALColorTable *poCT;
340 :
341 : int nHFADataType;
342 :
343 : int nOverviews;
344 : int nThisOverview;
345 : HFARasterBand **papoOverviewBands;
346 :
347 : CPLErr CleanOverviews();
348 :
349 : HFAHandle hHFA;
350 :
351 : int bMetadataDirty;
352 :
353 : GDALRasterAttributeTable *poDefaultRAT;
354 :
355 : void ReadAuxMetadata();
356 : void ReadHistogramMetadata();
357 : void EstablishOverviews();
358 :
359 : GDALRasterAttributeTable* ReadNamedRAT( const char *pszName );
360 :
361 : public:
362 :
363 : HFARasterBand( HFADataset *, int, int );
364 : virtual ~HFARasterBand();
365 :
366 : virtual CPLErr IReadBlock( int, int, void * );
367 : virtual CPLErr IWriteBlock( int, int, void * );
368 :
369 : virtual const char *GetDescription() const;
370 : virtual void SetDescription( const char * );
371 :
372 : virtual GDALColorInterp GetColorInterpretation();
373 : virtual GDALColorTable *GetColorTable();
374 : virtual CPLErr SetColorTable( GDALColorTable * );
375 : virtual int GetOverviewCount();
376 : virtual GDALRasterBand *GetOverview( int );
377 :
378 : virtual double GetMinimum( int *pbSuccess = NULL );
379 : virtual double GetMaximum(int *pbSuccess = NULL );
380 : virtual double GetNoDataValue( int *pbSuccess = NULL );
381 : virtual CPLErr SetNoDataValue( double dfValue );
382 :
383 : virtual CPLErr SetMetadata( char **, const char * = "" );
384 : virtual CPLErr SetMetadataItem( const char *, const char *, const char * = "" );
385 : virtual CPLErr BuildOverviews( const char *, int, int *,
386 : GDALProgressFunc, void * );
387 :
388 : virtual CPLErr GetDefaultHistogram( double *pdfMin, double *pdfMax,
389 : int *pnBuckets, int ** ppanHistogram,
390 : int bForce,
391 : GDALProgressFunc, void *pProgressData);
392 :
393 : virtual const GDALRasterAttributeTable *GetDefaultRAT();
394 : virtual CPLErr SetDefaultRAT( const GDALRasterAttributeTable * );
395 : };
396 :
397 : /************************************************************************/
398 : /* HFARasterBand() */
399 : /************************************************************************/
400 :
401 419 : HFARasterBand::HFARasterBand( HFADataset *poDS, int nBand, int iOverview )
402 :
403 : {
404 : int nCompression;
405 :
406 419 : if( iOverview == -1 )
407 403 : this->poDS = poDS;
408 : else
409 16 : this->poDS = NULL;
410 :
411 419 : this->hHFA = poDS->hHFA;
412 419 : this->nBand = nBand;
413 419 : this->poCT = NULL;
414 419 : this->nThisOverview = iOverview;
415 419 : this->papoOverviewBands = NULL;
416 419 : this->bMetadataDirty = FALSE;
417 419 : this->poDefaultRAT = NULL;
418 419 : this->nOverviews = -1;
419 :
420 : HFAGetBandInfo( hHFA, nBand, &nHFADataType,
421 419 : &nBlockXSize, &nBlockYSize, &nCompression );
422 :
423 419 : if( nCompression != 0 )
424 : GDALMajorObject::SetMetadataItem( "COMPRESSION", "RLE",
425 77 : "IMAGE_STRUCTURE" );
426 :
427 419 : switch( nHFADataType )
428 : {
429 : case EPT_u1:
430 : case EPT_u2:
431 : case EPT_u4:
432 : case EPT_u8:
433 : case EPT_s8:
434 200 : eDataType = GDT_Byte;
435 200 : break;
436 :
437 : case EPT_u16:
438 44 : eDataType = GDT_UInt16;
439 44 : break;
440 :
441 : case EPT_s16:
442 21 : eDataType = GDT_Int16;
443 21 : break;
444 :
445 : case EPT_u32:
446 21 : eDataType = GDT_UInt32;
447 21 : break;
448 :
449 : case EPT_s32:
450 26 : eDataType = GDT_Int32;
451 26 : break;
452 :
453 : case EPT_f32:
454 25 : eDataType = GDT_Float32;
455 25 : break;
456 :
457 : case EPT_f64:
458 42 : eDataType = GDT_Float64;
459 42 : break;
460 :
461 : case EPT_c64:
462 20 : eDataType = GDT_CFloat32;
463 20 : break;
464 :
465 : case EPT_c128:
466 20 : eDataType = GDT_CFloat64;
467 20 : break;
468 :
469 : default:
470 0 : eDataType = GDT_Byte;
471 : /* notdef: this should really report an error, but this isn't
472 : so easy from within constructors. */
473 : CPLDebug( "GDAL", "Unsupported pixel type in HFARasterBand: %d.",
474 0 : (int) nHFADataType );
475 : break;
476 : }
477 :
478 419 : if( HFAGetDataTypeBits( nHFADataType ) < 8 )
479 : {
480 : GDALMajorObject::SetMetadataItem(
481 : "NBITS",
482 6 : CPLString().Printf( "%d", HFAGetDataTypeBits( nHFADataType ) ), "IMAGE_STRUCTURE" );
483 : }
484 :
485 419 : if( nHFADataType == EPT_s8 )
486 : {
487 : GDALMajorObject::SetMetadataItem( "PIXELTYPE", "SIGNEDBYTE",
488 2 : "IMAGE_STRUCTURE" );
489 : }
490 :
491 : /* -------------------------------------------------------------------- */
492 : /* If this is an overview, we need to fetch the actual size, */
493 : /* and block size. */
494 : /* -------------------------------------------------------------------- */
495 419 : if( iOverview > -1 )
496 : {
497 : int nHFADataTypeO;
498 :
499 16 : nOverviews = 0;
500 : HFAGetOverviewInfo( hHFA, nBand, iOverview,
501 : &nRasterXSize, &nRasterYSize,
502 16 : &nBlockXSize, &nBlockYSize, &nHFADataTypeO );
503 :
504 : /* -------------------------------------------------------------------- */
505 : /* If we are an 8bit overview of a 1bit layer, we need to mark */
506 : /* ourselves as being "resample: average_bit2grayscale". */
507 : /* -------------------------------------------------------------------- */
508 16 : if( nHFADataType == EPT_u1 && nHFADataTypeO == EPT_u8 )
509 : {
510 : GDALMajorObject::SetMetadataItem( "RESAMPLING",
511 0 : "AVERAGE_BIT2GRAYSCALE" );
512 0 : GDALMajorObject::SetMetadataItem( "NBITS", "8" );
513 : }
514 : }
515 :
516 : /* -------------------------------------------------------------------- */
517 : /* Collect color table if present. */
518 : /* -------------------------------------------------------------------- */
519 : double *padfRed, *padfGreen, *padfBlue, *padfAlpha, *padfBins;
520 : int nColors;
521 :
522 419 : if( iOverview == -1
523 : && HFAGetPCT( hHFA, nBand, &nColors,
524 : &padfRed, &padfGreen, &padfBlue, &padfAlpha,
525 : &padfBins ) == CE_None
526 : && nColors > 0 )
527 : {
528 4 : poCT = new GDALColorTable();
529 485 : for( int iColor = 0; iColor < nColors; iColor++ )
530 : {
531 : GDALColorEntry sEntry;
532 :
533 : // The following mapping assigns "equal sized" section of
534 : // the [0...1] range to each possible output value and avoid
535 : // rounding issues for the "normal" values generated using n/255.
536 : // See bug #1732 for some discussion.
537 481 : sEntry.c1 = MIN(255,(short) (padfRed[iColor] * 256));
538 481 : sEntry.c2 = MIN(255,(short) (padfGreen[iColor] * 256));
539 481 : sEntry.c3 = MIN(255,(short) (padfBlue[iColor] * 256));
540 481 : sEntry.c4 = MIN(255,(short) (padfAlpha[iColor] * 256));
541 :
542 481 : if( padfBins != NULL )
543 225 : poCT->SetColorEntry( (int) padfBins[iColor], &sEntry );
544 : else
545 256 : poCT->SetColorEntry( iColor, &sEntry );
546 : }
547 : }
548 :
549 419 : poDefaultRAT = ReadNamedRAT( "Descriptor_Table" );
550 419 : }
551 :
552 : /************************************************************************/
553 : /* ~HFARasterBand() */
554 : /************************************************************************/
555 :
556 838 : HFARasterBand::~HFARasterBand()
557 :
558 : {
559 419 : FlushCache();
560 :
561 434 : for( int iOvIndex = 0; iOvIndex < nOverviews; iOvIndex++ )
562 : {
563 15 : delete papoOverviewBands[iOvIndex];
564 : }
565 419 : CPLFree( papoOverviewBands );
566 :
567 419 : if( poCT != NULL )
568 3 : delete poCT;
569 :
570 419 : if( poDefaultRAT )
571 40 : delete poDefaultRAT;
572 838 : }
573 :
574 : /************************************************************************/
575 : /* ReadAuxMetadata() */
576 : /************************************************************************/
577 :
578 403 : void HFARasterBand::ReadAuxMetadata()
579 :
580 : {
581 : int i;
582 403 : HFABand *poBand = hHFA->papoBand[nBand-1];
583 :
584 : // only load metadata for full resolution layer.
585 403 : if( nThisOverview != -1 )
586 0 : return;
587 :
588 403 : const char ** pszAuxMetaData = GetHFAAuxMetaDataList();
589 5642 : for( i = 0; pszAuxMetaData[i] != NULL; i += 4 )
590 : {
591 : HFAEntry *poEntry;
592 :
593 5239 : if( strlen(pszAuxMetaData[i]) > 0 )
594 4836 : poEntry = poBand->poNode->GetNamedChild( pszAuxMetaData[i] );
595 : else
596 403 : poEntry = poBand->poNode;
597 :
598 5239 : const char *pszFieldName = pszAuxMetaData[i+1] + 1;
599 5239 : CPLErr eErr = CE_None;
600 :
601 5239 : if( poEntry == NULL )
602 4464 : continue;
603 :
604 775 : switch( pszAuxMetaData[i+1][0] )
605 : {
606 : case 'd':
607 : {
608 : int nCount, iValue;
609 : double dfValue;
610 279 : CPLString osValueList;
611 :
612 279 : nCount = poEntry->GetFieldCount( pszFieldName, &eErr );
613 279 : for( iValue = 0; eErr == CE_None && iValue < nCount; iValue++ )
614 : {
615 272 : CPLString osSubFieldName;
616 272 : osSubFieldName.Printf( "%s[%d]", pszFieldName, iValue );
617 272 : dfValue = poEntry->GetDoubleField( osSubFieldName, &eErr );
618 272 : if( eErr != CE_None )
619 : break;
620 :
621 : char szValueAsString[100];
622 254 : sprintf( szValueAsString, "%.14g", dfValue );
623 :
624 254 : if( iValue > 0 )
625 2 : osValueList += ",";
626 254 : osValueList += szValueAsString;
627 : }
628 279 : if( eErr == CE_None )
629 261 : SetMetadataItem( pszAuxMetaData[i+2], osValueList );
630 : }
631 279 : break;
632 : case 'i':
633 : case 'l':
634 : {
635 : int nValue, nCount, iValue;
636 93 : CPLString osValueList;
637 :
638 93 : nCount = poEntry->GetFieldCount( pszFieldName, &eErr );
639 93 : for( iValue = 0; eErr == CE_None && iValue < nCount; iValue++ )
640 : {
641 93 : CPLString osSubFieldName;
642 93 : osSubFieldName.Printf( "%s[%d]", pszFieldName, iValue );
643 93 : nValue = poEntry->GetIntField( osSubFieldName, &eErr );
644 93 : if( eErr != CE_None )
645 : break;
646 :
647 : char szValueAsString[100];
648 84 : sprintf( szValueAsString, "%d", nValue );
649 :
650 84 : if( iValue > 0 )
651 0 : osValueList += ",";
652 84 : osValueList += szValueAsString;
653 : }
654 93 : if( eErr == CE_None )
655 84 : SetMetadataItem( pszAuxMetaData[i+2], osValueList );
656 : }
657 93 : break;
658 : case 's':
659 : case 'e':
660 : {
661 : const char *pszValue;
662 403 : pszValue = poEntry->GetStringField( pszFieldName, &eErr );
663 403 : if( eErr == CE_None )
664 403 : SetMetadataItem( pszAuxMetaData[i+2], pszValue );
665 : }
666 : break;
667 : default:
668 : CPLAssert( FALSE );
669 : }
670 : }
671 : }
672 :
673 : /************************************************************************/
674 : /* ReadHistogramMetadata() */
675 : /************************************************************************/
676 :
677 403 : void HFARasterBand::ReadHistogramMetadata()
678 :
679 : {
680 : int i;
681 403 : HFABand *poBand = hHFA->papoBand[nBand-1];
682 :
683 : // only load metadata for full resolution layer.
684 403 : if( nThisOverview != -1 )
685 0 : return;
686 :
687 : HFAEntry *poEntry =
688 403 : poBand->poNode->GetNamedChild( "Descriptor_Table.Histogram" );
689 403 : if ( poEntry == NULL )
690 370 : return;
691 :
692 33 : int nNumBins = poEntry->GetIntField( "numRows" );
693 :
694 : /* -------------------------------------------------------------------- */
695 : /* Fetch the histogram values. */
696 : /* -------------------------------------------------------------------- */
697 :
698 33 : int nOffset = poEntry->GetIntField( "columnDataPtr" );
699 33 : const char * pszType = poEntry->GetStringField( "dataType" );
700 33 : int nBinSize = 4;
701 :
702 33 : if( pszType != NULL && EQUALN( "real", pszType, 4 ) )
703 27 : nBinSize = 8;
704 :
705 33 : int *panHistValues = (int *) CPLMalloc(sizeof(int) * nNumBins);
706 33 : GByte *pabyWorkBuf = (GByte *) CPLMalloc(nBinSize * nNumBins);
707 :
708 33 : VSIFSeekL( hHFA->fp, nOffset, SEEK_SET );
709 :
710 33 : if( (int)VSIFReadL( pabyWorkBuf, nBinSize, nNumBins, hHFA->fp ) != nNumBins)
711 : {
712 : CPLError( CE_Failure, CPLE_FileIO,
713 0 : "Cannot read histogram values." );
714 0 : CPLFree( panHistValues );
715 0 : CPLFree( pabyWorkBuf );
716 0 : return;
717 : }
718 :
719 : // Swap into local order.
720 33 : for( i = 0; i < nNumBins; i++ )
721 : HFAStandard( nBinSize, pabyWorkBuf + i*nBinSize );
722 :
723 33 : if( nBinSize == 8 ) // source is doubles
724 : {
725 5013 : for( i = 0; i < nNumBins; i++ )
726 4986 : panHistValues[i] = (int) ((double *) pabyWorkBuf)[i];
727 : }
728 : else // source is 32bit integers
729 : {
730 6 : memcpy( panHistValues, pabyWorkBuf, 4 * nNumBins );
731 : }
732 :
733 33 : CPLFree( pabyWorkBuf );
734 33 : pabyWorkBuf = NULL;
735 :
736 : /* -------------------------------------------------------------------- */
737 : /* Do we have unique values for the bins? */
738 : /* -------------------------------------------------------------------- */
739 33 : double *padfBinValues = NULL;
740 33 : HFAEntry *poBinEntry = poBand->poNode->GetNamedChild( "Descriptor_Table.#Bin_Function840#" );
741 :
742 33 : if( poBinEntry != NULL
743 : && EQUAL(poBinEntry->GetType(),"Edsc_BinFunction840")
744 : && EQUAL(poBinEntry->GetStringField( "binFunction.type.string" ),
745 : "BFUnique") )
746 9 : padfBinValues = HFAReadBFUniqueBins( poBinEntry, nNumBins );
747 :
748 33 : if( padfBinValues )
749 : {
750 9 : int nMaxValue = 0;
751 9 : int nMinValue = 1000000;
752 9 : int bAllInteger = TRUE;
753 :
754 839 : for( i = 0; i < nNumBins; i++ )
755 : {
756 830 : if( padfBinValues[i] != floor(padfBinValues[i]) )
757 0 : bAllInteger = FALSE;
758 :
759 830 : nMaxValue = MAX(nMaxValue,(int)padfBinValues[i]);
760 830 : nMinValue = MIN(nMinValue,(int)padfBinValues[i]);
761 : }
762 :
763 9 : if( nMinValue < 0 || nMaxValue > 1000 || !bAllInteger )
764 : {
765 0 : CPLFree( padfBinValues );
766 0 : CPLFree( panHistValues );
767 0 : CPLDebug( "HFA", "Unable to offer histogram because unique values list is not convenient to reform as HISTOBINVALUES." );
768 0 : return;
769 : }
770 :
771 9 : int nNewBins = nMaxValue + 1;
772 9 : int *panNewHistValues = (int *) CPLCalloc(sizeof(int),nNewBins);
773 :
774 839 : for( i = 0; i < nNumBins; i++ )
775 830 : panNewHistValues[(int) padfBinValues[i]] = panHistValues[i];
776 :
777 9 : CPLFree( panHistValues );
778 9 : panHistValues = panNewHistValues;
779 9 : nNumBins = nNewBins;
780 :
781 9 : SetMetadataItem( "STATISTICS_HISTOMIN", "0" );
782 : SetMetadataItem( "STATISTICS_HISTOMAX",
783 9 : CPLString().Printf("%d", nMaxValue ) );
784 : SetMetadataItem( "STATISTICS_HISTONUMBINS",
785 18 : CPLString().Printf("%d", nMaxValue+1 ) );
786 :
787 9 : CPLFree(padfBinValues);
788 9 : padfBinValues = NULL;
789 : }
790 :
791 : /* -------------------------------------------------------------------- */
792 : /* Format into HISTOBINVALUES text format. */
793 : /* -------------------------------------------------------------------- */
794 33 : unsigned int nBufSize = 1024;
795 33 : char * pszBinValues = (char *)CPLMalloc( nBufSize );
796 33 : int nBinValuesLen = 0;
797 33 : pszBinValues[0] = 0;
798 :
799 7701 : for ( int nBin = 0; nBin < nNumBins; ++nBin )
800 : {
801 : char szBuf[32];
802 7668 : snprintf( szBuf, 31, "%d", panHistValues[nBin] );
803 7668 : if ( ( nBinValuesLen + strlen( szBuf ) + 2 ) > nBufSize )
804 : {
805 2 : nBufSize *= 2;
806 2 : char* pszNewBinValues = (char *)VSIRealloc( pszBinValues, nBufSize );
807 2 : if (pszNewBinValues == NULL)
808 : {
809 0 : CPLError(CE_Failure, CPLE_OutOfMemory, "Cannot allocate memory");
810 0 : break;
811 : }
812 2 : pszBinValues = pszNewBinValues;
813 : }
814 7668 : strcat( pszBinValues+nBinValuesLen, szBuf );
815 7668 : strcat( pszBinValues+nBinValuesLen, "|" );
816 7668 : nBinValuesLen += strlen(pszBinValues+nBinValuesLen);
817 : }
818 :
819 33 : SetMetadataItem( "STATISTICS_HISTOBINVALUES", pszBinValues );
820 33 : CPLFree( panHistValues );
821 33 : CPLFree( pszBinValues );
822 : }
823 :
824 : /************************************************************************/
825 : /* GetNoDataValue() */
826 : /************************************************************************/
827 :
828 142 : double HFARasterBand::GetNoDataValue( int *pbSuccess )
829 :
830 : {
831 : double dfNoData;
832 :
833 142 : if( HFAGetBandNoData( hHFA, nBand, &dfNoData ) )
834 : {
835 8 : if( pbSuccess )
836 7 : *pbSuccess = TRUE;
837 8 : return dfNoData;
838 : }
839 : else
840 134 : return GDALPamRasterBand::GetNoDataValue( pbSuccess );
841 : }
842 :
843 : /************************************************************************/
844 : /* SetNoDataValue() */
845 : /************************************************************************/
846 :
847 4 : CPLErr HFARasterBand::SetNoDataValue( double dfValue )
848 : {
849 4 : return HFASetBandNoData( hHFA, nBand, dfValue );
850 : }
851 :
852 : /************************************************************************/
853 : /* GetMinimum() */
854 : /************************************************************************/
855 :
856 1 : double HFARasterBand::GetMinimum( int *pbSuccess )
857 :
858 : {
859 1 : const char *pszValue = GetMetadataItem( "STATISTICS_MINIMUM" );
860 :
861 1 : if( pszValue != NULL )
862 : {
863 1 : if( pbSuccess )
864 1 : *pbSuccess = TRUE;
865 1 : return atof(pszValue);
866 : }
867 : else
868 : {
869 0 : return GDALRasterBand::GetMinimum( pbSuccess );
870 : }
871 : }
872 :
873 : /************************************************************************/
874 : /* GetMaximum() */
875 : /************************************************************************/
876 :
877 1 : double HFARasterBand::GetMaximum( int *pbSuccess )
878 :
879 : {
880 1 : const char *pszValue = GetMetadataItem( "STATISTICS_MAXIMUM" );
881 :
882 1 : if( pszValue != NULL )
883 : {
884 1 : if( pbSuccess )
885 1 : *pbSuccess = TRUE;
886 1 : return atof(pszValue);
887 : }
888 : else
889 : {
890 0 : return GDALRasterBand::GetMaximum( pbSuccess );
891 : }
892 : }
893 :
894 : /************************************************************************/
895 : /* EstablishOverviews() */
896 : /* */
897 : /* Delayed population of overview information. */
898 : /************************************************************************/
899 :
900 31 : void HFARasterBand::EstablishOverviews()
901 :
902 : {
903 31 : if( nOverviews != -1 )
904 18 : return;
905 :
906 13 : nOverviews = HFAGetOverviewCount( hHFA, nBand );
907 13 : if( nOverviews > 0 )
908 : {
909 : papoOverviewBands = (HFARasterBand **)
910 8 : CPLMalloc(sizeof(void*)*nOverviews);
911 :
912 38 : for( int iOvIndex = 0; iOvIndex < nOverviews; iOvIndex++ )
913 : {
914 11 : papoOverviewBands[iOvIndex] =
915 22 : new HFARasterBand( (HFADataset *) poDS, nBand, iOvIndex );
916 : }
917 : }
918 : }
919 :
920 : /************************************************************************/
921 : /* GetOverviewCount() */
922 : /************************************************************************/
923 :
924 15 : int HFARasterBand::GetOverviewCount()
925 :
926 : {
927 15 : EstablishOverviews();
928 :
929 15 : if( nOverviews == 0 )
930 3 : return GDALRasterBand::GetOverviewCount();
931 : else
932 12 : return nOverviews;
933 : }
934 :
935 : /************************************************************************/
936 : /* GetOverview() */
937 : /************************************************************************/
938 :
939 10 : GDALRasterBand *HFARasterBand::GetOverview( int i )
940 :
941 : {
942 10 : EstablishOverviews();
943 :
944 10 : if( nOverviews == 0 )
945 0 : return GDALRasterBand::GetOverview( i );
946 10 : else if( i < 0 || i >= nOverviews )
947 0 : return NULL;
948 : else
949 10 : return papoOverviewBands[i];
950 : }
951 :
952 : /************************************************************************/
953 : /* IReadBlock() */
954 : /************************************************************************/
955 :
956 335 : CPLErr HFARasterBand::IReadBlock( int nBlockXOff, int nBlockYOff,
957 : void * pImage )
958 :
959 : {
960 : CPLErr eErr;
961 335 : int nThisDataType = nHFADataType; // overview may differ.
962 :
963 335 : if( nThisOverview == -1 )
964 : eErr = HFAGetRasterBlockEx( hHFA, nBand, nBlockXOff, nBlockYOff,
965 : pImage,
966 323 : nBlockXSize * nBlockYSize * (GDALGetDataTypeSize(eDataType) / 8) );
967 : else
968 : {
969 : eErr = HFAGetOverviewRasterBlockEx( hHFA, nBand, nThisOverview,
970 : nBlockXOff, nBlockYOff,
971 : pImage,
972 12 : nBlockXSize * nBlockYSize * (GDALGetDataTypeSize(eDataType) / 8));
973 : nThisDataType =
974 12 : hHFA->papoBand[nBand-1]->papoOverviews[nThisOverview]->nDataType;
975 : }
976 :
977 335 : if( eErr == CE_None && nThisDataType == EPT_u4 )
978 : {
979 0 : GByte *pabyData = (GByte *) pImage;
980 :
981 0 : for( int ii = nBlockXSize * nBlockYSize - 2; ii >= 0; ii -= 2 )
982 : {
983 0 : int k = ii>>1;
984 0 : pabyData[ii+1] = (pabyData[k]>>4) & 0xf;
985 0 : pabyData[ii] = (pabyData[k]) & 0xf;
986 : }
987 : }
988 335 : if( eErr == CE_None && nThisDataType == EPT_u2 )
989 : {
990 0 : GByte *pabyData = (GByte *) pImage;
991 :
992 0 : for( int ii = nBlockXSize * nBlockYSize - 4; ii >= 0; ii -= 4 )
993 : {
994 0 : int k = ii>>2;
995 0 : pabyData[ii+3] = (pabyData[k]>>6) & 0x3;
996 0 : pabyData[ii+2] = (pabyData[k]>>4) & 0x3;
997 0 : pabyData[ii+1] = (pabyData[k]>>2) & 0x3;
998 0 : pabyData[ii] = (pabyData[k]) & 0x3;
999 : }
1000 : }
1001 335 : if( eErr == CE_None && nThisDataType == EPT_u1)
1002 : {
1003 4 : GByte *pabyData = (GByte *) pImage;
1004 :
1005 16388 : for( int ii = nBlockXSize * nBlockYSize - 1; ii >= 0; ii-- )
1006 : {
1007 16384 : if( (pabyData[ii>>3] & (1 << (ii & 0x7))) )
1008 504 : pabyData[ii] = 1;
1009 : else
1010 15880 : pabyData[ii] = 0;
1011 : }
1012 : }
1013 :
1014 335 : return eErr;
1015 : }
1016 :
1017 : /************************************************************************/
1018 : /* IWriteBlock() */
1019 : /************************************************************************/
1020 :
1021 119 : CPLErr HFARasterBand::IWriteBlock( int nBlockXOff, int nBlockYOff,
1022 : void * pImage )
1023 :
1024 : {
1025 119 : GByte *pabyOutBuf = (GByte *) pImage;
1026 :
1027 : /* -------------------------------------------------------------------- */
1028 : /* Do we need to pack 1/2/4 bit data? */
1029 : /* -------------------------------------------------------------------- */
1030 119 : if( nHFADataType == EPT_u1
1031 : || nHFADataType == EPT_u2
1032 : || nHFADataType == EPT_u4 )
1033 : {
1034 2 : int nPixCount = nBlockXSize * nBlockYSize;
1035 2 : pabyOutBuf = (GByte *) VSIMalloc2(nBlockXSize, nBlockYSize);
1036 2 : if (pabyOutBuf == NULL)
1037 0 : return CE_Failure;
1038 :
1039 2 : if( nHFADataType == EPT_u1 )
1040 : {
1041 1026 : for( int ii = 0; ii < nPixCount - 7; ii += 8 )
1042 : {
1043 1024 : int k = ii>>3;
1044 1024 : pabyOutBuf[k] =
1045 1024 : (((GByte *) pImage)[ii] & 0x1)
1046 1024 : | ((((GByte *) pImage)[ii+1]&0x1) << 1)
1047 1024 : | ((((GByte *) pImage)[ii+2]&0x1) << 2)
1048 1024 : | ((((GByte *) pImage)[ii+3]&0x1) << 3)
1049 1024 : | ((((GByte *) pImage)[ii+4]&0x1) << 4)
1050 1024 : | ((((GByte *) pImage)[ii+5]&0x1) << 5)
1051 1024 : | ((((GByte *) pImage)[ii+6]&0x1) << 6)
1052 7168 : | ((((GByte *) pImage)[ii+7]&0x1) << 7);
1053 : }
1054 : }
1055 0 : else if( nHFADataType == EPT_u2 )
1056 : {
1057 0 : for( int ii = 0; ii < nPixCount - 3; ii += 4 )
1058 : {
1059 0 : int k = ii>>2;
1060 0 : pabyOutBuf[k] =
1061 0 : (((GByte *) pImage)[ii] & 0x3)
1062 0 : | ((((GByte *) pImage)[ii+1]&0x3) << 2)
1063 0 : | ((((GByte *) pImage)[ii+2]&0x3) << 4)
1064 0 : | ((((GByte *) pImage)[ii+3]&0x3) << 6);
1065 : }
1066 : }
1067 0 : else if( nHFADataType == EPT_u4 )
1068 : {
1069 0 : for( int ii = 0; ii < nPixCount - 1; ii += 2 )
1070 : {
1071 0 : int k = ii>>1;
1072 0 : pabyOutBuf[k] =
1073 0 : (((GByte *) pImage)[ii] & 0xf)
1074 0 : | ((((GByte *) pImage)[ii+1]&0xf) << 4);
1075 : }
1076 : }
1077 : }
1078 :
1079 : /* -------------------------------------------------------------------- */
1080 : /* Actually write out. */
1081 : /* -------------------------------------------------------------------- */
1082 : CPLErr nRetCode;
1083 :
1084 119 : if( nThisOverview == -1 )
1085 : nRetCode = HFASetRasterBlock( hHFA, nBand, nBlockXOff, nBlockYOff,
1086 111 : pabyOutBuf );
1087 : else
1088 : nRetCode = HFASetOverviewRasterBlock( hHFA, nBand, nThisOverview,
1089 : nBlockXOff, nBlockYOff,
1090 8 : pabyOutBuf );
1091 :
1092 119 : if( pabyOutBuf != pImage )
1093 2 : CPLFree( pabyOutBuf );
1094 :
1095 119 : return nRetCode;
1096 : }
1097 :
1098 : /************************************************************************/
1099 : /* GetDescription() */
1100 : /************************************************************************/
1101 :
1102 83 : const char * HFARasterBand::GetDescription() const
1103 : {
1104 83 : const char *pszName = HFAGetBandName( hHFA, nBand );
1105 :
1106 83 : if( pszName == NULL )
1107 0 : return GDALPamRasterBand::GetDescription();
1108 : else
1109 83 : return pszName;
1110 : }
1111 :
1112 : /************************************************************************/
1113 : /* SetDescription() */
1114 : /************************************************************************/
1115 111 : void HFARasterBand::SetDescription( const char *pszName )
1116 : {
1117 111 : HFASetBandName( hHFA, nBand, pszName );
1118 111 : }
1119 :
1120 : /************************************************************************/
1121 : /* GetColorInterpretation() */
1122 : /************************************************************************/
1123 :
1124 43 : GDALColorInterp HFARasterBand::GetColorInterpretation()
1125 :
1126 : {
1127 43 : if( poCT != NULL )
1128 1 : return GCI_PaletteIndex;
1129 : else
1130 42 : return GCI_Undefined;
1131 : }
1132 :
1133 : /************************************************************************/
1134 : /* GetColorTable() */
1135 : /************************************************************************/
1136 :
1137 11 : GDALColorTable *HFARasterBand::GetColorTable()
1138 :
1139 : {
1140 11 : return poCT;
1141 : }
1142 :
1143 : /************************************************************************/
1144 : /* SetColorTable() */
1145 : /************************************************************************/
1146 :
1147 3 : CPLErr HFARasterBand::SetColorTable( GDALColorTable * poCTable )
1148 :
1149 : {
1150 3 : if( GetAccess() == GA_ReadOnly )
1151 : {
1152 : CPLError( CE_Failure, CPLE_NoWriteAccess,
1153 0 : "Unable to set color table on read-only file." );
1154 0 : return CE_Failure;
1155 : }
1156 :
1157 : /* -------------------------------------------------------------------- */
1158 : /* Special case if we are clearing the color table. */
1159 : /* -------------------------------------------------------------------- */
1160 3 : if( poCTable == NULL )
1161 : {
1162 2 : delete poCT;
1163 2 : poCT = NULL;
1164 :
1165 2 : HFASetPCT( hHFA, nBand, 0, NULL, NULL, NULL, NULL );
1166 :
1167 2 : return CE_None;
1168 : }
1169 :
1170 : /* -------------------------------------------------------------------- */
1171 : /* Write out the colortable, and update the configuration. */
1172 : /* -------------------------------------------------------------------- */
1173 1 : int nColors = poCTable->GetColorEntryCount();
1174 :
1175 : double *padfRed, *padfGreen, *padfBlue, *padfAlpha;
1176 :
1177 1 : padfRed = (double *) CPLMalloc(sizeof(double) * nColors);
1178 1 : padfGreen = (double *) CPLMalloc(sizeof(double) * nColors);
1179 1 : padfBlue = (double *) CPLMalloc(sizeof(double) * nColors);
1180 1 : padfAlpha = (double *) CPLMalloc(sizeof(double) * nColors);
1181 :
1182 257 : for( int iColor = 0; iColor < nColors; iColor++ )
1183 : {
1184 : GDALColorEntry sRGB;
1185 :
1186 256 : poCTable->GetColorEntryAsRGB( iColor, &sRGB );
1187 :
1188 256 : padfRed[iColor] = sRGB.c1 / 255.0;
1189 256 : padfGreen[iColor] = sRGB.c2 / 255.0;
1190 256 : padfBlue[iColor] = sRGB.c3 / 255.0;
1191 256 : padfAlpha[iColor] = sRGB.c4 / 255.0;
1192 : }
1193 :
1194 : HFASetPCT( hHFA, nBand, nColors,
1195 1 : padfRed, padfGreen, padfBlue, padfAlpha);
1196 :
1197 1 : CPLFree( padfRed );
1198 1 : CPLFree( padfGreen );
1199 1 : CPLFree( padfBlue );
1200 1 : CPLFree( padfAlpha );
1201 :
1202 1 : if( poCT )
1203 0 : delete poCT;
1204 :
1205 1 : poCT = poCTable->Clone();
1206 :
1207 1 : return CE_None;
1208 : }
1209 :
1210 : /************************************************************************/
1211 : /* SetMetadata() */
1212 : /************************************************************************/
1213 :
1214 52 : CPLErr HFARasterBand::SetMetadata( char **papszMDIn, const char *pszDomain )
1215 :
1216 : {
1217 52 : bMetadataDirty = TRUE;
1218 :
1219 52 : return GDALPamRasterBand::SetMetadata( papszMDIn, pszDomain );
1220 : }
1221 :
1222 : /************************************************************************/
1223 : /* SetMetadata() */
1224 : /************************************************************************/
1225 :
1226 816 : CPLErr HFARasterBand::SetMetadataItem( const char *pszTag, const char *pszValue,
1227 : const char *pszDomain )
1228 :
1229 : {
1230 816 : bMetadataDirty = TRUE;
1231 :
1232 816 : return GDALPamRasterBand::SetMetadataItem( pszTag, pszValue, pszDomain );
1233 : }
1234 :
1235 : /************************************************************************/
1236 : /* CleanOverviews() */
1237 : /************************************************************************/
1238 :
1239 1 : CPLErr HFARasterBand::CleanOverviews()
1240 :
1241 : {
1242 1 : if( nOverviews == 0 )
1243 0 : return CE_None;
1244 :
1245 : /* -------------------------------------------------------------------- */
1246 : /* Clear our reference to overviews as bands. */
1247 : /* -------------------------------------------------------------------- */
1248 : int iOverview;
1249 :
1250 2 : for( iOverview = 0; iOverview < nOverviews; iOverview++ )
1251 1 : delete papoOverviewBands[iOverview];
1252 :
1253 1 : CPLFree( papoOverviewBands );
1254 1 : papoOverviewBands = NULL;
1255 1 : nOverviews = 0;
1256 :
1257 : /* -------------------------------------------------------------------- */
1258 : /* Search for any RRDNamesList and destroy it. */
1259 : /* -------------------------------------------------------------------- */
1260 1 : HFABand *poBand = hHFA->papoBand[nBand-1];
1261 1 : HFAEntry *poEntry = poBand->poNode->GetNamedChild( "RRDNamesList" );
1262 1 : if( poEntry != NULL )
1263 : {
1264 1 : poEntry->RemoveAndDestroy();
1265 : }
1266 :
1267 : /* -------------------------------------------------------------------- */
1268 : /* Destroy and subsample layers under our band. */
1269 : /* -------------------------------------------------------------------- */
1270 : HFAEntry *poChild;
1271 6 : for( poChild = poBand->poNode->GetChild();
1272 : poChild != NULL; )
1273 : {
1274 4 : HFAEntry *poNext = poChild->GetNext();
1275 :
1276 4 : if( EQUAL(poChild->GetType(),"Eimg_Layer_SubSample") )
1277 0 : poChild->RemoveAndDestroy();
1278 :
1279 4 : poChild = poNext;
1280 : }
1281 :
1282 : /* -------------------------------------------------------------------- */
1283 : /* Clean up dependent file if we are the last band under the */
1284 : /* assumption there will be nothing else referencing it after */
1285 : /* this. */
1286 : /* -------------------------------------------------------------------- */
1287 1 : if( hHFA->psDependent != hHFA && hHFA->psDependent != NULL )
1288 : {
1289 : CPLString osFilename =
1290 : CPLFormFilename( hHFA->psDependent->pszPath,
1291 1 : hHFA->psDependent->pszFilename, NULL );
1292 :
1293 1 : HFAClose( hHFA->psDependent );
1294 1 : hHFA->psDependent = NULL;
1295 :
1296 1 : CPLDebug( "HFA", "Unlink(%s)", osFilename.c_str() );
1297 1 : VSIUnlink( osFilename );
1298 : }
1299 :
1300 1 : return CE_None;
1301 : }
1302 :
1303 : /************************************************************************/
1304 : /* BuildOverviews() */
1305 : /************************************************************************/
1306 :
1307 6 : CPLErr HFARasterBand::BuildOverviews( const char *pszResampling,
1308 : int nReqOverviews, int *panOverviewList,
1309 : GDALProgressFunc pfnProgress,
1310 : void *pProgressData )
1311 :
1312 : {
1313 : int iOverview;
1314 : GDALRasterBand **papoOvBands;
1315 6 : int bNoRegen = FALSE;
1316 :
1317 6 : EstablishOverviews();
1318 :
1319 6 : if( nThisOverview != -1 )
1320 : {
1321 : CPLError( CE_Failure, CPLE_AppDefined,
1322 0 : "Attempt to build overviews on an overview layer." );
1323 :
1324 0 : return CE_Failure;
1325 : }
1326 :
1327 6 : if( nReqOverviews == 0 )
1328 1 : return CleanOverviews();
1329 :
1330 5 : papoOvBands = (GDALRasterBand **) CPLCalloc(sizeof(void*),nReqOverviews);
1331 :
1332 5 : if( EQUALN(pszResampling,"NO_REGEN:",9) )
1333 : {
1334 2 : pszResampling += 9;
1335 2 : bNoRegen = TRUE;
1336 : }
1337 :
1338 : /* -------------------------------------------------------------------- */
1339 : /* Loop over overview levels requested. */
1340 : /* -------------------------------------------------------------------- */
1341 11 : for( iOverview = 0; iOverview < nReqOverviews; iOverview++ )
1342 : {
1343 : /* -------------------------------------------------------------------- */
1344 : /* Find this overview level. */
1345 : /* -------------------------------------------------------------------- */
1346 6 : int i, iResult = -1, nReqOvLevel;
1347 :
1348 : nReqOvLevel =
1349 6 : GDALOvLevelAdjust(panOverviewList[iOverview],nRasterXSize);
1350 :
1351 10 : for( i = 0; i < nOverviews && papoOvBands[iOverview] == NULL; i++ )
1352 : {
1353 : int nThisOvLevel;
1354 :
1355 : nThisOvLevel = (int) (0.5 + GetXSize()
1356 4 : / (double) papoOverviewBands[i]->GetXSize());
1357 :
1358 4 : if( nReqOvLevel == nThisOvLevel )
1359 1 : papoOvBands[iOverview] = papoOverviewBands[i];
1360 : }
1361 :
1362 : /* -------------------------------------------------------------------- */
1363 : /* If this overview level does not yet exist, create it now. */
1364 : /* -------------------------------------------------------------------- */
1365 6 : if( papoOvBands[iOverview] == NULL )
1366 : {
1367 5 : iResult = HFACreateOverview( hHFA, nBand, panOverviewList[iOverview] );
1368 5 : if( iResult < 0 )
1369 0 : return CE_Failure;
1370 :
1371 5 : nOverviews = iResult + 1;
1372 : papoOverviewBands = (HFARasterBand **)
1373 5 : CPLRealloc( papoOverviewBands, sizeof(void*) * nOverviews);
1374 5 : papoOverviewBands[iResult] = new HFARasterBand(
1375 10 : (HFADataset *) poDS, nBand, iResult );
1376 :
1377 5 : papoOvBands[iOverview] = papoOverviewBands[iResult];
1378 : }
1379 :
1380 : }
1381 :
1382 : /* -------------------------------------------------------------------- */
1383 : /* Regenerate the overviews. */
1384 : /* -------------------------------------------------------------------- */
1385 5 : CPLErr eErr = CE_None;
1386 :
1387 5 : if( !bNoRegen )
1388 : eErr = GDALRegenerateOverviews( (GDALRasterBandH) this,
1389 : nReqOverviews,
1390 : (GDALRasterBandH *) papoOvBands,
1391 : pszResampling,
1392 3 : pfnProgress, pProgressData );
1393 :
1394 5 : CPLFree( papoOvBands );
1395 :
1396 5 : return eErr;
1397 : }
1398 :
1399 : /************************************************************************/
1400 : /* GetDefaultHistogram() */
1401 : /************************************************************************/
1402 :
1403 : CPLErr
1404 3 : HFARasterBand::GetDefaultHistogram( double *pdfMin, double *pdfMax,
1405 : int *pnBuckets, int ** ppanHistogram,
1406 : int bForce,
1407 : GDALProgressFunc pfnProgress,
1408 : void *pProgressData)
1409 :
1410 : {
1411 5 : if( GetMetadataItem( "STATISTICS_HISTOBINVALUES" ) != NULL
1412 1 : && GetMetadataItem( "STATISTICS_HISTOMIN" ) != NULL
1413 1 : && GetMetadataItem( "STATISTICS_HISTOMAX" ) != NULL )
1414 : {
1415 : int i;
1416 : const char *pszNextBin;
1417 : const char *pszBinValues =
1418 1 : GetMetadataItem( "STATISTICS_HISTOBINVALUES" );
1419 :
1420 1 : *pdfMin = atof(GetMetadataItem("STATISTICS_HISTOMIN"));
1421 1 : *pdfMax = atof(GetMetadataItem("STATISTICS_HISTOMAX"));
1422 :
1423 1 : *pnBuckets = 0;
1424 573 : for( i = 0; pszBinValues[i] != '\0'; i++ )
1425 : {
1426 572 : if( pszBinValues[i] == '|' )
1427 220 : (*pnBuckets)++;
1428 : }
1429 :
1430 1 : *ppanHistogram = (int *) CPLCalloc(sizeof(int),*pnBuckets);
1431 :
1432 1 : pszNextBin = pszBinValues;
1433 221 : for( i = 0; i < *pnBuckets; i++ )
1434 : {
1435 220 : (*ppanHistogram)[i] = atoi(pszNextBin);
1436 :
1437 792 : while( *pszNextBin != '|' && *pszNextBin != '\0' )
1438 352 : pszNextBin++;
1439 220 : if( *pszNextBin == '|' )
1440 220 : pszNextBin++;
1441 : }
1442 :
1443 : // Adjust min/max to reflect outer edges of buckets.
1444 1 : double dfBucketWidth = (*pdfMax - *pdfMin) / (*pnBuckets-1);
1445 1 : *pdfMax += 0.5 * dfBucketWidth;
1446 1 : *pdfMin -= 0.5 * dfBucketWidth;
1447 :
1448 1 : return CE_None;
1449 : }
1450 : else
1451 : return GDALPamRasterBand::GetDefaultHistogram( pdfMin, pdfMax,
1452 : pnBuckets,ppanHistogram,
1453 : bForce,
1454 : pfnProgress,
1455 2 : pProgressData );
1456 : }
1457 :
1458 : /************************************************************************/
1459 : /* SetDefaultRAT() */
1460 : /************************************************************************/
1461 :
1462 0 : CPLErr HFARasterBand::SetDefaultRAT( const GDALRasterAttributeTable * poRAT )
1463 :
1464 : {
1465 0 : return GDALPamRasterBand::SetDefaultRAT( poRAT );
1466 : }
1467 :
1468 : /************************************************************************/
1469 : /* GetDefaultRAT() */
1470 : /************************************************************************/
1471 :
1472 11 : const GDALRasterAttributeTable *HFARasterBand::GetDefaultRAT()
1473 :
1474 : {
1475 11 : return poDefaultRAT;
1476 : }
1477 :
1478 : /************************************************************************/
1479 : /* ReadNamedRAT() */
1480 : /************************************************************************/
1481 :
1482 419 : GDALRasterAttributeTable *HFARasterBand::ReadNamedRAT( const char *pszName )
1483 :
1484 : {
1485 : /* -------------------------------------------------------------------- */
1486 : /* Find the requested table. */
1487 : /* -------------------------------------------------------------------- */
1488 419 : HFAEntry *poDT = hHFA->papoBand[nBand-1]->poNode->GetNamedChild(pszName);
1489 :
1490 419 : if( poDT == NULL )
1491 379 : return NULL;
1492 :
1493 : /* -------------------------------------------------------------------- */
1494 : /* Create a corresponding RAT. */
1495 : /* -------------------------------------------------------------------- */
1496 40 : GDALRasterAttributeTable *poRAT = NULL;
1497 40 : int nRowCount = poDT->GetIntField( "numRows" );
1498 :
1499 40 : poRAT = new GDALRasterAttributeTable();
1500 :
1501 : /* -------------------------------------------------------------------- */
1502 : /* Scan under table for columns. */
1503 : /* -------------------------------------------------------------------- */
1504 : HFAEntry *poDTChild;
1505 :
1506 127 : for( poDTChild = poDT->GetChild();
1507 : poDTChild != NULL;
1508 : poDTChild = poDTChild->GetNext() )
1509 : {
1510 : int i;
1511 :
1512 87 : if( EQUAL(poDTChild->GetType(),"Edsc_BinFunction") )
1513 : {
1514 26 : double dfMax = poDTChild->GetDoubleField( "maxLimit" );
1515 26 : double dfMin = poDTChild->GetDoubleField( "minLimit" );
1516 26 : int nBinCount = poDTChild->GetIntField( "numBins" );
1517 :
1518 26 : if( nBinCount == nRowCount
1519 : && dfMax != dfMin && nBinCount != 0 )
1520 : poRAT->SetLinearBinning( dfMin,
1521 26 : (dfMax-dfMin) / (nBinCount-1) );
1522 : }
1523 :
1524 87 : if( EQUAL(poDTChild->GetType(),"Edsc_BinFunction840")
1525 : && EQUAL(poDTChild->GetStringField( "binFunction.type.string" ),
1526 : "BFUnique") )
1527 : {
1528 9 : double *padfBinValues = HFAReadBFUniqueBins( poDTChild, nRowCount );
1529 :
1530 9 : if( padfBinValues != NULL )
1531 : {
1532 9 : poRAT->CreateColumn( "Value", GFT_Integer, GFU_MinMax );
1533 839 : for( i = 0; i < nRowCount; i++ )
1534 : poRAT->SetValue( i, poRAT->GetColumnCount()-1,
1535 830 : padfBinValues[i] );
1536 :
1537 9 : CPLFree( padfBinValues );
1538 : }
1539 : }
1540 :
1541 87 : if( !EQUAL(poDTChild->GetType(),"Edsc_Column") )
1542 35 : continue;
1543 :
1544 52 : int nOffset = poDTChild->GetIntField( "columnDataPtr" );
1545 52 : const char * pszType = poDTChild->GetStringField( "dataType" );
1546 52 : GDALRATFieldUsage eType = GFU_Generic;
1547 :
1548 52 : if( pszType == NULL || nOffset == 0 )
1549 0 : continue;
1550 :
1551 52 : if( EQUAL(poDTChild->GetName(),"Histogram") )
1552 34 : eType = GFU_Generic;
1553 18 : else if( EQUAL(poDTChild->GetName(),"Red") )
1554 4 : eType = GFU_Red;
1555 14 : else if( EQUAL(poDTChild->GetName(),"Green") )
1556 4 : eType = GFU_Green;
1557 10 : else if( EQUAL(poDTChild->GetName(),"Blue") )
1558 4 : eType = GFU_Blue;
1559 6 : else if( EQUAL(poDTChild->GetName(),"Alpha") )
1560 0 : eType = GFU_Alpha;
1561 6 : else if( EQUAL(poDTChild->GetName(),"Class_Names") )
1562 0 : eType = GFU_Name;
1563 :
1564 52 : if( EQUAL(pszType,"real") )
1565 : {
1566 44 : double *padfColData = (double*)VSIMalloc2(nRowCount, sizeof(double));
1567 44 : if (nRowCount != 0 && padfColData == NULL)
1568 : {
1569 : CPLError( CE_Failure, CPLE_OutOfMemory,
1570 0 : "HFARasterBand::ReadNamedRAT : Out of memory");
1571 0 : delete poRAT;
1572 0 : return NULL;
1573 : }
1574 :
1575 44 : VSIFSeekL( hHFA->fp, nOffset, SEEK_SET );
1576 44 : if ((int)VSIFReadL( padfColData, sizeof(double), nRowCount, hHFA->fp ) != nRowCount)
1577 : {
1578 : CPLError( CE_Failure, CPLE_AppDefined,
1579 0 : "HFARasterBand::ReadNamedRAT : Cannot read values");
1580 0 : CPLFree(padfColData);
1581 0 : delete poRAT;
1582 0 : return NULL;
1583 : }
1584 : #ifdef CPL_MSB
1585 : GDALSwapWords( padfColData, 8, nRowCount, 8 );
1586 : #endif
1587 44 : poRAT->CreateColumn( poDTChild->GetName(), GFT_Real, eType );
1588 7171 : for( i = 0; i < nRowCount; i++ )
1589 7127 : poRAT->SetValue( i, poRAT->GetColumnCount()-1, padfColData[i]);
1590 :
1591 44 : CPLFree( padfColData );
1592 : }
1593 8 : else if( EQUAL(pszType,"string") )
1594 : {
1595 1 : int nMaxNumChars = poDTChild->GetIntField( "maxNumChars" );
1596 1 : char *pachColData = (char*)VSICalloc(nRowCount+1,nMaxNumChars);
1597 1 : if (pachColData == NULL)
1598 : {
1599 : CPLError( CE_Failure, CPLE_OutOfMemory,
1600 0 : "HFARasterBand::ReadNamedRAT : Out of memory");
1601 0 : delete poRAT;
1602 0 : return NULL;
1603 : }
1604 :
1605 1 : VSIFSeekL( hHFA->fp, nOffset, SEEK_SET );
1606 1 : if ((int)VSIFReadL( pachColData, nMaxNumChars, nRowCount, hHFA->fp ) != nRowCount)
1607 : {
1608 : CPLError( CE_Failure, CPLE_AppDefined,
1609 0 : "HFARasterBand::ReadNamedRAT : Cannot read values");
1610 0 : CPLFree(pachColData);
1611 0 : delete poRAT;
1612 0 : return NULL;
1613 : }
1614 :
1615 1 : poRAT->CreateColumn(poDTChild->GetName(),GFT_String,eType);
1616 219 : for( i = 0; i < nRowCount; i++ )
1617 : {
1618 218 : CPLString oRowVal;
1619 :
1620 218 : oRowVal.assign( pachColData+nMaxNumChars*i, nMaxNumChars );
1621 : poRAT->SetValue( i, poRAT->GetColumnCount()-1,
1622 218 : oRowVal.c_str() );
1623 : }
1624 :
1625 1 : CPLFree( pachColData );
1626 : }
1627 7 : else if( EQUALN(pszType,"int",3) )
1628 : {
1629 7 : GInt32 *panColData = (GInt32*)VSIMalloc2(nRowCount, sizeof(GInt32));
1630 7 : if (nRowCount != 0 && panColData == NULL)
1631 : {
1632 : CPLError( CE_Failure, CPLE_OutOfMemory,
1633 0 : "HFARasterBand::ReadNamedRAT : Out of memory");
1634 0 : delete poRAT;
1635 0 : return NULL;
1636 : }
1637 :
1638 7 : VSIFSeekL( hHFA->fp, nOffset, SEEK_SET );
1639 7 : if ((int)VSIFReadL( panColData, sizeof(GInt32), nRowCount, hHFA->fp ) != nRowCount)
1640 : {
1641 : CPLError( CE_Failure, CPLE_AppDefined,
1642 0 : "HFARasterBand::ReadNamedRAT : Cannot read values");
1643 0 : CPLFree(panColData);
1644 0 : delete poRAT;
1645 0 : return NULL;
1646 : }
1647 : #ifdef CPL_MSB
1648 : GDALSwapWords( panColData, 4, nRowCount, 4 );
1649 : #endif
1650 7 : poRAT->CreateColumn(poDTChild->GetName(),GFT_Integer,eType);
1651 1739 : for( i = 0; i < nRowCount; i++ )
1652 1732 : poRAT->SetValue( i, poRAT->GetColumnCount()-1, panColData[i] );
1653 :
1654 7 : CPLFree( panColData );
1655 : }
1656 : }
1657 :
1658 40 : return poRAT;
1659 : }
1660 :
1661 : /************************************************************************/
1662 : /* ==================================================================== */
1663 : /* HFADataset */
1664 : /* ==================================================================== */
1665 : /************************************************************************/
1666 :
1667 : /************************************************************************/
1668 : /* HFADataset() */
1669 : /************************************************************************/
1670 :
1671 277 : HFADataset::HFADataset()
1672 :
1673 : {
1674 277 : hHFA = NULL;
1675 277 : bGeoDirty = FALSE;
1676 277 : pszProjection = CPLStrdup("");
1677 277 : bMetadataDirty = FALSE;
1678 277 : bIgnoreUTM = FALSE;
1679 277 : bForceToPEString = FALSE;
1680 :
1681 277 : nGCPCount = 0;
1682 277 : }
1683 :
1684 : /************************************************************************/
1685 : /* ~HFADataset() */
1686 : /************************************************************************/
1687 :
1688 554 : HFADataset::~HFADataset()
1689 :
1690 : {
1691 277 : FlushCache();
1692 :
1693 : /* -------------------------------------------------------------------- */
1694 : /* Destroy the raster bands if they exist. We forcably clean */
1695 : /* them up now to avoid any effort to write to them after the */
1696 : /* file is closed. */
1697 : /* -------------------------------------------------------------------- */
1698 : int i;
1699 :
1700 680 : for( i = 0; i < nBands && papoBands != NULL; i++ )
1701 : {
1702 403 : if( papoBands[i] != NULL )
1703 403 : delete papoBands[i];
1704 : }
1705 :
1706 277 : CPLFree( papoBands );
1707 277 : papoBands = NULL;
1708 :
1709 : /* -------------------------------------------------------------------- */
1710 : /* Close the file */
1711 : /* -------------------------------------------------------------------- */
1712 277 : if( hHFA != NULL )
1713 : {
1714 277 : HFAClose( hHFA );
1715 277 : hHFA = NULL;
1716 : }
1717 :
1718 : /* -------------------------------------------------------------------- */
1719 : /* Cleanup */
1720 : /* -------------------------------------------------------------------- */
1721 277 : CPLFree( pszProjection );
1722 277 : if( nGCPCount > 0 )
1723 1 : GDALDeinitGCPs( 36, asGCPList );
1724 554 : }
1725 :
1726 : /************************************************************************/
1727 : /* FlushCache() */
1728 : /************************************************************************/
1729 :
1730 281 : void HFADataset::FlushCache()
1731 :
1732 : {
1733 281 : GDALPamDataset::FlushCache();
1734 :
1735 281 : if( eAccess != GA_Update )
1736 177 : return;
1737 :
1738 104 : if( bGeoDirty )
1739 47 : WriteProjection();
1740 :
1741 104 : if( bMetadataDirty && GetMetadata() != NULL )
1742 : {
1743 42 : HFASetMetadata( hHFA, 0, GetMetadata() );
1744 42 : bMetadataDirty = FALSE;
1745 : }
1746 :
1747 258 : for( int iBand = 0; iBand < nBands; iBand++ )
1748 : {
1749 154 : HFARasterBand *poBand = (HFARasterBand *) GetRasterBand(iBand+1);
1750 154 : if( poBand->bMetadataDirty && poBand->GetMetadata() != NULL )
1751 : {
1752 4 : HFASetMetadata( hHFA, iBand+1, poBand->GetMetadata() );
1753 4 : poBand->bMetadataDirty = FALSE;
1754 : }
1755 : }
1756 :
1757 104 : if( nGCPCount > 0 )
1758 : {
1759 0 : GDALDeinitGCPs( nGCPCount, asGCPList );
1760 : }
1761 : }
1762 :
1763 : /************************************************************************/
1764 : /* WriteProjection() */
1765 : /************************************************************************/
1766 :
1767 47 : CPLErr HFADataset::WriteProjection()
1768 :
1769 : {
1770 : Eprj_Datum sDatum;
1771 : Eprj_ProParameters sPro;
1772 : Eprj_MapInfo sMapInfo;
1773 47 : OGRSpatialReference oSRS;
1774 47 : OGRSpatialReference *poGeogSRS = NULL;
1775 : int bHaveSRS;
1776 47 : char *pszP = pszProjection;
1777 47 : OGRBoolean peStrStored = FALSE;
1778 :
1779 47 : bGeoDirty = FALSE;
1780 :
1781 47 : if( pszProjection != NULL && strlen(pszProjection) > 0
1782 : && oSRS.importFromWkt( &pszP ) == OGRERR_NONE )
1783 36 : bHaveSRS = TRUE;
1784 : else
1785 11 : bHaveSRS = FALSE;
1786 :
1787 : /* -------------------------------------------------------------------- */
1788 : /* Initialize projection and datum. */
1789 : /* -------------------------------------------------------------------- */
1790 47 : memset( &sPro, 0, sizeof(sPro) );
1791 47 : memset( &sDatum, 0, sizeof(sDatum) );
1792 47 : memset( &sMapInfo, 0, sizeof(sMapInfo) );
1793 :
1794 : /* -------------------------------------------------------------------- */
1795 : /* Collect datum information. */
1796 : /* -------------------------------------------------------------------- */
1797 47 : if( bHaveSRS )
1798 : {
1799 36 : poGeogSRS = oSRS.CloneGeogCS();
1800 : }
1801 :
1802 47 : if( poGeogSRS )
1803 : {
1804 : int i;
1805 :
1806 36 : sDatum.datumname = (char *) poGeogSRS->GetAttrValue( "GEOGCS|DATUM" );
1807 :
1808 : /* WKT to Imagine translation */
1809 66 : for( i = 0; apszDatumMap[i] != NULL; i += 2 )
1810 : {
1811 66 : if( EQUAL(sDatum.datumname,apszDatumMap[i+1]) )
1812 : {
1813 36 : sDatum.datumname = (char *) apszDatumMap[i];
1814 36 : break;
1815 : }
1816 : }
1817 :
1818 : /* Map some EPSG datum codes directly to Imagine names */
1819 36 : int nGCS = poGeogSRS->GetEPSGGeogCS();
1820 :
1821 36 : if( nGCS == 4326 )
1822 14 : sDatum.datumname = (char*) "WGS 84";
1823 36 : if( nGCS == 4322 )
1824 0 : sDatum.datumname = (char*) "WGS 1972";
1825 36 : if( nGCS == 4267 )
1826 20 : sDatum.datumname = (char*) "NAD27";
1827 36 : if( nGCS == 4269 )
1828 2 : sDatum.datumname = (char*) "NAD83";
1829 :
1830 36 : if( poGeogSRS->GetTOWGS84( sDatum.params ) == OGRERR_NONE )
1831 1 : sDatum.type = EPRJ_DATUM_PARAMETRIC;
1832 35 : else if( EQUAL(sDatum.datumname,"NAD27") )
1833 : {
1834 19 : sDatum.type = EPRJ_DATUM_GRID;
1835 19 : sDatum.gridname = (char*) "nadcon.dat";
1836 : }
1837 : else
1838 : {
1839 : /* we will default to this (effectively WGS84) for now */
1840 16 : sDatum.type = EPRJ_DATUM_PARAMETRIC;
1841 : }
1842 :
1843 : /* Verify if we need to write a ESRI PE string */
1844 36 : peStrStored = WritePeStringIfNeeded(&oSRS, hHFA);
1845 :
1846 : sPro.proSpheroid.sphereName = (char *)
1847 36 : poGeogSRS->GetAttrValue( "GEOGCS|DATUM|SPHEROID" );
1848 36 : sPro.proSpheroid.a = poGeogSRS->GetSemiMajor();
1849 36 : sPro.proSpheroid.b = poGeogSRS->GetSemiMinor();
1850 36 : sPro.proSpheroid.radius = sPro.proSpheroid.a;
1851 :
1852 36 : double a2 = sPro.proSpheroid.a*sPro.proSpheroid.a;
1853 36 : double b2 = sPro.proSpheroid.b*sPro.proSpheroid.b;
1854 :
1855 36 : sPro.proSpheroid.eSquared = (a2-b2)/a2;
1856 : }
1857 :
1858 : /* -------------------------------------------------------------------- */
1859 : /* Recognise various projections. */
1860 : /* -------------------------------------------------------------------- */
1861 47 : const char * pszProjName = NULL;
1862 :
1863 47 : if( bHaveSRS )
1864 36 : pszProjName = oSRS.GetAttrValue( "PROJCS|PROJECTION" );
1865 :
1866 47 : if( bForceToPEString )
1867 : {
1868 0 : char *pszPEString = NULL;
1869 0 : oSRS.morphToESRI();
1870 0 : oSRS.exportToWkt( &pszPEString );
1871 : // need to transform this into ESRI format.
1872 0 : HFASetPEString( hHFA, pszPEString );
1873 0 : CPLFree( pszPEString );
1874 : }
1875 47 : else if( pszProjName == NULL )
1876 : {
1877 24 : if( bHaveSRS && oSRS.IsGeographic() )
1878 : {
1879 13 : sPro.proNumber = EPRJ_LATLONG;
1880 13 : sPro.proName = (char*) "Geographic (Lat/Lon)";
1881 : }
1882 : }
1883 :
1884 : /* FIXME/NOTDEF/TODO: Add State Plane */
1885 23 : else if( !bIgnoreUTM && oSRS.GetUTMZone( NULL ) != 0 )
1886 : {
1887 : int bNorth, nZone;
1888 :
1889 20 : nZone = oSRS.GetUTMZone( &bNorth );
1890 20 : sPro.proNumber = EPRJ_UTM;
1891 20 : sPro.proName = (char*) "UTM";
1892 20 : sPro.proZone = nZone;
1893 20 : if( bNorth )
1894 20 : sPro.proParams[3] = 1.0;
1895 : else
1896 0 : sPro.proParams[3] = -1.0;
1897 : }
1898 :
1899 3 : else if( EQUAL(pszProjName,SRS_PT_ALBERS_CONIC_EQUAL_AREA) )
1900 : {
1901 0 : sPro.proNumber = EPRJ_ALBERS_CONIC_EQUAL_AREA;
1902 0 : sPro.proName = (char*) "Albers Conical Equal Area";
1903 0 : sPro.proParams[2] = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_1)*D2R;
1904 0 : sPro.proParams[3] = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_2)*D2R;
1905 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_LONGITUDE_OF_CENTER)*D2R;
1906 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_CENTER)*D2R;
1907 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
1908 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
1909 : }
1910 3 : else if( EQUAL(pszProjName,SRS_PT_LAMBERT_CONFORMAL_CONIC_2SP) )
1911 : {
1912 0 : sPro.proNumber = EPRJ_LAMBERT_CONFORMAL_CONIC;
1913 0 : sPro.proName = (char*) "Lambert Conformal Conic";
1914 0 : sPro.proParams[2] = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_1)*D2R;
1915 0 : sPro.proParams[3] = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_2)*D2R;
1916 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
1917 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN)*D2R;
1918 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
1919 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
1920 : }
1921 3 : else if( EQUAL(pszProjName,SRS_PT_MERCATOR_1SP) )
1922 : {
1923 0 : sPro.proNumber = EPRJ_MERCATOR;
1924 0 : sPro.proName = (char*) "Mercator";
1925 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
1926 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN)*D2R;
1927 : /* hopefully the scale factor is 1.0! */
1928 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
1929 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
1930 : }
1931 3 : else if( EQUAL(pszProjName,SRS_PT_POLAR_STEREOGRAPHIC) )
1932 : {
1933 0 : sPro.proNumber = EPRJ_POLAR_STEREOGRAPHIC;
1934 0 : sPro.proName = (char*) "Polar Stereographic";
1935 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
1936 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN)*D2R;
1937 : /* hopefully the scale factor is 1.0! */
1938 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
1939 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
1940 : }
1941 3 : else if( EQUAL(pszProjName,SRS_PT_POLYCONIC) )
1942 : {
1943 0 : sPro.proNumber = EPRJ_POLYCONIC;
1944 0 : sPro.proName = (char*) "Polyconic";
1945 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
1946 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN)*D2R;
1947 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
1948 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
1949 : }
1950 3 : else if( EQUAL(pszProjName,SRS_PT_EQUIDISTANT_CONIC) )
1951 : {
1952 0 : sPro.proNumber = EPRJ_EQUIDISTANT_CONIC;
1953 0 : sPro.proName = (char*) "Equidistant Conic";
1954 0 : sPro.proParams[2] = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_1)*D2R;
1955 0 : sPro.proParams[3] = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_2)*D2R;
1956 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_LONGITUDE_OF_CENTER)*D2R;
1957 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_CENTER)*D2R;
1958 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
1959 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
1960 0 : sPro.proParams[8] = 1.0;
1961 : }
1962 3 : else if( EQUAL(pszProjName,SRS_PT_TRANSVERSE_MERCATOR) )
1963 : {
1964 0 : sPro.proNumber = EPRJ_TRANSVERSE_MERCATOR;
1965 0 : sPro.proName = (char*) "Transverse Mercator";
1966 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
1967 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN)*D2R;
1968 0 : sPro.proParams[2] = oSRS.GetProjParm(SRS_PP_SCALE_FACTOR,1.0);
1969 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
1970 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
1971 : }
1972 3 : else if( EQUAL(pszProjName,SRS_PT_STEREOGRAPHIC) )
1973 : {
1974 0 : sPro.proNumber = EPRJ_STEREOGRAPHIC_EXTENDED;
1975 0 : sPro.proName = (char*) "Stereographic (Extended)";
1976 0 : sPro.proParams[2] = oSRS.GetProjParm(SRS_PP_SCALE_FACTOR,1.0);
1977 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
1978 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN)*D2R;
1979 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
1980 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
1981 : }
1982 3 : else if( EQUAL(pszProjName,SRS_PT_LAMBERT_AZIMUTHAL_EQUAL_AREA) )
1983 : {
1984 0 : sPro.proNumber = EPRJ_LAMBERT_AZIMUTHAL_EQUAL_AREA;
1985 0 : sPro.proName = (char*) "Lambert Azimuthal Equal-area";
1986 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_LONGITUDE_OF_CENTER)*D2R;
1987 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_CENTER)*D2R;
1988 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
1989 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
1990 : }
1991 3 : else if( EQUAL(pszProjName,SRS_PT_AZIMUTHAL_EQUIDISTANT) )
1992 : {
1993 0 : sPro.proNumber = EPRJ_AZIMUTHAL_EQUIDISTANT;
1994 0 : sPro.proName = (char*) "Azimuthal Equidistant";
1995 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_LONGITUDE_OF_CENTER)*D2R;
1996 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_CENTER)*D2R;
1997 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
1998 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
1999 : }
2000 3 : else if( EQUAL(pszProjName,SRS_PT_GNOMONIC) )
2001 : {
2002 0 : sPro.proNumber = EPRJ_GNOMONIC;
2003 0 : sPro.proName = (char*) "Gnomonic";
2004 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2005 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN)*D2R;
2006 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2007 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2008 : }
2009 3 : else if( EQUAL(pszProjName,SRS_PT_ORTHOGRAPHIC) )
2010 : {
2011 0 : sPro.proNumber = EPRJ_ORTHOGRAPHIC;
2012 0 : sPro.proName = (char*) "Orthographic";
2013 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2014 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN)*D2R;
2015 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2016 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2017 : }
2018 3 : else if( EQUAL(pszProjName,SRS_PT_SINUSOIDAL) )
2019 : {
2020 0 : sPro.proNumber = EPRJ_SINUSOIDAL;
2021 0 : sPro.proName = (char*) "Sinusoidal";
2022 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_LONGITUDE_OF_CENTER)*D2R;
2023 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2024 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2025 : }
2026 3 : else if( EQUAL(pszProjName,SRS_PT_EQUIRECTANGULAR) )
2027 : {
2028 0 : sPro.proNumber = EPRJ_EQUIRECTANGULAR;
2029 0 : sPro.proName = (char*) "Equirectangular";
2030 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2031 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_1)*D2R;
2032 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2033 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2034 : }
2035 3 : else if( EQUAL(pszProjName,SRS_PT_MILLER_CYLINDRICAL) )
2036 : {
2037 0 : sPro.proNumber = EPRJ_MILLER_CYLINDRICAL;
2038 0 : sPro.proName = (char*) "Miller Cylindrical";
2039 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_LONGITUDE_OF_CENTER)*D2R;
2040 : /* hopefully the latitude is zero! */
2041 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2042 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2043 : }
2044 3 : else if( EQUAL(pszProjName,SRS_PT_VANDERGRINTEN) )
2045 : {
2046 0 : sPro.proNumber = EPRJ_VANDERGRINTEN;
2047 0 : sPro.proName = (char*) "Van der Grinten";
2048 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2049 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2050 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2051 : }
2052 3 : else if( EQUAL(pszProjName,SRS_PT_HOTINE_OBLIQUE_MERCATOR) )
2053 : {
2054 0 : sPro.proNumber = EPRJ_HOTINE_OBLIQUE_MERCATOR;
2055 0 : sPro.proName = (char*) "Oblique Mercator (Hotine)";
2056 0 : sPro.proParams[2] = oSRS.GetProjParm(SRS_PP_SCALE_FACTOR,1.0);
2057 0 : sPro.proParams[3] = oSRS.GetProjParm(SRS_PP_AZIMUTH)*D2R;
2058 : /* hopefully the rectified grid angle is zero */
2059 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_LONGITUDE_OF_CENTER)*D2R;
2060 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_CENTER)*D2R;
2061 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2062 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2063 0 : sPro.proParams[12] = 1.0;
2064 : }
2065 3 : else if( EQUAL(pszProjName,SRS_PT_ROBINSON) )
2066 : {
2067 0 : sPro.proNumber = EPRJ_ROBINSON;
2068 0 : sPro.proName = (char*) "Robinson";
2069 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_LONGITUDE_OF_CENTER)*D2R;
2070 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2071 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2072 : }
2073 3 : else if( EQUAL(pszProjName,SRS_PT_MOLLWEIDE) )
2074 : {
2075 0 : sPro.proNumber = EPRJ_MOLLWEIDE;
2076 0 : sPro.proName = (char*) "Mollweide";
2077 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2078 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2079 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2080 : }
2081 3 : else if( EQUAL(pszProjName,SRS_PT_ECKERT_I) )
2082 : {
2083 0 : sPro.proNumber = EPRJ_ECKERT_I;
2084 0 : sPro.proName = (char*) "Eckert I";
2085 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2086 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2087 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2088 : }
2089 3 : else if( EQUAL(pszProjName,SRS_PT_ECKERT_II) )
2090 : {
2091 0 : sPro.proNumber = EPRJ_ECKERT_II;
2092 0 : sPro.proName = (char*) "Eckert II";
2093 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2094 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2095 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2096 : }
2097 3 : else if( EQUAL(pszProjName,SRS_PT_ECKERT_III) )
2098 : {
2099 0 : sPro.proNumber = EPRJ_ECKERT_III;
2100 0 : sPro.proName = (char*) "Eckert III";
2101 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2102 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2103 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2104 : }
2105 3 : else if( EQUAL(pszProjName,SRS_PT_ECKERT_IV) )
2106 : {
2107 0 : sPro.proNumber = EPRJ_ECKERT_IV;
2108 0 : sPro.proName = (char*) "Eckert IV";
2109 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2110 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2111 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2112 : }
2113 3 : else if( EQUAL(pszProjName,SRS_PT_ECKERT_V) )
2114 : {
2115 0 : sPro.proNumber = EPRJ_ECKERT_V;
2116 0 : sPro.proName = (char*) "Eckert V";
2117 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2118 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2119 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2120 : }
2121 3 : else if( EQUAL(pszProjName,SRS_PT_ECKERT_VI) )
2122 : {
2123 0 : sPro.proNumber = EPRJ_ECKERT_VI;
2124 0 : sPro.proName = (char*) "Eckert VI";
2125 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2126 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2127 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2128 : }
2129 3 : else if( EQUAL(pszProjName,SRS_PT_GALL_STEREOGRAPHIC) )
2130 : {
2131 0 : sPro.proNumber = EPRJ_GALL_STEREOGRAPHIC;
2132 0 : sPro.proName = (char*) "Gall Stereographic";
2133 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2134 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2135 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2136 : }
2137 3 : else if( EQUAL(pszProjName,SRS_PT_CASSINI_SOLDNER) )
2138 : {
2139 0 : sPro.proNumber = EPRJ_CASSINI;
2140 0 : sPro.proName = (char*) "Cassini";
2141 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2142 0 : sPro.proParams[5] = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN)*D2R;
2143 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2144 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2145 : }
2146 3 : else if( EQUAL(pszProjName,SRS_PT_BONNE) )
2147 : {
2148 0 : sPro.proNumber = EPRJ_BONNE;
2149 0 : sPro.proName = (char*) "Bonne";
2150 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2151 0 : sPro.proParams[2] = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_1)*D2R;
2152 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2153 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2154 : }
2155 3 : else if( EQUAL(pszProjName,"Loximuthal") )
2156 : {
2157 0 : sPro.proNumber = EPRJ_LOXIMUTHAL;
2158 0 : sPro.proName = (char*) "Loximuthal";
2159 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2160 0 : sPro.proParams[5] = oSRS.GetProjParm("central_parallel")*D2R;
2161 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2162 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2163 : }
2164 3 : else if( EQUAL(pszProjName,"Quartic_Authalic") )
2165 : {
2166 0 : sPro.proNumber = EPRJ_QUARTIC_AUTHALIC;
2167 0 : sPro.proName = (char*) "Quartic Authalic";
2168 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2169 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2170 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2171 : }
2172 3 : else if( EQUAL(pszProjName,"Winkel_I") )
2173 : {
2174 0 : sPro.proNumber = EPRJ_WINKEL_I;
2175 0 : sPro.proName = (char*) "Winkel I";
2176 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2177 0 : sPro.proParams[2] = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_1)*D2R;
2178 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2179 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2180 : }
2181 3 : else if( EQUAL(pszProjName,"Winkel_II") )
2182 : {
2183 0 : sPro.proNumber = EPRJ_WINKEL_II;
2184 0 : sPro.proName = (char*) "Winkel II";
2185 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2186 0 : sPro.proParams[2] = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_1)*D2R;
2187 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2188 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2189 : }
2190 3 : else if( EQUAL(pszProjName,"Winkel_II") )
2191 : {
2192 0 : sPro.proNumber = EPRJ_WINKEL_II;
2193 0 : sPro.proName = (char*) "Winkel II";
2194 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2195 0 : sPro.proParams[2] = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_1)*D2R;
2196 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2197 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2198 : }
2199 3 : else if( EQUAL(pszProjName,"Behrmann") )
2200 : {
2201 0 : sPro.proNumber = EPRJ_BEHRMANN;
2202 0 : sPro.proName = (char*) "Behrmann";
2203 0 : sPro.proParams[4] = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN)*D2R;
2204 0 : sPro.proParams[6] = oSRS.GetProjParm(SRS_PP_FALSE_EASTING);
2205 0 : sPro.proParams[7] = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2206 : }
2207 : // Anything we can't map, we store as an ESRI PE_STRING
2208 3 : else if( oSRS.IsProjected() || oSRS.IsGeographic() )
2209 : {
2210 3 : if(!peStrStored)
2211 : {
2212 0 : char *pszPEString = NULL;
2213 0 : oSRS.morphToESRI();
2214 0 : oSRS.exportToWkt( &pszPEString );
2215 : // need to transform this into ESRI format.
2216 0 : HFASetPEString( hHFA, pszPEString );
2217 0 : CPLFree( pszPEString );
2218 0 : peStrStored = TRUE;
2219 : }
2220 : }
2221 : else
2222 : {
2223 : CPLError( CE_Warning, CPLE_NotSupported,
2224 : "Projection %s not supported for translation to Imagine.",
2225 0 : pszProjName );
2226 : }
2227 :
2228 : /* -------------------------------------------------------------------- */
2229 : /* MapInfo */
2230 : /* -------------------------------------------------------------------- */
2231 47 : const char *pszPROJCS = oSRS.GetAttrValue( "PROJCS" );
2232 :
2233 47 : if( pszPROJCS )
2234 23 : sMapInfo.proName = (char *) pszPROJCS;
2235 37 : else if( bHaveSRS && sPro.proName != NULL )
2236 13 : sMapInfo.proName = sPro.proName;
2237 : else
2238 11 : sMapInfo.proName = (char*) "Unknown";
2239 :
2240 : sMapInfo.upperLeftCenter.x =
2241 47 : adfGeoTransform[0] + adfGeoTransform[1]*0.5;
2242 : sMapInfo.upperLeftCenter.y =
2243 47 : adfGeoTransform[3] + adfGeoTransform[5]*0.5;
2244 :
2245 : sMapInfo.lowerRightCenter.x =
2246 47 : adfGeoTransform[0] + adfGeoTransform[1] * (GetRasterXSize()-0.5);
2247 : sMapInfo.lowerRightCenter.y =
2248 47 : adfGeoTransform[3] + adfGeoTransform[5] * (GetRasterYSize()-0.5);
2249 :
2250 47 : sMapInfo.pixelSize.width = ABS(adfGeoTransform[1]);
2251 47 : sMapInfo.pixelSize.height = ABS(adfGeoTransform[5]);
2252 :
2253 : /* -------------------------------------------------------------------- */
2254 : /* Handle units. Try to match up with a known name. */
2255 : /* -------------------------------------------------------------------- */
2256 47 : sMapInfo.units = (char*) "meters";
2257 :
2258 47 : if( bHaveSRS && oSRS.IsGeographic() )
2259 13 : sMapInfo.units = (char*) "dd";
2260 34 : else if( bHaveSRS && oSRS.GetLinearUnits() != 1.0 )
2261 : {
2262 2 : double dfClosestDiff = 100.0;
2263 2 : int iClosest=-1, iUnit;
2264 2 : char *pszUnitName = NULL;
2265 2 : double dfActualSize = oSRS.GetLinearUnits( &pszUnitName );
2266 :
2267 70 : for( iUnit = 0; apszUnitMap[iUnit] != NULL; iUnit += 2 )
2268 : {
2269 68 : if( fabs(atof(apszUnitMap[iUnit+1]) - dfActualSize) < dfClosestDiff )
2270 : {
2271 6 : iClosest = iUnit;
2272 6 : dfClosestDiff = fabs(atof(apszUnitMap[iUnit+1])-dfActualSize);
2273 : }
2274 : }
2275 :
2276 2 : if( iClosest == -1 || fabs(dfClosestDiff/dfActualSize) > 0.0001 )
2277 : {
2278 : CPLError( CE_Warning, CPLE_NotSupported,
2279 : "Unable to identify Erdas units matching %s/%gm,\n"
2280 : "output units will be wrong.",
2281 0 : pszUnitName, dfActualSize );
2282 : }
2283 : else
2284 2 : sMapInfo.units = (char *) apszUnitMap[iClosest];
2285 :
2286 : /* We need to convert false easting and northing to meters. */
2287 2 : sPro.proParams[6] *= dfActualSize;
2288 2 : sPro.proParams[7] *= dfActualSize;
2289 : }
2290 :
2291 : /* -------------------------------------------------------------------- */
2292 : /* Write out definitions. */
2293 : /* -------------------------------------------------------------------- */
2294 93 : if( adfGeoTransform[2] == 0.0 && adfGeoTransform[4] == 0.0 )
2295 : {
2296 46 : HFASetMapInfo( hHFA, &sMapInfo );
2297 : }
2298 : else
2299 : {
2300 : HFASetGeoTransform( hHFA,
2301 : sMapInfo.proName, sMapInfo.units,
2302 1 : adfGeoTransform );
2303 : }
2304 :
2305 80 : if( bHaveSRS && sPro.proName != NULL)
2306 : {
2307 33 : HFASetProParameters( hHFA, &sPro );
2308 33 : HFASetDatum( hHFA, &sDatum );
2309 : }
2310 14 : else if( !peStrStored )
2311 11 : ClearSR(hHFA);
2312 :
2313 : /* -------------------------------------------------------------------- */
2314 : /* Cleanup */
2315 : /* -------------------------------------------------------------------- */
2316 47 : if( poGeogSRS != NULL )
2317 36 : delete poGeogSRS;
2318 :
2319 47 : return CE_None;
2320 : }
2321 :
2322 : /************************************************************************/
2323 : /* WritePeStringIfNeeded() */
2324 : /************************************************************************/
2325 36 : OGRBoolean WritePeStringIfNeeded(OGRSpatialReference* poSRS, HFAHandle hHFA)
2326 : {
2327 36 : OGRBoolean ret = FALSE;
2328 36 : if(!poSRS || !hHFA)
2329 0 : return ret;
2330 :
2331 36 : const char *pszGEOGCS = poSRS->GetAttrValue( "GEOGCS" );
2332 36 : const char *pszDatum = poSRS->GetAttrValue( "DATUM" );
2333 36 : int gcsNameOffset = 0;
2334 36 : int datumNameOffset = 0;
2335 36 : if(strstr(pszGEOGCS, "GCS_"))
2336 3 : gcsNameOffset = strlen("GCS_");
2337 36 : if(strstr(pszDatum, "D_"))
2338 0 : datumNameOffset = strlen("D_");
2339 :
2340 36 : if(!EQUAL(pszGEOGCS+gcsNameOffset, pszDatum+datumNameOffset))
2341 35 : ret = TRUE;
2342 : else
2343 : {
2344 1 : const char* name = poSRS->GetAttrValue("PRIMEM");
2345 1 : if(name && !EQUAL(name,"Greenwich"))
2346 0 : ret = TRUE;
2347 1 : if(!ret)
2348 : {
2349 1 : OGR_SRSNode * poAUnits = poSRS->GetAttrNode( "GEOGCS|UNIT" );
2350 1 : name = poAUnits->GetChild(0)->GetValue();
2351 1 : if(name && !EQUAL(name,"Degree"))
2352 0 : ret = TRUE;
2353 : }
2354 1 : if(!ret)
2355 : {
2356 1 : name = poSRS->GetAttrValue("UNIT");
2357 1 : if(name)
2358 : {
2359 1 : ret = TRUE;
2360 35 : for(int i=0; apszUnitMap[i] != NULL; i+=2)
2361 34 : if(EQUAL(name, apszUnitMap[i]))
2362 1 : ret = FALSE;
2363 : }
2364 : }
2365 1 : if(!ret)
2366 : {
2367 1 : int nGCS = poSRS->GetEPSGGeogCS();
2368 1 : switch(nGCS)
2369 : {
2370 : case 4326:
2371 1 : if(!EQUAL(pszDatum+datumNameOffset, "WGS_84"))
2372 1 : ret = TRUE;
2373 1 : break;
2374 : case 4322:
2375 0 : if(!EQUAL(pszDatum+datumNameOffset, "WGS_72"))
2376 0 : ret = TRUE;
2377 0 : break;
2378 : case 4267:
2379 0 : if(!EQUAL(pszDatum+datumNameOffset, "North_America_1927"))
2380 0 : ret = TRUE;
2381 0 : break;
2382 : case 4269:
2383 0 : if(!EQUAL(pszDatum+datumNameOffset, "North_America_1983"))
2384 0 : ret = TRUE;
2385 : break;
2386 : }
2387 : }
2388 : }
2389 36 : if(ret)
2390 : {
2391 36 : char *pszPEString = NULL;
2392 36 : poSRS->morphToESRI();
2393 36 : poSRS->exportToWkt( &pszPEString );
2394 36 : HFASetPEString( hHFA, pszPEString );
2395 36 : CPLFree( pszPEString );
2396 : }
2397 :
2398 36 : return ret;
2399 : }
2400 :
2401 : /************************************************************************/
2402 : /* ClearSR() */
2403 : /************************************************************************/
2404 11 : void ClearSR(HFAHandle hHFA)
2405 : {
2406 22 : for( int iBand = 0; iBand < hHFA->nBands; iBand++ )
2407 : {
2408 : HFAEntry *poMIEntry;
2409 11 : if( hHFA->papoBand[iBand]->poNode && (poMIEntry = hHFA->papoBand[iBand]->poNode->GetNamedChild("Projection")) != NULL )
2410 : {
2411 0 : poMIEntry->MarkDirty();
2412 0 : poMIEntry->SetIntField( "proType", 0 );
2413 0 : poMIEntry->SetIntField( "proNumber", 0 );
2414 0 : poMIEntry->SetStringField( "proExeName", "");
2415 0 : poMIEntry->SetStringField( "proName", "");
2416 0 : poMIEntry->SetIntField( "proZone", 0 );
2417 0 : poMIEntry->SetDoubleField( "proParams[0]", 0.0 );
2418 0 : poMIEntry->SetDoubleField( "proParams[1]", 0.0 );
2419 0 : poMIEntry->SetDoubleField( "proParams[2]", 0.0 );
2420 0 : poMIEntry->SetDoubleField( "proParams[3]", 0.0 );
2421 0 : poMIEntry->SetDoubleField( "proParams[4]", 0.0 );
2422 0 : poMIEntry->SetDoubleField( "proParams[5]", 0.0 );
2423 0 : poMIEntry->SetDoubleField( "proParams[6]", 0.0 );
2424 0 : poMIEntry->SetDoubleField( "proParams[7]", 0.0 );
2425 0 : poMIEntry->SetDoubleField( "proParams[8]", 0.0 );
2426 0 : poMIEntry->SetDoubleField( "proParams[9]", 0.0 );
2427 0 : poMIEntry->SetDoubleField( "proParams[10]", 0.0 );
2428 0 : poMIEntry->SetDoubleField( "proParams[11]", 0.0 );
2429 0 : poMIEntry->SetDoubleField( "proParams[12]", 0.0 );
2430 0 : poMIEntry->SetDoubleField( "proParams[13]", 0.0 );
2431 0 : poMIEntry->SetDoubleField( "proParams[14]", 0.0 );
2432 0 : poMIEntry->SetStringField( "proSpheroid.sphereName", "" );
2433 0 : poMIEntry->SetDoubleField( "proSpheroid.a", 0.0 );
2434 0 : poMIEntry->SetDoubleField( "proSpheroid.b", 0.0 );
2435 0 : poMIEntry->SetDoubleField( "proSpheroid.eSquared", 0.0 );
2436 0 : poMIEntry->SetDoubleField( "proSpheroid.radius", 0.0 );
2437 0 : HFAEntry* poDatumEntry = poMIEntry->GetNamedChild("Datum");
2438 0 : if( poDatumEntry != NULL )
2439 : {
2440 0 : poDatumEntry->MarkDirty();
2441 0 : poDatumEntry->SetStringField( "datumname", "" );
2442 0 : poDatumEntry->SetIntField( "type", 0 );
2443 0 : poDatumEntry->SetDoubleField( "params[0]", 0.0 );
2444 0 : poDatumEntry->SetDoubleField( "params[1]", 0.0 );
2445 0 : poDatumEntry->SetDoubleField( "params[2]", 0.0 );
2446 0 : poDatumEntry->SetDoubleField( "params[3]", 0.0 );
2447 0 : poDatumEntry->SetDoubleField( "params[4]", 0.0 );
2448 0 : poDatumEntry->SetDoubleField( "params[5]", 0.0 );
2449 0 : poDatumEntry->SetDoubleField( "params[6]", 0.0 );
2450 0 : poDatumEntry->SetStringField( "gridname", "" );
2451 : }
2452 0 : poMIEntry->FlushToDisk();
2453 0 : char* peStr = HFAGetPEString( hHFA );
2454 0 : if( peStr != NULL && strlen(peStr) > 0 )
2455 0 : HFASetPEString( hHFA, "" );
2456 : }
2457 : }
2458 : return;
2459 : }
2460 :
2461 : /************************************************************************/
2462 : /* ESRIToUSGSZone() */
2463 : /* */
2464 : /* Convert ESRI style state plane zones to USGS style state */
2465 : /* plane zones. */
2466 : /************************************************************************/
2467 :
2468 3 : static int ESRIToUSGSZone( int nESRIZone )
2469 :
2470 : {
2471 3 : int nPairs = sizeof(anUsgsEsriZones) / (2*sizeof(int));
2472 : int i;
2473 :
2474 3 : if( nESRIZone < 0 )
2475 1 : return ABS(nESRIZone);
2476 :
2477 212 : for( i = 0; i < nPairs; i++ )
2478 : {
2479 212 : if( anUsgsEsriZones[i*2+1] == nESRIZone )
2480 2 : return anUsgsEsriZones[i*2];
2481 : }
2482 :
2483 0 : return 0;
2484 : }
2485 :
2486 : /************************************************************************/
2487 : /* PCSStructToWKT() */
2488 : /* */
2489 : /* Convert the datum, proparameters and mapinfo structures into */
2490 : /* WKT format. */
2491 : /************************************************************************/
2492 :
2493 : char *
2494 32 : HFAPCSStructToWKT( const Eprj_Datum *psDatum,
2495 : const Eprj_ProParameters *psPro,
2496 : const Eprj_MapInfo *psMapInfo,
2497 : HFAEntry *poMapInformation )
2498 :
2499 : {
2500 32 : OGRSpatialReference oSRS;
2501 32 : char *pszNewProj = NULL;
2502 :
2503 : /* -------------------------------------------------------------------- */
2504 : /* General case for Erdas style projections. */
2505 : /* */
2506 : /* We make a particular effort to adapt the mapinfo->proname as */
2507 : /* the PROJCS[] name per #2422. */
2508 : /* -------------------------------------------------------------------- */
2509 :
2510 32 : if( psPro == NULL && psMapInfo != NULL )
2511 : {
2512 0 : oSRS.SetLocalCS( psMapInfo->proName );
2513 : }
2514 :
2515 32 : else if( psPro == NULL )
2516 : {
2517 0 : return NULL;
2518 : }
2519 :
2520 32 : else if( psPro->proType == EPRJ_EXTERNAL )
2521 : {
2522 0 : oSRS.SetLocalCS( psPro->proName );
2523 : }
2524 :
2525 61 : else if( psPro->proNumber != EPRJ_LATLONG
2526 : && psMapInfo != NULL )
2527 : {
2528 29 : oSRS.SetProjCS( psMapInfo->proName );
2529 : }
2530 3 : else if( psPro->proNumber != EPRJ_LATLONG )
2531 : {
2532 3 : oSRS.SetProjCS( psPro->proName );
2533 : }
2534 :
2535 : /* -------------------------------------------------------------------- */
2536 : /* Handle units. It is important to deal with this first so */
2537 : /* that the projection Set methods will automatically do */
2538 : /* translation of linear values (like false easting) to PROJCS */
2539 : /* units from meters. Erdas linear projection values are */
2540 : /* always in meters. */
2541 : /* -------------------------------------------------------------------- */
2542 32 : int iUnitIndex = 0;
2543 :
2544 32 : if( oSRS.IsProjected() || oSRS.IsLocal() )
2545 : {
2546 32 : const char *pszUnits = NULL;
2547 :
2548 32 : if( psMapInfo )
2549 29 : pszUnits = psMapInfo->units;
2550 3 : else if( poMapInformation != NULL )
2551 3 : pszUnits = poMapInformation->GetStringField( "units.string" );
2552 :
2553 32 : if( pszUnits != NULL )
2554 : {
2555 136 : for( iUnitIndex = 0;
2556 68 : apszUnitMap[iUnitIndex] != NULL;
2557 : iUnitIndex += 2 )
2558 : {
2559 68 : if( EQUAL(apszUnitMap[iUnitIndex], pszUnits ) )
2560 32 : break;
2561 : }
2562 :
2563 32 : if( apszUnitMap[iUnitIndex] == NULL )
2564 0 : iUnitIndex = 0;
2565 :
2566 : oSRS.SetLinearUnits( pszUnits,
2567 32 : atof(apszUnitMap[iUnitIndex+1]) );
2568 : }
2569 : else
2570 0 : oSRS.SetLinearUnits( SRS_UL_METER, 1.0 );
2571 : }
2572 :
2573 32 : if( psPro == NULL )
2574 : {
2575 0 : if( oSRS.IsLocal() )
2576 : {
2577 0 : if( oSRS.exportToWkt( &pszNewProj ) == OGRERR_NONE )
2578 0 : return pszNewProj;
2579 : else
2580 : {
2581 0 : pszNewProj = NULL;
2582 0 : return NULL;
2583 : }
2584 : }
2585 : else
2586 0 : return NULL;
2587 : }
2588 :
2589 : /* -------------------------------------------------------------------- */
2590 : /* Try to work out ellipsoid and datum information. */
2591 : /* -------------------------------------------------------------------- */
2592 32 : const char *pszDatumName = psPro->proSpheroid.sphereName;
2593 32 : const char *pszEllipsoidName = psPro->proSpheroid.sphereName;
2594 : double dfInvFlattening;
2595 :
2596 32 : if( psDatum != NULL )
2597 : {
2598 : int i;
2599 :
2600 32 : pszDatumName = psDatum->datumname;
2601 :
2602 : /* Imagine to WKT translation */
2603 67 : for( i = 0; apszDatumMap[i] != NULL; i += 2 )
2604 : {
2605 59 : if( EQUAL(pszDatumName,apszDatumMap[i]) )
2606 : {
2607 24 : pszDatumName = apszDatumMap[i+1];
2608 24 : break;
2609 : }
2610 : }
2611 : }
2612 :
2613 32 : if( psPro->proSpheroid.a == 0.0 )
2614 0 : ((Eprj_ProParameters *) psPro)->proSpheroid.a = 6378137.0;
2615 32 : if( psPro->proSpheroid.b == 0.0 )
2616 0 : ((Eprj_ProParameters *) psPro)->proSpheroid.b = 6356752.3;
2617 :
2618 32 : if( fabs(psPro->proSpheroid.b - psPro->proSpheroid.a) < 0.001 )
2619 0 : dfInvFlattening = 0.0; /* special value for sphere. */
2620 : else
2621 32 : dfInvFlattening = 1.0/(1.0-psPro->proSpheroid.b/psPro->proSpheroid.a);
2622 :
2623 : /* -------------------------------------------------------------------- */
2624 : /* Handle different projection methods. */
2625 : /* -------------------------------------------------------------------- */
2626 32 : switch( psPro->proNumber )
2627 : {
2628 : case EPRJ_LATLONG:
2629 0 : break;
2630 :
2631 : case EPRJ_UTM:
2632 : // We change this to unnamed so that SetUTM will set the long
2633 : // UTM description.
2634 21 : oSRS.SetProjCS( "unnamed" );
2635 21 : oSRS.SetUTM( psPro->proZone, psPro->proParams[3] >= 0.0 );
2636 21 : break;
2637 :
2638 : case EPRJ_STATE_PLANE:
2639 : {
2640 3 : char *pszUnitsName = NULL;
2641 3 : double dfLinearUnits = oSRS.GetLinearUnits( &pszUnitsName );
2642 :
2643 3 : pszUnitsName = CPLStrdup( pszUnitsName );
2644 :
2645 : /* Set state plane zone. Set NAD83/27 on basis of spheroid */
2646 : oSRS.SetStatePlane( ESRIToUSGSZone(psPro->proZone),
2647 : fabs(psPro->proSpheroid.a - 6378137.0)< 1.0,
2648 3 : pszUnitsName, dfLinearUnits );
2649 :
2650 3 : CPLFree( pszUnitsName );
2651 : }
2652 3 : break;
2653 :
2654 : case EPRJ_ALBERS_CONIC_EQUAL_AREA:
2655 0 : oSRS.SetACEA( psPro->proParams[2]*R2D, psPro->proParams[3]*R2D,
2656 0 : psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2657 0 : psPro->proParams[6], psPro->proParams[7] );
2658 0 : break;
2659 :
2660 : case EPRJ_LAMBERT_CONFORMAL_CONIC:
2661 0 : oSRS.SetLCC( psPro->proParams[2]*R2D, psPro->proParams[3]*R2D,
2662 0 : psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2663 0 : psPro->proParams[6], psPro->proParams[7] );
2664 0 : break;
2665 :
2666 : case EPRJ_MERCATOR:
2667 0 : oSRS.SetMercator( psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2668 : 1.0,
2669 0 : psPro->proParams[6], psPro->proParams[7] );
2670 0 : break;
2671 :
2672 : case EPRJ_POLAR_STEREOGRAPHIC:
2673 0 : oSRS.SetPS( psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2674 : 1.0,
2675 0 : psPro->proParams[6], psPro->proParams[7] );
2676 0 : break;
2677 :
2678 : case EPRJ_POLYCONIC:
2679 0 : oSRS.SetPolyconic( psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2680 0 : psPro->proParams[6], psPro->proParams[7] );
2681 0 : break;
2682 :
2683 : case EPRJ_EQUIDISTANT_CONIC:
2684 : double dfStdParallel2;
2685 :
2686 0 : if( psPro->proParams[8] != 0.0 )
2687 0 : dfStdParallel2 = psPro->proParams[3]*R2D;
2688 : else
2689 0 : dfStdParallel2 = psPro->proParams[2]*R2D;
2690 0 : oSRS.SetEC( psPro->proParams[2]*R2D, dfStdParallel2,
2691 0 : psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2692 0 : psPro->proParams[6], psPro->proParams[7] );
2693 0 : break;
2694 :
2695 : case EPRJ_TRANSVERSE_MERCATOR:
2696 16 : oSRS.SetTM( psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2697 : psPro->proParams[2],
2698 24 : psPro->proParams[6], psPro->proParams[7] );
2699 8 : break;
2700 :
2701 : case EPRJ_STEREOGRAPHIC:
2702 0 : oSRS.SetStereographic( psPro->proParams[5]*R2D,psPro->proParams[4]*R2D,
2703 : 1.0,
2704 0 : psPro->proParams[6], psPro->proParams[7] );
2705 0 : break;
2706 :
2707 : case EPRJ_LAMBERT_AZIMUTHAL_EQUAL_AREA:
2708 0 : oSRS.SetLAEA( psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2709 0 : psPro->proParams[6], psPro->proParams[7] );
2710 0 : break;
2711 :
2712 : case EPRJ_AZIMUTHAL_EQUIDISTANT:
2713 0 : oSRS.SetAE( psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2714 0 : psPro->proParams[6], psPro->proParams[7] );
2715 0 : break;
2716 :
2717 : case EPRJ_GNOMONIC:
2718 0 : oSRS.SetGnomonic( psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2719 0 : psPro->proParams[6], psPro->proParams[7] );
2720 0 : break;
2721 :
2722 : case EPRJ_ORTHOGRAPHIC:
2723 0 : oSRS.SetOrthographic( psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2724 0 : psPro->proParams[6], psPro->proParams[7] );
2725 0 : break;
2726 :
2727 : case EPRJ_SINUSOIDAL:
2728 0 : oSRS.SetSinusoidal( psPro->proParams[4]*R2D,
2729 0 : psPro->proParams[6], psPro->proParams[7] );
2730 0 : break;
2731 :
2732 : case EPRJ_PLATE_CARREE:
2733 : case EPRJ_EQUIRECTANGULAR:
2734 : oSRS.SetEquirectangular2( 0.0,
2735 0 : psPro->proParams[4]*R2D,
2736 0 : psPro->proParams[5]*R2D,
2737 0 : psPro->proParams[6], psPro->proParams[7] );
2738 0 : break;
2739 :
2740 : case EPRJ_EQUIDISTANT_CYLINDRICAL:
2741 : oSRS.SetEquirectangular2( 0.0,
2742 0 : psPro->proParams[4]*R2D,
2743 0 : psPro->proParams[2]*R2D,
2744 0 : psPro->proParams[6], psPro->proParams[7] );
2745 0 : break;
2746 :
2747 : case EPRJ_MILLER_CYLINDRICAL:
2748 0 : oSRS.SetMC( 0.0, psPro->proParams[4]*R2D,
2749 0 : psPro->proParams[6], psPro->proParams[7] );
2750 0 : break;
2751 :
2752 : case EPRJ_VANDERGRINTEN:
2753 0 : oSRS.SetVDG( psPro->proParams[4]*R2D,
2754 0 : psPro->proParams[6], psPro->proParams[7] );
2755 0 : break;
2756 :
2757 : case EPRJ_HOTINE_OBLIQUE_MERCATOR:
2758 0 : if( psPro->proParams[12] > 0.0 )
2759 0 : oSRS.SetHOM( psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2760 0 : psPro->proParams[3]*R2D, 0.0,
2761 : psPro->proParams[2],
2762 0 : psPro->proParams[6], psPro->proParams[7] );
2763 0 : break;
2764 :
2765 : case EPRJ_ROBINSON:
2766 0 : oSRS.SetRobinson( psPro->proParams[4]*R2D,
2767 0 : psPro->proParams[6], psPro->proParams[7] );
2768 0 : break;
2769 :
2770 : case EPRJ_MOLLWEIDE:
2771 0 : oSRS.SetMollweide( psPro->proParams[4]*R2D,
2772 0 : psPro->proParams[6], psPro->proParams[7] );
2773 0 : break;
2774 :
2775 : case EPRJ_GALL_STEREOGRAPHIC:
2776 0 : oSRS.SetGS( psPro->proParams[4]*R2D,
2777 0 : psPro->proParams[6], psPro->proParams[7] );
2778 0 : break;
2779 :
2780 : case EPRJ_ECKERT_I:
2781 0 : oSRS.SetEckert( 1, psPro->proParams[4]*R2D,
2782 0 : psPro->proParams[6], psPro->proParams[7] );
2783 0 : break;
2784 :
2785 : case EPRJ_ECKERT_II:
2786 0 : oSRS.SetEckert( 2, psPro->proParams[4]*R2D,
2787 0 : psPro->proParams[6], psPro->proParams[7] );
2788 0 : break;
2789 :
2790 : case EPRJ_ECKERT_III:
2791 0 : oSRS.SetEckert( 3, psPro->proParams[4]*R2D,
2792 0 : psPro->proParams[6], psPro->proParams[7] );
2793 0 : break;
2794 :
2795 : case EPRJ_ECKERT_IV:
2796 0 : oSRS.SetEckert( 4, psPro->proParams[4]*R2D,
2797 0 : psPro->proParams[6], psPro->proParams[7] );
2798 0 : break;
2799 :
2800 : case EPRJ_ECKERT_V:
2801 0 : oSRS.SetEckert( 5, psPro->proParams[4]*R2D,
2802 0 : psPro->proParams[6], psPro->proParams[7] );
2803 0 : break;
2804 :
2805 : case EPRJ_ECKERT_VI:
2806 0 : oSRS.SetEckert( 6, psPro->proParams[4]*R2D,
2807 0 : psPro->proParams[6], psPro->proParams[7] );
2808 0 : break;
2809 :
2810 : case EPRJ_CASSINI:
2811 0 : oSRS.SetCS( psPro->proParams[5]*R2D, psPro->proParams[4]*R2D,
2812 0 : psPro->proParams[6], psPro->proParams[7] );
2813 0 : break;
2814 :
2815 : case EPRJ_STEREOGRAPHIC_EXTENDED:
2816 0 : oSRS.SetStereographic( psPro->proParams[5]*R2D,psPro->proParams[4]*R2D,
2817 : psPro->proParams[2],
2818 0 : psPro->proParams[6], psPro->proParams[7] );
2819 0 : break;
2820 :
2821 : case EPRJ_BONNE:
2822 0 : oSRS.SetBonne( psPro->proParams[2]*R2D, psPro->proParams[4]*R2D,
2823 0 : psPro->proParams[6], psPro->proParams[7] );
2824 0 : break;
2825 :
2826 : case EPRJ_LOXIMUTHAL:
2827 : {
2828 0 : oSRS.SetProjection( "Loximuthal" );
2829 : oSRS.SetNormProjParm( SRS_PP_CENTRAL_MERIDIAN,
2830 0 : psPro->proParams[4] * R2D );
2831 : oSRS.SetNormProjParm( "central_parallel",
2832 0 : psPro->proParams[5] * R2D );
2833 0 : oSRS.SetNormProjParm( SRS_PP_FALSE_EASTING, psPro->proParams[6] );
2834 0 : oSRS.SetNormProjParm( SRS_PP_FALSE_NORTHING, psPro->proParams[6] );
2835 : }
2836 0 : break;
2837 :
2838 : case EPRJ_QUARTIC_AUTHALIC:
2839 : {
2840 0 : oSRS.SetProjection( "Quartic_Authalic" );
2841 : oSRS.SetNormProjParm( SRS_PP_CENTRAL_MERIDIAN,
2842 0 : psPro->proParams[4] * R2D );
2843 0 : oSRS.SetNormProjParm( SRS_PP_FALSE_EASTING, psPro->proParams[6] );
2844 0 : oSRS.SetNormProjParm( SRS_PP_FALSE_NORTHING, psPro->proParams[6] );
2845 : }
2846 0 : break;
2847 :
2848 : case EPRJ_WINKEL_I:
2849 : {
2850 0 : oSRS.SetProjection( "Winkel_I" );
2851 : oSRS.SetNormProjParm( SRS_PP_CENTRAL_MERIDIAN,
2852 0 : psPro->proParams[4] * R2D );
2853 : oSRS.SetNormProjParm( SRS_PP_STANDARD_PARALLEL_1,
2854 0 : psPro->proParams[2] * R2D );
2855 0 : oSRS.SetNormProjParm( SRS_PP_FALSE_EASTING, psPro->proParams[6] );
2856 0 : oSRS.SetNormProjParm( SRS_PP_FALSE_NORTHING, psPro->proParams[6] );
2857 : }
2858 0 : break;
2859 :
2860 : case EPRJ_WINKEL_II:
2861 : {
2862 0 : oSRS.SetProjection( "Winkel_II" );
2863 : oSRS.SetNormProjParm( SRS_PP_CENTRAL_MERIDIAN,
2864 0 : psPro->proParams[4] * R2D );
2865 : oSRS.SetNormProjParm( SRS_PP_STANDARD_PARALLEL_1,
2866 0 : psPro->proParams[2] * R2D );
2867 0 : oSRS.SetNormProjParm( SRS_PP_FALSE_EASTING, psPro->proParams[6] );
2868 0 : oSRS.SetNormProjParm( SRS_PP_FALSE_NORTHING, psPro->proParams[6] );
2869 : }
2870 0 : break;
2871 :
2872 : case EPRJ_BEHRMANN:
2873 : {
2874 0 : oSRS.SetProjection( "Behrmann" );
2875 : oSRS.SetNormProjParm( SRS_PP_CENTRAL_MERIDIAN,
2876 0 : psPro->proParams[4] * R2D );
2877 0 : oSRS.SetNormProjParm( SRS_PP_FALSE_EASTING, psPro->proParams[6] );
2878 0 : oSRS.SetNormProjParm( SRS_PP_FALSE_NORTHING, psPro->proParams[6] );
2879 : }
2880 0 : break;
2881 :
2882 : default:
2883 0 : if( oSRS.IsProjected() )
2884 0 : oSRS.GetRoot()->SetValue( "LOCAL_CS" );
2885 : else
2886 0 : oSRS.SetLocalCS( psPro->proName );
2887 : break;
2888 : }
2889 :
2890 : /* -------------------------------------------------------------------- */
2891 : /* Try and set the GeogCS information. */
2892 : /* -------------------------------------------------------------------- */
2893 32 : if( oSRS.GetAttrNode("GEOGCS") == NULL
2894 : && oSRS.GetAttrNode("LOCAL_CS") == NULL )
2895 : {
2896 29 : if( EQUAL(pszDatumName,"WGS 84")
2897 : || EQUAL(pszDatumName,"WGS_1984") )
2898 0 : oSRS.SetWellKnownGeogCS( "WGS84" );
2899 50 : else if( strstr(pszDatumName,"NAD27") != NULL
2900 : || EQUAL(pszDatumName,"North_American_Datum_1927") )
2901 21 : oSRS.SetWellKnownGeogCS( "NAD27" );
2902 8 : else if( strstr(pszDatumName,"NAD83") != NULL
2903 : || EQUAL(pszDatumName,"North_American_Datum_1983"))
2904 0 : oSRS.SetWellKnownGeogCS( "NAD83" );
2905 : else
2906 : oSRS.SetGeogCS( pszDatumName, pszDatumName, pszEllipsoidName,
2907 8 : psPro->proSpheroid.a, dfInvFlattening );
2908 :
2909 29 : if( psDatum != NULL && psDatum->type == EPRJ_DATUM_PARAMETRIC )
2910 : {
2911 : oSRS.SetTOWGS84( psDatum->params[0],
2912 : psDatum->params[1],
2913 : psDatum->params[2],
2914 : psDatum->params[3],
2915 : psDatum->params[4],
2916 : psDatum->params[5],
2917 8 : psDatum->params[6] );
2918 : }
2919 : }
2920 :
2921 : /* -------------------------------------------------------------------- */
2922 : /* Try to insert authority information if possible. Fixup any */
2923 : /* ordering oddities. */
2924 : /* -------------------------------------------------------------------- */
2925 32 : oSRS.AutoIdentifyEPSG();
2926 32 : oSRS.Fixup();
2927 :
2928 : /* -------------------------------------------------------------------- */
2929 : /* Get the WKT representation of the coordinate system. */
2930 : /* -------------------------------------------------------------------- */
2931 32 : if( oSRS.exportToWkt( &pszNewProj ) == OGRERR_NONE )
2932 32 : return pszNewProj;
2933 : else
2934 : {
2935 0 : return NULL;
2936 0 : }
2937 : }
2938 :
2939 : /************************************************************************/
2940 : /* ReadProjection() */
2941 : /************************************************************************/
2942 :
2943 275 : CPLErr HFADataset::ReadProjection()
2944 :
2945 : {
2946 : const Eprj_Datum *psDatum;
2947 : const Eprj_ProParameters *psPro;
2948 : const Eprj_MapInfo *psMapInfo;
2949 275 : OGRSpatialReference oSRS;
2950 : char *pszPE_COORDSYS;
2951 :
2952 : /* -------------------------------------------------------------------- */
2953 : /* Special logic for PE string in ProjectionX node. */
2954 : /* -------------------------------------------------------------------- */
2955 275 : pszPE_COORDSYS = HFAGetPEString( hHFA );
2956 275 : if( pszPE_COORDSYS != NULL
2957 : && oSRS.SetFromUserInput( pszPE_COORDSYS ) == OGRERR_NONE )
2958 : {
2959 46 : CPLFree( pszPE_COORDSYS );
2960 :
2961 46 : oSRS.morphFromESRI();
2962 46 : oSRS.Fixup();
2963 :
2964 46 : CPLFree( pszProjection );
2965 46 : pszProjection = NULL;
2966 46 : oSRS.exportToWkt( &pszProjection );
2967 :
2968 46 : return CE_None;
2969 : }
2970 :
2971 : /* -------------------------------------------------------------------- */
2972 : /* General case for Erdas style projections. */
2973 : /* */
2974 : /* We make a particular effort to adapt the mapinfo->proname as */
2975 : /* the PROJCS[] name per #2422. */
2976 : /* -------------------------------------------------------------------- */
2977 229 : psDatum = HFAGetDatum( hHFA );
2978 229 : psPro = HFAGetProParameters( hHFA );
2979 229 : psMapInfo = HFAGetMapInfo( hHFA );
2980 :
2981 229 : HFAEntry *poMapInformation = NULL;
2982 229 : if( psMapInfo == NULL )
2983 : {
2984 177 : HFABand *poBand = hHFA->papoBand[0];
2985 177 : poMapInformation = poBand->poNode->GetNamedChild("MapInformation");
2986 : }
2987 :
2988 229 : CPLFree( pszProjection );
2989 :
2990 229 : if( !psDatum || !psPro ||
2991 : (psMapInfo == NULL && poMapInformation == NULL) ||
2992 : ((strlen(psDatum->datumname) == 0 || EQUAL(psDatum->datumname, "Unknown")) &&
2993 : (strlen(psPro->proName) == 0 || EQUAL(psPro->proName, "Unknown")) &&
2994 : (psMapInfo && (strlen(psMapInfo->proName) == 0 || EQUAL(psMapInfo->proName, "Unknown"))) &&
2995 : psPro->proZone == 0) )
2996 : {
2997 197 : pszProjection = CPLStrdup("");
2998 197 : return CE_None;
2999 : }
3000 :
3001 : pszProjection = HFAPCSStructToWKT( psDatum, psPro, psMapInfo,
3002 32 : poMapInformation );
3003 :
3004 32 : if( pszProjection != NULL )
3005 32 : return CE_None;
3006 : else
3007 : {
3008 0 : pszProjection = CPLStrdup("");
3009 0 : return CE_Failure;
3010 0 : }
3011 : }
3012 :
3013 : /************************************************************************/
3014 : /* IBuildOverviews() */
3015 : /************************************************************************/
3016 :
3017 6 : CPLErr HFADataset::IBuildOverviews( const char *pszResampling,
3018 : int nOverviews, int *panOverviewList,
3019 : int nListBands, int *panBandList,
3020 : GDALProgressFunc pfnProgress,
3021 : void * pProgressData )
3022 :
3023 : {
3024 : int i;
3025 :
3026 6 : if( GetAccess() == GA_ReadOnly )
3027 : return GDALDataset::IBuildOverviews( pszResampling,
3028 : nOverviews, panOverviewList,
3029 : nListBands, panBandList,
3030 0 : pfnProgress, pProgressData );
3031 :
3032 12 : for( i = 0; i < nListBands; i++ )
3033 : {
3034 : CPLErr eErr;
3035 : GDALRasterBand *poBand;
3036 :
3037 : void* pScaledProgressData = GDALCreateScaledProgress(
3038 : i * 1.0 / nListBands, (i + 1) * 1.0 / nListBands,
3039 6 : pfnProgress, pProgressData);
3040 :
3041 6 : poBand = GetRasterBand( panBandList[i] );
3042 : eErr =
3043 : poBand->BuildOverviews( pszResampling, nOverviews, panOverviewList,
3044 6 : GDALScaledProgress, pScaledProgressData );
3045 :
3046 6 : GDALDestroyScaledProgress(pScaledProgressData);
3047 :
3048 6 : if( eErr != CE_None )
3049 0 : return eErr;
3050 : }
3051 :
3052 6 : return CE_None;
3053 : }
3054 :
3055 : /************************************************************************/
3056 : /* Identify() */
3057 : /************************************************************************/
3058 :
3059 10222 : int HFADataset::Identify( GDALOpenInfo * poOpenInfo )
3060 :
3061 : {
3062 : /* -------------------------------------------------------------------- */
3063 : /* Verify that this is a HFA file. */
3064 : /* -------------------------------------------------------------------- */
3065 10222 : if( poOpenInfo->nHeaderBytes < 15
3066 : || !EQUALN((char *) poOpenInfo->pabyHeader,"EHFA_HEADER_TAG",15) )
3067 9944 : return FALSE;
3068 : else
3069 278 : return TRUE;
3070 : }
3071 :
3072 : /************************************************************************/
3073 : /* Open() */
3074 : /************************************************************************/
3075 :
3076 2491 : GDALDataset *HFADataset::Open( GDALOpenInfo * poOpenInfo )
3077 :
3078 : {
3079 : HFAHandle hHFA;
3080 : int i;
3081 :
3082 : /* -------------------------------------------------------------------- */
3083 : /* Verify that this is a HFA file. */
3084 : /* -------------------------------------------------------------------- */
3085 2491 : if( !Identify( poOpenInfo ) )
3086 2214 : return NULL;
3087 :
3088 : /* -------------------------------------------------------------------- */
3089 : /* Open the file. */
3090 : /* -------------------------------------------------------------------- */
3091 277 : if( poOpenInfo->eAccess == GA_Update )
3092 102 : hHFA = HFAOpen( poOpenInfo->pszFilename, "r+" );
3093 : else
3094 175 : hHFA = HFAOpen( poOpenInfo->pszFilename, "r" );
3095 :
3096 277 : if( hHFA == NULL )
3097 0 : return NULL;
3098 :
3099 : /* -------------------------------------------------------------------- */
3100 : /* Create a corresponding GDALDataset. */
3101 : /* -------------------------------------------------------------------- */
3102 : HFADataset *poDS;
3103 :
3104 277 : poDS = new HFADataset();
3105 :
3106 277 : poDS->hHFA = hHFA;
3107 277 : poDS->eAccess = poOpenInfo->eAccess;
3108 :
3109 : /* -------------------------------------------------------------------- */
3110 : /* Establish raster info. */
3111 : /* -------------------------------------------------------------------- */
3112 : HFAGetRasterInfo( hHFA, &poDS->nRasterXSize, &poDS->nRasterYSize,
3113 277 : &poDS->nBands );
3114 :
3115 277 : if( poDS->nBands == 0 )
3116 : {
3117 2 : delete poDS;
3118 : CPLError( CE_Failure, CPLE_AppDefined,
3119 : "Unable to open %s, it has zero usable bands.",
3120 2 : poOpenInfo->pszFilename );
3121 2 : return NULL;
3122 : }
3123 :
3124 275 : if( poDS->nRasterXSize == 0 || poDS->nRasterYSize == 0 )
3125 : {
3126 0 : delete poDS;
3127 : CPLError( CE_Failure, CPLE_AppDefined,
3128 : "Unable to open %s, it has no pixels.",
3129 0 : poOpenInfo->pszFilename );
3130 0 : return NULL;
3131 : }
3132 :
3133 : /* -------------------------------------------------------------------- */
3134 : /* Get geotransform, or if that fails, try to find XForms to */
3135 : /* build gcps, and metadata. */
3136 : /* -------------------------------------------------------------------- */
3137 275 : if( !HFAGetGeoTransform( hHFA, poDS->adfGeoTransform ) )
3138 : {
3139 175 : Efga_Polynomial *pasPolyListForward = NULL;
3140 175 : Efga_Polynomial *pasPolyListReverse = NULL;
3141 : int nStepCount =
3142 : HFAReadXFormStack( hHFA, &pasPolyListForward,
3143 175 : &pasPolyListReverse );
3144 :
3145 175 : if( nStepCount > 0 )
3146 : {
3147 : poDS->UseXFormStack( nStepCount,
3148 : pasPolyListForward,
3149 1 : pasPolyListReverse );
3150 1 : CPLFree( pasPolyListForward );
3151 1 : CPLFree( pasPolyListReverse );
3152 : }
3153 : }
3154 :
3155 : /* -------------------------------------------------------------------- */
3156 : /* Get the projection. */
3157 : /* -------------------------------------------------------------------- */
3158 275 : poDS->ReadProjection();
3159 :
3160 : /* -------------------------------------------------------------------- */
3161 : /* Read the camera model as metadata, if present. */
3162 : /* -------------------------------------------------------------------- */
3163 275 : char **papszCM = HFAReadCameraModel( hHFA );
3164 :
3165 275 : if( papszCM != NULL )
3166 : {
3167 0 : poDS->SetMetadata( papszCM, "CAMERA_MODEL" );
3168 0 : CSLDestroy( papszCM );
3169 : }
3170 :
3171 : /* -------------------------------------------------------------------- */
3172 : /* Create band information objects. */
3173 : /* -------------------------------------------------------------------- */
3174 1356 : for( i = 0; i < poDS->nBands; i++ )
3175 : {
3176 403 : poDS->SetBand( i+1, new HFARasterBand( poDS, i+1, -1 ) );
3177 : }
3178 :
3179 : /* -------------------------------------------------------------------- */
3180 : /* Collect GDAL custom Metadata, and "auxilary" metadata from */
3181 : /* well known HFA structures for the bands. We defer this till */
3182 : /* now to ensure that the bands are properly setup before */
3183 : /* interacting with PAM. */
3184 : /* -------------------------------------------------------------------- */
3185 678 : for( i = 0; i < poDS->nBands; i++ )
3186 : {
3187 403 : HFARasterBand *poBand = (HFARasterBand *) poDS->GetRasterBand( i+1 );
3188 :
3189 403 : char **papszMD = HFAGetMetadata( hHFA, i+1 );
3190 403 : if( papszMD != NULL )
3191 : {
3192 4 : poBand->SetMetadata( papszMD );
3193 4 : CSLDestroy( papszMD );
3194 : }
3195 :
3196 403 : poBand->ReadAuxMetadata();
3197 403 : poBand->ReadHistogramMetadata();
3198 : }
3199 :
3200 : /* -------------------------------------------------------------------- */
3201 : /* Check for GDAL style metadata. */
3202 : /* -------------------------------------------------------------------- */
3203 275 : char **papszMD = HFAGetMetadata( hHFA, 0 );
3204 275 : if( papszMD != NULL )
3205 : {
3206 66 : poDS->SetMetadata( papszMD );
3207 66 : CSLDestroy( papszMD );
3208 : }
3209 :
3210 : /* -------------------------------------------------------------------- */
3211 : /* Check for dependent dataset value. */
3212 : /* -------------------------------------------------------------------- */
3213 275 : HFAInfo_t *psInfo = (HFAInfo_t *) hHFA;
3214 275 : HFAEntry *poEntry = psInfo->poRoot->GetNamedChild("DependentFile");
3215 275 : if( poEntry != NULL )
3216 : {
3217 : poDS->SetMetadataItem( "HFA_DEPENDENT_FILE",
3218 : poEntry->GetStringField( "dependent.string" ),
3219 16 : "HFA" );
3220 : }
3221 :
3222 : /* -------------------------------------------------------------------- */
3223 : /* Initialize any PAM information. */
3224 : /* -------------------------------------------------------------------- */
3225 275 : poDS->SetDescription( poOpenInfo->pszFilename );
3226 275 : poDS->TryLoadXML();
3227 :
3228 : /* -------------------------------------------------------------------- */
3229 : /* Check for external overviews. */
3230 : /* -------------------------------------------------------------------- */
3231 275 : poDS->oOvManager.Initialize( poDS, poOpenInfo->pszFilename );
3232 :
3233 : /* -------------------------------------------------------------------- */
3234 : /* Clear dirty metadata flags. */
3235 : /* -------------------------------------------------------------------- */
3236 678 : for( i = 0; i < poDS->nBands; i++ )
3237 : {
3238 403 : HFARasterBand *poBand = (HFARasterBand *) poDS->GetRasterBand( i+1 );
3239 403 : poBand->bMetadataDirty = FALSE;
3240 : }
3241 275 : poDS->bMetadataDirty = FALSE;
3242 :
3243 275 : return( poDS );
3244 : }
3245 :
3246 : /************************************************************************/
3247 : /* GetProjectionRef() */
3248 : /************************************************************************/
3249 :
3250 85 : const char *HFADataset::GetProjectionRef()
3251 :
3252 : {
3253 85 : return pszProjection;
3254 : }
3255 :
3256 : /************************************************************************/
3257 : /* SetProjection() */
3258 : /************************************************************************/
3259 :
3260 36 : CPLErr HFADataset::SetProjection( const char * pszNewProjection )
3261 :
3262 : {
3263 36 : CPLFree( pszProjection );
3264 36 : pszProjection = CPLStrdup( pszNewProjection );
3265 36 : bGeoDirty = TRUE;
3266 :
3267 36 : return CE_None;
3268 : }
3269 :
3270 : /************************************************************************/
3271 : /* SetMetadata() */
3272 : /************************************************************************/
3273 :
3274 108 : CPLErr HFADataset::SetMetadata( char **papszMDIn, const char *pszDomain )
3275 :
3276 : {
3277 108 : bMetadataDirty = TRUE;
3278 :
3279 108 : return GDALPamDataset::SetMetadata( papszMDIn, pszDomain );
3280 : }
3281 :
3282 : /************************************************************************/
3283 : /* SetMetadata() */
3284 : /************************************************************************/
3285 :
3286 16 : CPLErr HFADataset::SetMetadataItem( const char *pszTag, const char *pszValue,
3287 : const char *pszDomain )
3288 :
3289 : {
3290 16 : bMetadataDirty = TRUE;
3291 :
3292 16 : return GDALPamDataset::SetMetadataItem( pszTag, pszValue, pszDomain );
3293 : }
3294 :
3295 : /************************************************************************/
3296 : /* GetGeoTransform() */
3297 : /************************************************************************/
3298 :
3299 62 : CPLErr HFADataset::GetGeoTransform( double * padfTransform )
3300 :
3301 : {
3302 77 : if( adfGeoTransform[0] != 0.0
3303 3 : || adfGeoTransform[1] != 1.0
3304 3 : || adfGeoTransform[2] != 0.0
3305 3 : || adfGeoTransform[3] != 0.0
3306 3 : || adfGeoTransform[4] != 0.0
3307 3 : || adfGeoTransform[5] != 1.0 )
3308 : {
3309 59 : memcpy( padfTransform, adfGeoTransform, sizeof(double)*6 );
3310 59 : return CE_None;
3311 : }
3312 : else
3313 3 : return GDALPamDataset::GetGeoTransform( padfTransform );
3314 : }
3315 :
3316 : /************************************************************************/
3317 : /* SetGeoTransform() */
3318 : /************************************************************************/
3319 :
3320 47 : CPLErr HFADataset::SetGeoTransform( double * padfTransform )
3321 :
3322 : {
3323 47 : memcpy( adfGeoTransform, padfTransform, sizeof(double)*6 );
3324 47 : bGeoDirty = TRUE;
3325 :
3326 47 : return CE_None;
3327 : }
3328 :
3329 : /************************************************************************/
3330 : /* IRasterIO() */
3331 : /* */
3332 : /* Multi-band raster io handler. Here we ensure that the block */
3333 : /* based loading is used for spill file rasters. That is */
3334 : /* because they are effectively pixel interleaved, so */
3335 : /* processing all bands for a given block together avoid extra */
3336 : /* seeks. */
3337 : /************************************************************************/
3338 :
3339 45 : CPLErr HFADataset::IRasterIO( GDALRWFlag eRWFlag,
3340 : int nXOff, int nYOff, int nXSize, int nYSize,
3341 : void *pData, int nBufXSize, int nBufYSize,
3342 : GDALDataType eBufType,
3343 : int nBandCount, int *panBandMap,
3344 : int nPixelSpace, int nLineSpace, int nBandSpace )
3345 :
3346 : {
3347 45 : if( hHFA->papoBand[panBandMap[0]-1]->fpExternal != NULL
3348 : && nBandCount > 1 )
3349 : return GDALDataset::BlockBasedRasterIO(
3350 : eRWFlag, nXOff, nYOff, nXSize, nYSize,
3351 : pData, nBufXSize, nBufYSize, eBufType,
3352 0 : nBandCount, panBandMap, nPixelSpace, nLineSpace, nBandSpace );
3353 : else
3354 : return
3355 : GDALDataset::IRasterIO( eRWFlag, nXOff, nYOff, nXSize, nYSize,
3356 : pData, nBufXSize, nBufYSize, eBufType,
3357 : nBandCount, panBandMap,
3358 45 : nPixelSpace, nLineSpace, nBandSpace );
3359 : }
3360 :
3361 : /************************************************************************/
3362 : /* UseXFormStack() */
3363 : /************************************************************************/
3364 :
3365 1 : void HFADataset::UseXFormStack( int nStepCount,
3366 : Efga_Polynomial *pasPLForward,
3367 : Efga_Polynomial *pasPLReverse )
3368 :
3369 : {
3370 : /* -------------------------------------------------------------------- */
3371 : /* Generate GCPs using the transform. */
3372 : /* -------------------------------------------------------------------- */
3373 : double dfXRatio, dfYRatio;
3374 :
3375 1 : nGCPCount = 0;
3376 1 : GDALInitGCPs( 36, asGCPList );
3377 :
3378 7 : for( dfYRatio = 0.0; dfYRatio < 1.001; dfYRatio += 0.2 )
3379 : {
3380 42 : for( dfXRatio = 0.0; dfXRatio < 1.001; dfXRatio += 0.2 )
3381 : {
3382 36 : double dfLine = 0.5 + (GetRasterYSize()-1) * dfYRatio;
3383 36 : double dfPixel = 0.5 + (GetRasterXSize()-1) * dfXRatio;
3384 36 : int iGCP = nGCPCount;
3385 :
3386 36 : asGCPList[iGCP].dfGCPPixel = dfPixel;
3387 36 : asGCPList[iGCP].dfGCPLine = dfLine;
3388 :
3389 36 : asGCPList[iGCP].dfGCPX = dfPixel;
3390 36 : asGCPList[iGCP].dfGCPY = dfLine;
3391 36 : asGCPList[iGCP].dfGCPZ = 0.0;
3392 :
3393 36 : if( HFAEvaluateXFormStack( nStepCount, FALSE, pasPLReverse,
3394 : &(asGCPList[iGCP].dfGCPX),
3395 : &(asGCPList[iGCP].dfGCPY) ) )
3396 36 : nGCPCount++;
3397 : }
3398 : }
3399 :
3400 : /* -------------------------------------------------------------------- */
3401 : /* Store the transform as metadata. */
3402 : /* -------------------------------------------------------------------- */
3403 : int iStep, i;
3404 :
3405 : GDALMajorObject::SetMetadataItem(
3406 : "XFORM_STEPS",
3407 : CPLString().Printf("%d",nStepCount),
3408 1 : "XFORMS" );
3409 :
3410 3 : for( iStep = 0; iStep < nStepCount; iStep++ )
3411 : {
3412 : GDALMajorObject::SetMetadataItem(
3413 : CPLString().Printf("XFORM%d_ORDER", iStep),
3414 2 : CPLString().Printf("%d",pasPLForward[iStep].order),
3415 2 : "XFORMS" );
3416 :
3417 2 : if( pasPLForward[iStep].order == 1 )
3418 : {
3419 10 : for( i = 0; i < 4; i++ )
3420 : GDALMajorObject::SetMetadataItem(
3421 : CPLString().Printf("XFORM%d_POLYCOEFMTX[%d]", iStep, i),
3422 : CPLString().Printf("%.15g",
3423 4 : pasPLForward[iStep].polycoefmtx[i]),
3424 4 : "XFORMS" );
3425 :
3426 6 : for( i = 0; i < 2; i++ )
3427 : GDALMajorObject::SetMetadataItem(
3428 : CPLString().Printf("XFORM%d_POLYCOEFVECTOR[%d]", iStep, i),
3429 : CPLString().Printf("%.15g",
3430 2 : pasPLForward[iStep].polycoefvector[i]),
3431 2 : "XFORMS" );
3432 :
3433 1 : continue;
3434 : }
3435 :
3436 : int nCoefCount;
3437 :
3438 1 : if( pasPLForward[iStep].order == 2 )
3439 0 : nCoefCount = 10;
3440 : else
3441 : {
3442 : CPLAssert( pasPLForward[iStep].order == 3 );
3443 1 : nCoefCount = 18;
3444 : }
3445 :
3446 38 : for( i = 0; i < nCoefCount; i++ )
3447 : GDALMajorObject::SetMetadataItem(
3448 : CPLString().Printf("XFORM%d_FWD_POLYCOEFMTX[%d]", iStep, i),
3449 : CPLString().Printf("%.15g",
3450 18 : pasPLForward[iStep].polycoefmtx[i]),
3451 18 : "XFORMS" );
3452 :
3453 6 : for( i = 0; i < 2; i++ )
3454 : GDALMajorObject::SetMetadataItem(
3455 : CPLString().Printf("XFORM%d_FWD_POLYCOEFVECTOR[%d]", iStep, i),
3456 : CPLString().Printf("%.15g",
3457 2 : pasPLForward[iStep].polycoefvector[i]),
3458 2 : "XFORMS" );
3459 :
3460 38 : for( i = 0; i < nCoefCount; i++ )
3461 : GDALMajorObject::SetMetadataItem(
3462 : CPLString().Printf("XFORM%d_REV_POLYCOEFMTX[%d]", iStep, i),
3463 : CPLString().Printf("%.15g",
3464 18 : pasPLReverse[iStep].polycoefmtx[i]),
3465 18 : "XFORMS" );
3466 :
3467 6 : for( i = 0; i < 2; i++ )
3468 : GDALMajorObject::SetMetadataItem(
3469 : CPLString().Printf("XFORM%d_REV_POLYCOEFVECTOR[%d]", iStep, i),
3470 : CPLString().Printf("%.15g",
3471 2 : pasPLReverse[iStep].polycoefvector[i]),
3472 2 : "XFORMS" );
3473 : }
3474 1 : }
3475 :
3476 : /************************************************************************/
3477 : /* GetGCPCount() */
3478 : /************************************************************************/
3479 :
3480 7 : int HFADataset::GetGCPCount()
3481 :
3482 : {
3483 7 : return nGCPCount;
3484 : }
3485 :
3486 : /************************************************************************/
3487 : /* GetGCPProjection() */
3488 : /************************************************************************/
3489 :
3490 0 : const char *HFADataset::GetGCPProjection()
3491 :
3492 : {
3493 0 : if( nGCPCount > 0 )
3494 0 : return pszProjection;
3495 : else
3496 0 : return "";
3497 : }
3498 :
3499 : /************************************************************************/
3500 : /* GetGCPs() */
3501 : /************************************************************************/
3502 :
3503 1 : const GDAL_GCP *HFADataset::GetGCPs()
3504 :
3505 : {
3506 1 : return asGCPList;
3507 : }
3508 :
3509 : /************************************************************************/
3510 : /* GetFileList() */
3511 : /************************************************************************/
3512 :
3513 73 : char **HFADataset::GetFileList()
3514 :
3515 : {
3516 73 : char **papszFileList = GDALPamDataset::GetFileList();
3517 :
3518 73 : if( hHFA->pszIGEFilename != NULL )
3519 : {
3520 : papszFileList = CSLAddString( papszFileList,
3521 : CPLFormFilename( hHFA->pszPath,
3522 : hHFA->pszIGEFilename,
3523 0 : NULL ));
3524 : }
3525 :
3526 73 : if( hHFA->psDependent != NULL )
3527 : {
3528 0 : HFAInfo_t *psDep = hHFA->psDependent;
3529 :
3530 : papszFileList =
3531 : CSLAddString( papszFileList,
3532 : CPLFormFilename( psDep->pszPath,
3533 0 : psDep->pszFilename, NULL ));
3534 0 : if( psDep->pszIGEFilename != NULL )
3535 : {
3536 : papszFileList =
3537 : CSLAddString( papszFileList,
3538 : CPLFormFilename( psDep->pszPath,
3539 0 : psDep->pszIGEFilename, NULL ));
3540 : }
3541 : }
3542 :
3543 73 : return papszFileList;
3544 : }
3545 :
3546 : /************************************************************************/
3547 : /* Create() */
3548 : /************************************************************************/
3549 :
3550 99 : GDALDataset *HFADataset::Create( const char * pszFilenameIn,
3551 : int nXSize, int nYSize, int nBands,
3552 : GDALDataType eType,
3553 : char ** papszParmList )
3554 :
3555 : {
3556 : int nHfaDataType;
3557 99 : int nBits = 0;
3558 : const char *pszPixelType;
3559 :
3560 :
3561 99 : if( CSLFetchNameValue( papszParmList, "NBITS" ) != NULL )
3562 2 : nBits = atoi(CSLFetchNameValue(papszParmList,"NBITS"));
3563 :
3564 : pszPixelType =
3565 99 : CSLFetchNameValue( papszParmList, "PIXELTYPE" );
3566 99 : if( pszPixelType == NULL )
3567 99 : pszPixelType = "";
3568 :
3569 : /* -------------------------------------------------------------------- */
3570 : /* Translate the data type. */
3571 : /* -------------------------------------------------------------------- */
3572 99 : switch( eType )
3573 : {
3574 : case GDT_Byte:
3575 40 : if( nBits == 1 )
3576 2 : nHfaDataType = EPT_u1;
3577 38 : else if( nBits == 2 )
3578 0 : nHfaDataType = EPT_u2;
3579 38 : else if( nBits == 4 )
3580 0 : nHfaDataType = EPT_u4;
3581 38 : else if( EQUAL(pszPixelType,"SIGNEDBYTE") )
3582 0 : nHfaDataType = EPT_s8;
3583 : else
3584 38 : nHfaDataType = EPT_u8;
3585 40 : break;
3586 :
3587 : case GDT_UInt16:
3588 10 : nHfaDataType = EPT_u16;
3589 10 : break;
3590 :
3591 : case GDT_Int16:
3592 6 : nHfaDataType = EPT_s16;
3593 6 : break;
3594 :
3595 : case GDT_Int32:
3596 6 : nHfaDataType = EPT_s32;
3597 6 : break;
3598 :
3599 : case GDT_UInt32:
3600 6 : nHfaDataType = EPT_u32;
3601 6 : break;
3602 :
3603 : case GDT_Float32:
3604 6 : nHfaDataType = EPT_f32;
3605 6 : break;
3606 :
3607 : case GDT_Float64:
3608 9 : nHfaDataType = EPT_f64;
3609 9 : break;
3610 :
3611 : case GDT_CFloat32:
3612 6 : nHfaDataType = EPT_c64;
3613 6 : break;
3614 :
3615 : case GDT_CFloat64:
3616 6 : nHfaDataType = EPT_c128;
3617 6 : break;
3618 :
3619 : default:
3620 : CPLError( CE_Failure, CPLE_NotSupported,
3621 : "Data type %s not supported by Erdas Imagine (HFA) format.\n",
3622 4 : GDALGetDataTypeName( eType ) );
3623 4 : return NULL;
3624 :
3625 : }
3626 :
3627 : /* -------------------------------------------------------------------- */
3628 : /* Create the new file. */
3629 : /* -------------------------------------------------------------------- */
3630 : HFAHandle hHFA;
3631 :
3632 : hHFA = HFACreate( pszFilenameIn, nXSize, nYSize, nBands,
3633 95 : nHfaDataType, papszParmList );
3634 95 : if( hHFA == NULL )
3635 0 : return NULL;
3636 :
3637 95 : HFAClose( hHFA );
3638 :
3639 : /* -------------------------------------------------------------------- */
3640 : /* Open the dataset normally. */
3641 : /* -------------------------------------------------------------------- */
3642 95 : HFADataset *poDS = (HFADataset *) GDALOpen( pszFilenameIn, GA_Update );
3643 :
3644 : /* -------------------------------------------------------------------- */
3645 : /* Special creation option to disable checking for UTM */
3646 : /* parameters when writing the projection. This is a special */
3647 : /* hack for sam.gillingham@nrm.qld.gov.au. */
3648 : /* -------------------------------------------------------------------- */
3649 95 : if( poDS != NULL )
3650 : {
3651 : poDS->bIgnoreUTM = CSLFetchBoolean( papszParmList, "IGNOREUTM",
3652 93 : FALSE );
3653 : }
3654 :
3655 : /* -------------------------------------------------------------------- */
3656 : /* Sometimes we can improve ArcGIS compatability by forcing */
3657 : /* generation of a PEString instead of traditional Imagine */
3658 : /* coordinate system descriptions. */
3659 : /* -------------------------------------------------------------------- */
3660 95 : if( poDS != NULL )
3661 : {
3662 : poDS->bForceToPEString =
3663 93 : CSLFetchBoolean( papszParmList, "FORCETOPESTRING", FALSE );
3664 : }
3665 :
3666 95 : return poDS;
3667 :
3668 : }
3669 :
3670 : /************************************************************************/
3671 : /* CreateCopy() */
3672 : /************************************************************************/
3673 :
3674 : GDALDataset *
3675 39 : HFADataset::CreateCopy( const char * pszFilename, GDALDataset *poSrcDS,
3676 : int bStrict, char ** papszOptions,
3677 : GDALProgressFunc pfnProgress, void * pProgressData )
3678 :
3679 : {
3680 : HFADataset *poDS;
3681 39 : GDALDataType eType = GDT_Byte;
3682 : int iBand;
3683 39 : int nBandCount = poSrcDS->GetRasterCount();
3684 39 : char **papszModOptions = CSLDuplicate( papszOptions );
3685 :
3686 : /* -------------------------------------------------------------------- */
3687 : /* Do we really just want to create an .aux file? */
3688 : /* -------------------------------------------------------------------- */
3689 39 : int bCreateAux = CSLFetchBoolean( papszOptions, "AUX", FALSE );
3690 :
3691 : /* -------------------------------------------------------------------- */
3692 : /* Establish a representative data type to use. */
3693 : /* -------------------------------------------------------------------- */
3694 39 : if( !pfnProgress( 0.0, NULL, pProgressData ) )
3695 0 : return NULL;
3696 :
3697 87 : for( iBand = 0; iBand < nBandCount; iBand++ )
3698 : {
3699 48 : GDALRasterBand *poBand = poSrcDS->GetRasterBand( iBand+1 );
3700 48 : eType = GDALDataTypeUnion( eType, poBand->GetRasterDataType() );
3701 : }
3702 :
3703 : /* -------------------------------------------------------------------- */
3704 : /* If we have PIXELTYPE metadadata in the source, pass it */
3705 : /* through as a creation option. */
3706 : /* -------------------------------------------------------------------- */
3707 57 : if( CSLFetchNameValue( papszOptions, "PIXELTYPE" ) == NULL
3708 : && nBandCount > 0
3709 : && eType == GDT_Byte
3710 : && poSrcDS->GetRasterBand(1)->GetMetadataItem( "PIXELTYPE",
3711 18 : "IMAGE_STRUCTURE" ) )
3712 : {
3713 : papszModOptions =
3714 : CSLSetNameValue( papszModOptions, "PIXELTYPE",
3715 : poSrcDS->GetRasterBand(1)->GetMetadataItem(
3716 0 : "PIXELTYPE", "IMAGE_STRUCTURE" ) );
3717 : }
3718 :
3719 : /* -------------------------------------------------------------------- */
3720 : /* Create the file. */
3721 : /* -------------------------------------------------------------------- */
3722 : poDS = (HFADataset *) Create( pszFilename,
3723 : poSrcDS->GetRasterXSize(),
3724 : poSrcDS->GetRasterYSize(),
3725 : nBandCount,
3726 39 : eType, papszModOptions );
3727 :
3728 39 : CSLDestroy( papszModOptions );
3729 :
3730 39 : if( poDS == NULL )
3731 3 : return NULL;
3732 :
3733 : /* -------------------------------------------------------------------- */
3734 : /* Does the source have a PCT for any of the bands? If so, */
3735 : /* copy it over. */
3736 : /* -------------------------------------------------------------------- */
3737 82 : for( iBand = 0; iBand < nBandCount; iBand++ )
3738 : {
3739 46 : GDALRasterBand *poBand = poSrcDS->GetRasterBand( iBand+1 );
3740 : GDALColorTable *poCT;
3741 :
3742 46 : poCT = poBand->GetColorTable();
3743 46 : if( poCT != NULL )
3744 : {
3745 1 : poDS->GetRasterBand(iBand+1)->SetColorTable(poCT);
3746 : }
3747 : }
3748 :
3749 : /* -------------------------------------------------------------------- */
3750 : /* Do we have metadata for any of the bands or the dataset as a */
3751 : /* whole? */
3752 : /* -------------------------------------------------------------------- */
3753 36 : if( poSrcDS->GetMetadata() != NULL )
3754 28 : poDS->SetMetadata( poSrcDS->GetMetadata() );
3755 :
3756 82 : for( iBand = 0; iBand < nBandCount; iBand++ )
3757 : {
3758 46 : GDALRasterBand *poSrcBand = poSrcDS->GetRasterBand( iBand+1 );
3759 46 : poDS->GetRasterBand(iBand+1)->SetMetadata( poSrcBand->GetMetadata() );
3760 46 : poDS->GetRasterBand(iBand+1)->SetDescription( poSrcBand->GetDescription() );
3761 : }
3762 :
3763 : /* -------------------------------------------------------------------- */
3764 : /* Copy projection information. */
3765 : /* -------------------------------------------------------------------- */
3766 : double adfGeoTransform[6];
3767 : const char *pszProj;
3768 :
3769 72 : if( poSrcDS->GetGeoTransform( adfGeoTransform ) == CE_None
3770 36 : && (adfGeoTransform[0] != 0.0 || adfGeoTransform[1] != 1.0
3771 0 : || adfGeoTransform[2] != 0.0 || adfGeoTransform[3] != 0.0
3772 0 : || adfGeoTransform[4] != 0.0 || fabs(adfGeoTransform[5]) != 1.0))
3773 36 : poDS->SetGeoTransform( adfGeoTransform );
3774 :
3775 36 : pszProj = poSrcDS->GetProjectionRef();
3776 36 : if( pszProj != NULL && strlen(pszProj) > 0 )
3777 35 : poDS->SetProjection( pszProj );
3778 :
3779 : /* -------------------------------------------------------------------- */
3780 : /* Copy the imagery. */
3781 : /* -------------------------------------------------------------------- */
3782 36 : if( !bCreateAux )
3783 : {
3784 : CPLErr eErr;
3785 :
3786 : eErr = GDALDatasetCopyWholeRaster( (GDALDatasetH) poSrcDS,
3787 : (GDALDatasetH) poDS,
3788 36 : NULL, pfnProgress, pProgressData );
3789 :
3790 36 : if( eErr != CE_None )
3791 0 : return NULL;
3792 : }
3793 :
3794 : /* -------------------------------------------------------------------- */
3795 : /* Do we want to generate statistics and a histogram? */
3796 : /* -------------------------------------------------------------------- */
3797 36 : if( CSLFetchBoolean( papszOptions, "STATISTICS", FALSE ) )
3798 : {
3799 0 : for( iBand = 0; iBand < nBandCount; iBand++ )
3800 : {
3801 0 : GDALRasterBand *poSrcBand = poSrcDS->GetRasterBand( iBand+1 );
3802 : double dfMin, dfMax, dfMean, dfStdDev;
3803 0 : char **papszStatsMD = NULL;
3804 :
3805 : // -----------------------------------------------------------
3806 : // Statistics
3807 : // -----------------------------------------------------------
3808 :
3809 0 : if( poSrcBand->GetStatistics( TRUE, FALSE, &dfMin, &dfMax,
3810 0 : &dfMean, &dfStdDev ) == CE_None
3811 : || poSrcBand->ComputeStatistics( TRUE, &dfMin, &dfMax,
3812 : &dfMean, &dfStdDev,
3813 0 : pfnProgress, pProgressData )
3814 : == CE_None )
3815 : {
3816 0 : CPLString osValue;
3817 :
3818 : papszStatsMD =
3819 : CSLSetNameValue( papszStatsMD, "STATISTICS_MINIMUM",
3820 0 : osValue.Printf( "%.15g", dfMin ) );
3821 : papszStatsMD =
3822 : CSLSetNameValue( papszStatsMD, "STATISTICS_MAXIMUM",
3823 0 : osValue.Printf( "%.15g", dfMax ) );
3824 : papszStatsMD =
3825 : CSLSetNameValue( papszStatsMD, "STATISTICS_MEAN",
3826 0 : osValue.Printf( "%.15g", dfMean ) );
3827 : papszStatsMD =
3828 : CSLSetNameValue( papszStatsMD, "STATISTICS_STDDEV",
3829 0 : osValue.Printf( "%.15g", dfStdDev ) );
3830 : }
3831 :
3832 : // -----------------------------------------------------------
3833 : // Histogram
3834 : // -----------------------------------------------------------
3835 :
3836 0 : int nBuckets, *panHistogram = NULL;
3837 :
3838 0 : if( poSrcBand->GetDefaultHistogram( &dfMin, &dfMax,
3839 : &nBuckets, &panHistogram,
3840 : TRUE,
3841 0 : pfnProgress, pProgressData )
3842 : == CE_None )
3843 : {
3844 0 : CPLString osValue;
3845 0 : char *pszBinValues = (char *) CPLCalloc(12,nBuckets+1);
3846 0 : int iBin, nBinValuesLen = 0;
3847 0 : double dfBinWidth = (dfMax - dfMin) / nBuckets;
3848 :
3849 : papszStatsMD = CSLSetNameValue(
3850 : papszStatsMD, "STATISTICS_HISTOMIN",
3851 0 : osValue.Printf( "%.15g", dfMin+dfBinWidth*0.5 ) );
3852 : papszStatsMD = CSLSetNameValue(
3853 : papszStatsMD, "STATISTICS_HISTOMAX",
3854 0 : osValue.Printf( "%.15g", dfMax-dfBinWidth*0.5 ) );
3855 : papszStatsMD =
3856 : CSLSetNameValue( papszStatsMD, "STATISTICS_HISTONUMBINS",
3857 0 : osValue.Printf( "%d", nBuckets ) );
3858 :
3859 0 : for( iBin = 0; iBin < nBuckets; iBin++ )
3860 : {
3861 :
3862 : strcat( pszBinValues+nBinValuesLen,
3863 0 : osValue.Printf( "%d", panHistogram[iBin]) );
3864 0 : strcat( pszBinValues+nBinValuesLen, "|" );
3865 0 : nBinValuesLen += strlen(pszBinValues+nBinValuesLen);
3866 : }
3867 : papszStatsMD =
3868 : CSLSetNameValue( papszStatsMD, "STATISTICS_HISTOBINVALUES",
3869 0 : pszBinValues );
3870 0 : CPLFree( pszBinValues );
3871 : }
3872 :
3873 0 : if( CSLCount(papszStatsMD) > 0 )
3874 0 : HFASetMetadata( poDS->hHFA, iBand+1, papszStatsMD );
3875 :
3876 0 : CSLDestroy( papszStatsMD );
3877 : }
3878 : }
3879 :
3880 : /* -------------------------------------------------------------------- */
3881 : /* All report completion. */
3882 : /* -------------------------------------------------------------------- */
3883 36 : if( !pfnProgress( 1.0, NULL, pProgressData ) )
3884 : {
3885 : CPLError( CE_Failure, CPLE_UserInterrupt,
3886 0 : "User terminated" );
3887 0 : delete poDS;
3888 :
3889 : GDALDriver *poHFADriver =
3890 0 : (GDALDriver *) GDALGetDriverByName( "HFA" );
3891 0 : poHFADriver->Delete( pszFilename );
3892 0 : return NULL;
3893 : }
3894 :
3895 36 : poDS->CloneInfo( poSrcDS, GCIF_PAM_DEFAULT );
3896 :
3897 36 : return poDS;
3898 : }
3899 :
3900 : /************************************************************************/
3901 : /* GDALRegister_HFA() */
3902 : /************************************************************************/
3903 :
3904 338 : void GDALRegister_HFA()
3905 :
3906 : {
3907 : GDALDriver *poDriver;
3908 :
3909 338 : if( GDALGetDriverByName( "HFA" ) == NULL )
3910 : {
3911 336 : poDriver = new GDALDriver();
3912 :
3913 336 : poDriver->SetDescription( "HFA" );
3914 : poDriver->SetMetadataItem( GDAL_DMD_LONGNAME,
3915 336 : "Erdas Imagine Images (.img)" );
3916 : poDriver->SetMetadataItem( GDAL_DMD_HELPTOPIC,
3917 336 : "frmt_hfa.html" );
3918 336 : poDriver->SetMetadataItem( GDAL_DMD_EXTENSION, "img" );
3919 : poDriver->SetMetadataItem( GDAL_DMD_CREATIONDATATYPES,
3920 336 : "Byte Int16 UInt16 Int32 UInt32 Float32 Float64 CFloat32 CFloat64" );
3921 :
3922 : poDriver->SetMetadataItem( GDAL_DMD_CREATIONOPTIONLIST,
3923 : "<CreationOptionList>"
3924 : " <Option name='BLOCKSIZE' type='integer' description='tile width/height (32-2048)' default='64'/>"
3925 : " <Option name='USE_SPILL' type='boolean' description='Force use of spill file'/>"
3926 : " <Option name='COMPRESSED' alias='COMPRESS' type='boolean' description='compress blocks'/>"
3927 : " <Option name='PIXELTYPE' type='string' description='By setting this to SIGNEDBYTE, a new Byte file can be forced to be written as signed byte'/>"
3928 : " <Option name='AUX' type='boolean' description='Create an .aux file'/>"
3929 : " <Option name='IGNOREUTM' type='boolean' description='Ignore UTM when selecting coordinate system - will use Transverse Mercator. Only used for Create() method'/>"
3930 : " <Option name='NBITS' type='integer' description='Create file with special sub-byte data type (1/2/4)'/>"
3931 : " <Option name='STATISTICS' type='boolean' description='Generate statistics and a histogram'/>"
3932 : " <Option name='DEPENDENT_FILE' type='string' description='Name of dependent file (must not have absolute path)'/>"
3933 : " <Option name='FORCETOPESTRING' type='boolean' description='Force use of ArcGIS PE String in file instead of Imagine coordinate system format'/>"
3934 336 : "</CreationOptionList>" );
3935 :
3936 336 : poDriver->SetMetadataItem( GDAL_DCAP_VIRTUALIO, "YES" );
3937 :
3938 336 : poDriver->pfnOpen = HFADataset::Open;
3939 336 : poDriver->pfnCreate = HFADataset::Create;
3940 336 : poDriver->pfnCreateCopy = HFADataset::CreateCopy;
3941 336 : poDriver->pfnIdentify = HFADataset::Identify;
3942 :
3943 336 : GetGDALDriverManager()->RegisterDriver( poDriver );
3944 : }
3945 338 : }
3946 :
|
