Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GXF Reader
4 : * Purpose: GDAL binding for GXF reader.
5 : * Author: Frank Warmerdam, warmerdam@pobox.com
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 1998, Frank Warmerdam
9 : * Copyright (c) 2008-2012, Even Rouault <even dot rouault at spatialys.com>
10 : *
11 : * SPDX-License-Identifier: MIT
12 : ****************************************************************************/
13 :
14 : #include "gdal_frmts.h"
15 : #include "gdal_pam.h"
16 : #include "gxfopen.h"
17 :
18 : /************************************************************************/
19 : /* ==================================================================== */
20 : /* GXFDataset */
21 : /* ==================================================================== */
22 : /************************************************************************/
23 :
24 : class GXFRasterBand;
25 :
26 : class GXFDataset final : public GDALPamDataset
27 : {
28 : friend class GXFRasterBand;
29 :
30 : GXFHandle hGXF;
31 :
32 : OGRSpatialReference m_oSRS{};
33 : double dfNoDataValue;
34 : GDALDataType eDataType;
35 :
36 : public:
37 : GXFDataset();
38 : ~GXFDataset();
39 :
40 : static GDALDataset *Open(GDALOpenInfo *);
41 :
42 : CPLErr GetGeoTransform(GDALGeoTransform >) const override;
43 : const OGRSpatialReference *GetSpatialRef() const override;
44 : };
45 :
46 : /************************************************************************/
47 : /* ==================================================================== */
48 : /* GXFRasterBand */
49 : /* ==================================================================== */
50 : /************************************************************************/
51 :
52 : class GXFRasterBand final : public GDALPamRasterBand
53 : {
54 : friend class GXFDataset;
55 :
56 : public:
57 : GXFRasterBand(GXFDataset *, int);
58 : double GetNoDataValue(int *bGotNoDataValue) override;
59 :
60 : virtual CPLErr IReadBlock(int, int, void *) override;
61 : };
62 :
63 : /************************************************************************/
64 : /* GXFRasterBand() */
65 : /************************************************************************/
66 :
67 4 : GXFRasterBand::GXFRasterBand(GXFDataset *poDSIn, int nBandIn)
68 :
69 : {
70 4 : poDS = poDSIn;
71 4 : nBand = nBandIn;
72 :
73 4 : eDataType = poDSIn->eDataType;
74 :
75 4 : nBlockXSize = poDS->GetRasterXSize();
76 4 : nBlockYSize = 1;
77 4 : }
78 :
79 : /************************************************************************/
80 : /* GetNoDataValue() */
81 : /************************************************************************/
82 :
83 0 : double GXFRasterBand::GetNoDataValue(int *bGotNoDataValue)
84 :
85 : {
86 0 : GXFDataset *poGXF_DS = cpl::down_cast<GXFDataset *>(poDS);
87 0 : if (bGotNoDataValue)
88 0 : *bGotNoDataValue = (fabs(poGXF_DS->dfNoDataValue - -1e12) > .1);
89 0 : if (eDataType == GDT_Float32)
90 0 : return (double)(float)poGXF_DS->dfNoDataValue;
91 :
92 0 : return poGXF_DS->dfNoDataValue;
93 : }
94 :
95 : /************************************************************************/
96 : /* IReadBlock() */
97 : /************************************************************************/
98 :
99 11 : CPLErr GXFRasterBand::IReadBlock(CPL_UNUSED int nBlockXOff, int nBlockYOff,
100 : void *pImage)
101 : {
102 11 : GXFDataset *const poGXF_DS = cpl::down_cast<GXFDataset *>(poDS);
103 :
104 11 : if (eDataType == GDT_Float32)
105 : {
106 11 : double *padfBuffer = (double *)VSIMalloc2(sizeof(double), nBlockXSize);
107 11 : if (padfBuffer == nullptr)
108 0 : return CE_Failure;
109 : const CPLErr eErr =
110 11 : GXFGetScanline(poGXF_DS->hGXF, nBlockYOff, padfBuffer);
111 :
112 11 : float *pafBuffer = (float *)pImage;
113 103 : for (int i = 0; i < nBlockXSize; i++)
114 92 : pafBuffer[i] = (float)padfBuffer[i];
115 :
116 11 : CPLFree(padfBuffer);
117 :
118 11 : return eErr;
119 : }
120 :
121 : const CPLErr eErr =
122 0 : eDataType == GDT_Float64
123 0 : ? GXFGetScanline(poGXF_DS->hGXF, nBlockYOff, (double *)pImage)
124 0 : : CE_Failure;
125 :
126 0 : return eErr;
127 : }
128 :
129 : /************************************************************************/
130 : /* ==================================================================== */
131 : /* GXFDataset */
132 : /* ==================================================================== */
133 : /************************************************************************/
134 :
135 : /************************************************************************/
136 : /* GXFDataset() */
137 : /************************************************************************/
138 :
139 4 : GXFDataset::GXFDataset()
140 4 : : hGXF(nullptr), dfNoDataValue(0), eDataType(GDT_Float32)
141 : {
142 4 : m_oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
143 4 : }
144 :
145 : /************************************************************************/
146 : /* ~GXFDataset() */
147 : /************************************************************************/
148 :
149 8 : GXFDataset::~GXFDataset()
150 :
151 : {
152 4 : FlushCache(true);
153 4 : if (hGXF != nullptr)
154 4 : GXFClose(hGXF);
155 8 : }
156 :
157 : /************************************************************************/
158 : /* GetGeoTransform() */
159 : /************************************************************************/
160 :
161 0 : CPLErr GXFDataset::GetGeoTransform(GDALGeoTransform >) const
162 :
163 : {
164 0 : double dfXOrigin = 0.0;
165 0 : double dfYOrigin = 0.0;
166 0 : double dfXSize = 0.0;
167 0 : double dfYSize = 0.0;
168 0 : double dfRotation = 0.0;
169 :
170 0 : const CPLErr eErr = GXFGetPosition(hGXF, &dfXOrigin, &dfYOrigin, &dfXSize,
171 : &dfYSize, &dfRotation);
172 :
173 0 : if (eErr != CE_None)
174 0 : return eErr;
175 :
176 : // Transform to radians.
177 0 : dfRotation = (dfRotation / 360.0) * 2.0 * M_PI;
178 :
179 0 : gt[1] = dfXSize * cos(dfRotation);
180 0 : gt[2] = dfYSize * sin(dfRotation);
181 0 : gt[4] = dfXSize * sin(dfRotation);
182 0 : gt[5] = -1 * dfYSize * cos(dfRotation);
183 :
184 : // take into account that GXF is point or center of pixel oriented.
185 0 : gt[0] = dfXOrigin - 0.5 * gt[1] - 0.5 * gt[2];
186 0 : gt[3] = dfYOrigin - 0.5 * gt[4] - 0.5 * gt[5];
187 :
188 0 : return CE_None;
189 : }
190 :
191 : /************************************************************************/
192 : /* GetSpatialRef() */
193 : /************************************************************************/
194 :
195 1 : const OGRSpatialReference *GXFDataset::GetSpatialRef() const
196 : {
197 1 : return m_oSRS.IsEmpty() ? nullptr : &m_oSRS;
198 : }
199 :
200 : /************************************************************************/
201 : /* Open() */
202 : /************************************************************************/
203 :
204 33216 : GDALDataset *GXFDataset::Open(GDALOpenInfo *poOpenInfo)
205 :
206 : {
207 : /* -------------------------------------------------------------------- */
208 : /* Before trying GXFOpen() we first verify that there is at */
209 : /* least one "\n#keyword" type signature in the first chunk of */
210 : /* the file. */
211 : /* -------------------------------------------------------------------- */
212 33216 : if (poOpenInfo->nHeaderBytes < 50 || poOpenInfo->fpL == nullptr)
213 30278 : return nullptr;
214 :
215 2938 : bool bFoundKeyword = false;
216 2938 : bool bFoundIllegal = false;
217 819861 : for (int i = 0; i < poOpenInfo->nHeaderBytes - 1; i++)
218 : {
219 819023 : if ((poOpenInfo->pabyHeader[i] == 10 ||
220 798008 : poOpenInfo->pabyHeader[i] == 13) &&
221 22170 : poOpenInfo->pabyHeader[i + 1] == '#')
222 : {
223 183 : if (STARTS_WITH((const char *)poOpenInfo->pabyHeader + i + 2,
224 : "include"))
225 0 : return nullptr;
226 183 : if (STARTS_WITH((const char *)poOpenInfo->pabyHeader + i + 2,
227 : "define"))
228 0 : return nullptr;
229 183 : if (STARTS_WITH((const char *)poOpenInfo->pabyHeader + i + 2,
230 : "ifdef"))
231 0 : return nullptr;
232 183 : bFoundKeyword = true;
233 : }
234 819023 : if (poOpenInfo->pabyHeader[i] == 0)
235 : {
236 2100 : bFoundIllegal = true;
237 2100 : break;
238 : }
239 : }
240 :
241 2938 : if (!bFoundKeyword || bFoundIllegal)
242 2895 : return nullptr;
243 :
244 : /* -------------------------------------------------------------------- */
245 : /* At this point it is plausible that this is a GXF file, but */
246 : /* we also now verify that there is a #GRID keyword before */
247 : /* passing it off to GXFOpen(). We check in the first 50K. */
248 : /* -------------------------------------------------------------------- */
249 43 : CPL_IGNORE_RET_VAL(poOpenInfo->TryToIngest(50000));
250 43 : bool bGotGrid = false;
251 :
252 43 : const char *pszBigBuf = (const char *)poOpenInfo->pabyHeader;
253 13085 : for (int i = 0; i < poOpenInfo->nHeaderBytes - 5 && !bGotGrid; i++)
254 : {
255 13042 : if (pszBigBuf[i] == '#' && STARTS_WITH_CI(pszBigBuf + i + 1, "GRID"))
256 4 : bGotGrid = true;
257 : }
258 :
259 43 : if (!bGotGrid)
260 39 : return nullptr;
261 :
262 4 : VSIFCloseL(poOpenInfo->fpL);
263 4 : poOpenInfo->fpL = nullptr;
264 :
265 : /* -------------------------------------------------------------------- */
266 : /* Try opening the dataset. */
267 : /* -------------------------------------------------------------------- */
268 :
269 4 : GXFHandle l_hGXF = GXFOpen(poOpenInfo->pszFilename);
270 :
271 4 : if (l_hGXF == nullptr)
272 0 : return nullptr;
273 :
274 : /* -------------------------------------------------------------------- */
275 : /* Confirm the requested access is supported. */
276 : /* -------------------------------------------------------------------- */
277 4 : if (poOpenInfo->eAccess == GA_Update)
278 : {
279 0 : GXFClose(l_hGXF);
280 0 : ReportUpdateNotSupportedByDriver("GXF");
281 0 : return nullptr;
282 : }
283 :
284 : /* -------------------------------------------------------------------- */
285 : /* Create a corresponding GDALDataset. */
286 : /* -------------------------------------------------------------------- */
287 4 : GXFDataset *poDS = new GXFDataset();
288 :
289 4 : const char *pszGXFDataType = CPLGetConfigOption("GXF_DATATYPE", "Float32");
290 4 : GDALDataType eDT = GDALGetDataTypeByName(pszGXFDataType);
291 4 : if (!(eDT == GDT_Float32 || eDT == GDT_Float64))
292 : {
293 0 : CPLError(CE_Warning, CPLE_NotSupported,
294 : "Unsupported value for GXF_DATATYPE : %s", pszGXFDataType);
295 0 : eDT = GDT_Float32;
296 : }
297 :
298 4 : poDS->hGXF = l_hGXF;
299 4 : poDS->eDataType = eDT;
300 :
301 : /* -------------------------------------------------------------------- */
302 : /* Establish the projection. */
303 : /* -------------------------------------------------------------------- */
304 4 : char *pszProjection = GXFGetMapProjectionAsOGCWKT(l_hGXF);
305 4 : if (pszProjection && pszProjection[0] != '\0')
306 2 : poDS->m_oSRS.importFromWkt(pszProjection);
307 4 : CPLFree(pszProjection);
308 :
309 : /* -------------------------------------------------------------------- */
310 : /* Capture some information from the file that is of interest. */
311 : /* -------------------------------------------------------------------- */
312 4 : GXFGetRawInfo(l_hGXF, &(poDS->nRasterXSize), &(poDS->nRasterYSize), nullptr,
313 : nullptr, nullptr, &(poDS->dfNoDataValue));
314 :
315 4 : if (poDS->nRasterXSize <= 0 || poDS->nRasterYSize <= 0)
316 : {
317 0 : CPLError(CE_Failure, CPLE_AppDefined, "Invalid dimensions : %d x %d",
318 : poDS->nRasterXSize, poDS->nRasterYSize);
319 0 : delete poDS;
320 0 : return nullptr;
321 : }
322 :
323 : /* -------------------------------------------------------------------- */
324 : /* Create band information objects. */
325 : /* -------------------------------------------------------------------- */
326 4 : poDS->nBands = 1;
327 4 : poDS->SetBand(1, new GXFRasterBand(poDS, 1));
328 :
329 : /* -------------------------------------------------------------------- */
330 : /* Initialize any PAM information. */
331 : /* -------------------------------------------------------------------- */
332 4 : poDS->SetDescription(poOpenInfo->pszFilename);
333 4 : poDS->TryLoadXML();
334 :
335 : /* -------------------------------------------------------------------- */
336 : /* Check for external overviews. */
337 : /* -------------------------------------------------------------------- */
338 8 : poDS->oOvManager.Initialize(poDS, poOpenInfo->pszFilename,
339 4 : poOpenInfo->GetSiblingFiles());
340 :
341 4 : return poDS;
342 : }
343 :
344 : /************************************************************************/
345 : /* GDALRegister_GXF() */
346 : /************************************************************************/
347 :
348 1928 : void GDALRegister_GXF()
349 :
350 : {
351 1928 : if (GDALGetDriverByName("GXF") != nullptr)
352 282 : return;
353 :
354 1646 : GDALDriver *poDriver = new GDALDriver();
355 :
356 1646 : poDriver->SetDescription("GXF");
357 1646 : poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
358 1646 : poDriver->SetMetadataItem(GDAL_DMD_LONGNAME,
359 1646 : "GeoSoft Grid Exchange Format");
360 1646 : poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/gxf.html");
361 1646 : poDriver->SetMetadataItem(GDAL_DMD_EXTENSION, "gxf");
362 1646 : poDriver->SetMetadataItem(GDAL_DCAP_VIRTUALIO, "YES");
363 :
364 1646 : poDriver->pfnOpen = GXFDataset::Open;
365 :
366 1646 : GetGDALDriverManager()->RegisterDriver(poDriver);
367 : }
|
