/src/gdal/frmts/pds/isis2dataset.cpp
Line | Count | Source (jump to first uncovered line) |
1 | | /****************************************************************************** |
2 | | * |
3 | | * Project: ISIS Version 2 Driver |
4 | | * Purpose: Implementation of ISIS2Dataset |
5 | | * Author: Trent Hare (thare at usgs.gov), |
6 | | * Robert Soricone (rsoricone at usgs.gov) |
7 | | * Ludovic Mercier (ludovic.mercier at gmail.com) |
8 | | * Frank Warmerdam (warmerdam at pobox.com) |
9 | | * |
10 | | * NOTE: Original code authored by Trent and Robert and placed in the public |
11 | | * domain as per US government policy. I have (within my rights) appropriated |
12 | | * it and placed it under the following license. This is not intended to |
13 | | * diminish Trent and Roberts contribution. |
14 | | ****************************************************************************** |
15 | | * Copyright (c) 2006, Frank Warmerdam <warmerdam at pobox.com> |
16 | | * Copyright (c) 2008-2011, Even Rouault <even dot rouault at spatialys.com> |
17 | | * |
18 | | * SPDX-License-Identifier: MIT |
19 | | ****************************************************************************/ |
20 | | |
21 | | constexpr int NULL1 = 0; |
22 | | constexpr int NULL2 = -32768; |
23 | | constexpr double NULL3 = -3.4028226550889044521e+38; |
24 | | |
25 | | #include "cpl_string.h" |
26 | | #include "gdal_frmts.h" |
27 | | #include "nasakeywordhandler.h" |
28 | | #include "ogr_spatialref.h" |
29 | | #include "rawdataset.h" |
30 | | #include "pdsdrivercore.h" |
31 | | |
32 | | /************************************************************************/ |
33 | | /* ==================================================================== */ |
34 | | /* ISISDataset version2 */ |
35 | | /* ==================================================================== */ |
36 | | /************************************************************************/ |
37 | | |
38 | | class ISIS2Dataset final : public RawDataset |
39 | | { |
40 | | VSILFILE *fpImage{}; // image data file. |
41 | | CPLString osExternalCube{}; |
42 | | |
43 | | NASAKeywordHandler oKeywords{}; |
44 | | |
45 | | bool bGotTransform{}; |
46 | | GDALGeoTransform m_gt{}; |
47 | | |
48 | | OGRSpatialReference m_oSRS{}; |
49 | | |
50 | | int parse_label(const char *file, char *keyword, char *value); |
51 | | int strstrip(char instr[], char outstr[], int position); |
52 | | |
53 | | CPLString oTempResult{}; |
54 | | |
55 | | static void CleanString(CPLString &osInput); |
56 | | |
57 | | const char *GetKeyword(const char *pszPath, const char *pszDefault = ""); |
58 | | const char *GetKeywordSub(const char *pszPath, int iSubscript, |
59 | | const char *pszDefault = ""); |
60 | | |
61 | | CPLErr Close() override; |
62 | | |
63 | | CPL_DISALLOW_COPY_ASSIGN(ISIS2Dataset) |
64 | | |
65 | | public: |
66 | | ISIS2Dataset(); |
67 | | virtual ~ISIS2Dataset(); |
68 | | |
69 | | virtual CPLErr GetGeoTransform(GDALGeoTransform >) const override; |
70 | | const OGRSpatialReference *GetSpatialRef() const override; |
71 | | |
72 | | virtual char **GetFileList() override; |
73 | | |
74 | | static GDALDataset *Open(GDALOpenInfo *); |
75 | | }; |
76 | | |
77 | | /************************************************************************/ |
78 | | /* ISIS2Dataset() */ |
79 | | /************************************************************************/ |
80 | | |
81 | | ISIS2Dataset::ISIS2Dataset() |
82 | 3.69k | { |
83 | 3.69k | m_oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER); |
84 | 3.69k | } |
85 | | |
86 | | /************************************************************************/ |
87 | | /* ~ISIS2Dataset() */ |
88 | | /************************************************************************/ |
89 | | |
90 | | ISIS2Dataset::~ISIS2Dataset() |
91 | | |
92 | 3.69k | { |
93 | 3.69k | ISIS2Dataset::Close(); |
94 | 3.69k | } |
95 | | |
96 | | /************************************************************************/ |
97 | | /* Close() */ |
98 | | /************************************************************************/ |
99 | | |
100 | | CPLErr ISIS2Dataset::Close() |
101 | 3.69k | { |
102 | 3.69k | CPLErr eErr = CE_None; |
103 | 3.69k | if (nOpenFlags != OPEN_FLAGS_CLOSED) |
104 | 3.69k | { |
105 | 3.69k | if (ISIS2Dataset::FlushCache(true) != CE_None) |
106 | 0 | eErr = CE_Failure; |
107 | 3.69k | if (fpImage != nullptr) |
108 | 4 | { |
109 | 4 | if (VSIFCloseL(fpImage) != 0) |
110 | 0 | eErr = CE_Failure; |
111 | 4 | } |
112 | 3.69k | if (GDALPamDataset::Close() != CE_None) |
113 | 0 | eErr = CE_Failure; |
114 | 3.69k | } |
115 | 3.69k | return eErr; |
116 | 3.69k | } |
117 | | |
118 | | /************************************************************************/ |
119 | | /* GetFileList() */ |
120 | | /************************************************************************/ |
121 | | |
122 | | char **ISIS2Dataset::GetFileList() |
123 | | |
124 | 8 | { |
125 | 8 | char **papszFileList = GDALPamDataset::GetFileList(); |
126 | | |
127 | 8 | if (!osExternalCube.empty()) |
128 | 0 | papszFileList = CSLAddString(papszFileList, osExternalCube); |
129 | | |
130 | 8 | return papszFileList; |
131 | 8 | } |
132 | | |
133 | | /************************************************************************/ |
134 | | /* GetSpatialRef() */ |
135 | | /************************************************************************/ |
136 | | |
137 | | const OGRSpatialReference *ISIS2Dataset::GetSpatialRef() const |
138 | 4 | { |
139 | 4 | if (!m_oSRS.IsEmpty()) |
140 | 2 | return &m_oSRS; |
141 | 2 | return GDALPamDataset::GetSpatialRef(); |
142 | 4 | } |
143 | | |
144 | | /************************************************************************/ |
145 | | /* GetGeoTransform() */ |
146 | | /************************************************************************/ |
147 | | |
148 | | CPLErr ISIS2Dataset::GetGeoTransform(GDALGeoTransform >) const |
149 | | |
150 | 4 | { |
151 | 4 | if (bGotTransform) |
152 | 3 | { |
153 | 3 | gt = m_gt; |
154 | 3 | return CE_None; |
155 | 3 | } |
156 | | |
157 | 1 | return GDALPamDataset::GetGeoTransform(gt); |
158 | 4 | } |
159 | | |
160 | | /************************************************************************/ |
161 | | /* Open() */ |
162 | | /************************************************************************/ |
163 | | |
164 | | GDALDataset *ISIS2Dataset::Open(GDALOpenInfo *poOpenInfo) |
165 | 3.69k | { |
166 | | /* -------------------------------------------------------------------- */ |
167 | | /* Does this look like a CUBE or an IMAGE Primary Data Object? */ |
168 | | /* -------------------------------------------------------------------- */ |
169 | 3.69k | if (!ISIS2DriverIdentify(poOpenInfo) || poOpenInfo->fpL == nullptr) |
170 | 0 | return nullptr; |
171 | | |
172 | 3.69k | VSILFILE *fpQube = poOpenInfo->fpL; |
173 | 3.69k | poOpenInfo->fpL = nullptr; |
174 | | |
175 | 3.69k | auto poDS = std::make_unique<ISIS2Dataset>(); |
176 | | |
177 | 3.69k | if (!poDS->oKeywords.Ingest(fpQube, 0)) |
178 | 3.27k | { |
179 | 3.27k | VSIFCloseL(fpQube); |
180 | 3.27k | return nullptr; |
181 | 3.27k | } |
182 | | |
183 | 421 | VSIFCloseL(fpQube); |
184 | | |
185 | | /* -------------------------------------------------------------------- */ |
186 | | /* We assume the user is pointing to the label (i.e. .lab) file. */ |
187 | | /* -------------------------------------------------------------------- */ |
188 | | // QUBE can be inline or detached and point to an image name |
189 | | // ^QUBE = 76 |
190 | | // ^QUBE = ("ui31s015.img",6441<BYTES>) - has another label on the image |
191 | | // ^QUBE = "ui31s015.img" - which implies no label or skip value |
192 | | |
193 | 421 | const char *pszQube = poDS->GetKeyword("^QUBE"); |
194 | 421 | int nQube = 0; |
195 | 421 | int bByteLocation = FALSE; |
196 | 421 | CPLString osTargetFile = poOpenInfo->pszFilename; |
197 | | |
198 | 421 | if (pszQube[0] == '"') |
199 | 31 | { |
200 | 31 | const CPLString osTPath = CPLGetPathSafe(poOpenInfo->pszFilename); |
201 | 31 | CPLString osFilename = pszQube; |
202 | 31 | poDS->CleanString(osFilename); |
203 | 31 | osTargetFile = CPLFormCIFilenameSafe(osTPath, osFilename, nullptr); |
204 | 31 | poDS->osExternalCube = osTargetFile; |
205 | 31 | } |
206 | 390 | else if (pszQube[0] == '(') |
207 | 0 | { |
208 | 0 | const CPLString osTPath = CPLGetPathSafe(poOpenInfo->pszFilename); |
209 | 0 | CPLString osFilename = poDS->GetKeywordSub("^QUBE", 1, ""); |
210 | 0 | poDS->CleanString(osFilename); |
211 | 0 | osTargetFile = CPLFormCIFilenameSafe(osTPath, osFilename, nullptr); |
212 | 0 | poDS->osExternalCube = osTargetFile; |
213 | |
|
214 | 0 | nQube = atoi(poDS->GetKeywordSub("^QUBE", 2, "1")); |
215 | 0 | if (strstr(poDS->GetKeywordSub("^QUBE", 2, "1"), "<BYTES>") != nullptr) |
216 | 0 | bByteLocation = true; |
217 | 0 | } |
218 | 390 | else |
219 | 390 | { |
220 | 390 | nQube = atoi(pszQube); |
221 | 390 | if (strstr(pszQube, "<BYTES>") != nullptr) |
222 | 0 | bByteLocation = true; |
223 | 390 | } |
224 | | |
225 | | /* -------------------------------------------------------------------- */ |
226 | | /* Check if file an ISIS2 header file? Read a few lines of text */ |
227 | | /* searching for something starting with nrows or ncols. */ |
228 | | /* -------------------------------------------------------------------- */ |
229 | | |
230 | | /* -------------------------------------------------------------------- */ |
231 | | /* Checks to see if this is valid ISIS2 cube */ |
232 | | /* SUFFIX_ITEM tag in .cub file should be (0,0,0); no side-planes */ |
233 | | /* -------------------------------------------------------------------- */ |
234 | 421 | const int s_ix = atoi(poDS->GetKeywordSub("QUBE.SUFFIX_ITEMS", 1)); |
235 | 421 | const int s_iy = atoi(poDS->GetKeywordSub("QUBE.SUFFIX_ITEMS", 2)); |
236 | 421 | const int s_iz = atoi(poDS->GetKeywordSub("QUBE.SUFFIX_ITEMS", 3)); |
237 | | |
238 | 421 | if (s_ix != 0 || s_iy != 0 || s_iz != 0) |
239 | 0 | { |
240 | 0 | CPLError(CE_Failure, CPLE_OpenFailed, |
241 | 0 | "*** ISIS 2 cube file has invalid SUFFIX_ITEMS parameters:\n" |
242 | 0 | "*** gdal isis2 driver requires (0, 0, 0), thus no sideplanes " |
243 | 0 | "or backplanes\n" |
244 | 0 | "found: (%i, %i, %i)\n\n", |
245 | 0 | s_ix, s_iy, s_iz); |
246 | 0 | return nullptr; |
247 | 0 | } |
248 | | |
249 | | /**************** end SUFFIX_ITEM check ***********************/ |
250 | | |
251 | | /*********** Grab layout type (BSQ, BIP, BIL) ************/ |
252 | | // AXIS_NAME = (SAMPLE,LINE,BAND) |
253 | | /***********************************************************/ |
254 | | |
255 | 421 | char szLayout[10] = "BSQ"; // default to band seq. |
256 | 421 | const char *value = poDS->GetKeyword("QUBE.AXIS_NAME", ""); |
257 | 421 | if (EQUAL(value, "(SAMPLE,LINE,BAND)")) |
258 | 5 | strcpy(szLayout, "BSQ"); |
259 | 416 | else if (EQUAL(value, "(BAND,LINE,SAMPLE)")) |
260 | 0 | strcpy(szLayout, "BIP"); |
261 | 416 | else if (EQUAL(value, "(SAMPLE,BAND,LINE)") || EQUAL(value, "")) |
262 | 416 | strcpy(szLayout, "BSQ"); |
263 | 0 | else |
264 | 0 | { |
265 | 0 | CPLError(CE_Failure, CPLE_OpenFailed, |
266 | 0 | "%s layout not supported. Abort\n\n", value); |
267 | 0 | return nullptr; |
268 | 0 | } |
269 | | |
270 | | /*********** Grab samples lines band ************/ |
271 | 421 | const int nCols = atoi(poDS->GetKeywordSub("QUBE.CORE_ITEMS", 1)); |
272 | 421 | const int nRows = atoi(poDS->GetKeywordSub("QUBE.CORE_ITEMS", 2)); |
273 | 421 | const int nBands = atoi(poDS->GetKeywordSub("QUBE.CORE_ITEMS", 3)); |
274 | | |
275 | | /*********** Grab Qube record bytes **********/ |
276 | 421 | const int record_bytes = atoi(poDS->GetKeyword("RECORD_BYTES")); |
277 | 421 | if (record_bytes < 0) |
278 | 0 | { |
279 | 0 | return nullptr; |
280 | 0 | } |
281 | | |
282 | 421 | GUIntBig nSkipBytes = 0; |
283 | 421 | if (nQube > 0 && bByteLocation) |
284 | 0 | nSkipBytes = (nQube - 1); |
285 | 421 | else if (nQube > 0) |
286 | 4 | nSkipBytes = static_cast<GUIntBig>(nQube - 1) * record_bytes; |
287 | 417 | else |
288 | 417 | nSkipBytes = 0; |
289 | | |
290 | | /*********** Grab samples lines band ************/ |
291 | 421 | char chByteOrder = 'M'; // default to MSB |
292 | 421 | CPLString osCoreItemType = poDS->GetKeyword("QUBE.CORE_ITEM_TYPE"); |
293 | 421 | if ((EQUAL(osCoreItemType, "PC_INTEGER")) || |
294 | 421 | (EQUAL(osCoreItemType, "PC_UNSIGNED_INTEGER")) || |
295 | 421 | (EQUAL(osCoreItemType, "PC_REAL"))) |
296 | 0 | { |
297 | 0 | chByteOrder = 'I'; |
298 | 0 | } |
299 | | |
300 | | /******** Grab format type - isis2 only supports 8,16,32 *******/ |
301 | 421 | GDALDataType eDataType = GDT_Byte; |
302 | 421 | bool bNoDataSet = false; |
303 | 421 | double dfNoData = 0.0; |
304 | | |
305 | 421 | int itype = atoi(poDS->GetKeyword("QUBE.CORE_ITEM_BYTES", "")); |
306 | 421 | switch (itype) |
307 | 421 | { |
308 | 0 | case 1: |
309 | 0 | eDataType = GDT_Byte; |
310 | 0 | dfNoData = NULL1; |
311 | 0 | bNoDataSet = true; |
312 | 0 | break; |
313 | 0 | case 2: |
314 | 0 | if (strstr(osCoreItemType, "UNSIGNED") != nullptr) |
315 | 0 | { |
316 | 0 | dfNoData = 0; |
317 | 0 | eDataType = GDT_UInt16; |
318 | 0 | } |
319 | 0 | else |
320 | 0 | { |
321 | 0 | dfNoData = NULL2; |
322 | 0 | eDataType = GDT_Int16; |
323 | 0 | } |
324 | 0 | bNoDataSet = true; |
325 | 0 | break; |
326 | 5 | case 4: |
327 | 5 | eDataType = GDT_Float32; |
328 | 5 | dfNoData = NULL3; |
329 | 5 | bNoDataSet = true; |
330 | 5 | break; |
331 | 0 | case 8: |
332 | 0 | eDataType = GDT_Float64; |
333 | 0 | dfNoData = NULL3; |
334 | 0 | bNoDataSet = true; |
335 | 0 | break; |
336 | 416 | default: |
337 | 416 | CPLError(CE_Failure, CPLE_AppDefined, |
338 | 416 | "Itype of %d is not supported in ISIS 2.", itype); |
339 | 416 | return nullptr; |
340 | 421 | } |
341 | | |
342 | | /*********** Grab Cellsize ************/ |
343 | 5 | double dfXDim = 1.0; |
344 | 5 | double dfYDim = 1.0; |
345 | | |
346 | 5 | value = poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.MAP_SCALE"); |
347 | 5 | if (strlen(value) > 0) |
348 | 3 | { |
349 | | // Convert km to m |
350 | 3 | dfXDim = static_cast<float>(CPLAtof(value) * 1000.0); |
351 | 3 | dfYDim = static_cast<float>(CPLAtof(value) * 1000.0 * -1); |
352 | 3 | } |
353 | | |
354 | | /*********** Grab LINE_PROJECTION_OFFSET ************/ |
355 | 5 | double dfULYMap = 0.5; |
356 | | |
357 | 5 | value = |
358 | 5 | poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.LINE_PROJECTION_OFFSET"); |
359 | 5 | if (strlen(value) > 0) |
360 | 4 | { |
361 | 4 | const double yulcenter = static_cast<float>(CPLAtof(value)) * dfYDim; |
362 | 4 | dfULYMap = yulcenter - (dfYDim / 2); |
363 | 4 | } |
364 | | |
365 | | /*********** Grab SAMPLE_PROJECTION_OFFSET ************/ |
366 | 5 | double dfULXMap = 0.5; |
367 | | |
368 | 5 | value = |
369 | 5 | poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.SAMPLE_PROJECTION_OFFSET"); |
370 | 5 | if (strlen(value) > 0) |
371 | 3 | { |
372 | 3 | const double xulcenter = static_cast<float>(CPLAtof(value)) * dfXDim; |
373 | 3 | dfULXMap = xulcenter - (dfXDim / 2); |
374 | 3 | } |
375 | | |
376 | | /*********** Grab TARGET_NAME ************/ |
377 | | /**** This is the planets name i.e. MARS ***/ |
378 | 5 | CPLString target_name = poDS->GetKeyword("QUBE.TARGET_NAME"); |
379 | | |
380 | | /*********** Grab MAP_PROJECTION_TYPE ************/ |
381 | 5 | CPLString map_proj_name = |
382 | 5 | poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.MAP_PROJECTION_TYPE"); |
383 | 5 | poDS->CleanString(map_proj_name); |
384 | | |
385 | | /*********** Grab SEMI-MAJOR ************/ |
386 | 5 | const double semi_major = |
387 | 5 | CPLAtof(poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.A_AXIS_RADIUS")) * |
388 | 5 | 1000.0; |
389 | | |
390 | | /*********** Grab semi-minor ************/ |
391 | 5 | const double semi_minor = |
392 | 5 | CPLAtof(poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.C_AXIS_RADIUS")) * |
393 | 5 | 1000.0; |
394 | | |
395 | | /*********** Grab CENTER_LAT ************/ |
396 | 5 | const double center_lat = |
397 | 5 | CPLAtof(poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.CENTER_LATITUDE")); |
398 | | |
399 | | /*********** Grab CENTER_LON ************/ |
400 | 5 | const double center_lon = |
401 | 5 | CPLAtof(poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.CENTER_LONGITUDE")); |
402 | | |
403 | | /*********** Grab 1st std parallel ************/ |
404 | 5 | const double first_std_parallel = CPLAtof( |
405 | 5 | poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.FIRST_STANDARD_PARALLEL")); |
406 | | |
407 | | /*********** Grab 2nd std parallel ************/ |
408 | 5 | const double second_std_parallel = CPLAtof( |
409 | 5 | poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.SECOND_STANDARD_PARALLEL")); |
410 | | |
411 | | /*** grab PROJECTION_LATITUDE_TYPE = "PLANETOCENTRIC" ****/ |
412 | | // Need to further study how ocentric/ographic will effect the gdal library. |
413 | | // So far we will use this fact to define a sphere or ellipse for some |
414 | | // projections Frank - may need to talk this over |
415 | 5 | bool bIsGeographic = true; |
416 | 5 | value = |
417 | 5 | poDS->GetKeyword("CUBE.IMAGE_MAP_PROJECTION.PROJECTION_LATITUDE_TYPE"); |
418 | 5 | if (EQUAL(value, "\"PLANETOCENTRIC\"")) |
419 | 0 | bIsGeographic = false; |
420 | | |
421 | 5 | CPLDebug("ISIS2", "using projection %s", map_proj_name.c_str()); |
422 | | |
423 | 5 | OGRSpatialReference oSRS; |
424 | 5 | oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER); |
425 | 5 | bool bProjectionSet = true; |
426 | | |
427 | | // Set oSRS projection and parameters |
428 | 5 | if ((EQUAL(map_proj_name, "EQUIRECTANGULAR_CYLINDRICAL")) || |
429 | 5 | (EQUAL(map_proj_name, "EQUIRECTANGULAR")) || |
430 | 5 | (EQUAL(map_proj_name, "SIMPLE_CYLINDRICAL"))) |
431 | 3 | { |
432 | 3 | oSRS.OGRSpatialReference::SetEquirectangular2(0.0, center_lon, |
433 | 3 | center_lat, 0, 0); |
434 | 3 | } |
435 | 2 | else if (EQUAL(map_proj_name, "ORTHOGRAPHIC")) |
436 | 0 | { |
437 | 0 | oSRS.OGRSpatialReference::SetOrthographic(center_lat, center_lon, 0, 0); |
438 | 0 | } |
439 | 2 | else if ((EQUAL(map_proj_name, "SINUSOIDAL")) || |
440 | 2 | (EQUAL(map_proj_name, "SINUSOIDAL_EQUAL-AREA"))) |
441 | 0 | { |
442 | 0 | oSRS.OGRSpatialReference::SetSinusoidal(center_lon, 0, 0); |
443 | 0 | } |
444 | 2 | else if (EQUAL(map_proj_name, "MERCATOR")) |
445 | 0 | { |
446 | 0 | oSRS.OGRSpatialReference::SetMercator(center_lat, center_lon, 1, 0, 0); |
447 | 0 | } |
448 | 2 | else if (EQUAL(map_proj_name, "POLAR_STEREOGRAPHIC")) |
449 | 0 | { |
450 | 0 | oSRS.OGRSpatialReference::SetPS(center_lat, center_lon, 1, 0, 0); |
451 | 0 | } |
452 | 2 | else if (EQUAL(map_proj_name, "TRANSVERSE_MERCATOR")) |
453 | 0 | { |
454 | 0 | oSRS.OGRSpatialReference::SetTM(center_lat, center_lon, 1, 0, 0); |
455 | 0 | } |
456 | 2 | else if (EQUAL(map_proj_name, "LAMBERT_CONFORMAL_CONIC")) |
457 | 0 | { |
458 | 0 | oSRS.OGRSpatialReference::SetLCC(first_std_parallel, |
459 | 0 | second_std_parallel, center_lat, |
460 | 0 | center_lon, 0, 0); |
461 | 0 | } |
462 | 2 | else if (EQUAL(map_proj_name, "")) |
463 | 1 | { |
464 | | /* no projection */ |
465 | 1 | bProjectionSet = false; |
466 | 1 | } |
467 | 1 | else |
468 | 1 | { |
469 | 1 | CPLDebug("ISIS2", |
470 | 1 | "Dataset projection %s is not supported. Continuing...", |
471 | 1 | map_proj_name.c_str()); |
472 | 1 | bProjectionSet = false; |
473 | 1 | } |
474 | | |
475 | 5 | if (bProjectionSet) |
476 | 3 | { |
477 | | // Create projection name, i.e. MERCATOR MARS and set as ProjCS keyword |
478 | 3 | const CPLString proj_target_name = map_proj_name + " " + target_name; |
479 | 3 | oSRS.SetProjCS(proj_target_name); // set ProjCS keyword |
480 | | |
481 | | // The geographic/geocentric name will be the same basic name as the |
482 | | // body name 'GCS' = Geographic/Geocentric Coordinate System |
483 | 3 | const CPLString geog_name = "GCS_" + target_name; |
484 | | |
485 | | // The datum and sphere names will be the same basic name aas the planet |
486 | 3 | const CPLString datum_name = "D_" + target_name; |
487 | | |
488 | 3 | CPLString sphere_name = std::move(target_name); |
489 | | |
490 | | // calculate inverse flattening from major and minor axis: 1/f = a/(a-b) |
491 | 3 | double iflattening = 0.0; |
492 | 3 | if ((semi_major - semi_minor) < 0.0000001) |
493 | 3 | iflattening = 0; |
494 | 0 | else |
495 | 0 | iflattening = semi_major / (semi_major - semi_minor); |
496 | | |
497 | | // Set the body size but take into consideration which proj is being |
498 | | // used to help w/ proj4 compatibility The use of a Sphere, polar radius |
499 | | // or ellipse here is based on how ISIS does it internally |
500 | 3 | if (((EQUAL(map_proj_name, "STEREOGRAPHIC") && |
501 | 3 | (fabs(center_lat) == 90))) || |
502 | 3 | (EQUAL(map_proj_name, "POLAR_STEREOGRAPHIC"))) |
503 | 0 | { |
504 | 0 | if (bIsGeographic) |
505 | 0 | { |
506 | | // Geograpraphic, so set an ellipse |
507 | 0 | oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_major, |
508 | 0 | iflattening, "Reference_Meridian", 0.0); |
509 | 0 | } |
510 | 0 | else |
511 | 0 | { |
512 | | // Geocentric, so force a sphere using the semi-minor axis. I |
513 | | // hope... |
514 | 0 | sphere_name += "_polarRadius"; |
515 | 0 | oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_minor, |
516 | 0 | 0.0, "Reference_Meridian", 0.0); |
517 | 0 | } |
518 | 0 | } |
519 | 3 | else if ((EQUAL(map_proj_name, "SIMPLE_CYLINDRICAL")) || |
520 | 3 | (EQUAL(map_proj_name, "ORTHOGRAPHIC")) || |
521 | 3 | (EQUAL(map_proj_name, "STEREOGRAPHIC")) || |
522 | 3 | (EQUAL(map_proj_name, "SINUSOIDAL_EQUAL-AREA")) || |
523 | 3 | (EQUAL(map_proj_name, "SINUSOIDAL"))) |
524 | 3 | { |
525 | | // ISIS uses the spherical equation for these projections so force |
526 | | // a sphere. |
527 | 3 | oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_major, 0.0, |
528 | 3 | "Reference_Meridian", 0.0); |
529 | 3 | } |
530 | 0 | else if ((EQUAL(map_proj_name, "EQUIRECTANGULAR_CYLINDRICAL")) || |
531 | 0 | (EQUAL(map_proj_name, "EQUIRECTANGULAR"))) |
532 | 0 | { |
533 | | // Calculate localRadius using ISIS3 simple elliptical method |
534 | | // not the more standard Radius of Curvature method |
535 | | // PI = 4 * atan(1); |
536 | 0 | if (center_lon == 0) |
537 | 0 | { // No need to calculate local radius |
538 | 0 | oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_major, |
539 | 0 | 0.0, "Reference_Meridian", 0.0); |
540 | 0 | } |
541 | 0 | else |
542 | 0 | { |
543 | 0 | const double radLat = center_lat * M_PI / 180; // in radians |
544 | 0 | const double localRadius = |
545 | 0 | semi_major * semi_minor / |
546 | 0 | sqrt(pow(semi_minor * cos(radLat), 2) + |
547 | 0 | pow(semi_major * sin(radLat), 2)); |
548 | 0 | sphere_name += "_localRadius"; |
549 | 0 | oSRS.SetGeogCS(geog_name, datum_name, sphere_name, localRadius, |
550 | 0 | 0.0, "Reference_Meridian", 0.0); |
551 | 0 | CPLDebug("ISIS2", "local radius: %f", localRadius); |
552 | 0 | } |
553 | 0 | } |
554 | 0 | else |
555 | 0 | { |
556 | | // All other projections: Mercator, Transverse Mercator, Lambert |
557 | | // Conformal, etc. Geographic, so set an ellipse |
558 | 0 | if (bIsGeographic) |
559 | 0 | { |
560 | 0 | oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_major, |
561 | 0 | iflattening, "Reference_Meridian", 0.0); |
562 | 0 | } |
563 | 0 | else |
564 | 0 | { |
565 | | // Geocentric, so force a sphere. I hope... |
566 | 0 | oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_major, |
567 | 0 | 0.0, "Reference_Meridian", 0.0); |
568 | 0 | } |
569 | 0 | } |
570 | | |
571 | | // translate back into a projection string. |
572 | 3 | poDS->m_oSRS = std::move(oSRS); |
573 | 3 | } |
574 | | |
575 | | /* END ISIS2 Label Read */ |
576 | | /*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*/ |
577 | | |
578 | | /* -------------------------------------------------------------------- */ |
579 | | /* Did we get the required keywords? If not we return with */ |
580 | | /* this never having been considered to be a match. This isn't */ |
581 | | /* an error! */ |
582 | | /* -------------------------------------------------------------------- */ |
583 | 5 | if (!GDALCheckDatasetDimensions(nCols, nRows) || |
584 | 5 | !GDALCheckBandCount(nBands, false)) |
585 | 1 | { |
586 | 1 | return nullptr; |
587 | 1 | } |
588 | | |
589 | | /* -------------------------------------------------------------------- */ |
590 | | /* Capture some information from the file that is of interest. */ |
591 | | /* -------------------------------------------------------------------- */ |
592 | 4 | poDS->nRasterXSize = nCols; |
593 | 4 | poDS->nRasterYSize = nRows; |
594 | | |
595 | | /* -------------------------------------------------------------------- */ |
596 | | /* Open target binary file. */ |
597 | | /* -------------------------------------------------------------------- */ |
598 | | |
599 | 4 | if (poOpenInfo->eAccess == GA_ReadOnly) |
600 | 4 | poDS->fpImage = VSIFOpenL(osTargetFile, "rb"); |
601 | 0 | else |
602 | 0 | poDS->fpImage = VSIFOpenL(osTargetFile, "r+b"); |
603 | | |
604 | 4 | if (poDS->fpImage == nullptr) |
605 | 0 | { |
606 | 0 | CPLError(CE_Failure, CPLE_OpenFailed, |
607 | 0 | "Failed to open %s with write permission.\n%s", |
608 | 0 | osTargetFile.c_str(), VSIStrerror(errno)); |
609 | 0 | return nullptr; |
610 | 0 | } |
611 | | |
612 | 4 | poDS->eAccess = poOpenInfo->eAccess; |
613 | | |
614 | | /* -------------------------------------------------------------------- */ |
615 | | /* Compute the line offset. */ |
616 | | /* -------------------------------------------------------------------- */ |
617 | 4 | int nItemSize = GDALGetDataTypeSizeBytes(eDataType); |
618 | 4 | int nLineOffset, nPixelOffset; |
619 | 4 | vsi_l_offset nBandOffset; |
620 | | |
621 | 4 | if (EQUAL(szLayout, "BIP")) |
622 | 0 | { |
623 | 0 | nPixelOffset = nItemSize * nBands; |
624 | 0 | if (nPixelOffset > INT_MAX / nBands) |
625 | 0 | { |
626 | 0 | return nullptr; |
627 | 0 | } |
628 | 0 | nLineOffset = nPixelOffset * nCols; |
629 | 0 | nBandOffset = nItemSize; |
630 | 0 | } |
631 | 4 | else if (EQUAL(szLayout, "BSQ")) |
632 | 4 | { |
633 | 4 | nPixelOffset = nItemSize; |
634 | 4 | if (nPixelOffset > INT_MAX / nCols) |
635 | 0 | { |
636 | 0 | return nullptr; |
637 | 0 | } |
638 | 4 | nLineOffset = nPixelOffset * nCols; |
639 | 4 | nBandOffset = static_cast<vsi_l_offset>(nLineOffset) * nRows; |
640 | 4 | } |
641 | 0 | else /* assume BIL */ |
642 | 0 | { |
643 | 0 | nPixelOffset = nItemSize; |
644 | 0 | if (nPixelOffset > INT_MAX / nBands || |
645 | 0 | nPixelOffset * nBands > INT_MAX / nCols) |
646 | 0 | { |
647 | 0 | return nullptr; |
648 | 0 | } |
649 | 0 | nLineOffset = nItemSize * nBands * nCols; |
650 | 0 | nBandOffset = static_cast<vsi_l_offset>(nItemSize) * nCols; |
651 | 0 | } |
652 | | |
653 | | /* -------------------------------------------------------------------- */ |
654 | | /* Create band information objects. */ |
655 | | /* -------------------------------------------------------------------- */ |
656 | 8 | for (int i = 0; i < nBands; i++) |
657 | 4 | { |
658 | 4 | auto poBand = RawRasterBand::Create( |
659 | 4 | poDS.get(), i + 1, poDS->fpImage, nSkipBytes + nBandOffset * i, |
660 | 4 | nPixelOffset, nLineOffset, eDataType, |
661 | 4 | chByteOrder == 'I' || chByteOrder == 'L' |
662 | 4 | ? RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN |
663 | 4 | : RawRasterBand::ByteOrder::ORDER_BIG_ENDIAN, |
664 | 4 | RawRasterBand::OwnFP::NO); |
665 | 4 | if (!poBand) |
666 | 0 | return nullptr; |
667 | | |
668 | 4 | if (bNoDataSet) |
669 | 4 | poBand->SetNoDataValue(dfNoData); |
670 | | |
671 | | // Set offset/scale values at the PAM level. |
672 | 4 | poBand->SetOffset(CPLAtofM(poDS->GetKeyword("QUBE.CORE_BASE", "0.0"))); |
673 | 4 | poBand->SetScale( |
674 | 4 | CPLAtofM(poDS->GetKeyword("QUBE.CORE_MULTIPLIER", "1.0"))); |
675 | | |
676 | 4 | poDS->SetBand(i + 1, std::move(poBand)); |
677 | 4 | } |
678 | | |
679 | | /* -------------------------------------------------------------------- */ |
680 | | /* Check for a .prj file. For isis2 I would like to keep this in */ |
681 | | /* -------------------------------------------------------------------- */ |
682 | 4 | const CPLString osPath = CPLGetPathSafe(poOpenInfo->pszFilename); |
683 | 4 | const CPLString osName = CPLGetBasenameSafe(poOpenInfo->pszFilename); |
684 | 4 | const std::string osPrjFile = CPLFormCIFilenameSafe(osPath, osName, "prj"); |
685 | | |
686 | 4 | VSILFILE *fp = VSIFOpenL(osPrjFile.c_str(), "r"); |
687 | 4 | if (fp != nullptr) |
688 | 0 | { |
689 | 0 | VSIFCloseL(fp); |
690 | |
|
691 | 0 | char **papszLines = CSLLoad(osPrjFile.c_str()); |
692 | |
|
693 | 0 | poDS->m_oSRS.importFromESRI(papszLines); |
694 | |
|
695 | 0 | CSLDestroy(papszLines); |
696 | 0 | } |
697 | | |
698 | 4 | if (dfULXMap != 0.5 || dfULYMap != 0.5 || dfXDim != 1.0 || dfYDim != 1.0) |
699 | 3 | { |
700 | 3 | poDS->bGotTransform = TRUE; |
701 | 3 | poDS->m_gt[0] = dfULXMap; |
702 | 3 | poDS->m_gt[1] = dfXDim; |
703 | 3 | poDS->m_gt[2] = 0.0; |
704 | 3 | poDS->m_gt[3] = dfULYMap; |
705 | 3 | poDS->m_gt[4] = 0.0; |
706 | 3 | poDS->m_gt[5] = dfYDim; |
707 | 3 | } |
708 | | |
709 | 4 | if (!poDS->bGotTransform) |
710 | 1 | poDS->bGotTransform = GDALReadWorldFile(poOpenInfo->pszFilename, "cbw", |
711 | 1 | poDS->m_gt.data()); |
712 | | |
713 | 4 | if (!poDS->bGotTransform) |
714 | 1 | poDS->bGotTransform = GDALReadWorldFile(poOpenInfo->pszFilename, "wld", |
715 | 1 | poDS->m_gt.data()); |
716 | | |
717 | | /* -------------------------------------------------------------------- */ |
718 | | /* Initialize any PAM information. */ |
719 | | /* -------------------------------------------------------------------- */ |
720 | 4 | poDS->SetDescription(poOpenInfo->pszFilename); |
721 | 4 | poDS->TryLoadXML(); |
722 | | |
723 | | /* -------------------------------------------------------------------- */ |
724 | | /* Check for overviews. */ |
725 | | /* -------------------------------------------------------------------- */ |
726 | 4 | poDS->oOvManager.Initialize(poDS.get(), poOpenInfo->pszFilename); |
727 | | |
728 | 4 | return poDS.release(); |
729 | 4 | } |
730 | | |
731 | | /************************************************************************/ |
732 | | /* GetKeyword() */ |
733 | | /************************************************************************/ |
734 | | |
735 | | const char *ISIS2Dataset::GetKeyword(const char *pszPath, |
736 | | const char *pszDefault) |
737 | | |
738 | 2.17k | { |
739 | 2.17k | return oKeywords.GetKeyword(pszPath, pszDefault); |
740 | 2.17k | } |
741 | | |
742 | | /************************************************************************/ |
743 | | /* GetKeywordSub() */ |
744 | | /************************************************************************/ |
745 | | |
746 | | const char *ISIS2Dataset::GetKeywordSub(const char *pszPath, int iSubscript, |
747 | | const char *pszDefault) |
748 | | |
749 | 2.52k | { |
750 | 2.52k | const char *pszResult = oKeywords.GetKeyword(pszPath, nullptr); |
751 | | |
752 | 2.52k | if (pszResult == nullptr) |
753 | 2.50k | return pszDefault; |
754 | | |
755 | 24 | if (pszResult[0] != '(') |
756 | 0 | return pszDefault; |
757 | | |
758 | 24 | char **papszTokens = |
759 | 24 | CSLTokenizeString2(pszResult, "(,)", CSLT_HONOURSTRINGS); |
760 | | |
761 | 24 | if (iSubscript <= CSLCount(papszTokens)) |
762 | 24 | { |
763 | 24 | oTempResult = papszTokens[iSubscript - 1]; |
764 | 24 | CSLDestroy(papszTokens); |
765 | 24 | return oTempResult.c_str(); |
766 | 24 | } |
767 | | |
768 | 0 | CSLDestroy(papszTokens); |
769 | 0 | return pszDefault; |
770 | 24 | } |
771 | | |
772 | | /************************************************************************/ |
773 | | /* CleanString() */ |
774 | | /* */ |
775 | | /* Removes single or double quotes, and converts spaces to underscores. */ |
776 | | /* The change is made in-place to CPLString. */ |
777 | | /************************************************************************/ |
778 | | |
779 | | void ISIS2Dataset::CleanString(CPLString &osInput) |
780 | | |
781 | 36 | { |
782 | 36 | if ((osInput.size() < 2) || |
783 | 36 | ((osInput.at(0) != '"' || osInput.back() != '"') && |
784 | 35 | (osInput.at(0) != '\'' || osInput.back() != '\''))) |
785 | 5 | return; |
786 | | |
787 | 31 | char *pszWrk = CPLStrdup(osInput.c_str() + 1); |
788 | | |
789 | 31 | pszWrk[strlen(pszWrk) - 1] = '\0'; |
790 | | |
791 | 4.25M | for (int i = 0; pszWrk[i] != '\0'; i++) |
792 | 4.25M | { |
793 | 4.25M | if (pszWrk[i] == ' ') |
794 | 1.12k | pszWrk[i] = '_'; |
795 | 4.25M | } |
796 | | |
797 | 31 | osInput = pszWrk; |
798 | 31 | CPLFree(pszWrk); |
799 | 31 | } |
800 | | |
801 | | /************************************************************************/ |
802 | | /* GDALRegister_ISIS2() */ |
803 | | /************************************************************************/ |
804 | | |
805 | | void GDALRegister_ISIS2() |
806 | | |
807 | 22 | { |
808 | 22 | if (GDALGetDriverByName(ISIS2_DRIVER_NAME) != nullptr) |
809 | 0 | return; |
810 | | |
811 | 22 | GDALDriver *poDriver = new GDALDriver(); |
812 | 22 | ISIS2DriverSetCommonMetadata(poDriver); |
813 | | |
814 | 22 | poDriver->pfnOpen = ISIS2Dataset::Open; |
815 | | |
816 | 22 | GetGDALDriverManager()->RegisterDriver(poDriver); |
817 | 22 | } |