1 : /******************************************************************************
2 : * $Id: ehdrdataset.cpp 25701 2013-03-07 14:28:02Z rouault $
3 : *
4 : * Project: ESRI .hdr Driver
5 : * Purpose: Implementation of EHdrDataset
6 : * Author: Frank Warmerdam, warmerdam@pobox.com
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 1999, Frank Warmerdam <warmerdam@pobox.com>
10 : *
11 : * Permission is hereby granted, free of charge, to any person obtaining a
12 : * copy of this software and associated documentation files (the "Software"),
13 : * to deal in the Software without restriction, including without limitation
14 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
15 : * and/or sell copies of the Software, and to permit persons to whom the
16 : * Software is furnished to do so, subject to the following conditions:
17 : *
18 : * The above copyright notice and this permission notice shall be included
19 : * in all copies or substantial portions of the Software.
20 : *
21 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
22 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
23 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
24 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
25 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
26 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
27 : * DEALINGS IN THE SOFTWARE.
28 : ****************************************************************************/
29 :
30 : #include "rawdataset.h"
31 : #include "ogr_spatialref.h"
32 : #include "cpl_string.h"
33 :
34 : CPL_CVSID("$Id: ehdrdataset.cpp 25701 2013-03-07 14:28:02Z rouault $");
35 :
36 : CPL_C_START
37 : void GDALRegister_EHdr(void);
38 : CPL_C_END
39 :
40 : #define HAS_MIN_FLAG 0x1
41 : #define HAS_MAX_FLAG 0x2
42 : #define HAS_MEAN_FLAG 0x4
43 : #define HAS_STDDEV_FLAG 0x8
44 : #define HAS_ALL_FLAGS (HAS_MIN_FLAG | HAS_MAX_FLAG | HAS_MEAN_FLAG | HAS_STDDEV_FLAG)
45 :
46 : /************************************************************************/
47 : /* ==================================================================== */
48 : /* EHdrDataset */
49 : /* ==================================================================== */
50 : /************************************************************************/
51 :
52 : class EHdrRasterBand;
53 :
54 : class EHdrDataset : public RawDataset
55 : {
56 : friend class EHdrRasterBand;
57 :
58 : VSILFILE *fpImage; // image data file.
59 :
60 : CPLString osHeaderExt;
61 :
62 : int bGotTransform;
63 : double adfGeoTransform[6];
64 : char *pszProjection;
65 :
66 : int bHDRDirty;
67 : char **papszHDR;
68 :
69 : int bCLRDirty;
70 :
71 : CPLErr ReadSTX();
72 : CPLErr RewriteSTX();
73 : CPLErr RewriteHDR();
74 : void ResetKeyValue( const char *pszKey, const char *pszValue );
75 : const char *GetKeyValue( const char *pszKey, const char *pszDefault = "" );
76 : void RewriteColorTable( GDALColorTable * );
77 :
78 : public:
79 : EHdrDataset();
80 : ~EHdrDataset();
81 :
82 : virtual CPLErr GetGeoTransform( double * padfTransform );
83 : virtual CPLErr SetGeoTransform( double *padfTransform );
84 : virtual const char *GetProjectionRef(void);
85 : virtual CPLErr SetProjection( const char * );
86 :
87 : virtual char **GetFileList();
88 :
89 : static GDALDataset *Open( GDALOpenInfo * );
90 : static GDALDataset *Create( const char * pszFilename,
91 : int nXSize, int nYSize, int nBands,
92 : GDALDataType eType, char ** papszParmList );
93 : static GDALDataset *CreateCopy( const char * pszFilename,
94 : GDALDataset * poSrcDS,
95 : int bStrict, char ** papszOptions,
96 : GDALProgressFunc pfnProgress,
97 : void * pProgressData );
98 : static CPLString GetImageRepFilename(const char* pszFilename);
99 : };
100 :
101 : /************************************************************************/
102 : /* ==================================================================== */
103 : /* EHdrRasterBand */
104 : /* ==================================================================== */
105 : /************************************************************************/
106 :
107 : class EHdrRasterBand : public RawRasterBand
108 126 : {
109 : friend class EHdrDataset;
110 :
111 : int nBits;
112 : long nStartBit;
113 : int nPixelOffsetBits;
114 : int nLineOffsetBits;
115 :
116 : int bNoDataSet;
117 : double dfNoData;
118 : double dfMin;
119 : double dfMax;
120 : double dfMean;
121 : double dfStdDev;
122 :
123 : int minmaxmeanstddev;
124 :
125 : virtual CPLErr IRasterIO( GDALRWFlag, int, int, int, int,
126 : void *, int, int, GDALDataType,
127 : int, int );
128 :
129 : public:
130 : EHdrRasterBand( GDALDataset *poDS, int nBand, VSILFILE * fpRaw,
131 : vsi_l_offset nImgOffset, int nPixelOffset,
132 : int nLineOffset,
133 : GDALDataType eDataType, int bNativeOrder,
134 : int nBits);
135 :
136 : virtual CPLErr IReadBlock( int, int, void * );
137 : virtual CPLErr IWriteBlock( int, int, void * );
138 :
139 : virtual double GetNoDataValue( int *pbSuccess = NULL );
140 : virtual double GetMinimum( int *pbSuccess = NULL );
141 : virtual double GetMaximum(int *pbSuccess = NULL );
142 : virtual CPLErr GetStatistics( int bApproxOK, int bForce,
143 : double *pdfMin, double *pdfMax,
144 : double *pdfMean, double *pdfStdDev );
145 : virtual CPLErr SetStatistics( double dfMin, double dfMax,
146 : double dfMean, double dfStdDev );
147 : virtual CPLErr SetColorTable( GDALColorTable *poNewCT );
148 :
149 : };
150 :
151 : /************************************************************************/
152 : /* EHdrRasterBand() */
153 : /************************************************************************/
154 :
155 126 : EHdrRasterBand::EHdrRasterBand( GDALDataset *poDS,
156 : int nBand, VSILFILE * fpRaw,
157 : vsi_l_offset nImgOffset, int nPixelOffset,
158 : int nLineOffset,
159 : GDALDataType eDataType, int bNativeOrder,
160 : int nBits)
161 : : RawRasterBand( poDS, nBand, fpRaw, nImgOffset, nPixelOffset, nLineOffset,
162 : eDataType, bNativeOrder, TRUE ),
163 : nBits(nBits),
164 : bNoDataSet(FALSE),
165 : dfNoData(0),
166 : dfMin(0),
167 : dfMax(0),
168 126 : minmaxmeanstddev(0)
169 : {
170 126 : EHdrDataset* poEDS = (EHdrDataset*)poDS;
171 :
172 126 : if (nBits < 8)
173 : {
174 3 : nStartBit = atoi(poEDS->GetKeyValue("SKIPBYTES")) * 8;
175 3 : if (nBand >= 2)
176 : {
177 0 : long nRowBytes = atoi(poEDS->GetKeyValue("BANDROWBYTES"));
178 0 : if (nRowBytes == 0)
179 0 : nRowBytes = (nBits * poDS->GetRasterXSize() + 7) / 8;
180 :
181 0 : nStartBit += nRowBytes * (nBand-1) * 8;
182 : }
183 :
184 3 : nPixelOffsetBits = nBits;
185 3 : nLineOffsetBits = atoi(poEDS->GetKeyValue("TOTALROWBYTES")) * 8;
186 3 : if( nLineOffsetBits == 0 )
187 0 : nLineOffsetBits = nPixelOffsetBits * poDS->GetRasterXSize();
188 :
189 3 : nBlockXSize = poDS->GetRasterXSize();
190 3 : nBlockYSize = 1;
191 :
192 : SetMetadataItem( "NBITS",
193 : CPLString().Printf( "%d", nBits ),
194 3 : "IMAGE_STRUCTURE" );
195 : }
196 :
197 126 : if( eDataType == GDT_Byte
198 : && EQUAL(poEDS->GetKeyValue("PIXELTYPE",""),"SIGNEDINT") )
199 : SetMetadataItem( "PIXELTYPE", "SIGNEDBYTE",
200 2 : "IMAGE_STRUCTURE" );
201 126 : }
202 :
203 : /************************************************************************/
204 : /* IReadBlock() */
205 : /************************************************************************/
206 :
207 256 : CPLErr EHdrRasterBand::IReadBlock( int nBlockXOff, int nBlockYOff,
208 : void * pImage )
209 :
210 : {
211 256 : if (nBits >= 8)
212 224 : return RawRasterBand::IReadBlock(nBlockXOff, nBlockYOff, pImage);
213 :
214 : vsi_l_offset nLineStart;
215 : unsigned int nLineBytes;
216 : int iBitOffset;
217 : GByte *pabyBuffer;
218 :
219 : /* -------------------------------------------------------------------- */
220 : /* Establish desired position. */
221 : /* -------------------------------------------------------------------- */
222 32 : nLineBytes = (nPixelOffsetBits*nBlockXSize + 7)/8;
223 32 : nLineStart = (nStartBit + ((vsi_l_offset)nLineOffsetBits) * nBlockYOff) / 8;
224 : iBitOffset = (int)
225 32 : ((nStartBit + ((vsi_l_offset)nLineOffsetBits) * nBlockYOff) % 8);
226 :
227 : /* -------------------------------------------------------------------- */
228 : /* Read data into buffer. */
229 : /* -------------------------------------------------------------------- */
230 32 : pabyBuffer = (GByte *) CPLCalloc(nLineBytes,1);
231 :
232 32 : if( VSIFSeekL( GetFPL(), nLineStart, SEEK_SET ) != 0
233 : || VSIFReadL( pabyBuffer, 1, nLineBytes, GetFPL() ) != nLineBytes )
234 : {
235 : CPLError( CE_Failure, CPLE_FileIO,
236 : "Failed to read %u bytes at offset %lu.\n%s",
237 : nLineBytes, (unsigned long)nLineStart,
238 0 : VSIStrerror( errno ) );
239 0 : return CE_Failure;
240 : }
241 :
242 : /* -------------------------------------------------------------------- */
243 : /* Copy data, promoting to 8bit. */
244 : /* -------------------------------------------------------------------- */
245 32 : int iPixel = 0, iX;
246 :
247 1056 : for( iX = 0; iX < nBlockXSize; iX++ )
248 : {
249 1024 : int nOutWord = 0, iBit;
250 :
251 2048 : for( iBit = 0; iBit < nBits; iBit++ )
252 : {
253 1024 : if( pabyBuffer[iBitOffset>>3] & (0x80 >>(iBitOffset & 7)) )
254 200 : nOutWord |= (1 << (nBits - 1 - iBit));
255 1024 : iBitOffset++;
256 : }
257 :
258 1024 : iBitOffset = iBitOffset + nPixelOffsetBits - nBits;
259 :
260 1024 : ((GByte *) pImage)[iPixel++] = (GByte) nOutWord;
261 : }
262 :
263 32 : CPLFree( pabyBuffer );
264 :
265 32 : return CE_None;
266 : }
267 :
268 : /************************************************************************/
269 : /* IWriteBlock() */
270 : /************************************************************************/
271 :
272 472 : CPLErr EHdrRasterBand::IWriteBlock( int nBlockXOff, int nBlockYOff,
273 : void * pImage )
274 :
275 : {
276 472 : if (nBits >= 8)
277 440 : return RawRasterBand::IWriteBlock(nBlockXOff, nBlockYOff, pImage);
278 :
279 : vsi_l_offset nLineStart;
280 : unsigned int nLineBytes;
281 : int iBitOffset;
282 : GByte *pabyBuffer;
283 :
284 : /* -------------------------------------------------------------------- */
285 : /* Establish desired position. */
286 : /* -------------------------------------------------------------------- */
287 32 : nLineBytes = (nPixelOffsetBits*nBlockXSize + 7)/8;
288 32 : nLineStart = (nStartBit + ((vsi_l_offset)nLineOffsetBits) * nBlockYOff) / 8;
289 : iBitOffset = (int)
290 32 : ((nStartBit + ((vsi_l_offset)nLineOffsetBits) * nBlockYOff) % 8);
291 :
292 : /* -------------------------------------------------------------------- */
293 : /* Read data into buffer. */
294 : /* -------------------------------------------------------------------- */
295 32 : pabyBuffer = (GByte *) CPLCalloc(nLineBytes,1);
296 :
297 32 : if( VSIFSeekL( GetFPL(), nLineStart, SEEK_SET ) != 0 )
298 : {
299 : CPLError( CE_Failure, CPLE_FileIO,
300 : "Failed to read %u bytes at offset %lu.\n%s",
301 : nLineBytes, (unsigned long)nLineStart,
302 0 : VSIStrerror( errno ) );
303 0 : return CE_Failure;
304 : }
305 :
306 32 : VSIFReadL( pabyBuffer, 1, nLineBytes, GetFPL() );
307 :
308 : /* -------------------------------------------------------------------- */
309 : /* Copy data, promoting to 8bit. */
310 : /* -------------------------------------------------------------------- */
311 32 : int iPixel = 0, iX;
312 :
313 1056 : for( iX = 0; iX < nBlockXSize; iX++ )
314 : {
315 : int iBit;
316 1024 : int nOutWord = ((GByte *) pImage)[iPixel++];
317 :
318 2048 : for( iBit = 0; iBit < nBits; iBit++ )
319 : {
320 1024 : if( nOutWord & (1 << (nBits - 1 - iBit)) )
321 200 : pabyBuffer[iBitOffset>>3] |= (0x80 >>(iBitOffset & 7));
322 : else
323 824 : pabyBuffer[iBitOffset>>3] &= ~((0x80 >>(iBitOffset & 7)));
324 :
325 1024 : iBitOffset++;
326 : }
327 :
328 1024 : iBitOffset = iBitOffset + nPixelOffsetBits - nBits;
329 : }
330 :
331 : /* -------------------------------------------------------------------- */
332 : /* Write the data back out. */
333 : /* -------------------------------------------------------------------- */
334 32 : if( VSIFSeekL( GetFPL(), nLineStart, SEEK_SET ) != 0
335 : || VSIFWriteL( pabyBuffer, 1, nLineBytes, GetFPL() ) != nLineBytes )
336 : {
337 : CPLError( CE_Failure, CPLE_FileIO,
338 : "Failed to write %u bytes at offset %lu.\n%s",
339 : nLineBytes, (unsigned long)nLineStart,
340 0 : VSIStrerror( errno ) );
341 0 : return CE_Failure;
342 : }
343 :
344 32 : CPLFree( pabyBuffer );
345 :
346 32 : return CE_None;
347 : }
348 :
349 : /************************************************************************/
350 : /* IRasterIO() */
351 : /************************************************************************/
352 :
353 408 : CPLErr EHdrRasterBand::IRasterIO( GDALRWFlag eRWFlag,
354 : int nXOff, int nYOff, int nXSize, int nYSize,
355 : void * pData, int nBufXSize, int nBufYSize,
356 : GDALDataType eBufType,
357 : int nPixelSpace, int nLineSpace )
358 :
359 : {
360 : // Defer to RawRasterBand
361 408 : if (nBits >= 8)
362 : return RawRasterBand::IRasterIO( eRWFlag,
363 : nXOff, nYOff, nXSize, nYSize,
364 : pData, nBufXSize, nBufYSize,
365 375 : eBufType, nPixelSpace, nLineSpace );
366 :
367 : // Force use of IReadBlock() and IWriteBlock()
368 : else
369 : return GDALRasterBand::IRasterIO( eRWFlag,
370 : nXOff, nYOff, nXSize, nYSize,
371 : pData, nBufXSize, nBufYSize,
372 33 : eBufType, nPixelSpace, nLineSpace );
373 : }
374 :
375 : /************************************************************************/
376 : /* OSR_GDS() */
377 : /************************************************************************/
378 :
379 16 : static const char*OSR_GDS( char* pszResult, int nResultLen,
380 : char **papszNV, const char * pszField,
381 : const char *pszDefaultValue )
382 :
383 : {
384 : int iLine;
385 :
386 16 : if( papszNV == NULL || papszNV[0] == NULL )
387 0 : return pszDefaultValue;
388 :
389 125 : for( iLine = 0;
390 62 : papszNV[iLine] != NULL &&
391 47 : !EQUALN(papszNV[iLine],pszField,strlen(pszField));
392 : iLine++ ) {}
393 :
394 16 : if( papszNV[iLine] == NULL )
395 15 : return pszDefaultValue;
396 : else
397 : {
398 : char **papszTokens;
399 :
400 1 : papszTokens = CSLTokenizeString(papszNV[iLine]);
401 :
402 1 : if( CSLCount(papszTokens) > 1 )
403 1 : strncpy( pszResult, papszTokens[1], nResultLen);
404 : else
405 0 : strncpy( pszResult, pszDefaultValue, nResultLen);
406 1 : pszResult[nResultLen-1] = '\0';
407 :
408 1 : CSLDestroy( papszTokens );
409 1 : return pszResult;
410 : }
411 : }
412 :
413 :
414 : /************************************************************************/
415 : /* ==================================================================== */
416 : /* EHdrDataset */
417 : /* ==================================================================== */
418 : /************************************************************************/
419 :
420 : /************************************************************************/
421 : /* EHdrDataset() */
422 : /************************************************************************/
423 :
424 76 : EHdrDataset::EHdrDataset()
425 : {
426 76 : fpImage = NULL;
427 76 : pszProjection = CPLStrdup("");
428 76 : bGotTransform = FALSE;
429 76 : adfGeoTransform[0] = 0.0;
430 76 : adfGeoTransform[1] = 1.0;
431 76 : adfGeoTransform[2] = 0.0;
432 76 : adfGeoTransform[3] = 0.0;
433 76 : adfGeoTransform[4] = 0.0;
434 76 : adfGeoTransform[5] = 1.0;
435 76 : papszHDR = NULL;
436 76 : bHDRDirty = FALSE;
437 76 : bCLRDirty = FALSE;
438 76 : osHeaderExt = "hdr";
439 76 : }
440 :
441 : /************************************************************************/
442 : /* ~EHdrDataset() */
443 : /************************************************************************/
444 :
445 76 : EHdrDataset::~EHdrDataset()
446 :
447 : {
448 76 : FlushCache();
449 :
450 76 : if( nBands > 0 && GetAccess() == GA_Update )
451 : {
452 : int bNoDataSet;
453 : double dfNoData;
454 34 : RawRasterBand *poBand = (RawRasterBand *) GetRasterBand( 1 );
455 :
456 34 : dfNoData = poBand->GetNoDataValue(&bNoDataSet);
457 34 : if( bNoDataSet )
458 : {
459 : ResetKeyValue( "NODATA",
460 1 : CPLString().Printf( "%.8g", dfNoData ) );
461 : }
462 :
463 34 : if( bCLRDirty )
464 2 : RewriteColorTable( poBand->GetColorTable() );
465 :
466 34 : if( bHDRDirty )
467 33 : RewriteHDR();
468 : }
469 :
470 76 : if( fpImage != NULL )
471 76 : VSIFCloseL( fpImage );
472 :
473 76 : CPLFree( pszProjection );
474 76 : CSLDestroy( papszHDR );
475 76 : }
476 :
477 : /************************************************************************/
478 : /* GetKeyValue() */
479 : /************************************************************************/
480 :
481 66 : const char *EHdrDataset::GetKeyValue( const char *pszKey,
482 : const char *pszDefault )
483 :
484 : {
485 : int i;
486 :
487 602 : for( i = 0; papszHDR[i] != NULL; i++ )
488 : {
489 662 : if( EQUALN(pszKey,papszHDR[i],strlen(pszKey))
490 63 : && isspace((unsigned char)papszHDR[i][strlen(pszKey)]) )
491 : {
492 63 : const char *pszValue = papszHDR[i] + strlen(pszKey);
493 492 : while( isspace((unsigned char)*pszValue) )
494 366 : pszValue++;
495 :
496 63 : return pszValue;
497 : }
498 : }
499 :
500 3 : return pszDefault;
501 : }
502 :
503 : /************************************************************************/
504 : /* ResetKeyValue() */
505 : /* */
506 : /* Replace or add the keyword with the indicated value in the */
507 : /* papszHDR list. */
508 : /************************************************************************/
509 :
510 129 : void EHdrDataset::ResetKeyValue( const char *pszKey, const char *pszValue )
511 :
512 : {
513 : int i;
514 : char szNewLine[82];
515 :
516 129 : if( strlen(pszValue) > 65 )
517 : {
518 0 : CPLAssert( strlen(pszValue) <= 65 );
519 0 : return;
520 : }
521 :
522 129 : sprintf( szNewLine, "%-15s%s", pszKey, pszValue );
523 :
524 1482 : for( i = CSLCount(papszHDR)-1; i >= 0; i-- )
525 : {
526 1353 : if( EQUALN(papszHDR[i],szNewLine,strlen(pszKey)+1 ) )
527 : {
528 0 : if( strcmp(papszHDR[i],szNewLine) != 0 )
529 : {
530 0 : CPLFree( papszHDR[i] );
531 0 : papszHDR[i] = CPLStrdup( szNewLine );
532 0 : bHDRDirty = TRUE;
533 : }
534 0 : return;
535 : }
536 : }
537 :
538 129 : bHDRDirty = TRUE;
539 129 : papszHDR = CSLAddString( papszHDR, szNewLine );
540 : }
541 :
542 : /************************************************************************/
543 : /* RewriteColorTable() */
544 : /************************************************************************/
545 :
546 2 : void EHdrDataset::RewriteColorTable( GDALColorTable *poTable )
547 :
548 : {
549 2 : CPLString osCLRFilename = CPLResetExtension( GetDescription(), "clr" );
550 2 : if( poTable )
551 : {
552 2 : VSILFILE *fp = VSIFOpenL( osCLRFilename, "wt" );
553 2 : if( fp != NULL )
554 : {
555 8 : for( int iColor = 0; iColor < poTable->GetColorEntryCount(); iColor++ )
556 : {
557 6 : CPLString oLine;
558 : GDALColorEntry sEntry;
559 :
560 6 : poTable->GetColorEntryAsRGB( iColor, &sEntry );
561 :
562 : // I wish we had a way to mark transparency.
563 : oLine.Printf( "%3d %3d %3d %3d\n",
564 6 : iColor, sEntry.c1, sEntry.c2, sEntry.c3 );
565 6 : VSIFWriteL( (void *) oLine.c_str(), 1, strlen(oLine), fp );
566 : }
567 2 : VSIFCloseL( fp );
568 : }
569 : else
570 : {
571 : CPLError( CE_Failure, CPLE_OpenFailed,
572 : "Unable to create color file %s.",
573 0 : osCLRFilename.c_str() );
574 : }
575 : }
576 : else
577 0 : VSIUnlink( osCLRFilename );
578 2 : }
579 :
580 : /************************************************************************/
581 : /* GetProjectionRef() */
582 : /************************************************************************/
583 :
584 16 : const char *EHdrDataset::GetProjectionRef()
585 :
586 : {
587 16 : if (pszProjection && strlen(pszProjection) > 0)
588 16 : return pszProjection;
589 :
590 0 : return GDALPamDataset::GetProjectionRef();
591 : }
592 :
593 : /************************************************************************/
594 : /* SetProjection() */
595 : /************************************************************************/
596 :
597 30 : CPLErr EHdrDataset::SetProjection( const char *pszSRS )
598 :
599 : {
600 : /* -------------------------------------------------------------------- */
601 : /* Reset coordinate system on the dataset. */
602 : /* -------------------------------------------------------------------- */
603 30 : CPLFree( pszProjection );
604 30 : pszProjection = CPLStrdup( pszSRS );
605 :
606 30 : if( strlen(pszSRS) == 0 )
607 0 : return CE_None;
608 :
609 : /* -------------------------------------------------------------------- */
610 : /* Convert to ESRI WKT. */
611 : /* -------------------------------------------------------------------- */
612 30 : OGRSpatialReference oSRS( pszSRS );
613 30 : char *pszESRI_SRS = NULL;
614 :
615 30 : oSRS.morphToESRI();
616 30 : oSRS.exportToWkt( &pszESRI_SRS );
617 :
618 : /* -------------------------------------------------------------------- */
619 : /* Write to .prj file. */
620 : /* -------------------------------------------------------------------- */
621 30 : CPLString osPrjFilename = CPLResetExtension( GetDescription(), "prj" );
622 : VSILFILE *fp;
623 :
624 30 : fp = VSIFOpenL( osPrjFilename.c_str(), "wt" );
625 30 : if( fp != NULL )
626 : {
627 30 : VSIFWriteL( pszESRI_SRS, 1, strlen(pszESRI_SRS), fp );
628 30 : VSIFWriteL( (void *) "\n", 1, 1, fp );
629 30 : VSIFCloseL( fp );
630 : }
631 :
632 30 : CPLFree( pszESRI_SRS );
633 :
634 30 : return CE_None;
635 : }
636 :
637 : /************************************************************************/
638 : /* GetGeoTransform() */
639 : /************************************************************************/
640 :
641 24 : CPLErr EHdrDataset::GetGeoTransform( double * padfTransform )
642 :
643 : {
644 24 : if( bGotTransform )
645 : {
646 24 : memcpy( padfTransform, adfGeoTransform, sizeof(double) * 6 );
647 24 : return CE_None;
648 : }
649 : else
650 : {
651 0 : return GDALPamDataset::GetGeoTransform( padfTransform );
652 : }
653 : }
654 :
655 : /************************************************************************/
656 : /* SetGeoTransform() */
657 : /************************************************************************/
658 :
659 32 : CPLErr EHdrDataset::SetGeoTransform( double *padfGeoTransform )
660 :
661 : {
662 : /* -------------------------------------------------------------------- */
663 : /* We only support non-rotated images with info in the .HDR file. */
664 : /* -------------------------------------------------------------------- */
665 64 : if( padfGeoTransform[2] != 0.0
666 32 : || padfGeoTransform[4] != 0.0 )
667 : {
668 0 : return GDALPamDataset::SetGeoTransform( padfGeoTransform );
669 : }
670 :
671 : /* -------------------------------------------------------------------- */
672 : /* Record new geotransform. */
673 : /* -------------------------------------------------------------------- */
674 32 : bGotTransform = TRUE;
675 32 : memcpy( adfGeoTransform, padfGeoTransform, sizeof(double) * 6 );
676 :
677 : /* -------------------------------------------------------------------- */
678 : /* Strip out all old geotransform keywords from HDR records. */
679 : /* -------------------------------------------------------------------- */
680 : int i;
681 320 : for( i = CSLCount(papszHDR)-1; i >= 0; i-- )
682 : {
683 1152 : if( EQUALN(papszHDR[i],"ul",2)
684 288 : || EQUALN(papszHDR[i]+1,"ll",2)
685 288 : || EQUALN(papszHDR[i],"cell",4)
686 288 : || EQUALN(papszHDR[i]+1,"dim",3) )
687 : {
688 0 : papszHDR = CSLRemoveStrings( papszHDR, i, 1, NULL );
689 : }
690 : }
691 :
692 : /* -------------------------------------------------------------------- */
693 : /* Set the transformation information. */
694 : /* -------------------------------------------------------------------- */
695 32 : CPLString oValue;
696 :
697 32 : oValue.Printf( "%.15g", adfGeoTransform[0] + adfGeoTransform[1] * 0.5 );
698 32 : ResetKeyValue( "ULXMAP", oValue );
699 :
700 32 : oValue.Printf( "%.15g", adfGeoTransform[3] + adfGeoTransform[5] * 0.5 );
701 32 : ResetKeyValue( "ULYMAP", oValue );
702 :
703 32 : oValue.Printf( "%.15g", adfGeoTransform[1] );
704 32 : ResetKeyValue( "XDIM", oValue );
705 :
706 32 : oValue.Printf( "%.15g", fabs(adfGeoTransform[5]) );
707 32 : ResetKeyValue( "YDIM", oValue );
708 :
709 32 : return CE_None;
710 : }
711 :
712 : /************************************************************************/
713 : /* RewriteHDR() */
714 : /************************************************************************/
715 :
716 33 : CPLErr EHdrDataset::RewriteHDR()
717 :
718 : {
719 33 : CPLString osPath = CPLGetPath( GetDescription() );
720 33 : CPLString osName = CPLGetBasename( GetDescription() );
721 33 : CPLString osHDRFilename = CPLFormCIFilename( osPath, osName, osHeaderExt );
722 :
723 : /* -------------------------------------------------------------------- */
724 : /* Write .hdr file. */
725 : /* -------------------------------------------------------------------- */
726 : VSILFILE *fp;
727 : int i;
728 :
729 33 : fp = VSIFOpenL( osHDRFilename, "wt" );
730 :
731 33 : if( fp == NULL )
732 : {
733 : CPLError( CE_Failure, CPLE_OpenFailed,
734 : "Failed to rewrite .hdr file %s.",
735 0 : osHDRFilename.c_str() );
736 0 : return CE_Failure;
737 : }
738 :
739 459 : for( i = 0; papszHDR[i] != NULL; i++ )
740 : {
741 426 : VSIFWriteL( papszHDR[i], 1, strlen(papszHDR[i]), fp );
742 426 : VSIFWriteL( (void *) "\n", 1, 1, fp );
743 : }
744 :
745 33 : VSIFCloseL( fp );
746 :
747 33 : bHDRDirty = FALSE;
748 :
749 33 : return CE_None;
750 : }
751 :
752 : /************************************************************************/
753 : /* RewriteSTX() */
754 : /************************************************************************/
755 :
756 2 : CPLErr EHdrDataset::RewriteSTX()
757 : {
758 2 : CPLString osPath = CPLGetPath( GetDescription() );
759 2 : CPLString osName = CPLGetBasename( GetDescription() );
760 2 : CPLString osSTXFilename = CPLFormCIFilename( osPath, osName, "stx" );
761 :
762 : /* -------------------------------------------------------------------- */
763 : /* Write .stx file. */
764 : /* -------------------------------------------------------------------- */
765 : VSILFILE *fp;
766 2 : fp = VSIFOpenL( osSTXFilename, "wt" );
767 2 : if( fp == NULL )
768 : {
769 : CPLDebug( "EHDR", "Failed to rewrite .stx file %s.",
770 0 : osSTXFilename.c_str() );
771 0 : return CE_Failure;
772 : }
773 :
774 4 : for (int i = 0; i < nBands; ++i)
775 : {
776 2 : EHdrRasterBand* poBand = (EHdrRasterBand*)papoBands[i];
777 2 : VSIFPrintfL( fp, "%d %.10f %.10f ", i+1, poBand->dfMin, poBand->dfMax );
778 2 : if ( poBand->minmaxmeanstddev & HAS_MEAN_FLAG )
779 2 : VSIFPrintfL( fp, "%.10f ", poBand->dfMean);
780 : else
781 0 : VSIFPrintfL( fp, "# ");
782 :
783 2 : if ( poBand->minmaxmeanstddev & HAS_STDDEV_FLAG )
784 2 : VSIFPrintfL( fp, "%.10f\n", poBand->dfStdDev);
785 : else
786 0 : VSIFPrintfL( fp, "#\n");
787 : }
788 :
789 2 : VSIFCloseL( fp );
790 :
791 2 : return CE_None;
792 : }
793 :
794 : /************************************************************************/
795 : /* ReadSTX() */
796 : /************************************************************************/
797 :
798 76 : CPLErr EHdrDataset::ReadSTX()
799 : {
800 76 : CPLString osPath = CPLGetPath( GetDescription() );
801 76 : CPLString osName = CPLGetBasename( GetDescription() );
802 76 : CPLString osSTXFilename = CPLFormCIFilename( osPath, osName, "stx" );
803 :
804 : /* -------------------------------------------------------------------- */
805 : /* Read .stx file. */
806 : /* -------------------------------------------------------------------- */
807 : VSILFILE *fp;
808 76 : if ((fp = VSIFOpenL( osSTXFilename, "rt" )) != NULL)
809 : {
810 : const char * pszLine;
811 6 : while( (pszLine = CPLReadLineL( fp )) != NULL )
812 : {
813 : char **papszTokens;
814 2 : papszTokens = CSLTokenizeStringComplex( pszLine, " \t", TRUE, FALSE );
815 2 : int nTokens = CSLCount( papszTokens );
816 2 : if( nTokens >= 5 )
817 : {
818 2 : int i = atoi(papszTokens[0]);
819 2 : if (i > 0 && i <= nBands)
820 : {
821 2 : EHdrRasterBand* poBand = (EHdrRasterBand*)papoBands[i-1];
822 2 : poBand->dfMin = atof(papszTokens[1]);
823 2 : poBand->dfMax = atof(papszTokens[2]);
824 :
825 2 : int bNoDataSet = FALSE;
826 2 : double dfNoData = poBand->GetNoDataValue(&bNoDataSet);
827 2 : if (bNoDataSet && dfNoData == poBand->dfMin)
828 : {
829 : /* Triggered by /vsicurl/http://eros.usgs.gov/archive/nslrsda/GeoTowns/HongKong/srtm/n22e113.zip/n22e113.bil */
830 : CPLDebug("EHDr", "Ignoring .stx file where min == nodata. "
831 : "The nodata value shouldn't be taken into account "
832 0 : "in minimum value computation.");
833 0 : CSLDestroy( papszTokens );
834 0 : papszTokens = NULL;
835 0 : break;
836 : }
837 :
838 2 : poBand->minmaxmeanstddev = HAS_MIN_FLAG | HAS_MAX_FLAG;
839 : // reads optional mean and stddev
840 2 : if ( !EQUAL(papszTokens[3], "#") )
841 : {
842 2 : poBand->dfMean = atof(papszTokens[3]);
843 2 : poBand->minmaxmeanstddev |= HAS_MEAN_FLAG;
844 : }
845 2 : if ( !EQUAL(papszTokens[4], "#") )
846 : {
847 2 : poBand->dfStdDev = atof(papszTokens[4]);
848 2 : poBand->minmaxmeanstddev |= HAS_STDDEV_FLAG;
849 : }
850 :
851 2 : if( nTokens >= 6 && !EQUAL(papszTokens[5], "#") )
852 0 : poBand->SetMetadataItem( "STRETCHMIN", papszTokens[5], "RENDERING_HINTS" );
853 :
854 2 : if( nTokens >= 7 && !EQUAL(papszTokens[6], "#") )
855 0 : poBand->SetMetadataItem( "STRETCHMAX", papszTokens[6], "RENDERING_HINTS" );
856 : }
857 : }
858 :
859 2 : CSLDestroy( papszTokens );
860 : }
861 :
862 2 : VSIFCloseL( fp );
863 : }
864 :
865 76 : return CE_None;
866 : }
867 :
868 :
869 : /************************************************************************/
870 : /* GetImageRepFilename() */
871 : /************************************************************************/
872 :
873 : /* -------------------------------------------------------------------- */
874 : /* Check for IMAGE.REP (Spatiocarte Defense 1.0) or name_of_image.rep */
875 : /* if it's a GIS-GeoSPOT image */
876 : /* For the specification of SPDF (in French), */
877 : /* see http://eden.ign.fr/download/pub/doc/emabgi/spdf10.pdf/download */
878 : /* -------------------------------------------------------------------- */
879 :
880 53 : CPLString EHdrDataset::GetImageRepFilename(const char* pszFilename)
881 : {
882 : VSIStatBufL sStatBuf;
883 :
884 53 : CPLString osPath = CPLGetPath( pszFilename );
885 53 : CPLString osName = CPLGetBasename( pszFilename );
886 : CPLString osREPFilename =
887 53 : CPLFormCIFilename( osPath, osName, "rep" );
888 53 : if( VSIStatExL( osREPFilename.c_str(), &sStatBuf, VSI_STAT_EXISTS_FLAG ) == 0 )
889 0 : return osREPFilename;
890 :
891 53 : if (EQUAL(CPLGetFilename(pszFilename), "imspatio.bil") ||
892 : EQUAL(CPLGetFilename(pszFilename), "haspatio.bil"))
893 : {
894 0 : CPLString osImageRepFilename(CPLFormCIFilename( osPath, "image", "rep" ));
895 0 : if( VSIStatExL( osImageRepFilename.c_str(), &sStatBuf, VSI_STAT_EXISTS_FLAG ) == 0 )
896 0 : return osImageRepFilename;
897 :
898 : /* Try in the upper directories if not found in the BIL image directory */
899 0 : CPLString dirName(CPLGetDirname(osPath));
900 0 : if (CPLIsFilenameRelative(osPath.c_str()))
901 : {
902 0 : char* cwd = CPLGetCurrentDir();
903 0 : if (cwd)
904 : {
905 0 : dirName = CPLFormFilename(cwd, dirName.c_str(), NULL);
906 0 : CPLFree(cwd);
907 : }
908 : }
909 0 : while (dirName[0] != 0 && EQUAL(dirName, ".") == FALSE && EQUAL(dirName, "/") == FALSE)
910 : {
911 0 : osImageRepFilename = CPLFormCIFilename( dirName.c_str(), "image", "rep" );
912 0 : if( VSIStatExL( osImageRepFilename.c_str(), &sStatBuf, VSI_STAT_EXISTS_FLAG ) == 0 )
913 0 : return osImageRepFilename;
914 :
915 : /* Don't try to recurse above the 'image' subdirectory */
916 0 : if (EQUAL(dirName, "image"))
917 : {
918 0 : break;
919 : }
920 0 : dirName = CPLString(CPLGetDirname(dirName));
921 0 : }
922 : }
923 53 : return CPLString();
924 : }
925 :
926 : /************************************************************************/
927 : /* GetFileList() */
928 : /************************************************************************/
929 :
930 9 : char **EHdrDataset::GetFileList()
931 :
932 : {
933 : VSIStatBufL sStatBuf;
934 9 : CPLString osPath = CPLGetPath( GetDescription() );
935 9 : CPLString osName = CPLGetBasename( GetDescription() );
936 9 : char **papszFileList = NULL;
937 :
938 : // Main data file, etc.
939 9 : papszFileList = GDALPamDataset::GetFileList();
940 :
941 : // Header file.
942 9 : CPLString osFilename = CPLFormCIFilename( osPath, osName, osHeaderExt );
943 9 : papszFileList = CSLAddString( papszFileList, osFilename );
944 :
945 : // Statistics file
946 9 : osFilename = CPLFormCIFilename( osPath, osName, "stx" );
947 9 : if( VSIStatExL( osFilename, &sStatBuf, VSI_STAT_EXISTS_FLAG ) == 0 )
948 1 : papszFileList = CSLAddString( papszFileList, osFilename );
949 :
950 : // color table file.
951 9 : osFilename = CPLFormCIFilename( osPath, osName, "clr" );
952 9 : if( VSIStatExL( osFilename, &sStatBuf, VSI_STAT_EXISTS_FLAG ) == 0 )
953 2 : papszFileList = CSLAddString( papszFileList, osFilename );
954 :
955 : // projections file.
956 9 : osFilename = CPLFormCIFilename( osPath, osName, "prj" );
957 9 : if( VSIStatExL( osFilename, &sStatBuf, VSI_STAT_EXISTS_FLAG ) == 0 )
958 5 : papszFileList = CSLAddString( papszFileList, osFilename );
959 :
960 9 : CPLString imageRepFilename = GetImageRepFilename( GetDescription() );
961 9 : if (!imageRepFilename.empty())
962 0 : papszFileList = CSLAddString( papszFileList, imageRepFilename.c_str() );
963 :
964 9 : return papszFileList;
965 : }
966 :
967 : /************************************************************************/
968 : /* Open() */
969 : /************************************************************************/
970 :
971 12268 : GDALDataset *EHdrDataset::Open( GDALOpenInfo * poOpenInfo )
972 :
973 : {
974 : int i, bSelectedHDR;
975 :
976 : /* -------------------------------------------------------------------- */
977 : /* We assume the user is pointing to the binary (ie. .bil) file. */
978 : /* -------------------------------------------------------------------- */
979 12268 : if( poOpenInfo->nHeaderBytes < 2 )
980 11361 : return NULL;
981 :
982 : /* -------------------------------------------------------------------- */
983 : /* Now we need to tear apart the filename to form a .HDR */
984 : /* filename. */
985 : /* -------------------------------------------------------------------- */
986 907 : CPLString osPath = CPLGetPath( poOpenInfo->pszFilename );
987 907 : CPLString osName = CPLGetBasename( poOpenInfo->pszFilename );
988 907 : CPLString osHDRFilename;
989 :
990 907 : const char* pszHeaderExt = "hdr";
991 907 : if( EQUAL( CPLGetExtension( poOpenInfo->pszFilename ), "SRC" ) &&
992 : osName.size() == 7 &&
993 : (osName[0] == 'e' || osName[0] == 'E' || osName[0] == 'w' || osName[0] == 'W') &&
994 : (osName[4] == 'n' || osName[4] == 'N' || osName[4] == 's' || osName[4] == 'S') )
995 : {
996 : /* It is a GTOPO30 or SRTM30 source file, whose header extension is .sch */
997 : /* see http://dds.cr.usgs.gov/srtm/version1/SRTM30/GTOPO30_Documentation */
998 0 : pszHeaderExt = "sch";
999 : }
1000 :
1001 907 : if( poOpenInfo->papszSiblingFiles )
1002 : {
1003 : int iFile = CSLFindString(poOpenInfo->papszSiblingFiles,
1004 903 : CPLFormFilename( NULL, osName, pszHeaderExt ) );
1005 903 : if( iFile < 0 ) // return if there is no corresponding .hdr file
1006 796 : return NULL;
1007 :
1008 : osHDRFilename =
1009 107 : CPLFormFilename( osPath, poOpenInfo->papszSiblingFiles[iFile],
1010 214 : NULL );
1011 : }
1012 : else
1013 : {
1014 4 : osHDRFilename = CPLFormCIFilename( osPath, osName, pszHeaderExt );
1015 : }
1016 :
1017 111 : bSelectedHDR = EQUAL( osHDRFilename, poOpenInfo->pszFilename );
1018 :
1019 : /* -------------------------------------------------------------------- */
1020 : /* Do we have a .hdr file? */
1021 : /* -------------------------------------------------------------------- */
1022 : VSILFILE *fp;
1023 :
1024 111 : fp = VSIFOpenL( osHDRFilename, "r" );
1025 :
1026 111 : if( fp == NULL )
1027 : {
1028 4 : return NULL;
1029 : }
1030 :
1031 : /* -------------------------------------------------------------------- */
1032 : /* Is this file an ESRI header file? Read a few lines of text */
1033 : /* searching for something starting with nrows or ncols. */
1034 : /* -------------------------------------------------------------------- */
1035 : const char * pszLine;
1036 107 : int nRows = -1, nCols = -1, nBands = 1;
1037 107 : int nSkipBytes = 0;
1038 107 : double dfULXMap=0.5, dfULYMap = 0.5, dfYLLCorner = -123.456;
1039 107 : int bCenter = TRUE;
1040 107 : double dfXDim = 1.0, dfYDim = 1.0, dfNoData = 0.0;
1041 107 : int nLineCount = 0, bNoDataSet = FALSE;
1042 107 : GDALDataType eDataType = GDT_Byte;
1043 107 : int nBits = -1;
1044 107 : char chByteOrder = 'M';
1045 107 : char chPixelType = 'N'; // not defined
1046 107 : char szLayout[10] = "BIL";
1047 107 : char **papszHDR = NULL;
1048 107 : int bHasInternalProjection = FALSE;
1049 107 : int bHasMin = FALSE;
1050 107 : int bHasMax = FALSE;
1051 107 : double dfMin = 0, dfMax = 0;
1052 :
1053 1375 : while( (pszLine = CPLReadLineL( fp )) != NULL )
1054 : {
1055 : char **papszTokens;
1056 :
1057 1163 : nLineCount++;
1058 :
1059 1163 : if( nLineCount > 50 || strlen(pszLine) > 1000 )
1060 2 : break;
1061 :
1062 1161 : papszHDR = CSLAddString( papszHDR, pszLine );
1063 :
1064 1161 : papszTokens = CSLTokenizeStringComplex( pszLine, " \t", TRUE, FALSE );
1065 1161 : if( CSLCount( papszTokens ) < 2 )
1066 : {
1067 49 : CSLDestroy( papszTokens );
1068 49 : continue;
1069 : }
1070 :
1071 1112 : if( EQUAL(papszTokens[0],"ncols") )
1072 : {
1073 76 : nCols = atoi(papszTokens[1]);
1074 : }
1075 1036 : else if( EQUAL(papszTokens[0],"nrows") )
1076 : {
1077 76 : nRows = atoi(papszTokens[1]);
1078 : }
1079 960 : else if( EQUAL(papszTokens[0],"skipbytes") )
1080 : {
1081 0 : nSkipBytes = atoi(papszTokens[1]);
1082 : }
1083 2844 : else if( EQUAL(papszTokens[0],"ulxmap")
1084 924 : || EQUAL(papszTokens[0],"xllcorner")
1085 922 : || EQUAL(papszTokens[0],"xllcenter") )
1086 : {
1087 38 : dfULXMap = CPLAtofM(papszTokens[1]);
1088 38 : if( EQUAL(papszTokens[0],"xllcorner") )
1089 2 : bCenter = FALSE;
1090 : }
1091 922 : else if( EQUAL(papszTokens[0],"ulymap") )
1092 : {
1093 36 : dfULYMap = CPLAtofM(papszTokens[1]);
1094 : }
1095 1772 : else if( EQUAL(papszTokens[0],"yllcorner")
1096 884 : || EQUAL(papszTokens[0],"yllcenter") )
1097 : {
1098 2 : dfYLLCorner = CPLAtofM(papszTokens[1]);
1099 2 : if( EQUAL(papszTokens[0],"yllcorner") )
1100 2 : bCenter = FALSE;
1101 : }
1102 884 : else if( EQUAL(papszTokens[0],"xdim") )
1103 : {
1104 36 : dfXDim = CPLAtofM(papszTokens[1]);
1105 : }
1106 848 : else if( EQUAL(papszTokens[0],"ydim") )
1107 : {
1108 36 : dfYDim = CPLAtofM(papszTokens[1]);
1109 : }
1110 812 : else if( EQUAL(papszTokens[0],"cellsize") )
1111 : {
1112 2 : dfXDim = dfYDim = CPLAtofM(papszTokens[1]);
1113 : }
1114 810 : else if( EQUAL(papszTokens[0],"nbands") )
1115 : {
1116 74 : nBands = atoi(papszTokens[1]);
1117 : }
1118 736 : else if( EQUAL(papszTokens[0],"layout") )
1119 : {
1120 74 : strncpy( szLayout, papszTokens[1], sizeof(szLayout) );
1121 74 : szLayout[sizeof(szLayout)-1] = '\0';
1122 : }
1123 1328 : else if( EQUAL(papszTokens[0],"NODATA_value")
1124 662 : || EQUAL(papszTokens[0],"NODATA") )
1125 : {
1126 4 : dfNoData = CPLAtofM(papszTokens[1]);
1127 4 : bNoDataSet = TRUE;
1128 : }
1129 658 : else if( EQUAL(papszTokens[0],"NBITS") )
1130 : {
1131 75 : nBits = atoi(papszTokens[1]);
1132 : }
1133 583 : else if( EQUAL(papszTokens[0],"PIXELTYPE") )
1134 : {
1135 72 : chPixelType = (char) toupper(papszTokens[1][0]);
1136 : }
1137 511 : else if( EQUAL(papszTokens[0],"byteorder") )
1138 : {
1139 76 : chByteOrder = (char) toupper(papszTokens[1][0]);
1140 : }
1141 :
1142 : /* http://www.worldclim.org/futdown.htm have the projection extensions */
1143 435 : else if( EQUAL(papszTokens[0],"Projection") )
1144 : {
1145 3 : bHasInternalProjection = TRUE;
1146 : }
1147 864 : else if( EQUAL(papszTokens[0],"MinValue") ||
1148 431 : EQUAL(papszTokens[0],"MIN_VALUE") )
1149 : {
1150 1 : dfMin = CPLAtofM(papszTokens[1]);
1151 1 : bHasMin = TRUE;
1152 : }
1153 861 : else if( EQUAL(papszTokens[0],"MaxValue") ||
1154 430 : EQUAL(papszTokens[0],"MAX_VALUE") )
1155 : {
1156 1 : dfMax = CPLAtofM(papszTokens[1]);
1157 1 : bHasMax = TRUE;
1158 : }
1159 :
1160 1112 : CSLDestroy( papszTokens );
1161 : }
1162 :
1163 107 : VSIFCloseL( fp );
1164 :
1165 : /* -------------------------------------------------------------------- */
1166 : /* Did we get the required keywords? If not we return with */
1167 : /* this never having been considered to be a match. This isn't */
1168 : /* an error! */
1169 : /* -------------------------------------------------------------------- */
1170 107 : if( nRows == -1 || nCols == -1 )
1171 : {
1172 31 : CSLDestroy( papszHDR );
1173 31 : return NULL;
1174 : }
1175 :
1176 76 : if (!GDALCheckDatasetDimensions(nCols, nRows) ||
1177 : !GDALCheckBandCount(nBands, FALSE))
1178 : {
1179 0 : CSLDestroy( papszHDR );
1180 0 : return NULL;
1181 : }
1182 :
1183 : /* -------------------------------------------------------------------- */
1184 : /* Has the user selected the .hdr file to open? */
1185 : /* -------------------------------------------------------------------- */
1186 76 : if( bSelectedHDR )
1187 : {
1188 : CPLError( CE_Failure, CPLE_AppDefined,
1189 : "The selected file is an ESRI BIL header file, but to\n"
1190 : "open ESRI BIL datasets, the data file should be selected\n"
1191 : "instead of the .hdr file. Please try again selecting\n"
1192 : "the data file (often with the extension .bil) corresponding\n"
1193 : "to the header file: %s\n",
1194 0 : poOpenInfo->pszFilename );
1195 0 : CSLDestroy( papszHDR );
1196 0 : return NULL;
1197 : }
1198 :
1199 : /* -------------------------------------------------------------------- */
1200 : /* If we aren't sure of the file type, check the data file */
1201 : /* size. If it is 4 bytes or more per pixel then we assume it */
1202 : /* is floating point data. */
1203 : /* -------------------------------------------------------------------- */
1204 76 : if( nBits == -1 && chPixelType == 'N' )
1205 : {
1206 : VSIStatBufL sStatBuf;
1207 1 : if( VSIStatL( poOpenInfo->pszFilename, &sStatBuf ) == 0 )
1208 : {
1209 1 : size_t nBytes = (size_t) (sStatBuf.st_size/nCols/nRows/nBands);
1210 1 : if( nBytes > 0 && nBytes != 3 )
1211 1 : nBits = nBytes*8;
1212 :
1213 1 : if( nBytes == 4 )
1214 1 : chPixelType = 'F';
1215 : }
1216 : }
1217 :
1218 : /* -------------------------------------------------------------------- */
1219 : /* If the extension is FLT it is likely a floating point file. */
1220 : /* -------------------------------------------------------------------- */
1221 76 : if( chPixelType == 'N' )
1222 : {
1223 3 : if( EQUAL( CPLGetExtension( poOpenInfo->pszFilename ), "FLT" ) )
1224 1 : chPixelType = 'F';
1225 : }
1226 :
1227 : /* -------------------------------------------------------------------- */
1228 : /* If we have a negative nodata value, let's assume that the */
1229 : /* pixel type is signed. This is necessary for datasets from */
1230 : /* http://www.worldclim.org/futdown.htm */
1231 : /* -------------------------------------------------------------------- */
1232 76 : if( bNoDataSet && dfNoData < 0 && chPixelType == 'N' )
1233 : {
1234 2 : chPixelType = 'S';
1235 : }
1236 :
1237 : /* -------------------------------------------------------------------- */
1238 : /* Create a corresponding GDALDataset. */
1239 : /* -------------------------------------------------------------------- */
1240 : EHdrDataset *poDS;
1241 :
1242 76 : poDS = new EHdrDataset();
1243 :
1244 152 : poDS->osHeaderExt = pszHeaderExt;
1245 :
1246 : /* -------------------------------------------------------------------- */
1247 : /* Capture some information from the file that is of interest. */
1248 : /* -------------------------------------------------------------------- */
1249 76 : poDS->nRasterXSize = nCols;
1250 76 : poDS->nRasterYSize = nRows;
1251 76 : poDS->papszHDR = papszHDR;
1252 :
1253 : /* -------------------------------------------------------------------- */
1254 : /* Open target binary file. */
1255 : /* -------------------------------------------------------------------- */
1256 76 : if( poOpenInfo->eAccess == GA_ReadOnly )
1257 42 : poDS->fpImage = VSIFOpenL( poOpenInfo->pszFilename, "rb" );
1258 : else
1259 34 : poDS->fpImage = VSIFOpenL( poOpenInfo->pszFilename, "r+b" );
1260 :
1261 76 : if( poDS->fpImage == NULL )
1262 : {
1263 : CPLError( CE_Failure, CPLE_OpenFailed,
1264 : "Failed to open %s with write permission.\n%s",
1265 0 : osName.c_str(), VSIStrerror( errno ) );
1266 0 : delete poDS;
1267 0 : return NULL;
1268 : }
1269 :
1270 76 : poDS->eAccess = poOpenInfo->eAccess;
1271 :
1272 : /* -------------------------------------------------------------------- */
1273 : /* Figure out the data type. */
1274 : /* -------------------------------------------------------------------- */
1275 76 : if( nBits == 16 )
1276 : {
1277 15 : if ( chPixelType == 'S' )
1278 7 : eDataType = GDT_Int16;
1279 : else
1280 8 : eDataType = GDT_UInt16; // default
1281 : }
1282 61 : else if( nBits == 32 )
1283 : {
1284 31 : if( chPixelType == 'S' )
1285 8 : eDataType = GDT_Int32;
1286 23 : else if( chPixelType == 'F' )
1287 18 : eDataType = GDT_Float32;
1288 : else
1289 5 : eDataType = GDT_UInt32; // default
1290 : }
1291 30 : else if( nBits == 8 )
1292 : {
1293 27 : eDataType = GDT_Byte;
1294 27 : nBits = 8;
1295 : }
1296 6 : else if( nBits < 8 && nBits >= 1 )
1297 3 : eDataType = GDT_Byte;
1298 0 : else if( nBits == -1 )
1299 : {
1300 0 : if( chPixelType == 'F' )
1301 : {
1302 0 : eDataType = GDT_Float32;
1303 0 : nBits = 32;
1304 : }
1305 : else
1306 : {
1307 0 : eDataType = GDT_Byte;
1308 0 : nBits = 8;
1309 : }
1310 : }
1311 : else
1312 : {
1313 : CPLError( CE_Failure, CPLE_NotSupported,
1314 : "EHdr driver does not support %d NBITS value.",
1315 0 : nBits );
1316 0 : delete poDS;
1317 0 : return NULL;
1318 : }
1319 :
1320 : /* -------------------------------------------------------------------- */
1321 : /* Compute the line offset. */
1322 : /* -------------------------------------------------------------------- */
1323 76 : int nItemSize = GDALGetDataTypeSize(eDataType)/8;
1324 : int nPixelOffset, nLineOffset;
1325 : vsi_l_offset nBandOffset;
1326 :
1327 76 : if( EQUAL(szLayout,"BIP") )
1328 : {
1329 0 : nPixelOffset = nItemSize * nBands;
1330 0 : nLineOffset = nPixelOffset * nCols;
1331 0 : nBandOffset = (vsi_l_offset)nItemSize;
1332 : }
1333 76 : else if( EQUAL(szLayout,"BSQ") )
1334 : {
1335 0 : nPixelOffset = nItemSize;
1336 0 : nLineOffset = nPixelOffset * nCols;
1337 0 : nBandOffset = (vsi_l_offset)nLineOffset * nRows;
1338 : }
1339 : else /* assume BIL */
1340 : {
1341 76 : nPixelOffset = nItemSize;
1342 76 : nLineOffset = nItemSize * nBands * nCols;
1343 76 : nBandOffset = (vsi_l_offset)nItemSize * nCols;
1344 : }
1345 :
1346 76 : poDS->SetDescription( poOpenInfo->pszFilename );
1347 76 : poDS->PamInitialize();
1348 :
1349 : /* -------------------------------------------------------------------- */
1350 : /* Create band information objects. */
1351 : /* -------------------------------------------------------------------- */
1352 76 : poDS->nBands = nBands;
1353 202 : for( i = 0; i < poDS->nBands; i++ )
1354 : {
1355 : EHdrRasterBand *poBand;
1356 :
1357 : poBand =
1358 : new EHdrRasterBand( poDS, i+1, poDS->fpImage,
1359 : nSkipBytes + nBandOffset * i,
1360 : nPixelOffset, nLineOffset, eDataType,
1361 : #ifdef CPL_LSB
1362 : chByteOrder == 'I' || chByteOrder == 'L',
1363 : #else
1364 : chByteOrder == 'M',
1365 : #endif
1366 126 : nBits);
1367 :
1368 126 : poBand->bNoDataSet = bNoDataSet;
1369 126 : poBand->dfNoData = dfNoData;
1370 :
1371 127 : if( bHasMin && bHasMax )
1372 : {
1373 1 : poBand->dfMin = dfMin;
1374 1 : poBand->dfMax = dfMax;
1375 1 : poBand->minmaxmeanstddev = HAS_MIN_FLAG | HAS_MAX_FLAG;
1376 : }
1377 :
1378 126 : poDS->SetBand( i+1, poBand );
1379 : }
1380 :
1381 : /* -------------------------------------------------------------------- */
1382 : /* If we didn't get bounds in the .hdr, look for a worldfile. */
1383 : /* -------------------------------------------------------------------- */
1384 76 : if( dfYLLCorner != -123.456 )
1385 : {
1386 2 : if( bCenter )
1387 0 : dfULYMap = dfYLLCorner + (nRows-1) * dfYDim;
1388 : else
1389 2 : dfULYMap = dfYLLCorner + nRows * dfYDim;
1390 : }
1391 :
1392 76 : if( dfULYMap != 0.5 || dfULYMap != 0.5 || dfXDim != 1.0 || dfYDim != 1.0 )
1393 : {
1394 38 : poDS->bGotTransform = TRUE;
1395 :
1396 38 : if( bCenter )
1397 : {
1398 36 : poDS->adfGeoTransform[0] = dfULXMap - dfXDim * 0.5;
1399 36 : poDS->adfGeoTransform[1] = dfXDim;
1400 36 : poDS->adfGeoTransform[2] = 0.0;
1401 36 : poDS->adfGeoTransform[3] = dfULYMap + dfYDim * 0.5;
1402 36 : poDS->adfGeoTransform[4] = 0.0;
1403 36 : poDS->adfGeoTransform[5] = - dfYDim;
1404 : }
1405 : else
1406 : {
1407 2 : poDS->adfGeoTransform[0] = dfULXMap;
1408 2 : poDS->adfGeoTransform[1] = dfXDim;
1409 2 : poDS->adfGeoTransform[2] = 0.0;
1410 2 : poDS->adfGeoTransform[3] = dfULYMap;
1411 2 : poDS->adfGeoTransform[4] = 0.0;
1412 2 : poDS->adfGeoTransform[5] = - dfYDim;
1413 : }
1414 : }
1415 :
1416 76 : if( !poDS->bGotTransform )
1417 : poDS->bGotTransform =
1418 : GDALReadWorldFile( poOpenInfo->pszFilename, 0,
1419 38 : poDS->adfGeoTransform );
1420 :
1421 76 : if( !poDS->bGotTransform )
1422 : poDS->bGotTransform =
1423 : GDALReadWorldFile( poOpenInfo->pszFilename, "wld",
1424 38 : poDS->adfGeoTransform );
1425 :
1426 : /* -------------------------------------------------------------------- */
1427 : /* Check for a .prj file. */
1428 : /* -------------------------------------------------------------------- */
1429 76 : const char *pszPrjFilename = CPLFormCIFilename( osPath, osName, "prj" );
1430 :
1431 76 : fp = VSIFOpenL( pszPrjFilename, "r" );
1432 :
1433 : /* .hdr files from http://www.worldclim.org/futdown.htm have the projection */
1434 : /* info in the .hdr file itself ! */
1435 76 : if (fp == NULL && bHasInternalProjection)
1436 : {
1437 1 : pszPrjFilename = osHDRFilename;
1438 1 : fp = VSIFOpenL( pszPrjFilename, "r" );
1439 : }
1440 :
1441 76 : if( fp != NULL )
1442 : {
1443 : char **papszLines;
1444 32 : OGRSpatialReference oSRS;
1445 :
1446 32 : VSIFCloseL( fp );
1447 :
1448 32 : papszLines = CSLLoad( pszPrjFilename );
1449 :
1450 32 : if( oSRS.importFromESRI( papszLines ) == OGRERR_NONE )
1451 : {
1452 : // If geographic values are in seconds, we must transform.
1453 : // Is there a code for minutes too?
1454 : char szResult[80];
1455 32 : if( oSRS.IsGeographic()
1456 : && EQUAL(OSR_GDS( szResult, sizeof(szResult),
1457 : papszLines, "Units", ""), "DS") )
1458 : {
1459 0 : poDS->adfGeoTransform[0] /= 3600.0;
1460 0 : poDS->adfGeoTransform[1] /= 3600.0;
1461 0 : poDS->adfGeoTransform[2] /= 3600.0;
1462 0 : poDS->adfGeoTransform[3] /= 3600.0;
1463 0 : poDS->adfGeoTransform[4] /= 3600.0;
1464 0 : poDS->adfGeoTransform[5] /= 3600.0;
1465 : }
1466 :
1467 32 : CPLFree( poDS->pszProjection );
1468 32 : oSRS.exportToWkt( &(poDS->pszProjection) );
1469 : }
1470 :
1471 32 : CSLDestroy( papszLines );
1472 : }
1473 : else
1474 : {
1475 : /* -------------------------------------------------------------------- */
1476 : /* Check for IMAGE.REP (Spatiocarte Defense 1.0) or name_of_image.rep */
1477 : /* if it's a GIS-GeoSPOT image */
1478 : /* For the specification of SPDF (in French), */
1479 : /* see http://eden.ign.fr/download/pub/doc/emabgi/spdf10.pdf/download */
1480 : /* -------------------------------------------------------------------- */
1481 44 : CPLString szImageRepFilename = GetImageRepFilename(poOpenInfo->pszFilename );
1482 44 : if (!szImageRepFilename.empty())
1483 : {
1484 0 : fp = VSIFOpenL( szImageRepFilename.c_str(), "r" );
1485 : }
1486 44 : if (fp != NULL)
1487 : {
1488 : const char *pszLine;
1489 0 : int bUTM = FALSE;
1490 0 : int bWGS84 = FALSE;
1491 0 : int bNorth = FALSE;
1492 0 : int bSouth = FALSE;
1493 0 : int utmZone = 0;
1494 :
1495 0 : while( (pszLine = CPLReadLineL( fp )) != NULL )
1496 : {
1497 0 : if (strncmp(pszLine, "PROJ_ID", strlen("PROJ_ID")) == 0 &&
1498 : strstr(pszLine, "UTM"))
1499 : {
1500 0 : bUTM = TRUE;
1501 : }
1502 0 : else if (strncmp(pszLine, "PROJ_ZONE", strlen("PROJ_ZONE")) == 0)
1503 : {
1504 0 : const char* c = strchr(pszLine, '"');
1505 0 : if (c)
1506 : {
1507 0 : c++;
1508 0 : if (*c >= '0' && *c <= '9')
1509 : {
1510 0 : utmZone = atoi(c);
1511 0 : if (utmZone >= 1 && utmZone <= 60)
1512 : {
1513 0 : if (strstr(pszLine, "Nord") || strstr(pszLine, "NORD"))
1514 : {
1515 0 : bNorth = TRUE;
1516 : }
1517 0 : else if (strstr(pszLine, "Sud") || strstr(pszLine, "SUD"))
1518 : {
1519 0 : bSouth = TRUE;
1520 : }
1521 : }
1522 : }
1523 : }
1524 : }
1525 0 : else if (strncmp(pszLine, "PROJ_CODE", strlen("PROJ_CODE")) == 0 &&
1526 : strstr(pszLine, "FR-MINDEF"))
1527 : {
1528 0 : const char* c = strchr(pszLine, 'A');
1529 0 : if (c)
1530 : {
1531 0 : c++;
1532 0 : if (*c >= '0' && *c <= '9')
1533 : {
1534 0 : utmZone = atoi(c);
1535 0 : if (utmZone >= 1 && utmZone <= 60)
1536 : {
1537 0 : if (c[1] == 'N' ||
1538 0 : (c[1] != '\0' && c[2] == 'N'))
1539 : {
1540 0 : bNorth = TRUE;
1541 : }
1542 0 : else if (c[1] == 'S' ||
1543 0 : (c[1] != '\0' && c[2] == 'S'))
1544 : {
1545 0 : bSouth = TRUE;
1546 : }
1547 : }
1548 : }
1549 : }
1550 : }
1551 0 : else if (strncmp(pszLine, "HORIZ_DATUM", strlen("HORIZ_DATUM")) == 0 &&
1552 : (strstr(pszLine, "WGS 84") || strstr(pszLine, "WGS84")))
1553 : {
1554 0 : bWGS84 = TRUE;
1555 : }
1556 0 : else if (strncmp(pszLine, "MAP_NUMBER", strlen("MAP_NUMBER")) == 0)
1557 : {
1558 0 : const char* c = strchr(pszLine, '"');
1559 0 : if (c)
1560 : {
1561 0 : char* pszMapNumber = CPLStrdup(c+1);
1562 0 : char* c2 = strchr(pszMapNumber, '"');
1563 0 : if (c2) *c2 = 0;
1564 0 : poDS->SetMetadataItem("SPDF_MAP_NUMBER", pszMapNumber);
1565 0 : CPLFree(pszMapNumber);
1566 : }
1567 : }
1568 0 : else if (strncmp(pszLine, "PRODUCTION_DATE", strlen("PRODUCTION_DATE")) == 0)
1569 : {
1570 0 : const char* c = pszLine + strlen("PRODUCTION_DATE");
1571 0 : while(*c == ' ')
1572 0 : c++;
1573 0 : if (*c)
1574 : {
1575 0 : poDS->SetMetadataItem("SPDF_PRODUCTION_DATE", c );
1576 : }
1577 : }
1578 : }
1579 :
1580 0 : VSIFCloseL( fp );
1581 :
1582 0 : if (utmZone != 0 && bUTM && bWGS84 && (bNorth || bSouth))
1583 : {
1584 : char projCSStr[64];
1585 0 : OGRSpatialReference oSRS;
1586 :
1587 : sprintf(projCSStr, "WGS 84 / UTM zone %d%c",
1588 0 : utmZone, (bNorth) ? 'N' : 'S');
1589 0 : oSRS.SetProjCS(projCSStr);
1590 0 : oSRS.SetWellKnownGeogCS( "WGS84" );
1591 0 : oSRS.SetUTM(utmZone, bNorth);
1592 0 : oSRS.SetAuthority("PROJCS", "EPSG", ((bNorth) ? 32600 : 32700) + utmZone);
1593 0 : oSRS.AutoIdentifyEPSG();
1594 :
1595 0 : CPLFree( poDS->pszProjection );
1596 0 : oSRS.exportToWkt( &(poDS->pszProjection) );
1597 : }
1598 : else
1599 : {
1600 0 : CPLError( CE_Warning, CPLE_NotSupported, "Cannot retrive projection from IMAGE.REP");
1601 : }
1602 44 : }
1603 : }
1604 :
1605 : /* -------------------------------------------------------------------- */
1606 : /* Check for a color table. */
1607 : /* -------------------------------------------------------------------- */
1608 76 : const char *pszCLRFilename = CPLFormCIFilename( osPath, osName, "clr" );
1609 :
1610 : /* Only read the .clr for byte, int16 or uint16 bands */
1611 76 : if (nItemSize <= 2)
1612 45 : fp = VSIFOpenL( pszCLRFilename, "r" );
1613 : else
1614 31 : fp = NULL;
1615 :
1616 76 : if( fp != NULL )
1617 : {
1618 4 : GDALColorTable oColorTable;
1619 4 : int bHasWarned = FALSE;
1620 :
1621 12 : while(TRUE)
1622 : {
1623 16 : const char *pszLine = CPLReadLineL(fp);
1624 16 : if ( !pszLine )
1625 : break;
1626 :
1627 12 : if( *pszLine == '#' || *pszLine == '!' )
1628 0 : continue;
1629 :
1630 : char **papszValues = CSLTokenizeString2(pszLine, "\t ",
1631 12 : CSLT_HONOURSTRINGS);
1632 : GDALColorEntry oEntry;
1633 :
1634 12 : if ( CSLCount(papszValues) >= 4 )
1635 : {
1636 12 : int nIndex = atoi( papszValues[0] ); // Index
1637 24 : if (nIndex >= 0 && nIndex < 65536)
1638 : {
1639 12 : oEntry.c1 = (short) atoi( papszValues[1] ); // Red
1640 12 : oEntry.c2 = (short) atoi( papszValues[2] ); // Green
1641 12 : oEntry.c3 = (short) atoi( papszValues[3] ); // Blue
1642 12 : oEntry.c4 = 255;
1643 :
1644 12 : oColorTable.SetColorEntry( nIndex, &oEntry );
1645 : }
1646 : else
1647 : {
1648 : /* Negative values are valid. At least we can find use of */
1649 : /* them here : http://www.ngdc.noaa.gov/mgg/topo/elev/esri/clr/ */
1650 : /* but there's no way of representing them with GDAL color */
1651 : /* table model */
1652 0 : if (!bHasWarned)
1653 0 : CPLDebug("EHdr", "Ignoring color index : %d", nIndex);
1654 0 : bHasWarned = TRUE;
1655 : }
1656 : }
1657 :
1658 12 : CSLDestroy( papszValues );
1659 : }
1660 :
1661 4 : VSIFCloseL( fp );
1662 :
1663 8 : for( i = 1; i <= poDS->nBands; i++ )
1664 : {
1665 4 : GDALRasterBand *poBand = poDS->GetRasterBand( i );
1666 4 : poBand->SetColorTable( &oColorTable );
1667 4 : poBand->SetColorInterpretation( GCI_PaletteIndex );
1668 : }
1669 :
1670 4 : poDS->bCLRDirty = FALSE;
1671 : }
1672 :
1673 : /* -------------------------------------------------------------------- */
1674 : /* Read statistics (.STX) */
1675 : /* -------------------------------------------------------------------- */
1676 76 : poDS->ReadSTX();
1677 :
1678 : /* -------------------------------------------------------------------- */
1679 : /* Initialize any PAM information. */
1680 : /* -------------------------------------------------------------------- */
1681 76 : poDS->TryLoadXML();
1682 :
1683 : /* -------------------------------------------------------------------- */
1684 : /* Check for overviews. */
1685 : /* -------------------------------------------------------------------- */
1686 76 : poDS->oOvManager.Initialize( poDS, poOpenInfo->pszFilename );
1687 :
1688 76 : return( poDS );
1689 : }
1690 :
1691 : /************************************************************************/
1692 : /* Create() */
1693 : /************************************************************************/
1694 :
1695 62 : GDALDataset *EHdrDataset::Create( const char * pszFilename,
1696 : int nXSize, int nYSize, int nBands,
1697 : GDALDataType eType,
1698 : char **papszParmList )
1699 :
1700 : {
1701 : /* -------------------------------------------------------------------- */
1702 : /* Verify input options. */
1703 : /* -------------------------------------------------------------------- */
1704 62 : if (nBands <= 0)
1705 : {
1706 : CPLError( CE_Failure, CPLE_NotSupported,
1707 1 : "EHdr driver does not support %d bands.\n", nBands);
1708 1 : return NULL;
1709 : }
1710 :
1711 61 : if( eType != GDT_Byte && eType != GDT_Float32 && eType != GDT_UInt16
1712 : && eType != GDT_Int16 && eType != GDT_Int32 && eType != GDT_UInt32 )
1713 : {
1714 : CPLError( CE_Failure, CPLE_AppDefined,
1715 : "Attempt to create ESRI .hdr labelled dataset with an illegal\n"
1716 : "data type (%s).\n",
1717 15 : GDALGetDataTypeName(eType) );
1718 :
1719 15 : return NULL;
1720 : }
1721 :
1722 : /* -------------------------------------------------------------------- */
1723 : /* Try to create the file. */
1724 : /* -------------------------------------------------------------------- */
1725 : VSILFILE *fp;
1726 :
1727 46 : fp = VSIFOpenL( pszFilename, "wb" );
1728 :
1729 46 : if( fp == NULL )
1730 : {
1731 : CPLError( CE_Failure, CPLE_OpenFailed,
1732 : "Attempt to create file `%s' failed.\n",
1733 12 : pszFilename );
1734 12 : return NULL;
1735 : }
1736 :
1737 : /* -------------------------------------------------------------------- */
1738 : /* Just write out a couple of bytes to establish the binary */
1739 : /* file, and then close it. */
1740 : /* -------------------------------------------------------------------- */
1741 34 : VSIFWriteL( (void *) "\0\0", 2, 1, fp );
1742 34 : VSIFCloseL( fp );
1743 :
1744 : /* -------------------------------------------------------------------- */
1745 : /* Create the hdr filename. */
1746 : /* -------------------------------------------------------------------- */
1747 : char *pszHdrFilename;
1748 :
1749 : pszHdrFilename =
1750 34 : CPLStrdup( CPLResetExtension( pszFilename, "hdr" ) );
1751 :
1752 : /* -------------------------------------------------------------------- */
1753 : /* Open the file. */
1754 : /* -------------------------------------------------------------------- */
1755 34 : fp = VSIFOpenL( pszHdrFilename, "wt" );
1756 34 : if( fp == NULL )
1757 : {
1758 : CPLError( CE_Failure, CPLE_OpenFailed,
1759 : "Attempt to create file `%s' failed.\n",
1760 0 : pszHdrFilename );
1761 0 : CPLFree( pszHdrFilename );
1762 0 : return NULL;
1763 : }
1764 :
1765 : /* -------------------------------------------------------------------- */
1766 : /* Decide how many bits the file should have. */
1767 : /* -------------------------------------------------------------------- */
1768 : int nRowBytes;
1769 34 : int nBits = GDALGetDataTypeSize(eType);
1770 :
1771 34 : if( CSLFetchNameValue( papszParmList, "NBITS" ) != NULL )
1772 1 : nBits = atoi(CSLFetchNameValue( papszParmList, "NBITS" ));
1773 :
1774 34 : nRowBytes = (nBits * nXSize + 7) / 8;
1775 :
1776 : /* -------------------------------------------------------------------- */
1777 : /* Check for signed byte. */
1778 : /* -------------------------------------------------------------------- */
1779 34 : const char *pszPixelType = CSLFetchNameValue( papszParmList, "PIXELTYPE" );
1780 34 : if( pszPixelType == NULL )
1781 33 : pszPixelType = "";
1782 :
1783 : /* -------------------------------------------------------------------- */
1784 : /* Write out the raw definition for the dataset as a whole. */
1785 : /* -------------------------------------------------------------------- */
1786 34 : VSIFPrintfL( fp, "BYTEORDER I\n" );
1787 34 : VSIFPrintfL( fp, "LAYOUT BIL\n" );
1788 34 : VSIFPrintfL( fp, "NROWS %d\n", nYSize );
1789 34 : VSIFPrintfL( fp, "NCOLS %d\n", nXSize );
1790 34 : VSIFPrintfL( fp, "NBANDS %d\n", nBands );
1791 34 : VSIFPrintfL( fp, "NBITS %d\n", nBits );
1792 34 : VSIFPrintfL( fp, "BANDROWBYTES %d\n", nRowBytes );
1793 34 : VSIFPrintfL( fp, "TOTALROWBYTES %d\n", nRowBytes * nBands );
1794 :
1795 34 : if( eType == GDT_Float32 )
1796 5 : VSIFPrintfL( fp, "PIXELTYPE FLOAT\n");
1797 36 : else if( eType == GDT_Int16 || eType == GDT_Int32 )
1798 7 : VSIFPrintfL( fp, "PIXELTYPE SIGNEDINT\n");
1799 23 : else if( eType == GDT_Byte && EQUAL(pszPixelType,"SIGNEDBYTE") )
1800 1 : VSIFPrintfL( fp, "PIXELTYPE SIGNEDINT\n");
1801 : else
1802 21 : VSIFPrintfL( fp, "PIXELTYPE UNSIGNEDINT\n");
1803 :
1804 : /* -------------------------------------------------------------------- */
1805 : /* Cleanup */
1806 : /* -------------------------------------------------------------------- */
1807 34 : VSIFCloseL( fp );
1808 :
1809 34 : CPLFree( pszHdrFilename );
1810 :
1811 34 : return (GDALDataset *) GDALOpen( pszFilename, GA_Update );
1812 : }
1813 :
1814 : /************************************************************************/
1815 : /* CreateCopy() */
1816 : /************************************************************************/
1817 :
1818 35 : GDALDataset *EHdrDataset::CreateCopy( const char * pszFilename,
1819 : GDALDataset * poSrcDS,
1820 : int bStrict, char ** papszOptions,
1821 : GDALProgressFunc pfnProgress,
1822 : void * pProgressData )
1823 :
1824 : {
1825 35 : char **papszAdjustedOptions = CSLDuplicate( papszOptions );
1826 : GDALDataset *poOutDS;
1827 :
1828 35 : int nBands = poSrcDS->GetRasterCount();
1829 35 : if (nBands == 0)
1830 : {
1831 : CPLError( CE_Failure, CPLE_NotSupported,
1832 1 : "EHdr driver does not support source dataset with zero band.\n");
1833 1 : return NULL;
1834 : }
1835 :
1836 : /* -------------------------------------------------------------------- */
1837 : /* Ensure we pass on NBITS and PIXELTYPE structure information. */
1838 : /* -------------------------------------------------------------------- */
1839 68 : if( poSrcDS->GetRasterBand(1)->GetMetadataItem( "NBITS",
1840 34 : "IMAGE_STRUCTURE" ) !=NULL
1841 : && CSLFetchNameValue( papszOptions, "NBITS" ) == NULL )
1842 : {
1843 : papszAdjustedOptions =
1844 : CSLSetNameValue( papszAdjustedOptions,
1845 : "NBITS",
1846 0 : poSrcDS->GetRasterBand(1)->GetMetadataItem("NBITS","IMAGE_STRUCTURE") );
1847 : }
1848 :
1849 68 : if( poSrcDS->GetRasterBand(1)->GetMetadataItem( "PIXELTYPE",
1850 34 : "IMAGE_STRUCTURE" ) !=NULL
1851 : && CSLFetchNameValue( papszOptions, "PIXELTYPE" ) == NULL )
1852 : {
1853 : papszAdjustedOptions =
1854 : CSLSetNameValue( papszAdjustedOptions,
1855 : "PIXELTYPE",
1856 1 : poSrcDS->GetRasterBand(1)->GetMetadataItem("PIXELTYPE","IMAGE_STRUCTURE") );
1857 : }
1858 :
1859 : /* -------------------------------------------------------------------- */
1860 : /* Proceed with normal copying using the default createcopy */
1861 : /* operators. */
1862 : /* -------------------------------------------------------------------- */
1863 34 : GDALDriver *poDriver = (GDALDriver *) GDALGetDriverByName( "EHdr" );
1864 :
1865 : poOutDS = poDriver->DefaultCreateCopy( pszFilename, poSrcDS, bStrict,
1866 : papszAdjustedOptions,
1867 34 : pfnProgress, pProgressData );
1868 34 : CSLDestroy( papszAdjustedOptions );
1869 :
1870 34 : return poOutDS;
1871 : }
1872 :
1873 : /************************************************************************/
1874 : /* GetNoDataValue() */
1875 : /************************************************************************/
1876 :
1877 61 : double EHdrRasterBand::GetNoDataValue( int *pbSuccess )
1878 : {
1879 61 : if( pbSuccess )
1880 61 : *pbSuccess = bNoDataSet;
1881 :
1882 61 : if( bNoDataSet )
1883 1 : return dfNoData;
1884 :
1885 60 : return RawRasterBand::GetNoDataValue( pbSuccess );
1886 : }
1887 :
1888 : /************************************************************************/
1889 : /* GetMinimum() */
1890 : /************************************************************************/
1891 :
1892 3 : double EHdrRasterBand::GetMinimum( int *pbSuccess )
1893 : {
1894 3 : if( pbSuccess != NULL )
1895 3 : *pbSuccess = (minmaxmeanstddev & HAS_MIN_FLAG) != 0;
1896 :
1897 3 : if( minmaxmeanstddev & HAS_MIN_FLAG )
1898 2 : return dfMin;
1899 :
1900 1 : return RawRasterBand::GetMinimum( pbSuccess );
1901 : }
1902 :
1903 : /************************************************************************/
1904 : /* GetMaximum() */
1905 : /************************************************************************/
1906 :
1907 2 : double EHdrRasterBand::GetMaximum( int *pbSuccess )
1908 : {
1909 2 : if( pbSuccess != NULL )
1910 2 : *pbSuccess = (minmaxmeanstddev & HAS_MAX_FLAG) != 0;
1911 :
1912 2 : if( minmaxmeanstddev & HAS_MAX_FLAG )
1913 1 : return dfMax;
1914 :
1915 1 : return RawRasterBand::GetMaximum( pbSuccess );
1916 : }
1917 :
1918 : /************************************************************************/
1919 : /* GetStatistics() */
1920 : /************************************************************************/
1921 :
1922 2 : CPLErr EHdrRasterBand::GetStatistics( int bApproxOK, int bForce, double *pdfMin, double *pdfMax, double *pdfMean, double *pdfStdDev )
1923 : {
1924 2 : if( (minmaxmeanstddev & HAS_ALL_FLAGS) == HAS_ALL_FLAGS)
1925 : {
1926 1 : if ( pdfMin ) *pdfMin = dfMin;
1927 1 : if ( pdfMax ) *pdfMax = dfMax;
1928 1 : if ( pdfMean ) *pdfMean = dfMean;
1929 1 : if ( pdfStdDev ) *pdfStdDev = dfStdDev;
1930 1 : return CE_None;
1931 : }
1932 :
1933 : CPLErr eErr = RawRasterBand::GetStatistics( bApproxOK, bForce,
1934 : &dfMin, &dfMax,
1935 1 : &dfMean, &dfStdDev );
1936 :
1937 1 : if( eErr == CE_None )
1938 : {
1939 1 : EHdrDataset* poEDS = (EHdrDataset *) poDS;
1940 :
1941 1 : minmaxmeanstddev = HAS_ALL_FLAGS;
1942 :
1943 1 : if( poEDS->RewriteSTX() != CE_None )
1944 0 : RawRasterBand::SetStatistics( dfMin, dfMax, dfMean, dfStdDev );
1945 :
1946 1 : if( pdfMin )
1947 1 : *pdfMin = dfMin;
1948 1 : if( pdfMax )
1949 1 : *pdfMax = dfMax;
1950 1 : if( pdfMean )
1951 1 : *pdfMean = dfMean;
1952 1 : if( pdfStdDev )
1953 1 : *pdfStdDev = dfStdDev;
1954 : }
1955 :
1956 1 : return eErr;
1957 : }
1958 :
1959 : /************************************************************************/
1960 : /* SetStatistics() */
1961 : /************************************************************************/
1962 :
1963 1 : CPLErr EHdrRasterBand::SetStatistics( double dfMin, double dfMax, double dfMean, double dfStdDev )
1964 : {
1965 : // avoid churn if nothing is changing.
1966 1 : if( dfMin == this->dfMin
1967 : && dfMax == this->dfMax
1968 : && dfMean == this->dfMean
1969 : && dfStdDev == this->dfStdDev )
1970 0 : return CE_None;
1971 :
1972 1 : this->dfMin = dfMin;
1973 1 : this->dfMax = dfMax;
1974 1 : this->dfMean = dfMean;
1975 1 : this->dfStdDev = dfStdDev;
1976 :
1977 : // marks stats valid
1978 1 : minmaxmeanstddev = HAS_ALL_FLAGS;
1979 :
1980 1 : EHdrDataset* poEDS = (EHdrDataset *) poDS;
1981 :
1982 1 : if( poEDS->RewriteSTX() != CE_None )
1983 0 : return RawRasterBand::SetStatistics( dfMin, dfMax, dfMean, dfStdDev );
1984 : else
1985 1 : return CE_None;
1986 : }
1987 :
1988 : /************************************************************************/
1989 : /* SetColorTable() */
1990 : /************************************************************************/
1991 :
1992 6 : CPLErr EHdrRasterBand::SetColorTable( GDALColorTable *poNewCT )
1993 : {
1994 6 : CPLErr err = RawRasterBand::SetColorTable( poNewCT );
1995 6 : if( err != CE_None )
1996 0 : return err;
1997 :
1998 6 : ((EHdrDataset*)poDS)->bCLRDirty = TRUE;
1999 :
2000 6 : return CE_None;
2001 : }
2002 :
2003 : /************************************************************************/
2004 : /* GDALRegister_EHdr() */
2005 : /************************************************************************/
2006 :
2007 610 : void GDALRegister_EHdr()
2008 :
2009 : {
2010 : GDALDriver *poDriver;
2011 :
2012 610 : if( GDALGetDriverByName( "EHdr" ) == NULL )
2013 : {
2014 588 : poDriver = new GDALDriver();
2015 :
2016 588 : poDriver->SetDescription( "EHdr" );
2017 : poDriver->SetMetadataItem( GDAL_DMD_LONGNAME,
2018 588 : "ESRI .hdr Labelled" );
2019 : poDriver->SetMetadataItem( GDAL_DMD_HELPTOPIC,
2020 588 : "frmt_various.html#EHdr" );
2021 : poDriver->SetMetadataItem( GDAL_DMD_CREATIONDATATYPES,
2022 588 : "Byte Int16 UInt16 Int32 UInt32 Float32" );
2023 :
2024 : poDriver->SetMetadataItem( GDAL_DMD_CREATIONOPTIONLIST,
2025 : "<CreationOptionList>"
2026 : " <Option name='NBITS' type='int' description='Special pixel bits (1-7)'/>"
2027 : " <Option name='PIXELTYPE' type='string' description='By setting this to SIGNEDBYTE, a new Byte file can be forced to be written as signed byte'/>"
2028 588 : "</CreationOptionList>" );
2029 :
2030 588 : poDriver->SetMetadataItem( GDAL_DCAP_VIRTUALIO, "YES" );
2031 588 : poDriver->pfnOpen = EHdrDataset::Open;
2032 588 : poDriver->pfnCreate = EHdrDataset::Create;
2033 588 : poDriver->pfnCreateCopy = EHdrDataset::CreateCopy;
2034 :
2035 588 : GetGDALDriverManager()->RegisterDriver( poDriver );
2036 : }
2037 610 : }
|
