Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: eCognition
4 : * Purpose: Implementation of Erdas .LAN / .GIS format.
5 : * Author: Frank Warmerdam, warmerdam@pobox.com
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2004, Frank Warmerdam
9 : * Copyright (c) 2008-2011, Even Rouault <even dot rouault at spatialys.com>
10 : *
11 : * Permission is hereby granted, free of charge, to any person obtaining a
12 : * copy of this software and associated documentation files (the "Software"),
13 : * to deal in the Software without restriction, including without limitation
14 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
15 : * and/or sell copies of the Software, and to permit persons to whom the
16 : * Software is furnished to do so, subject to the following conditions:
17 : *
18 : * The above copyright notice and this permission notice shall be included
19 : * in all copies or substantial portions of the Software.
20 : *
21 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
22 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
23 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
24 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
25 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
26 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
27 : * DEALINGS IN THE SOFTWARE.
28 : ****************************************************************************/
29 :
30 : #include "cpl_string.h"
31 : #include "gdal_frmts.h"
32 : #include "ogr_spatialref.h"
33 : #include "rawdataset.h"
34 :
35 : #include <cmath>
36 :
37 : #include <algorithm>
38 :
39 : /**
40 :
41 : Erdas Header format: "HEAD74"
42 :
43 : Offset Size Type Description
44 : ------ ---- ---- -----------
45 : 0 6 char magic cookie / version (i.e. HEAD74).
46 : 6 2 Int16 Pixel type, 0=8bit, 1=4bit, 2=16bit
47 : 8 2 Int16 Number of Bands.
48 : 10 6 char Unknown.
49 : 16 4 Int32 Width
50 : 20 4 Int32 Height
51 : 24 4 Int32 X Start (offset in original file?)
52 : 28 4 Int32 Y Start (offset in original file?)
53 : 32 56 char Unknown.
54 : 88 2 Int16 0=LAT, 1=UTM, 2=StatePlane, 3- are projections?
55 : 90 2 Int16 Classes in coverage.
56 : 92 14 char Unknown.
57 : 106 2 Int16 Area Unit (0=none, 1=Acre, 2=Hectare, 3=Other)
58 : 108 4 Float32 Pixel area.
59 : 112 4 Float32 Upper Left corner X (center of pixel?)
60 : 116 4 Float32 Upper Left corner Y (center of pixel?)
61 : 120 4 Float32 Width of a pixel.
62 : 124 4 Float32 Height of a pixel.
63 :
64 : Erdas Header format: "HEADER"
65 :
66 : Offset Size Type Description
67 : ------ ---- ---- -----------
68 : 0 6 char magic cookie / version (i.e. HEAD74).
69 : 6 2 Int16 Pixel type, 0=8bit, 1=4bit, 2=16bit
70 : 8 2 Int16 Number of Bands.
71 : 10 6 char Unknown.
72 : 16 4 Float32 Width
73 : 20 4 Float32 Height
74 : 24 4 Int32 X Start (offset in original file?)
75 : 28 4 Int32 Y Start (offset in original file?)
76 : 32 56 char Unknown.
77 : 88 2 Int16 0=LAT, 1=UTM, 2=StatePlane, 3- are projections?
78 : 90 2 Int16 Classes in coverage.
79 : 92 14 char Unknown.
80 : 106 2 Int16 Area Unit (0=none, 1=Acre, 2=Hectare, 3=Other)
81 : 108 4 Float32 Pixel area.
82 : 112 4 Float32 Upper Left corner X (center of pixel?)
83 : 116 4 Float32 Upper Left corner Y (center of pixel?)
84 : 120 4 Float32 Width of a pixel.
85 : 124 4 Float32 Height of a pixel.
86 :
87 : All binary fields are in the same byte order but it may be big endian or
88 : little endian depending on what platform the file was written on. Usually
89 : this can be checked against the number of bands though this test won't work
90 : if there are more than 255 bands.
91 :
92 : There is also some information on .STA and .TRL files at:
93 :
94 : http://www.pcigeomatics.com/cgi-bin/pcihlp/ERDASWR%7CTRAILER+FORMAT
95 :
96 : **/
97 :
98 : constexpr int ERD_HEADER_SIZE = 128;
99 :
100 : /************************************************************************/
101 : /* ==================================================================== */
102 : /* LAN4BitRasterBand */
103 : /* ==================================================================== */
104 : /************************************************************************/
105 :
106 : class LANDataset;
107 :
108 : class LAN4BitRasterBand final : public GDALPamRasterBand
109 : {
110 : GDALColorTable *poCT;
111 : GDALColorInterp eInterp;
112 :
113 : CPL_DISALLOW_COPY_ASSIGN(LAN4BitRasterBand)
114 :
115 : public:
116 : LAN4BitRasterBand(LANDataset *, int);
117 : ~LAN4BitRasterBand() override;
118 :
119 : GDALColorTable *GetColorTable() override;
120 : GDALColorInterp GetColorInterpretation() override;
121 : CPLErr SetColorTable(GDALColorTable *) override;
122 : CPLErr SetColorInterpretation(GDALColorInterp) override;
123 :
124 : CPLErr IReadBlock(int, int, void *) override;
125 : };
126 :
127 : /************************************************************************/
128 : /* ==================================================================== */
129 : /* LANDataset */
130 : /* ==================================================================== */
131 : /************************************************************************/
132 :
133 : class LANDataset final : public RawDataset
134 : {
135 : CPL_DISALLOW_COPY_ASSIGN(LANDataset)
136 :
137 : public:
138 : VSILFILE *fpImage; // Image data file.
139 :
140 : char pachHeader[ERD_HEADER_SIZE];
141 :
142 : OGRSpatialReference *m_poSRS = nullptr;
143 :
144 : double adfGeoTransform[6];
145 :
146 : CPLString osSTAFilename{};
147 : void CheckForStatistics(void);
148 :
149 : char **GetFileList() override;
150 :
151 : CPLErr Close() override;
152 :
153 : public:
154 : LANDataset();
155 : ~LANDataset() override;
156 :
157 : CPLErr GetGeoTransform(double *padfTransform) override;
158 : CPLErr SetGeoTransform(double *padfTransform) override;
159 :
160 : const OGRSpatialReference *GetSpatialRef() const override;
161 : CPLErr SetSpatialRef(const OGRSpatialReference *poSRS) override;
162 :
163 : static GDALDataset *Open(GDALOpenInfo *);
164 : static GDALDataset *Create(const char *pszFilename, int nXSize, int nYSize,
165 : int nBandsIn, GDALDataType eType,
166 : char **papszOptions);
167 : };
168 :
169 : /************************************************************************/
170 : /* ==================================================================== */
171 : /* LAN4BitRasterBand */
172 : /* ==================================================================== */
173 : /************************************************************************/
174 :
175 : /************************************************************************/
176 : /* LAN4BitRasterBand() */
177 : /************************************************************************/
178 :
179 2 : LAN4BitRasterBand::LAN4BitRasterBand(LANDataset *poDSIn, int nBandIn)
180 2 : : poCT(nullptr), eInterp(GCI_Undefined)
181 : {
182 2 : poDS = poDSIn;
183 2 : nBand = nBandIn;
184 2 : eDataType = GDT_Byte;
185 :
186 2 : nBlockXSize = poDSIn->GetRasterXSize();
187 2 : nBlockYSize = 1;
188 2 : }
189 :
190 : /************************************************************************/
191 : /* ~LAN4BitRasterBand() */
192 : /************************************************************************/
193 :
194 4 : LAN4BitRasterBand::~LAN4BitRasterBand()
195 :
196 : {
197 2 : if (poCT)
198 0 : delete poCT;
199 4 : }
200 :
201 : /************************************************************************/
202 : /* IReadBlock() */
203 : /************************************************************************/
204 :
205 2 : CPLErr LAN4BitRasterBand::IReadBlock(CPL_UNUSED int nBlockXOff, int nBlockYOff,
206 : void *pImage)
207 :
208 : {
209 2 : LANDataset *poLAN_DS = reinterpret_cast<LANDataset *>(poDS);
210 2 : CPLAssert(nBlockXOff == 0);
211 :
212 : /* -------------------------------------------------------------------- */
213 : /* Seek to profile. */
214 : /* -------------------------------------------------------------------- */
215 : const vsi_l_offset nOffset =
216 : ERD_HEADER_SIZE +
217 4 : (static_cast<vsi_l_offset>(nBlockYOff) * nRasterXSize *
218 2 : poLAN_DS->GetRasterCount()) /
219 2 : 2 +
220 2 : (static_cast<vsi_l_offset>(nBand - 1) * nRasterXSize) / 2;
221 :
222 2 : if (VSIFSeekL(poLAN_DS->fpImage, nOffset, SEEK_SET) != 0)
223 : {
224 0 : CPLError(CE_Failure, CPLE_FileIO, "LAN Seek failed:%s",
225 0 : VSIStrerror(errno));
226 0 : return CE_Failure;
227 : }
228 :
229 : /* -------------------------------------------------------------------- */
230 : /* Read the profile. */
231 : /* -------------------------------------------------------------------- */
232 2 : if (VSIFReadL(pImage, 1, nRasterXSize / 2, poLAN_DS->fpImage) !=
233 2 : static_cast<size_t>(nRasterXSize) / 2)
234 : {
235 0 : CPLError(CE_Failure, CPLE_FileIO, "LAN Read failed:%s",
236 0 : VSIStrerror(errno));
237 0 : return CE_Failure;
238 : }
239 :
240 : /* -------------------------------------------------------------------- */
241 : /* Convert 4bit to 8bit. */
242 : /* -------------------------------------------------------------------- */
243 6 : for (int i = nRasterXSize - 1; i >= 0; i--)
244 : {
245 4 : if ((i & 0x01) != 0)
246 2 : reinterpret_cast<GByte *>(pImage)[i] =
247 2 : reinterpret_cast<GByte *>(pImage)[i / 2] & 0x0f;
248 : else
249 2 : reinterpret_cast<GByte *>(pImage)[i] =
250 2 : (reinterpret_cast<GByte *>(pImage)[i / 2] & 0xf0) / 16;
251 : }
252 :
253 2 : return CE_None;
254 : }
255 :
256 : /************************************************************************/
257 : /* SetColorTable() */
258 : /************************************************************************/
259 :
260 0 : CPLErr LAN4BitRasterBand::SetColorTable(GDALColorTable *poNewCT)
261 :
262 : {
263 0 : if (poCT)
264 0 : delete poCT;
265 0 : if (poNewCT == nullptr)
266 0 : poCT = nullptr;
267 : else
268 0 : poCT = poNewCT->Clone();
269 :
270 0 : return CE_None;
271 : }
272 :
273 : /************************************************************************/
274 : /* GetColorTable() */
275 : /************************************************************************/
276 :
277 0 : GDALColorTable *LAN4BitRasterBand::GetColorTable()
278 :
279 : {
280 0 : if (poCT != nullptr)
281 0 : return poCT;
282 :
283 0 : return GDALPamRasterBand::GetColorTable();
284 : }
285 :
286 : /************************************************************************/
287 : /* SetColorInterpretation() */
288 : /************************************************************************/
289 :
290 0 : CPLErr LAN4BitRasterBand::SetColorInterpretation(GDALColorInterp eNewInterp)
291 :
292 : {
293 0 : eInterp = eNewInterp;
294 :
295 0 : return CE_None;
296 : }
297 :
298 : /************************************************************************/
299 : /* GetColorInterpretation() */
300 : /************************************************************************/
301 :
302 0 : GDALColorInterp LAN4BitRasterBand::GetColorInterpretation()
303 :
304 : {
305 0 : return eInterp;
306 : }
307 :
308 : /************************************************************************/
309 : /* ==================================================================== */
310 : /* LANDataset */
311 : /* ==================================================================== */
312 : /************************************************************************/
313 :
314 : /************************************************************************/
315 : /* LANDataset() */
316 : /************************************************************************/
317 :
318 26 : LANDataset::LANDataset() : fpImage(nullptr)
319 : {
320 26 : memset(pachHeader, 0, sizeof(pachHeader));
321 26 : adfGeoTransform[0] = 0.0;
322 26 : adfGeoTransform[1] = 0.0; // TODO(schwehr): Should this be 1.0?
323 26 : adfGeoTransform[2] = 0.0;
324 26 : adfGeoTransform[3] = 0.0;
325 26 : adfGeoTransform[4] = 0.0;
326 26 : adfGeoTransform[5] = 0.0; // TODO(schwehr): Should this be 1.0?
327 26 : }
328 :
329 : /************************************************************************/
330 : /* ~LANDataset() */
331 : /************************************************************************/
332 :
333 52 : LANDataset::~LANDataset()
334 :
335 : {
336 26 : LANDataset::Close();
337 52 : }
338 :
339 : /************************************************************************/
340 : /* Close() */
341 : /************************************************************************/
342 :
343 50 : CPLErr LANDataset::Close()
344 : {
345 50 : CPLErr eErr = CE_None;
346 50 : if (nOpenFlags != OPEN_FLAGS_CLOSED)
347 : {
348 26 : if (LANDataset::FlushCache(true) != CE_None)
349 0 : eErr = CE_Failure;
350 :
351 26 : if (fpImage)
352 : {
353 26 : if (VSIFCloseL(fpImage) != 0)
354 : {
355 0 : CPLError(CE_Failure, CPLE_FileIO, "I/O error");
356 0 : eErr = CE_Failure;
357 : }
358 : }
359 :
360 26 : if (m_poSRS)
361 24 : m_poSRS->Release();
362 :
363 26 : if (GDALPamDataset::Close() != CE_None)
364 0 : eErr = CE_Failure;
365 : }
366 50 : return eErr;
367 : }
368 :
369 : /************************************************************************/
370 : /* Open() */
371 : /************************************************************************/
372 :
373 29185 : GDALDataset *LANDataset::Open(GDALOpenInfo *poOpenInfo)
374 :
375 : {
376 : /* -------------------------------------------------------------------- */
377 : /* We assume the user is pointing to the header (.pcb) file. */
378 : /* Does this appear to be a pcb file? */
379 : /* -------------------------------------------------------------------- */
380 29185 : if (poOpenInfo->nHeaderBytes < ERD_HEADER_SIZE ||
381 2930 : poOpenInfo->fpL == nullptr)
382 26293 : return nullptr;
383 :
384 2892 : if (!STARTS_WITH_CI(reinterpret_cast<char *>(poOpenInfo->pabyHeader),
385 2892 : "HEADER") &&
386 2892 : !STARTS_WITH_CI(reinterpret_cast<char *>(poOpenInfo->pabyHeader),
387 : "HEAD74"))
388 2866 : return nullptr;
389 :
390 26 : if (memcmp(poOpenInfo->pabyHeader + 16, "S LAT ", 8) == 0)
391 : {
392 : // NTV1 format
393 0 : return nullptr;
394 : }
395 :
396 : /* -------------------------------------------------------------------- */
397 : /* Create a corresponding GDALDataset. */
398 : /* -------------------------------------------------------------------- */
399 52 : auto poDS = std::make_unique<LANDataset>();
400 :
401 26 : poDS->eAccess = poOpenInfo->eAccess;
402 26 : std::swap(poDS->fpImage, poOpenInfo->fpL);
403 :
404 : /* -------------------------------------------------------------------- */
405 : /* Do we need to byte swap the headers to local machine order? */
406 : /* -------------------------------------------------------------------- */
407 26 : const RawRasterBand::ByteOrder eByteOrder =
408 26 : poOpenInfo->pabyHeader[8] == 0
409 26 : ? RawRasterBand::ByteOrder::ORDER_BIG_ENDIAN
410 : : RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN;
411 :
412 26 : memcpy(poDS->pachHeader, poOpenInfo->pabyHeader, ERD_HEADER_SIZE);
413 :
414 26 : if (eByteOrder != RawRasterBand::NATIVE_BYTE_ORDER)
415 : {
416 2 : CPL_SWAP16PTR(poDS->pachHeader + 6);
417 2 : CPL_SWAP16PTR(poDS->pachHeader + 8);
418 :
419 2 : CPL_SWAP32PTR(poDS->pachHeader + 16);
420 2 : CPL_SWAP32PTR(poDS->pachHeader + 20);
421 2 : CPL_SWAP32PTR(poDS->pachHeader + 24);
422 2 : CPL_SWAP32PTR(poDS->pachHeader + 28);
423 :
424 2 : CPL_SWAP16PTR(poDS->pachHeader + 88);
425 2 : CPL_SWAP16PTR(poDS->pachHeader + 90);
426 :
427 2 : CPL_SWAP16PTR(poDS->pachHeader + 106);
428 2 : CPL_SWAP32PTR(poDS->pachHeader + 108);
429 2 : CPL_SWAP32PTR(poDS->pachHeader + 112);
430 2 : CPL_SWAP32PTR(poDS->pachHeader + 116);
431 2 : CPL_SWAP32PTR(poDS->pachHeader + 120);
432 2 : CPL_SWAP32PTR(poDS->pachHeader + 124);
433 : }
434 :
435 : /* -------------------------------------------------------------------- */
436 : /* Capture some information from the file that is of interest. */
437 : /* -------------------------------------------------------------------- */
438 26 : if (STARTS_WITH_CI(poDS->pachHeader, "HEADER"))
439 : {
440 0 : float fTmp = 0.0;
441 0 : memcpy(&fTmp, poDS->pachHeader + 16, 4);
442 0 : poDS->nRasterXSize = static_cast<int>(fTmp);
443 0 : memcpy(&fTmp, poDS->pachHeader + 20, 4);
444 0 : poDS->nRasterYSize = static_cast<int>(fTmp);
445 : }
446 : else
447 : {
448 26 : GInt32 nTmp = 0;
449 26 : memcpy(&nTmp, poDS->pachHeader + 16, 4);
450 26 : poDS->nRasterXSize = nTmp;
451 26 : memcpy(&nTmp, poDS->pachHeader + 20, 4);
452 26 : poDS->nRasterYSize = nTmp;
453 : }
454 :
455 26 : GInt16 nTmp16 = 0;
456 26 : memcpy(&nTmp16, poDS->pachHeader + 6, 2);
457 :
458 26 : int nPixelOffset = 0;
459 26 : GDALDataType eDataType = GDT_Unknown;
460 26 : if (nTmp16 == 0)
461 : {
462 19 : eDataType = GDT_Byte;
463 19 : nPixelOffset = 1;
464 : }
465 7 : else if (nTmp16 == 1) // 4 bit
466 : {
467 2 : eDataType = GDT_Byte;
468 2 : nPixelOffset = -1;
469 : }
470 5 : else if (nTmp16 == 2)
471 : {
472 5 : nPixelOffset = 2;
473 5 : eDataType = GDT_Int16;
474 : }
475 : else
476 : {
477 0 : CPLError(CE_Failure, CPLE_AppDefined, "Unsupported pixel type (%d).",
478 : nTmp16);
479 0 : return nullptr;
480 : }
481 :
482 26 : memcpy(&nTmp16, poDS->pachHeader + 8, 2);
483 26 : const int nBandCount = nTmp16;
484 :
485 52 : if (!GDALCheckDatasetDimensions(poDS->nRasterXSize, poDS->nRasterYSize) ||
486 26 : !GDALCheckBandCount(nBandCount, FALSE))
487 : {
488 2 : return nullptr;
489 : }
490 :
491 : // cppcheck-suppress knownConditionTrueFalse
492 46 : if (nPixelOffset != -1 &&
493 22 : poDS->nRasterXSize > INT_MAX / (nPixelOffset * nBandCount))
494 : {
495 0 : CPLError(CE_Failure, CPLE_AppDefined, "Int overflow occurred.");
496 0 : return nullptr;
497 : }
498 :
499 : /* -------------------------------------------------------------------- */
500 : /* Create band information object. */
501 : /* -------------------------------------------------------------------- */
502 82 : for (int iBand = 1; iBand <= nBandCount; iBand++)
503 : {
504 58 : if (nPixelOffset == -1) /* 4 bit case */
505 2 : poDS->SetBand(iBand, new LAN4BitRasterBand(poDS.get(), iBand));
506 : else
507 : {
508 : auto poBand = RawRasterBand::Create(
509 112 : poDS.get(), iBand, poDS->fpImage,
510 56 : ERD_HEADER_SIZE +
511 112 : (iBand - 1) * nPixelOffset * poDS->nRasterXSize,
512 56 : nPixelOffset, poDS->nRasterXSize * nPixelOffset * nBandCount,
513 56 : eDataType, eByteOrder, RawRasterBand::OwnFP::NO);
514 56 : if (!poBand)
515 0 : return nullptr;
516 56 : poDS->SetBand(iBand, std::move(poBand));
517 : }
518 : }
519 :
520 : /* -------------------------------------------------------------------- */
521 : /* Initialize any PAM information. */
522 : /* -------------------------------------------------------------------- */
523 24 : poDS->SetDescription(poOpenInfo->pszFilename);
524 24 : poDS->CheckForStatistics();
525 24 : poDS->TryLoadXML();
526 :
527 : /* -------------------------------------------------------------------- */
528 : /* Check for overviews. */
529 : /* -------------------------------------------------------------------- */
530 24 : poDS->oOvManager.Initialize(poDS.get(), poOpenInfo->pszFilename);
531 :
532 : /* -------------------------------------------------------------------- */
533 : /* Try to interpret georeferencing. */
534 : /* -------------------------------------------------------------------- */
535 24 : float fTmp = 0.0;
536 :
537 24 : memcpy(&fTmp, poDS->pachHeader + 112, 4);
538 24 : poDS->adfGeoTransform[0] = fTmp;
539 24 : memcpy(&fTmp, poDS->pachHeader + 120, 4);
540 24 : poDS->adfGeoTransform[1] = fTmp;
541 24 : poDS->adfGeoTransform[2] = 0.0;
542 24 : memcpy(&fTmp, poDS->pachHeader + 116, 4);
543 24 : poDS->adfGeoTransform[3] = fTmp;
544 24 : poDS->adfGeoTransform[4] = 0.0;
545 24 : memcpy(&fTmp, poDS->pachHeader + 124, 4);
546 24 : poDS->adfGeoTransform[5] = -fTmp;
547 :
548 : // adjust for center of pixel vs. top left corner of pixel.
549 24 : poDS->adfGeoTransform[0] -= poDS->adfGeoTransform[1] * 0.5;
550 24 : poDS->adfGeoTransform[3] -= poDS->adfGeoTransform[5] * 0.5;
551 :
552 : /* -------------------------------------------------------------------- */
553 : /* If we didn't get any georeferencing, try for a worldfile. */
554 : /* -------------------------------------------------------------------- */
555 24 : if (poDS->adfGeoTransform[1] == 0.0 || poDS->adfGeoTransform[5] == 0.0)
556 : {
557 0 : if (!GDALReadWorldFile(poOpenInfo->pszFilename, nullptr,
558 0 : poDS->adfGeoTransform))
559 0 : GDALReadWorldFile(poOpenInfo->pszFilename, ".wld",
560 0 : poDS->adfGeoTransform);
561 : }
562 :
563 : /* -------------------------------------------------------------------- */
564 : /* Try to come up with something for the coordinate system. */
565 : /* -------------------------------------------------------------------- */
566 24 : memcpy(&nTmp16, poDS->pachHeader + 88, 2);
567 24 : int nCoordSys = nTmp16;
568 :
569 24 : poDS->m_poSRS = new OGRSpatialReference();
570 24 : poDS->m_poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
571 24 : if (nCoordSys == 0)
572 : {
573 24 : poDS->m_poSRS->SetFromUserInput(SRS_WKT_WGS84_LAT_LONG);
574 : }
575 0 : else if (nCoordSys == 1)
576 : {
577 0 : poDS->m_poSRS->SetFromUserInput(
578 : "LOCAL_CS[\"UTM - Zone Unknown\",UNIT[\"Meter\",1]]");
579 : }
580 0 : else if (nCoordSys == 2)
581 : {
582 0 : poDS->m_poSRS->SetFromUserInput(
583 : "LOCAL_CS[\"State Plane - Zone Unknown\","
584 : "UNIT[\"US survey foot\",0.3048006096012192]]");
585 : }
586 : else
587 : {
588 0 : poDS->m_poSRS->SetFromUserInput(
589 : "LOCAL_CS[\"Unknown\",UNIT[\"Meter\",1]]");
590 : }
591 :
592 : /* -------------------------------------------------------------------- */
593 : /* Check for a trailer file with a colormap in it. */
594 : /* -------------------------------------------------------------------- */
595 24 : char *pszPath = CPLStrdup(CPLGetPath(poOpenInfo->pszFilename));
596 24 : char *pszBasename = CPLStrdup(CPLGetBasename(poOpenInfo->pszFilename));
597 24 : const char *pszTRLFilename = CPLFormCIFilename(pszPath, pszBasename, "trl");
598 24 : VSILFILE *fpTRL = VSIFOpenL(pszTRLFilename, "rb");
599 24 : if (fpTRL != nullptr)
600 : {
601 0 : char szTRLData[896] = {'\0'};
602 :
603 0 : CPL_IGNORE_RET_VAL(VSIFReadL(szTRLData, 1, 896, fpTRL));
604 0 : CPL_IGNORE_RET_VAL(VSIFCloseL(fpTRL));
605 :
606 0 : GDALColorTable oCT;
607 0 : for (int iColor = 0; iColor < 256; iColor++)
608 : {
609 0 : GDALColorEntry sEntry = {0, 0, 0, 0};
610 :
611 0 : sEntry.c2 = reinterpret_cast<GByte *>(szTRLData)[iColor + 128];
612 0 : sEntry.c1 =
613 0 : reinterpret_cast<GByte *>(szTRLData)[iColor + 128 + 256];
614 0 : sEntry.c3 =
615 0 : reinterpret_cast<GByte *>(szTRLData)[iColor + 128 + 512];
616 0 : sEntry.c4 = 255;
617 0 : oCT.SetColorEntry(iColor, &sEntry);
618 :
619 : // Only 16 colors in 4bit files.
620 0 : if (nPixelOffset == -1 && iColor == 15)
621 0 : break;
622 : }
623 :
624 0 : poDS->GetRasterBand(1)->SetColorTable(&oCT);
625 0 : poDS->GetRasterBand(1)->SetColorInterpretation(GCI_PaletteIndex);
626 : }
627 :
628 24 : CPLFree(pszPath);
629 24 : CPLFree(pszBasename);
630 :
631 24 : return poDS.release();
632 : }
633 :
634 : /************************************************************************/
635 : /* GetGeoTransform() */
636 : /************************************************************************/
637 :
638 7 : CPLErr LANDataset::GetGeoTransform(double *padfTransform)
639 :
640 : {
641 7 : if (adfGeoTransform[1] != 0.0 && adfGeoTransform[5] != 0.0)
642 : {
643 7 : memcpy(padfTransform, adfGeoTransform, sizeof(double) * 6);
644 7 : return CE_None;
645 : }
646 :
647 0 : return GDALPamDataset::GetGeoTransform(padfTransform);
648 : }
649 :
650 : /************************************************************************/
651 : /* SetGeoTransform() */
652 : /************************************************************************/
653 :
654 13 : CPLErr LANDataset::SetGeoTransform(double *padfTransform)
655 :
656 : {
657 13 : unsigned char abyHeader[128] = {'\0'};
658 :
659 13 : memcpy(adfGeoTransform, padfTransform, sizeof(double) * 6);
660 :
661 13 : CPL_IGNORE_RET_VAL(VSIFSeekL(fpImage, 0, SEEK_SET));
662 13 : CPL_IGNORE_RET_VAL(VSIFReadL(abyHeader, 128, 1, fpImage));
663 :
664 : // Upper Left X.
665 13 : float f32Val =
666 13 : static_cast<float>(adfGeoTransform[0] + 0.5 * adfGeoTransform[1]);
667 13 : memcpy(abyHeader + 112, &f32Val, 4);
668 :
669 : // Upper Left Y.
670 13 : f32Val = static_cast<float>(adfGeoTransform[3] + 0.5 * adfGeoTransform[5]);
671 13 : memcpy(abyHeader + 116, &f32Val, 4);
672 :
673 : // Width of pixel.
674 13 : f32Val = static_cast<float>(adfGeoTransform[1]);
675 13 : memcpy(abyHeader + 120, &f32Val, 4);
676 :
677 : // Height of pixel.
678 13 : f32Val = static_cast<float>(std::abs(adfGeoTransform[5]));
679 13 : memcpy(abyHeader + 124, &f32Val, 4);
680 :
681 26 : if (VSIFSeekL(fpImage, 0, SEEK_SET) != 0 ||
682 13 : VSIFWriteL(abyHeader, 128, 1, fpImage) != 1)
683 : {
684 0 : CPLError(CE_Failure, CPLE_FileIO,
685 : "File IO Error writing header with new geotransform.");
686 0 : return CE_Failure;
687 : }
688 :
689 13 : return CE_None;
690 : }
691 :
692 : /************************************************************************/
693 : /* GetSpatialRef() */
694 : /* */
695 : /* Use PAM coordinate system if available in preference to the */
696 : /* generally poor value derived from the file itself. */
697 : /************************************************************************/
698 :
699 0 : const OGRSpatialReference *LANDataset::GetSpatialRef() const
700 :
701 : {
702 0 : const auto poSRS = GDALPamDataset::GetSpatialRef();
703 0 : if (poSRS)
704 0 : return poSRS;
705 :
706 0 : return m_poSRS;
707 : }
708 :
709 : /************************************************************************/
710 : /* SetSpatialRef() */
711 : /************************************************************************/
712 :
713 13 : CPLErr LANDataset::SetSpatialRef(const OGRSpatialReference *poSRS)
714 :
715 : {
716 13 : if (poSRS == nullptr)
717 0 : return GDALPamDataset::SetSpatialRef(poSRS);
718 :
719 13 : unsigned char abyHeader[128] = {'\0'};
720 :
721 13 : CPL_IGNORE_RET_VAL(VSIFSeekL(fpImage, 0, SEEK_SET));
722 13 : CPL_IGNORE_RET_VAL(VSIFReadL(abyHeader, 128, 1, fpImage));
723 :
724 13 : GUInt16 nProjCode = 0;
725 :
726 13 : if (poSRS->IsGeographic())
727 : {
728 13 : nProjCode = 0;
729 : }
730 0 : else if (poSRS->GetUTMZone() != 0)
731 : {
732 0 : nProjCode = 1;
733 : }
734 : // Too bad we have no way of recognising state plane projections.
735 : else
736 : {
737 0 : const char *l_pszProjection = poSRS->GetAttrValue("PROJECTION");
738 :
739 0 : if (l_pszProjection == nullptr)
740 : ;
741 0 : else if (EQUAL(l_pszProjection, SRS_PT_ALBERS_CONIC_EQUAL_AREA))
742 0 : nProjCode = 3;
743 0 : else if (EQUAL(l_pszProjection, SRS_PT_LAMBERT_CONFORMAL_CONIC_1SP))
744 0 : nProjCode = 4;
745 0 : else if (EQUAL(l_pszProjection, SRS_PT_MERCATOR_1SP))
746 0 : nProjCode = 5;
747 0 : else if (EQUAL(l_pszProjection, SRS_PT_POLAR_STEREOGRAPHIC))
748 0 : nProjCode = 6;
749 0 : else if (EQUAL(l_pszProjection, SRS_PT_POLYCONIC))
750 0 : nProjCode = 7;
751 0 : else if (EQUAL(l_pszProjection, SRS_PT_EQUIDISTANT_CONIC))
752 0 : nProjCode = 8;
753 0 : else if (EQUAL(l_pszProjection, SRS_PT_TRANSVERSE_MERCATOR))
754 0 : nProjCode = 9;
755 0 : else if (EQUAL(l_pszProjection, SRS_PT_STEREOGRAPHIC))
756 0 : nProjCode = 10;
757 0 : else if (EQUAL(l_pszProjection, SRS_PT_LAMBERT_AZIMUTHAL_EQUAL_AREA))
758 0 : nProjCode = 11;
759 0 : else if (EQUAL(l_pszProjection, SRS_PT_AZIMUTHAL_EQUIDISTANT))
760 0 : nProjCode = 12;
761 0 : else if (EQUAL(l_pszProjection, SRS_PT_GNOMONIC))
762 0 : nProjCode = 13;
763 0 : else if (EQUAL(l_pszProjection, SRS_PT_ORTHOGRAPHIC))
764 0 : nProjCode = 14;
765 : // We do not have GVNP.
766 0 : else if (EQUAL(l_pszProjection, SRS_PT_SINUSOIDAL))
767 0 : nProjCode = 16;
768 0 : else if (EQUAL(l_pszProjection, SRS_PT_EQUIRECTANGULAR))
769 0 : nProjCode = 17;
770 0 : else if (EQUAL(l_pszProjection, SRS_PT_MILLER_CYLINDRICAL))
771 0 : nProjCode = 18;
772 0 : else if (EQUAL(l_pszProjection, SRS_PT_VANDERGRINTEN))
773 0 : nProjCode = 19;
774 0 : else if (EQUAL(l_pszProjection, SRS_PT_HOTINE_OBLIQUE_MERCATOR))
775 0 : nProjCode = 20;
776 : }
777 :
778 13 : memcpy(abyHeader + 88, &nProjCode, 2);
779 :
780 13 : CPL_IGNORE_RET_VAL(VSIFSeekL(fpImage, 0, SEEK_SET));
781 13 : CPL_IGNORE_RET_VAL(VSIFWriteL(abyHeader, 128, 1, fpImage));
782 :
783 13 : return GDALPamDataset::SetSpatialRef(poSRS);
784 : }
785 :
786 : /************************************************************************/
787 : /* GetFileList() */
788 : /************************************************************************/
789 :
790 2 : char **LANDataset::GetFileList()
791 :
792 : {
793 : // Main data file, etc.
794 2 : char **papszFileList = GDALPamDataset::GetFileList();
795 :
796 2 : if (!osSTAFilename.empty())
797 0 : papszFileList = CSLAddString(papszFileList, osSTAFilename);
798 :
799 2 : return papszFileList;
800 : }
801 :
802 : /************************************************************************/
803 : /* CheckForStatistics() */
804 : /************************************************************************/
805 :
806 24 : void LANDataset::CheckForStatistics()
807 :
808 : {
809 : /* -------------------------------------------------------------------- */
810 : /* Do we have a statistics file? */
811 : /* -------------------------------------------------------------------- */
812 24 : osSTAFilename = CPLResetExtension(GetDescription(), "sta");
813 :
814 24 : VSILFILE *fpSTA = VSIFOpenL(osSTAFilename, "r");
815 :
816 24 : if (fpSTA == nullptr && VSIIsCaseSensitiveFS(osSTAFilename))
817 : {
818 24 : osSTAFilename = CPLResetExtension(GetDescription(), "STA");
819 24 : fpSTA = VSIFOpenL(osSTAFilename, "r");
820 : }
821 :
822 24 : if (fpSTA == nullptr)
823 : {
824 24 : osSTAFilename = "";
825 24 : return;
826 : }
827 :
828 : /* -------------------------------------------------------------------- */
829 : /* Read it one band at a time. */
830 : /* -------------------------------------------------------------------- */
831 0 : GByte abyBandInfo[1152] = {'\0'};
832 :
833 0 : for (int iBand = 0; iBand < nBands; iBand++)
834 : {
835 0 : if (VSIFReadL(abyBandInfo, 1152, 1, fpSTA) != 1)
836 0 : break;
837 :
838 0 : const int nBandNumber = abyBandInfo[7];
839 0 : GDALRasterBand *poBand = GetRasterBand(nBandNumber);
840 0 : if (poBand == nullptr)
841 0 : break;
842 :
843 0 : GInt16 nMin = 0;
844 0 : GInt16 nMax = 0;
845 :
846 0 : if (poBand->GetRasterDataType() != GDT_Byte)
847 : {
848 0 : memcpy(&nMin, abyBandInfo + 28, 2);
849 0 : memcpy(&nMax, abyBandInfo + 30, 2);
850 0 : CPL_LSBPTR16(&nMin);
851 0 : CPL_LSBPTR16(&nMax);
852 : }
853 : else
854 : {
855 0 : nMin = abyBandInfo[9];
856 0 : nMax = abyBandInfo[8];
857 : }
858 :
859 0 : float fMean = 0.0;
860 0 : float fStdDev = 0.0;
861 0 : memcpy(&fMean, abyBandInfo + 12, 4);
862 0 : memcpy(&fStdDev, abyBandInfo + 24, 4);
863 0 : CPL_LSBPTR32(&fMean);
864 0 : CPL_LSBPTR32(&fStdDev);
865 :
866 0 : poBand->SetStatistics(nMin, nMax, fMean, fStdDev);
867 : }
868 :
869 0 : CPL_IGNORE_RET_VAL(VSIFCloseL(fpSTA));
870 : }
871 :
872 : /************************************************************************/
873 : /* Create() */
874 : /************************************************************************/
875 :
876 60 : GDALDataset *LANDataset::Create(const char *pszFilename, int nXSize, int nYSize,
877 : int nBandsIn, GDALDataType eType,
878 : char ** /* papszOptions */)
879 : {
880 60 : if (eType != GDT_Byte && eType != GDT_Int16)
881 : {
882 33 : CPLError(CE_Failure, CPLE_AppDefined,
883 : "Attempt to create .GIS file with unsupported data type '%s'.",
884 : GDALGetDataTypeName(eType));
885 33 : return nullptr;
886 : }
887 :
888 : /* -------------------------------------------------------------------- */
889 : /* Try to create the file. */
890 : /* -------------------------------------------------------------------- */
891 27 : VSILFILE *fp = VSIFOpenL(pszFilename, "wb");
892 27 : if (fp == nullptr)
893 : {
894 3 : CPLError(CE_Failure, CPLE_OpenFailed,
895 : "Attempt to create file `%s' failed.\n", pszFilename);
896 3 : return nullptr;
897 : }
898 :
899 : /* -------------------------------------------------------------------- */
900 : /* Write out the header. */
901 : /* -------------------------------------------------------------------- */
902 24 : unsigned char abyHeader[128] = {'\0'};
903 :
904 24 : memset(abyHeader, 0, sizeof(abyHeader));
905 :
906 24 : memcpy(abyHeader + 0, "HEAD74", 6);
907 :
908 : // Pixel type.
909 24 : GInt16 n16Val = 0;
910 24 : if (eType == GDT_Byte) // Do we want 4bit?
911 21 : n16Val = 0;
912 : else
913 3 : n16Val = 2;
914 24 : memcpy(abyHeader + 6, &n16Val, 2);
915 :
916 : // Number of Bands.
917 24 : n16Val = static_cast<GInt16>(nBandsIn);
918 24 : memcpy(abyHeader + 8, &n16Val, 2);
919 :
920 : // Unknown (6).
921 :
922 : // Width.
923 24 : GInt32 n32Val = nXSize;
924 24 : memcpy(abyHeader + 16, &n32Val, 4);
925 :
926 : // Height.
927 24 : n32Val = nYSize;
928 24 : memcpy(abyHeader + 20, &n32Val, 4);
929 :
930 : // X Start (4).
931 : // Y Start (4).
932 :
933 : // Unknown (56).
934 :
935 : // Coordinate System.
936 24 : n16Val = 0;
937 24 : memcpy(abyHeader + 88, &n16Val, 2);
938 :
939 : // Classes in coverage.
940 24 : n16Val = 0;
941 24 : memcpy(abyHeader + 90, &n16Val, 2);
942 :
943 : // Unknown (14).
944 :
945 : // Area Unit.
946 24 : n16Val = 0;
947 24 : memcpy(abyHeader + 106, &n16Val, 2);
948 :
949 : // Pixel Area.
950 24 : float f32Val = 0.0f;
951 24 : memcpy(abyHeader + 108, &f32Val, 4);
952 :
953 : // Upper Left X.
954 24 : f32Val = 0.5f;
955 24 : memcpy(abyHeader + 112, &f32Val, 4);
956 :
957 : // Upper Left Y
958 24 : f32Val = static_cast<float>(nYSize - 0.5);
959 24 : memcpy(abyHeader + 116, &f32Val, 4);
960 :
961 : // Width of pixel.
962 24 : f32Val = 1.0f;
963 24 : memcpy(abyHeader + 120, &f32Val, 4);
964 :
965 : // Height of pixel.
966 24 : f32Val = 1.0f;
967 24 : memcpy(abyHeader + 124, &f32Val, 4);
968 :
969 24 : CPL_IGNORE_RET_VAL(VSIFWriteL(abyHeader, sizeof(abyHeader), 1, fp));
970 :
971 : /* -------------------------------------------------------------------- */
972 : /* Extend the file to the target size. */
973 : /* -------------------------------------------------------------------- */
974 24 : vsi_l_offset nImageBytes = 0;
975 :
976 24 : if (eType != GDT_Byte)
977 3 : nImageBytes = nXSize * static_cast<vsi_l_offset>(nYSize) * 2;
978 : else
979 21 : nImageBytes = nXSize * static_cast<vsi_l_offset>(nYSize);
980 :
981 24 : memset(abyHeader, 0, sizeof(abyHeader));
982 :
983 820 : while (nImageBytes > 0)
984 : {
985 : const vsi_l_offset nWriteThisTime =
986 805 : std::min(static_cast<size_t>(nImageBytes), sizeof(abyHeader));
987 :
988 805 : if (VSIFWriteL(abyHeader, 1, static_cast<size_t>(nWriteThisTime), fp) !=
989 : nWriteThisTime)
990 : {
991 9 : CPL_IGNORE_RET_VAL(VSIFCloseL(fp));
992 9 : CPLError(CE_Failure, CPLE_FileIO,
993 : "Failed to write whole Istar file.");
994 9 : return nullptr;
995 : }
996 796 : nImageBytes -= nWriteThisTime;
997 : }
998 :
999 15 : if (VSIFCloseL(fp) != 0)
1000 : {
1001 0 : CPLError(CE_Failure, CPLE_FileIO, "Failed to write whole Istar file.");
1002 0 : return nullptr;
1003 : }
1004 :
1005 15 : return static_cast<GDALDataset *>(GDALOpen(pszFilename, GA_Update));
1006 : }
1007 :
1008 : /************************************************************************/
1009 : /* GDALRegister_LAN() */
1010 : /************************************************************************/
1011 :
1012 1520 : void GDALRegister_LAN()
1013 :
1014 : {
1015 1520 : if (GDALGetDriverByName("LAN") != nullptr)
1016 301 : return;
1017 :
1018 1219 : GDALDriver *poDriver = new GDALDriver();
1019 :
1020 1219 : poDriver->SetDescription("LAN");
1021 1219 : poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
1022 1219 : poDriver->SetMetadataItem(GDAL_DMD_LONGNAME, "Erdas .LAN/.GIS");
1023 1219 : poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/lan.html");
1024 1219 : poDriver->SetMetadataItem(GDAL_DCAP_VIRTUALIO, "YES");
1025 1219 : poDriver->SetMetadataItem(GDAL_DMD_CREATIONDATATYPES, "Byte Int16");
1026 :
1027 1219 : poDriver->pfnOpen = LANDataset::Open;
1028 1219 : poDriver->pfnCreate = LANDataset::Create;
1029 :
1030 1219 : GetGDALDriverManager()->RegisterDriver(poDriver);
1031 : }
|
