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