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