Line data Source code
1 : /****************************************************************************
2 : *
3 : * Project: GDAL
4 : * Purpose: Implements the Golden Software Surfer 7 Binary Grid Format.
5 : * Author: Adam Guernsey, adam@ctech.com
6 : * (Based almost entirely on gsbgdataset.cpp by Kevin Locke)
7 : * Create functions added by Russell Jurgensen.
8 : *
9 : ****************************************************************************
10 : * Copyright (c) 2007, Adam Guernsey <adam@ctech.com>
11 : * Copyright (c) 2009-2011, Even Rouault <even dot rouault at spatialys.com>
12 : *
13 : * SPDX-License-Identifier: MIT
14 : ****************************************************************************/
15 :
16 : #include <cassert>
17 : #include <cfloat>
18 : #include <climits>
19 : #include <cmath>
20 : #include <limits>
21 :
22 : #include "gdal_frmts.h"
23 : #include "gdal_pam.h"
24 :
25 : /************************************************************************/
26 : /* ==================================================================== */
27 : /* GS7BGDataset */
28 : /* ==================================================================== */
29 : /************************************************************************/
30 :
31 : class GS7BGRasterBand;
32 :
33 : constexpr double dfDefaultNoDataValue = 1.701410009187828e+38f;
34 :
35 : class GS7BGDataset final : public GDALPamDataset
36 : {
37 : friend class GS7BGRasterBand;
38 :
39 : double dfNoData_Value;
40 : static const size_t nHEADER_SIZE;
41 : size_t nData_Position;
42 :
43 : static CPLErr WriteHeader(VSILFILE *fp, GInt32 nXSize, GInt32 nYSize,
44 : double dfMinX, double dfMaxX, double dfMinY,
45 : double dfMaxY, double dfMinZ, double dfMaxZ);
46 :
47 : VSILFILE *fp;
48 :
49 : public:
50 33 : GS7BGDataset()
51 33 : : /* NOTE: This is not mentioned in the spec, but Surfer 8 uses this
52 : value */
53 : /* 0x7effffee (Little Endian: eeffff7e) */
54 33 : dfNoData_Value(dfDefaultNoDataValue), nData_Position(0), fp(nullptr)
55 : {
56 33 : }
57 :
58 : ~GS7BGDataset();
59 :
60 : static int Identify(GDALOpenInfo *);
61 : static GDALDataset *Open(GDALOpenInfo *);
62 : static GDALDataset *Create(const char *pszFilename, int nXSize, int nYSize,
63 : int nBandsIn, GDALDataType eType,
64 : char **papszParamList);
65 : static GDALDataset *CreateCopy(const char *pszFilename,
66 : GDALDataset *poSrcDS, int bStrict,
67 : char **papszOptions,
68 : GDALProgressFunc pfnProgress,
69 : void *pProgressData);
70 :
71 : CPLErr GetGeoTransform(GDALGeoTransform >) const override;
72 : CPLErr SetGeoTransform(const GDALGeoTransform >) override;
73 : };
74 :
75 : const size_t GS7BGDataset::nHEADER_SIZE = 100;
76 :
77 : constexpr long nHEADER_TAG = 0x42525344;
78 : constexpr long nGRID_TAG = 0x44495247;
79 : constexpr long nDATA_TAG = 0x41544144;
80 : #if 0 /* Unused */
81 : const long nFAULT_TAG = 0x49544c46;
82 : #endif
83 :
84 : /************************************************************************/
85 : /* ==================================================================== */
86 : /* GS7BGRasterBand */
87 : /* ==================================================================== */
88 : /************************************************************************/
89 :
90 : class GS7BGRasterBand final : public GDALPamRasterBand
91 : {
92 : friend class GS7BGDataset;
93 :
94 : double dfMinX;
95 : double dfMaxX;
96 : double dfMinY;
97 : double dfMaxY;
98 : double dfMinZ;
99 : double dfMaxZ;
100 :
101 : double *pafRowMinZ;
102 : double *pafRowMaxZ;
103 : int nMinZRow;
104 : int nMaxZRow;
105 :
106 : CPLErr ScanForMinMaxZ();
107 :
108 : public:
109 : GS7BGRasterBand(GS7BGDataset *, int);
110 : ~GS7BGRasterBand();
111 :
112 : CPLErr IReadBlock(int, int, void *) override;
113 : CPLErr IWriteBlock(int, int, void *) override;
114 : double GetMinimum(int *pbSuccess = nullptr) override;
115 : double GetMaximum(int *pbSuccess = nullptr) override;
116 :
117 : double GetNoDataValue(int *pbSuccess = nullptr) override;
118 : };
119 :
120 : /************************************************************************/
121 : /* GS7BGRasterBand() */
122 : /************************************************************************/
123 :
124 33 : GS7BGRasterBand::GS7BGRasterBand(GS7BGDataset *poDSIn, int nBandIn)
125 : : dfMinX(0.0), dfMaxX(0.0), dfMinY(0.0), dfMaxY(0.0), dfMinZ(0.0),
126 : dfMaxZ(0.0), pafRowMinZ(nullptr), pafRowMaxZ(nullptr), nMinZRow(-1),
127 33 : nMaxZRow(-1)
128 :
129 : {
130 33 : poDS = poDSIn;
131 33 : nBand = nBandIn;
132 :
133 33 : eDataType = GDT_Float64;
134 :
135 33 : nBlockXSize = poDS->GetRasterXSize();
136 33 : nBlockYSize = 1;
137 33 : }
138 :
139 : /************************************************************************/
140 : /* ~GSBGRasterBand() */
141 : /************************************************************************/
142 :
143 66 : GS7BGRasterBand::~GS7BGRasterBand()
144 :
145 : {
146 33 : CPLFree(pafRowMinZ);
147 33 : CPLFree(pafRowMaxZ);
148 66 : }
149 :
150 : /************************************************************************/
151 : /* ScanForMinMaxZ() */
152 : /************************************************************************/
153 :
154 1 : CPLErr GS7BGRasterBand::ScanForMinMaxZ()
155 :
156 : {
157 1 : GS7BGDataset *poGDS = reinterpret_cast<GS7BGDataset *>(poDS);
158 : double *pafRowVals =
159 1 : (double *)VSI_MALLOC2_VERBOSE(nRasterXSize, sizeof(double));
160 :
161 1 : if (pafRowVals == nullptr)
162 : {
163 0 : return CE_Failure;
164 : }
165 :
166 1 : double dfNewMinZ = std::numeric_limits<double>::max();
167 1 : double dfNewMaxZ = std::numeric_limits<double>::lowest();
168 1 : int nNewMinZRow = 0;
169 1 : int nNewMaxZRow = 0;
170 :
171 : /* Since we have to scan, lets calc. statistics too */
172 1 : double dfSum = 0.0;
173 1 : double dfSum2 = 0.0;
174 1 : unsigned long nValuesRead = 0;
175 21 : for (int iRow = 0; iRow < nRasterYSize; iRow++)
176 : {
177 20 : CPLErr eErr = IReadBlock(0, iRow, pafRowVals);
178 20 : if (eErr != CE_None)
179 : {
180 0 : VSIFree(pafRowVals);
181 0 : return CE_Failure;
182 : }
183 :
184 20 : pafRowMinZ[iRow] = std::numeric_limits<float>::max();
185 20 : pafRowMaxZ[iRow] = std::numeric_limits<float>::lowest();
186 420 : for (int iCol = 0; iCol < nRasterXSize; iCol++)
187 : {
188 400 : if (pafRowVals[iCol] == poGDS->dfNoData_Value)
189 400 : continue;
190 :
191 0 : if (pafRowVals[iCol] < pafRowMinZ[iRow])
192 0 : pafRowMinZ[iRow] = pafRowVals[iCol];
193 :
194 0 : if (pafRowVals[iCol] > pafRowMinZ[iRow])
195 0 : pafRowMaxZ[iRow] = pafRowVals[iCol];
196 :
197 0 : dfSum += pafRowVals[iCol];
198 0 : dfSum2 += pafRowVals[iCol] * pafRowVals[iCol];
199 0 : nValuesRead++;
200 : }
201 :
202 20 : if (pafRowMinZ[iRow] < dfNewMinZ)
203 : {
204 1 : dfNewMinZ = pafRowMinZ[iRow];
205 1 : nNewMinZRow = iRow;
206 : }
207 :
208 20 : if (pafRowMaxZ[iRow] > dfNewMaxZ)
209 : {
210 1 : dfNewMaxZ = pafRowMaxZ[iRow];
211 1 : nNewMaxZRow = iRow;
212 : }
213 : }
214 :
215 1 : VSIFree(pafRowVals);
216 :
217 1 : if (nValuesRead == 0)
218 : {
219 1 : dfMinZ = 0.0;
220 1 : dfMaxZ = 0.0;
221 1 : nMinZRow = 0;
222 1 : nMaxZRow = 0;
223 1 : return CE_None;
224 : }
225 :
226 0 : dfMinZ = dfNewMinZ;
227 0 : dfMaxZ = dfNewMaxZ;
228 0 : nMinZRow = nNewMinZRow;
229 0 : nMaxZRow = nNewMaxZRow;
230 :
231 0 : double dfMean = dfSum / nValuesRead;
232 0 : double dfStdDev = sqrt((dfSum2 / nValuesRead) - (dfMean * dfMean));
233 0 : SetStatistics(dfMinZ, dfMaxZ, dfMean, dfStdDev);
234 :
235 0 : return CE_None;
236 : }
237 :
238 : /************************************************************************/
239 : /* IReadBlock() */
240 : /************************************************************************/
241 :
242 140 : CPLErr GS7BGRasterBand::IReadBlock(int nBlockXOff, int nBlockYOff, void *pImage)
243 :
244 : {
245 140 : if (nBlockYOff < 0 || nBlockYOff > nRasterYSize - 1 || nBlockXOff != 0)
246 0 : return CE_Failure;
247 :
248 140 : GS7BGDataset *poGDS = cpl::down_cast<GS7BGDataset *>(poDS);
249 :
250 280 : if (VSIFSeekL(poGDS->fp,
251 140 : (poGDS->nData_Position +
252 140 : sizeof(double) * static_cast<vsi_l_offset>(nRasterXSize) *
253 140 : (nRasterYSize - nBlockYOff - 1)),
254 140 : SEEK_SET) != 0)
255 : {
256 0 : CPLError(CE_Failure, CPLE_FileIO,
257 : "Unable to seek to beginning of grid row.\n");
258 0 : return CE_Failure;
259 : }
260 :
261 140 : if (VSIFReadL(pImage, sizeof(double), nBlockXSize, poGDS->fp) !=
262 140 : static_cast<unsigned>(nBlockXSize))
263 : {
264 0 : CPLError(CE_Failure, CPLE_FileIO,
265 : "Unable to read block from grid file.\n");
266 0 : return CE_Failure;
267 : }
268 :
269 : #ifdef CPL_MSB
270 : double *pfImage = (double *)pImage;
271 : for (int iPixel = 0; iPixel < nBlockXSize; iPixel++)
272 : CPL_LSBPTR64(pfImage + iPixel);
273 : #endif
274 :
275 140 : return CE_None;
276 : }
277 :
278 : /************************************************************************/
279 : /* IWriteBlock() */
280 : /************************************************************************/
281 :
282 20 : CPLErr GS7BGRasterBand::IWriteBlock(int nBlockXOff, int nBlockYOff,
283 : void *pImage)
284 :
285 : {
286 20 : if (eAccess == GA_ReadOnly)
287 : {
288 0 : CPLError(CE_Failure, CPLE_NoWriteAccess,
289 : "Unable to write block, dataset opened read only.\n");
290 0 : return CE_Failure;
291 : }
292 :
293 20 : if (nBlockYOff < 0 || nBlockYOff > nRasterYSize - 1 || nBlockXOff != 0)
294 0 : return CE_Failure;
295 :
296 20 : GS7BGDataset *poGDS = cpl::down_cast<GS7BGDataset *>(poDS);
297 :
298 20 : if (pafRowMinZ == nullptr || pafRowMaxZ == nullptr || nMinZRow < 0 ||
299 19 : nMaxZRow < 0)
300 : {
301 1 : pafRowMinZ =
302 1 : (double *)VSI_MALLOC2_VERBOSE(nRasterYSize, sizeof(double));
303 1 : if (pafRowMinZ == nullptr)
304 : {
305 0 : return CE_Failure;
306 : }
307 :
308 1 : pafRowMaxZ =
309 1 : (double *)VSI_MALLOC2_VERBOSE(nRasterYSize, sizeof(double));
310 1 : if (pafRowMaxZ == nullptr)
311 : {
312 0 : VSIFree(pafRowMinZ);
313 0 : pafRowMinZ = nullptr;
314 0 : return CE_Failure;
315 : }
316 :
317 1 : CPLErr eErr = ScanForMinMaxZ();
318 1 : if (eErr != CE_None)
319 0 : return eErr;
320 : }
321 :
322 40 : if (VSIFSeekL(poGDS->fp,
323 20 : GS7BGDataset::nHEADER_SIZE +
324 20 : sizeof(double) * nRasterXSize *
325 20 : (nRasterYSize - nBlockYOff - 1),
326 20 : SEEK_SET) != 0)
327 : {
328 0 : CPLError(CE_Failure, CPLE_FileIO,
329 : "Unable to seek to beginning of grid row.\n");
330 0 : return CE_Failure;
331 : }
332 :
333 20 : double *pdfImage = (double *)pImage;
334 20 : pafRowMinZ[nBlockYOff] = std::numeric_limits<double>::max();
335 20 : pafRowMaxZ[nBlockYOff] = std::numeric_limits<double>::lowest();
336 420 : for (int iPixel = 0; iPixel < nBlockXSize; iPixel++)
337 : {
338 400 : if (pdfImage[iPixel] != poGDS->dfNoData_Value)
339 : {
340 400 : if (pdfImage[iPixel] < pafRowMinZ[nBlockYOff])
341 78 : pafRowMinZ[nBlockYOff] = pdfImage[iPixel];
342 :
343 400 : if (pdfImage[iPixel] > pafRowMaxZ[nBlockYOff])
344 43 : pafRowMaxZ[nBlockYOff] = pdfImage[iPixel];
345 : }
346 :
347 400 : CPL_LSBPTR64(pdfImage + iPixel);
348 : }
349 :
350 20 : if (VSIFWriteL(pImage, sizeof(double), nBlockXSize, poGDS->fp) !=
351 20 : static_cast<unsigned>(nBlockXSize))
352 : {
353 0 : CPLError(CE_Failure, CPLE_FileIO,
354 : "Unable to write block to grid file.\n");
355 0 : return CE_Failure;
356 : }
357 :
358 : /* Update min/max Z values as appropriate */
359 20 : bool bHeaderNeedsUpdate = false;
360 20 : if (nMinZRow == nBlockYOff && pafRowMinZ[nBlockYOff] > dfMinZ)
361 : {
362 1 : double dfNewMinZ = std::numeric_limits<double>::max();
363 21 : for (int iRow = 0; iRow < nRasterYSize; iRow++)
364 : {
365 20 : if (pafRowMinZ[iRow] < dfNewMinZ)
366 : {
367 1 : dfNewMinZ = pafRowMinZ[iRow];
368 1 : nMinZRow = iRow;
369 : }
370 : }
371 :
372 1 : if (dfNewMinZ != dfMinZ)
373 : {
374 1 : dfMinZ = dfNewMinZ;
375 1 : bHeaderNeedsUpdate = true;
376 : }
377 : }
378 :
379 20 : if (nMaxZRow == nBlockYOff && pafRowMaxZ[nBlockYOff] < dfMaxZ)
380 : {
381 0 : double dfNewMaxZ = std::numeric_limits<double>::lowest();
382 0 : for (int iRow = 0; iRow < nRasterYSize; iRow++)
383 : {
384 0 : if (pafRowMaxZ[iRow] > dfNewMaxZ)
385 : {
386 0 : dfNewMaxZ = pafRowMaxZ[iRow];
387 0 : nMaxZRow = iRow;
388 : }
389 : }
390 :
391 0 : if (dfNewMaxZ != dfMaxZ)
392 : {
393 0 : dfMaxZ = dfNewMaxZ;
394 0 : bHeaderNeedsUpdate = true;
395 : }
396 : }
397 :
398 20 : if (pafRowMinZ[nBlockYOff] < dfMinZ || pafRowMaxZ[nBlockYOff] > dfMaxZ)
399 : {
400 6 : if (pafRowMinZ[nBlockYOff] < dfMinZ)
401 : {
402 3 : dfMinZ = pafRowMinZ[nBlockYOff];
403 3 : nMinZRow = nBlockYOff;
404 : }
405 :
406 6 : if (pafRowMaxZ[nBlockYOff] > dfMaxZ)
407 : {
408 4 : dfMaxZ = pafRowMaxZ[nBlockYOff];
409 4 : nMaxZRow = nBlockYOff;
410 : }
411 :
412 6 : bHeaderNeedsUpdate = true;
413 : }
414 :
415 20 : if (bHeaderNeedsUpdate && dfMaxZ > dfMinZ)
416 : {
417 : CPLErr eErr =
418 6 : poGDS->WriteHeader(poGDS->fp, nRasterXSize, nRasterYSize, dfMinX,
419 : dfMaxX, dfMinY, dfMaxY, dfMinZ, dfMaxZ);
420 6 : return eErr;
421 : }
422 :
423 14 : return CE_None;
424 : }
425 :
426 : /************************************************************************/
427 : /* GetNoDataValue() */
428 : /************************************************************************/
429 :
430 11 : double GS7BGRasterBand::GetNoDataValue(int *pbSuccess)
431 : {
432 11 : GS7BGDataset *poGDS = reinterpret_cast<GS7BGDataset *>(poDS);
433 11 : if (pbSuccess)
434 10 : *pbSuccess = TRUE;
435 :
436 11 : return poGDS->dfNoData_Value;
437 : }
438 :
439 : /************************************************************************/
440 : /* GetMinimum() */
441 : /************************************************************************/
442 :
443 0 : double GS7BGRasterBand::GetMinimum(int *pbSuccess)
444 : {
445 0 : if (pbSuccess)
446 0 : *pbSuccess = TRUE;
447 :
448 0 : return dfMinZ;
449 : }
450 :
451 : /************************************************************************/
452 : /* GetMaximum() */
453 : /************************************************************************/
454 :
455 0 : double GS7BGRasterBand::GetMaximum(int *pbSuccess)
456 : {
457 0 : if (pbSuccess)
458 0 : *pbSuccess = TRUE;
459 :
460 0 : return dfMaxZ;
461 : }
462 :
463 : /************************************************************************/
464 : /* ==================================================================== */
465 : /* GS7BGDataset */
466 : /* ==================================================================== */
467 : /************************************************************************/
468 :
469 66 : GS7BGDataset::~GS7BGDataset()
470 :
471 : {
472 33 : FlushCache(true);
473 33 : if (fp != nullptr)
474 33 : VSIFCloseL(fp);
475 66 : }
476 :
477 : /************************************************************************/
478 : /* Identify() */
479 : /************************************************************************/
480 :
481 58870 : int GS7BGDataset::Identify(GDALOpenInfo *poOpenInfo)
482 :
483 : {
484 : /* Check for signature - for GS7BG the signature is the */
485 : /* nHEADER_TAG with reverse byte order. */
486 58870 : if (poOpenInfo->nHeaderBytes < 4 ||
487 4568 : !STARTS_WITH_CI((const char *)poOpenInfo->pabyHeader, "DSRB"))
488 : {
489 58804 : return FALSE;
490 : }
491 :
492 66 : return TRUE;
493 : }
494 :
495 : /************************************************************************/
496 : /* Open() */
497 : /************************************************************************/
498 :
499 33 : GDALDataset *GS7BGDataset::Open(GDALOpenInfo *poOpenInfo)
500 :
501 : {
502 33 : if (!Identify(poOpenInfo) || poOpenInfo->fpL == nullptr)
503 : {
504 0 : return nullptr;
505 : }
506 :
507 : /* ------------------------------------------------------------------- */
508 : /* Create a corresponding GDALDataset. */
509 : /* ------------------------------------------------------------------- */
510 33 : GS7BGDataset *poDS = new GS7BGDataset();
511 33 : poDS->eAccess = poOpenInfo->eAccess;
512 33 : poDS->fp = poOpenInfo->fpL;
513 33 : poOpenInfo->fpL = nullptr;
514 :
515 : /* ------------------------------------------------------------------- */
516 : /* Read the header. The Header section must be the first section */
517 : /* in the file. */
518 : /* ------------------------------------------------------------------- */
519 33 : if (VSIFSeekL(poDS->fp, 0, SEEK_SET) != 0)
520 : {
521 0 : delete poDS;
522 0 : CPLError(CE_Failure, CPLE_FileIO,
523 : "Unable to seek to start of grid file header.\n");
524 0 : return nullptr;
525 : }
526 :
527 : GInt32 nTag;
528 33 : if (VSIFReadL((void *)&nTag, sizeof(GInt32), 1, poDS->fp) != 1)
529 : {
530 0 : delete poDS;
531 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to read Tag.\n");
532 0 : return nullptr;
533 : }
534 :
535 33 : CPL_LSBPTR32(&nTag);
536 :
537 33 : if (nTag != nHEADER_TAG)
538 : {
539 0 : delete poDS;
540 0 : CPLError(CE_Failure, CPLE_FileIO, "Header tag not found.\n");
541 0 : return nullptr;
542 : }
543 :
544 : GUInt32 nSize;
545 33 : if (VSIFReadL((void *)&nSize, sizeof(GUInt32), 1, poDS->fp) != 1)
546 : {
547 0 : delete poDS;
548 0 : CPLError(CE_Failure, CPLE_FileIO,
549 : "Unable to read file section size.\n");
550 0 : return nullptr;
551 : }
552 :
553 33 : CPL_LSBPTR32(&nSize);
554 :
555 : GInt32 nVersion;
556 33 : if (VSIFReadL((void *)&nVersion, sizeof(GInt32), 1, poDS->fp) != 1)
557 : {
558 0 : delete poDS;
559 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to read file version.\n");
560 0 : return nullptr;
561 : }
562 :
563 33 : CPL_LSBPTR32(&nVersion);
564 :
565 33 : if (nVersion != 1 && nVersion != 2)
566 : {
567 0 : delete poDS;
568 0 : CPLError(CE_Failure, CPLE_FileIO, "Incorrect file version (%d).",
569 : nVersion);
570 0 : return nullptr;
571 : }
572 :
573 : // advance until the grid tag is found
574 66 : while (nTag != nGRID_TAG)
575 : {
576 33 : if (VSIFReadL((void *)&nTag, sizeof(GInt32), 1, poDS->fp) != 1)
577 : {
578 0 : delete poDS;
579 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to read Tag.\n");
580 0 : return nullptr;
581 : }
582 :
583 33 : CPL_LSBPTR32(&nTag);
584 :
585 33 : if (VSIFReadL((void *)&nSize, sizeof(GUInt32), 1, poDS->fp) != 1)
586 : {
587 0 : delete poDS;
588 0 : CPLError(CE_Failure, CPLE_FileIO,
589 : "Unable to read file section size.\n");
590 0 : return nullptr;
591 : }
592 :
593 33 : CPL_LSBPTR32(&nSize);
594 :
595 33 : if (nTag != nGRID_TAG)
596 : {
597 0 : if (VSIFSeekL(poDS->fp, nSize, SEEK_CUR) != 0)
598 : {
599 0 : delete poDS;
600 0 : CPLError(CE_Failure, CPLE_FileIO,
601 : "Unable to seek to end of file section.\n");
602 0 : return nullptr;
603 : }
604 : }
605 : }
606 :
607 : /* --------------------------------------------------------------------*/
608 : /* Read the grid. */
609 : /* --------------------------------------------------------------------*/
610 : /* Parse number of Y axis grid rows */
611 : GInt32 nRows;
612 33 : if (VSIFReadL((void *)&nRows, sizeof(GInt32), 1, poDS->fp) != 1)
613 : {
614 0 : delete poDS;
615 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to read raster Y size.\n");
616 0 : return nullptr;
617 : }
618 33 : CPL_LSBPTR32(&nRows);
619 33 : poDS->nRasterYSize = nRows;
620 :
621 : /* Parse number of X axis grid columns */
622 : GInt32 nCols;
623 33 : if (VSIFReadL((void *)&nCols, sizeof(GInt32), 1, poDS->fp) != 1)
624 : {
625 0 : delete poDS;
626 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to read raster X size.\n");
627 0 : return nullptr;
628 : }
629 33 : CPL_LSBPTR32(&nCols);
630 33 : poDS->nRasterXSize = nCols;
631 :
632 33 : if (!GDALCheckDatasetDimensions(poDS->nRasterXSize, poDS->nRasterYSize))
633 : {
634 0 : delete poDS;
635 0 : return nullptr;
636 : }
637 :
638 : /* --------------------------------------------------------------------*/
639 : /* Create band information objects. */
640 : /* --------------------------------------------------------------------*/
641 33 : GS7BGRasterBand *poBand = new GS7BGRasterBand(poDS, 1);
642 33 : poDS->SetBand(1, poBand);
643 :
644 : // find the min X Value of the grid
645 : double dfTemp;
646 33 : if (VSIFReadL((void *)&dfTemp, sizeof(double), 1, poDS->fp) != 1)
647 : {
648 0 : delete poDS;
649 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to read minimum X value.\n");
650 0 : return nullptr;
651 : }
652 33 : CPL_LSBPTR64(&dfTemp);
653 33 : poBand->dfMinX = dfTemp;
654 :
655 : // find the min Y value of the grid
656 33 : if (VSIFReadL((void *)&dfTemp, sizeof(double), 1, poDS->fp) != 1)
657 : {
658 0 : delete poDS;
659 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to read minimum X value.\n");
660 0 : return nullptr;
661 : }
662 33 : CPL_LSBPTR64(&dfTemp);
663 33 : poBand->dfMinY = dfTemp;
664 :
665 : // find the spacing between adjacent nodes in the X direction
666 : // (between columns)
667 33 : if (VSIFReadL((void *)&dfTemp, sizeof(double), 1, poDS->fp) != 1)
668 : {
669 0 : delete poDS;
670 0 : CPLError(CE_Failure, CPLE_FileIO,
671 : "Unable to read spacing in X value.\n");
672 0 : return nullptr;
673 : }
674 33 : CPL_LSBPTR64(&dfTemp);
675 33 : poBand->dfMaxX = poBand->dfMinX + (dfTemp * (nCols - 1));
676 :
677 : // find the spacing between adjacent nodes in the Y direction
678 : // (between rows)
679 33 : if (VSIFReadL((void *)&dfTemp, sizeof(double), 1, poDS->fp) != 1)
680 : {
681 0 : delete poDS;
682 0 : CPLError(CE_Failure, CPLE_FileIO,
683 : "Unable to read spacing in Y value.\n");
684 0 : return nullptr;
685 : }
686 33 : CPL_LSBPTR64(&dfTemp);
687 33 : poBand->dfMaxY = poBand->dfMinY + (dfTemp * (nRows - 1));
688 :
689 : // set the z min
690 33 : if (VSIFReadL((void *)&dfTemp, sizeof(double), 1, poDS->fp) != 1)
691 : {
692 0 : delete poDS;
693 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to read Z min value.\n");
694 0 : return nullptr;
695 : }
696 33 : CPL_LSBPTR64(&dfTemp);
697 33 : poBand->dfMinZ = dfTemp;
698 :
699 : // set the z max
700 33 : if (VSIFReadL((void *)&dfTemp, sizeof(double), 1, poDS->fp) != 1)
701 : {
702 0 : delete poDS;
703 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to read Z max value.\n");
704 0 : return nullptr;
705 : }
706 33 : CPL_LSBPTR64(&dfTemp);
707 33 : poBand->dfMaxZ = dfTemp;
708 :
709 : // read and ignore the rotation value
710 : //(This is not used in the current version).
711 33 : if (VSIFReadL((void *)&dfTemp, sizeof(double), 1, poDS->fp) != 1)
712 : {
713 0 : delete poDS;
714 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to read rotation value.\n");
715 0 : return nullptr;
716 : }
717 :
718 : // read and set the cell blank value
719 33 : if (VSIFReadL((void *)&dfTemp, sizeof(double), 1, poDS->fp) != 1)
720 : {
721 0 : delete poDS;
722 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to Blank value.\n");
723 0 : return nullptr;
724 : }
725 33 : CPL_LSBPTR64(&dfTemp);
726 33 : poDS->dfNoData_Value = dfTemp;
727 :
728 : /* --------------------------------------------------------------------*/
729 : /* Set the current offset of the grid data. */
730 : /* --------------------------------------------------------------------*/
731 33 : if (VSIFReadL((void *)&nTag, sizeof(GInt32), 1, poDS->fp) != 1)
732 : {
733 0 : delete poDS;
734 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to read Tag.\n");
735 0 : return nullptr;
736 : }
737 :
738 33 : CPL_LSBPTR32(&nTag);
739 33 : if (nTag != nDATA_TAG)
740 : {
741 0 : delete poDS;
742 0 : CPLError(CE_Failure, CPLE_FileIO, "Data tag not found.\n");
743 0 : return nullptr;
744 : }
745 :
746 33 : if (VSIFReadL((void *)&nSize, sizeof(GInt32), 1, poDS->fp) != 1)
747 : {
748 0 : delete poDS;
749 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to data section size.\n");
750 0 : return nullptr;
751 : }
752 :
753 33 : poDS->nData_Position = (size_t)VSIFTellL(poDS->fp);
754 :
755 : /* --------------------------------------------------------------------*/
756 : /* Initialize any PAM information. */
757 : /* --------------------------------------------------------------------*/
758 33 : poDS->SetDescription(poOpenInfo->pszFilename);
759 33 : poDS->TryLoadXML();
760 :
761 : /* -------------------------------------------------------------------- */
762 : /* Check for external overviews. */
763 : /* -------------------------------------------------------------------- */
764 66 : poDS->oOvManager.Initialize(poDS, poOpenInfo->pszFilename,
765 33 : poOpenInfo->GetSiblingFiles());
766 :
767 33 : return poDS;
768 : }
769 :
770 : /************************************************************************/
771 : /* GetGeoTransform() */
772 : /************************************************************************/
773 :
774 23 : CPLErr GS7BGDataset::GetGeoTransform(GDALGeoTransform >) const
775 : {
776 : const GS7BGRasterBand *poGRB =
777 23 : cpl::down_cast<const GS7BGRasterBand *>(GetRasterBand(1));
778 :
779 23 : if (poGRB == nullptr)
780 : {
781 0 : gt = GDALGeoTransform();
782 0 : return CE_Failure;
783 : }
784 :
785 : /* check if we have a PAM GeoTransform stored */
786 23 : CPLPushErrorHandler(CPLQuietErrorHandler);
787 23 : CPLErr eErr = GDALPamDataset::GetGeoTransform(gt);
788 23 : CPLPopErrorHandler();
789 :
790 23 : if (eErr == CE_None)
791 0 : return CE_None;
792 :
793 23 : if (nRasterXSize == 1 || nRasterYSize == 1)
794 0 : return CE_Failure;
795 :
796 : /* calculate pixel size first */
797 23 : gt[1] = (poGRB->dfMaxX - poGRB->dfMinX) / (nRasterXSize - 1);
798 23 : gt[5] = (poGRB->dfMinY - poGRB->dfMaxY) / (nRasterYSize - 1);
799 :
800 : /* then calculate image origin */
801 23 : gt[0] = poGRB->dfMinX - gt[1] / 2;
802 23 : gt[3] = poGRB->dfMaxY - gt[5] / 2;
803 :
804 : /* tilt/rotation does not supported by the GS grids */
805 23 : gt[4] = 0.0;
806 23 : gt[2] = 0.0;
807 :
808 23 : return CE_None;
809 : }
810 :
811 : /************************************************************************/
812 : /* SetGeoTransform() */
813 : /************************************************************************/
814 :
815 6 : CPLErr GS7BGDataset::SetGeoTransform(const GDALGeoTransform >)
816 : {
817 6 : if (eAccess == GA_ReadOnly)
818 : {
819 0 : CPLError(CE_Failure, CPLE_NoWriteAccess,
820 : "Unable to set GeoTransform, dataset opened read only.\n");
821 0 : return CE_Failure;
822 : }
823 :
824 : GS7BGRasterBand *poGRB =
825 6 : cpl::down_cast<GS7BGRasterBand *>(GetRasterBand(1));
826 :
827 : /* non-zero transform 2 or 4 or negative 1 or 5 not supported natively */
828 : /*if( gt[2] != 0.0 || gt[4] != 0.0
829 : || gt[1] < 0.0 || gt[5] < 0.0 )
830 : eErr = GDALPamDataset::SetGeoTransform( gt );
831 :
832 : if( eErr != CE_None )
833 : return eErr;*/
834 :
835 6 : double dfMinX = gt[0] + gt[1] / 2;
836 6 : double dfMaxX = gt[1] * (nRasterXSize - 0.5) + gt[0];
837 6 : double dfMinY = gt[5] * (nRasterYSize - 0.5) + gt[3];
838 6 : double dfMaxY = gt[3] + gt[5] / 2;
839 :
840 : CPLErr eErr =
841 6 : WriteHeader(fp, poGRB->nRasterXSize, poGRB->nRasterYSize, dfMinX,
842 : dfMaxX, dfMinY, dfMaxY, poGRB->dfMinZ, poGRB->dfMaxZ);
843 :
844 6 : if (eErr == CE_None)
845 : {
846 6 : poGRB->dfMinX = dfMinX;
847 6 : poGRB->dfMaxX = dfMaxX;
848 6 : poGRB->dfMinY = dfMinY;
849 6 : poGRB->dfMaxY = dfMaxY;
850 : }
851 :
852 6 : return eErr;
853 : }
854 :
855 : /************************************************************************/
856 : /* WriteHeader() */
857 : /************************************************************************/
858 :
859 54 : CPLErr GS7BGDataset::WriteHeader(VSILFILE *fp, GInt32 nXSize, GInt32 nYSize,
860 : double dfMinX, double dfMaxX, double dfMinY,
861 : double dfMaxY, double dfMinZ, double dfMaxZ)
862 :
863 : {
864 54 : if (VSIFSeekL(fp, 0, SEEK_SET) != 0)
865 : {
866 0 : CPLError(CE_Failure, CPLE_FileIO,
867 : "Unable to seek to start of grid file.\n");
868 0 : return CE_Failure;
869 : }
870 :
871 54 : GInt32 nTemp = CPL_LSBWORD32(nHEADER_TAG);
872 54 : if (VSIFWriteL((void *)&nTemp, sizeof(GInt32), 1, fp) != 1)
873 : {
874 1 : CPLError(CE_Failure, CPLE_FileIO,
875 : "Unable to write header tag to grid file.\n");
876 1 : return CE_Failure;
877 : }
878 :
879 53 : nTemp = CPL_LSBWORD32(sizeof(GInt32)); // Size of version section.
880 53 : if (VSIFWriteL((void *)&nTemp, sizeof(GInt32), 1, fp) != 1)
881 : {
882 0 : CPLError(CE_Failure, CPLE_FileIO,
883 : "Unable to write size to grid file.\n");
884 0 : return CE_Failure;
885 : }
886 :
887 53 : nTemp = CPL_LSBWORD32(1); // Version
888 53 : if (VSIFWriteL((void *)&nTemp, sizeof(GInt32), 1, fp) != 1)
889 : {
890 0 : CPLError(CE_Failure, CPLE_FileIO,
891 : "Unable to write size to grid file.\n");
892 0 : return CE_Failure;
893 : }
894 :
895 53 : nTemp = CPL_LSBWORD32(nGRID_TAG); // Mark start of grid
896 53 : if (VSIFWriteL((void *)&nTemp, sizeof(GInt32), 1, fp) != 1)
897 : {
898 0 : CPLError(CE_Failure, CPLE_FileIO,
899 : "Unable to write size to grid file.\n");
900 0 : return CE_Failure;
901 : }
902 :
903 53 : nTemp = CPL_LSBWORD32(72); // Grid info size (the remainder of the header)
904 53 : if (VSIFWriteL((void *)&nTemp, sizeof(GInt32), 1, fp) != 1)
905 : {
906 0 : CPLError(CE_Failure, CPLE_FileIO,
907 : "Unable to write size to grid file.\n");
908 0 : return CE_Failure;
909 : }
910 :
911 53 : nTemp = CPL_LSBWORD32(nYSize);
912 53 : if (VSIFWriteL((void *)&nTemp, sizeof(GInt32), 1, fp) != 1)
913 : {
914 0 : CPLError(CE_Failure, CPLE_FileIO,
915 : "Unable to write Y size to grid file.\n");
916 0 : return CE_Failure;
917 : }
918 :
919 53 : nTemp = CPL_LSBWORD32(nXSize);
920 53 : if (VSIFWriteL((void *)&nTemp, sizeof(GInt32), 1, fp) != 1)
921 : {
922 0 : CPLError(CE_Failure, CPLE_FileIO,
923 : "Unable to write X size to grid file.\n");
924 0 : return CE_Failure;
925 : }
926 :
927 53 : double dfTemp = dfMinX;
928 53 : CPL_LSBPTR64(&dfTemp);
929 53 : if (VSIFWriteL((void *)&dfTemp, sizeof(double), 1, fp) != 1)
930 : {
931 0 : CPLError(CE_Failure, CPLE_FileIO,
932 : "Unable to write minimum X value to grid file.\n");
933 0 : return CE_Failure;
934 : }
935 :
936 53 : dfTemp = dfMinY;
937 53 : CPL_LSBPTR64(&dfTemp);
938 53 : if (VSIFWriteL((void *)&dfTemp, sizeof(double), 1, fp) != 1)
939 : {
940 0 : CPLError(CE_Failure, CPLE_FileIO,
941 : "Unable to write minimum Y value to grid file.\n");
942 0 : return CE_Failure;
943 : }
944 :
945 : // Write node spacing in x direction
946 53 : dfTemp = (dfMaxX - dfMinX) / (nXSize - 1);
947 53 : CPL_LSBPTR64(&dfTemp);
948 53 : if (VSIFWriteL((void *)&dfTemp, sizeof(double), 1, fp) != 1)
949 : {
950 0 : CPLError(CE_Failure, CPLE_FileIO,
951 : "Unable to write spacing in X value.\n");
952 0 : return CE_Failure;
953 : }
954 :
955 : // Write node spacing in y direction
956 53 : dfTemp = (dfMaxY - dfMinY) / (nYSize - 1);
957 53 : CPL_LSBPTR64(&dfTemp);
958 53 : if (VSIFWriteL((void *)&dfTemp, sizeof(double), 1, fp) != 1)
959 : {
960 0 : CPLError(CE_Failure, CPLE_FileIO,
961 : "Unable to write spacing in Y value.\n");
962 0 : return CE_Failure;
963 : }
964 :
965 53 : dfTemp = dfMinZ;
966 53 : CPL_LSBPTR64(&dfTemp);
967 53 : if (VSIFWriteL((void *)&dfTemp, sizeof(double), 1, fp) != 1)
968 : {
969 0 : CPLError(CE_Failure, CPLE_FileIO,
970 : "Unable to write minimum Z value to grid file.\n");
971 0 : return CE_Failure;
972 : }
973 :
974 53 : dfTemp = dfMaxZ;
975 53 : CPL_LSBPTR64(&dfTemp);
976 53 : if (VSIFWriteL((void *)&dfTemp, sizeof(double), 1, fp) != 1)
977 : {
978 0 : CPLError(CE_Failure, CPLE_FileIO,
979 : "Unable to write maximum Z value to grid file.\n");
980 0 : return CE_Failure;
981 : }
982 :
983 53 : dfTemp = 0; // Rotation value is zero
984 53 : CPL_LSBPTR64(&dfTemp);
985 53 : if (VSIFWriteL((void *)&dfTemp, sizeof(double), 1, fp) != 1)
986 : {
987 0 : CPLError(CE_Failure, CPLE_FileIO,
988 : "Unable to write rotation value to grid file.\n");
989 0 : return CE_Failure;
990 : }
991 :
992 53 : dfTemp = dfDefaultNoDataValue;
993 53 : CPL_LSBPTR64(&dfTemp);
994 53 : if (VSIFWriteL((void *)&dfTemp, sizeof(double), 1, fp) != 1)
995 : {
996 1 : CPLError(CE_Failure, CPLE_FileIO,
997 : "Unable to write cell blank value to grid file.\n");
998 1 : return CE_Failure;
999 : }
1000 :
1001 : // Only supports 1 band so go ahead and write band info here
1002 52 : nTemp = CPL_LSBWORD32(nDATA_TAG); // Mark start of data
1003 52 : if (VSIFWriteL((void *)&nTemp, sizeof(GInt32), 1, fp) != 1)
1004 : {
1005 0 : CPLError(CE_Failure, CPLE_FileIO, "Unable to data tag to grid file.\n");
1006 0 : return CE_Failure;
1007 : }
1008 :
1009 52 : int nSize = nXSize * nYSize * (int)sizeof(double);
1010 52 : nTemp = CPL_LSBWORD32(nSize); // Mark size of data
1011 52 : if (VSIFWriteL((void *)&nTemp, sizeof(GInt32), 1, fp) != 1)
1012 : {
1013 0 : CPLError(CE_Failure, CPLE_FileIO,
1014 : "Unable to write data size to grid file.\n");
1015 0 : return CE_Failure;
1016 : }
1017 :
1018 52 : return CE_None;
1019 : }
1020 :
1021 : /************************************************************************/
1022 : /* GS7BGCreateCheckDims() */
1023 : /************************************************************************/
1024 :
1025 64 : static bool GS7BGCreateCheckDims(int nXSize, int nYSize)
1026 : {
1027 64 : if (nXSize <= 1 || nYSize <= 1)
1028 : {
1029 3 : CPLError(CE_Failure, CPLE_IllegalArg,
1030 : "Unable to create grid, both X and Y size must be "
1031 : "larger or equal to 2.");
1032 3 : return false;
1033 : }
1034 61 : if (nXSize > INT_MAX / nYSize / static_cast<int>(sizeof(double)))
1035 : {
1036 2 : CPLError(CE_Failure, CPLE_IllegalArg,
1037 : "Unable to create grid, too large X and Y size.");
1038 2 : return false;
1039 : }
1040 59 : return true;
1041 : }
1042 :
1043 : /************************************************************************/
1044 : /* Create() */
1045 : /************************************************************************/
1046 :
1047 37 : GDALDataset *GS7BGDataset::Create(const char *pszFilename, int nXSize,
1048 : int nYSize, int nBandsIn, GDALDataType eType,
1049 : char ** /* papszParamList*/)
1050 :
1051 : {
1052 37 : if (!GS7BGCreateCheckDims(nXSize, nYSize))
1053 : {
1054 4 : return nullptr;
1055 : }
1056 :
1057 33 : if (eType != GDT_Byte && eType != GDT_Float32 && eType != GDT_UInt16 &&
1058 21 : eType != GDT_Int16 && eType != GDT_Float64)
1059 : {
1060 18 : CPLError(
1061 : CE_Failure, CPLE_AppDefined,
1062 : "GS7BG Grid only supports Byte, Int16, "
1063 : "Uint16, Float32, and Float64 datatypes. Unable to create with "
1064 : "type %s.\n",
1065 : GDALGetDataTypeName(eType));
1066 :
1067 18 : return nullptr;
1068 : }
1069 :
1070 15 : if (nBandsIn > 1)
1071 : {
1072 8 : CPLError(CE_Failure, CPLE_NotSupported,
1073 : "Unable to create copy, "
1074 : "format only supports one raster band.\n");
1075 8 : return nullptr;
1076 : }
1077 :
1078 7 : VSILFILE *fp = VSIFOpenL(pszFilename, "w+b");
1079 :
1080 7 : if (fp == nullptr)
1081 : {
1082 0 : CPLError(CE_Failure, CPLE_OpenFailed,
1083 : "Attempt to create file '%s' failed.\n", pszFilename);
1084 0 : return nullptr;
1085 : }
1086 :
1087 : CPLErr eErr =
1088 7 : WriteHeader(fp, nXSize, nYSize, 0.0, nXSize, 0.0, nYSize, 0.0, 0.0);
1089 7 : if (eErr != CE_None)
1090 : {
1091 0 : VSIFCloseL(fp);
1092 0 : return nullptr;
1093 : }
1094 :
1095 7 : double dfVal = dfDefaultNoDataValue;
1096 7 : CPL_LSBPTR64(&dfVal);
1097 627 : for (int iRow = 0; iRow < nYSize; iRow++)
1098 : {
1099 61020 : for (int iCol = 0; iCol < nXSize; iCol++)
1100 : {
1101 60400 : if (VSIFWriteL((void *)&dfVal, sizeof(double), 1, fp) != 1)
1102 : {
1103 0 : VSIFCloseL(fp);
1104 0 : CPLError(CE_Failure, CPLE_FileIO,
1105 : "Unable to write grid cell. Disk full?\n");
1106 0 : return nullptr;
1107 : }
1108 : }
1109 : }
1110 :
1111 7 : VSIFCloseL(fp);
1112 :
1113 7 : return (GDALDataset *)GDALOpen(pszFilename, GA_Update);
1114 : }
1115 :
1116 : /************************************************************************/
1117 : /* CreateCopy() */
1118 : /************************************************************************/
1119 :
1120 32 : GDALDataset *GS7BGDataset::CreateCopy(const char *pszFilename,
1121 : GDALDataset *poSrcDS, int bStrict,
1122 : char ** /*papszOptions*/,
1123 : GDALProgressFunc pfnProgress,
1124 : void *pProgressData)
1125 : {
1126 32 : if (pfnProgress == nullptr)
1127 0 : pfnProgress = GDALDummyProgress;
1128 :
1129 32 : int nBands = poSrcDS->GetRasterCount();
1130 32 : if (nBands == 0)
1131 : {
1132 1 : CPLError(CE_Failure, CPLE_NotSupported,
1133 : "Driver does not support source dataset with zero band.\n");
1134 1 : return nullptr;
1135 : }
1136 31 : else if (nBands > 1)
1137 : {
1138 4 : if (bStrict)
1139 : {
1140 4 : CPLError(CE_Failure, CPLE_NotSupported,
1141 : "Unable to create copy, "
1142 : "format only supports one raster band.\n");
1143 4 : return nullptr;
1144 : }
1145 : else
1146 0 : CPLError(CE_Warning, CPLE_NotSupported,
1147 : "Format only supports one "
1148 : "raster band, first band will be copied.\n");
1149 : }
1150 :
1151 27 : const int nXSize = poSrcDS->GetRasterXSize();
1152 27 : const int nYSize = poSrcDS->GetRasterXSize();
1153 27 : if (!GS7BGCreateCheckDims(nXSize, nYSize))
1154 : {
1155 1 : return nullptr;
1156 : }
1157 :
1158 26 : GDALRasterBand *poSrcBand = poSrcDS->GetRasterBand(1);
1159 :
1160 26 : if (!pfnProgress(0.0, nullptr, pProgressData))
1161 : {
1162 0 : CPLError(CE_Failure, CPLE_UserInterrupt, "User terminated\n");
1163 0 : return nullptr;
1164 : }
1165 :
1166 26 : VSILFILE *fp = VSIFOpenL(pszFilename, "w+b");
1167 :
1168 26 : if (fp == nullptr)
1169 : {
1170 3 : CPLError(CE_Failure, CPLE_OpenFailed,
1171 : "Attempt to create file '%s' failed.\n", pszFilename);
1172 3 : return nullptr;
1173 : }
1174 :
1175 23 : GDALGeoTransform gt;
1176 23 : poSrcDS->GetGeoTransform(gt);
1177 :
1178 23 : double dfMinX = gt[0] + gt[1] / 2;
1179 23 : double dfMaxX = gt[1] * (nXSize - 0.5) + gt[0];
1180 23 : double dfMinY = gt[5] * (nYSize - 0.5) + gt[3];
1181 23 : double dfMaxY = gt[3] + gt[5] / 2;
1182 23 : CPLErr eErr = WriteHeader(fp, nXSize, nYSize, dfMinX, dfMaxX, dfMinY,
1183 : dfMaxY, 0.0, 0.0);
1184 :
1185 23 : if (eErr != CE_None)
1186 : {
1187 2 : VSIFCloseL(fp);
1188 2 : return nullptr;
1189 : }
1190 :
1191 : /* -------------------------------------------------------------------- */
1192 : /* Copy band data. */
1193 : /* -------------------------------------------------------------------- */
1194 21 : double *pfData = (double *)VSI_MALLOC2_VERBOSE(nXSize, sizeof(double));
1195 21 : if (pfData == nullptr)
1196 : {
1197 0 : VSIFCloseL(fp);
1198 0 : return nullptr;
1199 : }
1200 :
1201 : int bSrcHasNDValue;
1202 21 : double dfSrcNoDataValue = poSrcBand->GetNoDataValue(&bSrcHasNDValue);
1203 21 : double dfMinZ = std::numeric_limits<double>::max();
1204 21 : double dfMaxZ = std::numeric_limits<double>::lowest();
1205 187 : for (GInt32 iRow = nYSize - 1; iRow >= 0; iRow--)
1206 : {
1207 175 : eErr = poSrcBand->RasterIO(GF_Read, 0, iRow, nXSize, 1, pfData, nXSize,
1208 : 1, GDT_Float64, 0, 0, nullptr);
1209 :
1210 175 : if (eErr != CE_None)
1211 : {
1212 1 : VSIFCloseL(fp);
1213 1 : VSIFree(pfData);
1214 1 : return nullptr;
1215 : }
1216 :
1217 2114 : for (int iCol = 0; iCol < nXSize; iCol++)
1218 : {
1219 1940 : if (bSrcHasNDValue && pfData[iCol] == dfSrcNoDataValue)
1220 : {
1221 0 : pfData[iCol] = dfDefaultNoDataValue;
1222 : }
1223 : else
1224 : {
1225 1940 : if (pfData[iCol] > dfMaxZ)
1226 22 : dfMaxZ = pfData[iCol];
1227 :
1228 1940 : if (pfData[iCol] < dfMinZ)
1229 27 : dfMinZ = pfData[iCol];
1230 : }
1231 :
1232 1940 : CPL_LSBPTR64(pfData + iCol);
1233 : }
1234 :
1235 174 : if (VSIFWriteL((void *)pfData, sizeof(double), nXSize, fp) !=
1236 174 : static_cast<unsigned>(nXSize))
1237 : {
1238 8 : VSIFCloseL(fp);
1239 8 : VSIFree(pfData);
1240 8 : CPLError(CE_Failure, CPLE_FileIO,
1241 : "Unable to write grid row. Disk full?\n");
1242 8 : return nullptr;
1243 : }
1244 :
1245 166 : if (!pfnProgress(static_cast<double>(nYSize - iRow) / nYSize, nullptr,
1246 : pProgressData))
1247 : {
1248 0 : VSIFCloseL(fp);
1249 0 : VSIFree(pfData);
1250 0 : CPLError(CE_Failure, CPLE_UserInterrupt, "User terminated");
1251 0 : return nullptr;
1252 : }
1253 : }
1254 :
1255 12 : VSIFree(pfData);
1256 :
1257 : /* write out the min and max values */
1258 12 : eErr = WriteHeader(fp, nXSize, nYSize, dfMinX, dfMaxX, dfMinY, dfMaxY,
1259 : dfMinZ, dfMaxZ);
1260 :
1261 12 : if (eErr != CE_None)
1262 : {
1263 0 : VSIFCloseL(fp);
1264 0 : return nullptr;
1265 : }
1266 :
1267 12 : VSIFCloseL(fp);
1268 :
1269 12 : GDALPamDataset *poDS = (GDALPamDataset *)GDALOpen(pszFilename, GA_Update);
1270 12 : if (poDS)
1271 : {
1272 12 : poDS->CloneInfo(poSrcDS, GCIF_PAM_DEFAULT);
1273 : }
1274 :
1275 12 : return poDS;
1276 : }
1277 :
1278 : /************************************************************************/
1279 : /* GDALRegister_GS7BG() */
1280 : /************************************************************************/
1281 1935 : void GDALRegister_GS7BG()
1282 :
1283 : {
1284 1935 : if (GDALGetDriverByName("GS7BG") != nullptr)
1285 282 : return;
1286 :
1287 1653 : GDALDriver *poDriver = new GDALDriver();
1288 :
1289 1653 : poDriver->SetDescription("GS7BG");
1290 1653 : poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
1291 1653 : poDriver->SetMetadataItem(GDAL_DMD_LONGNAME,
1292 1653 : "Golden Software 7 Binary Grid (.grd)");
1293 1653 : poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/gs7bg.html");
1294 1653 : poDriver->SetMetadataItem(GDAL_DMD_EXTENSION, "grd");
1295 1653 : poDriver->SetMetadataItem(GDAL_DMD_CREATIONDATATYPES,
1296 1653 : "Byte Int16 UInt16 Float32 Float64");
1297 1653 : poDriver->SetMetadataItem(GDAL_DCAP_VIRTUALIO, "YES");
1298 :
1299 1653 : poDriver->pfnIdentify = GS7BGDataset::Identify;
1300 1653 : poDriver->pfnOpen = GS7BGDataset::Open;
1301 1653 : poDriver->pfnCreate = GS7BGDataset::Create;
1302 1653 : poDriver->pfnCreateCopy = GS7BGDataset::CreateCopy;
1303 :
1304 1653 : GetGDALDriverManager()->RegisterDriver(poDriver);
1305 : }
|