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