Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: NDF Driver
4 : * Purpose: Implementation of NLAPS Data Format read support.
5 : * Author: Frank Warmerdam, warmerdam@pobox.com
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2005, Frank Warmerdam
9 : * Copyright (c) 2008-2011, Even Rouault <even dot rouault at spatialys.com>
10 : *
11 : * SPDX-License-Identifier: MIT
12 : ****************************************************************************/
13 :
14 : #include "cpl_string.h"
15 : #include "gdal_frmts.h"
16 : #include "gdal_priv.h"
17 : #include "ogr_spatialref.h"
18 : #include "rawdataset.h"
19 :
20 : /************************************************************************/
21 : /* ==================================================================== */
22 : /* NDFDataset */
23 : /* ==================================================================== */
24 : /************************************************************************/
25 :
26 : class NDFDataset final : public RawDataset
27 : {
28 : GDALGeoTransform m_gt{};
29 :
30 : OGRSpatialReference m_oSRS{};
31 : char **papszExtraFiles;
32 :
33 : char **papszHeader;
34 : const char *Get(const char *pszKey, const char *pszDefault);
35 :
36 : CPL_DISALLOW_COPY_ASSIGN(NDFDataset)
37 :
38 : CPLErr Close() override;
39 :
40 : public:
41 : NDFDataset();
42 : ~NDFDataset() override;
43 :
44 : CPLErr GetGeoTransform(GDALGeoTransform >) const override;
45 :
46 1 : const OGRSpatialReference *GetSpatialRef() const override
47 : {
48 1 : return m_oSRS.IsEmpty() ? nullptr : &m_oSRS;
49 : }
50 :
51 : char **GetFileList(void) override;
52 :
53 : static GDALDataset *Open(GDALOpenInfo *);
54 : static int Identify(GDALOpenInfo *);
55 : };
56 :
57 : /************************************************************************/
58 : /* NDFDataset() */
59 : /************************************************************************/
60 :
61 2 : NDFDataset::NDFDataset() : papszExtraFiles(nullptr), papszHeader(nullptr)
62 : {
63 2 : m_oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
64 2 : }
65 :
66 : /************************************************************************/
67 : /* ~NDFDataset() */
68 : /************************************************************************/
69 :
70 4 : NDFDataset::~NDFDataset()
71 :
72 : {
73 2 : NDFDataset::Close();
74 4 : }
75 :
76 : /************************************************************************/
77 : /* Close() */
78 : /************************************************************************/
79 :
80 4 : CPLErr NDFDataset::Close()
81 : {
82 4 : CPLErr eErr = CE_None;
83 4 : if (nOpenFlags != OPEN_FLAGS_CLOSED)
84 : {
85 2 : if (NDFDataset::FlushCache(true) != CE_None)
86 0 : eErr = CE_Failure;
87 :
88 2 : CSLDestroy(papszHeader);
89 2 : CSLDestroy(papszExtraFiles);
90 :
91 2 : if (GDALPamDataset::Close() != CE_None)
92 0 : eErr = CE_Failure;
93 : }
94 4 : return eErr;
95 : }
96 :
97 : /************************************************************************/
98 : /* GetGeoTransform() */
99 : /************************************************************************/
100 :
101 1 : CPLErr NDFDataset::GetGeoTransform(GDALGeoTransform >) const
102 :
103 : {
104 1 : gt = m_gt;
105 1 : return CE_None;
106 : }
107 :
108 : /************************************************************************/
109 : /* Get() */
110 : /* */
111 : /* Fetch a value from the header by keyword. */
112 : /************************************************************************/
113 :
114 26 : const char *NDFDataset::Get(const char *pszKey, const char *pszDefault)
115 :
116 : {
117 26 : const char *pszResult = CSLFetchNameValue(papszHeader, pszKey);
118 :
119 26 : if (pszResult == nullptr)
120 0 : return pszDefault;
121 :
122 26 : return pszResult;
123 : }
124 :
125 : /************************************************************************/
126 : /* GetFileList() */
127 : /************************************************************************/
128 :
129 1 : char **NDFDataset::GetFileList()
130 :
131 : {
132 : // Main data file, etc.
133 1 : char **papszFileList = GDALPamDataset::GetFileList();
134 :
135 : // Header file.
136 1 : papszFileList = CSLInsertStrings(papszFileList, -1, papszExtraFiles);
137 :
138 1 : return papszFileList;
139 : }
140 :
141 : /************************************************************************/
142 : /* Identify() */
143 : /************************************************************************/
144 :
145 59137 : int NDFDataset::Identify(GDALOpenInfo *poOpenInfo)
146 :
147 : {
148 : /* -------------------------------------------------------------------- */
149 : /* The user must select the header file (i.e. .H1). */
150 : /* -------------------------------------------------------------------- */
151 62911 : return poOpenInfo->nHeaderBytes >= 50 &&
152 3774 : (STARTS_WITH_CI(
153 : reinterpret_cast<const char *>(poOpenInfo->pabyHeader),
154 3770 : "NDF_REVISION=2") ||
155 3770 : STARTS_WITH_CI(
156 : reinterpret_cast<const char *>(poOpenInfo->pabyHeader),
157 59137 : "NDF_REVISION=0"));
158 : }
159 :
160 : /************************************************************************/
161 : /* Open() */
162 : /************************************************************************/
163 :
164 2 : GDALDataset *NDFDataset::Open(GDALOpenInfo *poOpenInfo)
165 :
166 : {
167 : /* -------------------------------------------------------------------- */
168 : /* The user must select the header file (i.e. .H1). */
169 : /* -------------------------------------------------------------------- */
170 2 : if (!Identify(poOpenInfo) || poOpenInfo->fpL == nullptr)
171 0 : return nullptr;
172 :
173 : /* -------------------------------------------------------------------- */
174 : /* Confirm the requested access is supported. */
175 : /* -------------------------------------------------------------------- */
176 2 : if (poOpenInfo->eAccess == GA_Update)
177 : {
178 0 : ReportUpdateNotSupportedByDriver("NDF");
179 0 : return nullptr;
180 : }
181 : /* -------------------------------------------------------------------- */
182 : /* Read and process the header into a local name/value */
183 : /* stringlist. We just take off the trailing semicolon. The */
184 : /* keyword is already separated from the value by an equal */
185 : /* sign. */
186 : /* -------------------------------------------------------------------- */
187 :
188 2 : const char *pszLine = nullptr;
189 2 : const int nHeaderMax = 1000;
190 2 : int nHeaderLines = 0;
191 : char **papszHeader =
192 2 : static_cast<char **>(CPLMalloc(sizeof(char *) * (nHeaderMax + 1)));
193 :
194 210 : while (nHeaderLines < nHeaderMax &&
195 212 : (pszLine = CPLReadLineL(poOpenInfo->fpL)) != nullptr &&
196 106 : !EQUAL(pszLine, "END_OF_HDR;"))
197 : {
198 104 : if (strstr(pszLine, "=") == nullptr)
199 0 : break;
200 :
201 104 : char *pszFixed = CPLStrdup(pszLine);
202 104 : if (pszFixed[strlen(pszFixed) - 1] == ';')
203 104 : pszFixed[strlen(pszFixed) - 1] = '\0';
204 :
205 104 : papszHeader[nHeaderLines++] = pszFixed;
206 104 : papszHeader[nHeaderLines] = nullptr;
207 : }
208 2 : CPL_IGNORE_RET_VAL(VSIFCloseL(poOpenInfo->fpL));
209 2 : poOpenInfo->fpL = nullptr;
210 :
211 2 : if (CSLFetchNameValue(papszHeader, "PIXELS_PER_LINE") == nullptr ||
212 2 : CSLFetchNameValue(papszHeader, "LINES_PER_DATA_FILE") == nullptr ||
213 6 : CSLFetchNameValue(papszHeader, "BITS_PER_PIXEL") == nullptr ||
214 2 : CSLFetchNameValue(papszHeader, "PIXEL_FORMAT") == nullptr)
215 : {
216 0 : CPLError(CE_Failure, CPLE_AppDefined,
217 : "Dataset appears to be NDF but is missing a required field.");
218 0 : CSLDestroy(papszHeader);
219 0 : return nullptr;
220 : }
221 :
222 4 : if (!EQUAL(CSLFetchNameValue(papszHeader, "PIXEL_FORMAT"), "BYTE") ||
223 2 : !EQUAL(CSLFetchNameValue(papszHeader, "BITS_PER_PIXEL"), "8"))
224 : {
225 0 : CPLError(CE_Failure, CPLE_AppDefined,
226 : "Currently NDF driver supports only 8bit BYTE format.");
227 0 : CSLDestroy(papszHeader);
228 0 : return nullptr;
229 : }
230 :
231 : /* -------------------------------------------------------------------- */
232 : /* Confirm the requested access is supported. */
233 : /* -------------------------------------------------------------------- */
234 2 : if (poOpenInfo->eAccess == GA_Update)
235 : {
236 0 : CSLDestroy(papszHeader);
237 0 : ReportUpdateNotSupportedByDriver("NDF");
238 0 : return nullptr;
239 : }
240 :
241 : /* -------------------------------------------------------------------- */
242 : /* Create a corresponding GDALDataset. */
243 : /* -------------------------------------------------------------------- */
244 4 : auto poDS = std::make_unique<NDFDataset>();
245 2 : poDS->papszHeader = papszHeader;
246 :
247 2 : poDS->nRasterXSize = atoi(poDS->Get("PIXELS_PER_LINE", ""));
248 2 : poDS->nRasterYSize = atoi(poDS->Get("LINES_PER_DATA_FILE", ""));
249 :
250 : /* -------------------------------------------------------------------- */
251 : /* Create a raw raster band for each file. */
252 : /* -------------------------------------------------------------------- */
253 : const char *pszBand =
254 2 : CSLFetchNameValue(papszHeader, "NUMBER_OF_BANDS_IN_VOLUME");
255 2 : if (pszBand == nullptr)
256 : {
257 0 : CPLError(CE_Failure, CPLE_AppDefined, "Cannot find band count");
258 0 : return nullptr;
259 : }
260 2 : const int nBands = atoi(pszBand);
261 :
262 4 : if (!GDALCheckDatasetDimensions(poDS->nRasterXSize, poDS->nRasterYSize) ||
263 2 : !GDALCheckBandCount(nBands, FALSE))
264 : {
265 0 : return nullptr;
266 : }
267 :
268 4 : for (int iBand = 0; iBand < nBands; iBand++)
269 : {
270 : char szKey[100];
271 2 : snprintf(szKey, sizeof(szKey), "BAND%d_FILENAME", iBand + 1);
272 2 : CPLString osFilename = poDS->Get(szKey, "");
273 :
274 : // NDF1 file do not include the band filenames.
275 2 : if (osFilename.empty())
276 : {
277 : char szBandExtension[15];
278 0 : snprintf(szBandExtension, sizeof(szBandExtension), "I%d",
279 : iBand + 1);
280 : osFilename =
281 0 : CPLResetExtensionSafe(poOpenInfo->pszFilename, szBandExtension);
282 : }
283 : else
284 : {
285 2 : if (CPLHasPathTraversal(osFilename.c_str()))
286 : {
287 0 : CPLError(CE_Failure, CPLE_AppDefined,
288 : "Path traversal detected in %s", osFilename.c_str());
289 0 : return nullptr;
290 : }
291 2 : CPLString osBasePath = CPLGetPathSafe(poOpenInfo->pszFilename);
292 2 : osFilename = CPLFormFilenameSafe(osBasePath, osFilename, nullptr);
293 : }
294 :
295 2 : VSILFILE *fpRaw = VSIFOpenL(osFilename, "rb");
296 2 : if (fpRaw == nullptr)
297 : {
298 0 : CPLError(CE_Failure, CPLE_AppDefined,
299 : "Failed to open band file: %s", osFilename.c_str());
300 0 : return nullptr;
301 : }
302 2 : poDS->papszExtraFiles = CSLAddString(poDS->papszExtraFiles, osFilename);
303 :
304 : auto poBand = RawRasterBand::Create(
305 2 : poDS.get(), iBand + 1, fpRaw, 0, 1, poDS->nRasterXSize, GDT_Byte,
306 : RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN,
307 2 : RawRasterBand::OwnFP::YES);
308 2 : if (!poBand)
309 0 : return nullptr;
310 :
311 2 : snprintf(szKey, sizeof(szKey), "BAND%d_NAME", iBand + 1);
312 2 : poBand->SetDescription(poDS->Get(szKey, ""));
313 :
314 2 : snprintf(szKey, sizeof(szKey), "BAND%d_WAVELENGTHS", iBand + 1);
315 2 : poBand->SetMetadataItem("WAVELENGTHS", poDS->Get(szKey, ""));
316 :
317 2 : snprintf(szKey, sizeof(szKey), "BAND%d_RADIOMETRIC_GAINS/BIAS",
318 : iBand + 1);
319 2 : poBand->SetMetadataItem("RADIOMETRIC_GAINS_BIAS", poDS->Get(szKey, ""));
320 :
321 2 : poDS->SetBand(iBand + 1, std::move(poBand));
322 : }
323 :
324 : /* -------------------------------------------------------------------- */
325 : /* Fetch and parse USGS projection parameters. */
326 : /* -------------------------------------------------------------------- */
327 2 : double adfUSGSParams[15] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
328 : const CPLStringList aosParamTokens(CSLTokenizeStringComplex(
329 4 : poDS->Get("USGS_PROJECTION_NUMBER", ""), ",", FALSE, TRUE));
330 :
331 2 : if (aosParamTokens.size() >= 15)
332 : {
333 0 : for (int i = 0; i < 15; i++)
334 0 : adfUSGSParams[i] = CPLAtof(aosParamTokens[i]);
335 : }
336 :
337 : /* -------------------------------------------------------------------- */
338 : /* Minimal georef support ... should add full USGS style */
339 : /* support at some point. */
340 : /* -------------------------------------------------------------------- */
341 2 : const int nUSGSProjection = atoi(poDS->Get("USGS_PROJECTION_NUMBER", ""));
342 2 : const int nZone = atoi(poDS->Get("USGS_MAP_ZONE", "0"));
343 :
344 4 : OGRSpatialReference oSRS;
345 2 : oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
346 2 : oSRS.importFromUSGS(nUSGSProjection, nZone, adfUSGSParams, 12);
347 :
348 4 : const CPLString osDatum = poDS->Get("HORIZONTAL_DATUM", "");
349 2 : if (EQUAL(osDatum, "WGS84") || EQUAL(osDatum, "NAD83") ||
350 0 : EQUAL(osDatum, "NAD27"))
351 : {
352 2 : oSRS.SetWellKnownGeogCS(osDatum);
353 : }
354 0 : else if (STARTS_WITH_CI(osDatum, "NAD27"))
355 : {
356 0 : oSRS.SetWellKnownGeogCS("NAD27");
357 : }
358 : else
359 : {
360 0 : CPLError(CE_Warning, CPLE_AppDefined,
361 : "Unrecognized datum name in NLAPS/NDF file:%s, "
362 : "assuming WGS84.",
363 : osDatum.c_str());
364 0 : oSRS.SetWellKnownGeogCS("WGS84");
365 : }
366 :
367 2 : if (oSRS.GetRoot() != nullptr)
368 : {
369 2 : poDS->m_oSRS = std::move(oSRS);
370 : }
371 :
372 : /* -------------------------------------------------------------------- */
373 : /* Get geotransform. */
374 : /* -------------------------------------------------------------------- */
375 : char **papszUL =
376 2 : CSLTokenizeString2(poDS->Get("UPPER_LEFT_CORNER", ""), ",", 0);
377 : char **papszUR =
378 2 : CSLTokenizeString2(poDS->Get("UPPER_RIGHT_CORNER", ""), ",", 0);
379 : char **papszLL =
380 2 : CSLTokenizeString2(poDS->Get("LOWER_LEFT_CORNER", ""), ",", 0);
381 :
382 4 : if (CSLCount(papszUL) == 4 && CSLCount(papszUR) == 4 &&
383 2 : CSLCount(papszLL) == 4)
384 : {
385 2 : poDS->m_gt[0] = CPLAtof(papszUL[2]);
386 4 : poDS->m_gt[1] = (CPLAtof(papszUR[2]) - CPLAtof(papszUL[2])) /
387 2 : (poDS->nRasterXSize - 1);
388 4 : poDS->m_gt[2] = (CPLAtof(papszUR[3]) - CPLAtof(papszUL[3])) /
389 2 : (poDS->nRasterXSize - 1);
390 :
391 2 : poDS->m_gt[3] = CPLAtof(papszUL[3]);
392 4 : poDS->m_gt[4] = (CPLAtof(papszLL[2]) - CPLAtof(papszUL[2])) /
393 2 : (poDS->nRasterYSize - 1);
394 4 : poDS->m_gt[5] = (CPLAtof(papszLL[3]) - CPLAtof(papszUL[3])) /
395 2 : (poDS->nRasterYSize - 1);
396 :
397 : // Move origin up-left half a pixel.
398 2 : poDS->m_gt[0] -= poDS->m_gt[1] * 0.5;
399 2 : poDS->m_gt[0] -= poDS->m_gt[4] * 0.5;
400 2 : poDS->m_gt[3] -= poDS->m_gt[2] * 0.5;
401 2 : poDS->m_gt[3] -= poDS->m_gt[5] * 0.5;
402 : }
403 :
404 2 : CSLDestroy(papszUL);
405 2 : CSLDestroy(papszLL);
406 2 : CSLDestroy(papszUR);
407 :
408 : /* -------------------------------------------------------------------- */
409 : /* Initialize any PAM information. */
410 : /* -------------------------------------------------------------------- */
411 2 : poDS->SetDescription(poOpenInfo->pszFilename);
412 2 : poDS->TryLoadXML();
413 :
414 : /* -------------------------------------------------------------------- */
415 : /* Check for overviews. */
416 : /* -------------------------------------------------------------------- */
417 2 : poDS->oOvManager.Initialize(poDS.get(), poOpenInfo->pszFilename);
418 :
419 2 : return poDS.release();
420 : }
421 :
422 : /************************************************************************/
423 : /* GDALRegister_NDF() */
424 : /************************************************************************/
425 :
426 2038 : void GDALRegister_NDF()
427 :
428 : {
429 2038 : if (GDALGetDriverByName("NDF") != nullptr)
430 283 : return;
431 :
432 1755 : GDALDriver *poDriver = new GDALDriver();
433 :
434 1755 : poDriver->SetDescription("NDF");
435 1755 : poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
436 1755 : poDriver->SetMetadataItem(GDAL_DMD_LONGNAME, "NLAPS Data Format");
437 1755 : poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/ndf.html");
438 1755 : poDriver->SetMetadataItem(GDAL_DCAP_VIRTUALIO, "YES");
439 :
440 1755 : poDriver->pfnIdentify = NDFDataset::Identify;
441 1755 : poDriver->pfnOpen = NDFDataset::Open;
442 :
443 1755 : GetGDALDriverManager()->RegisterDriver(poDriver);
444 : }
|