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