Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GDAL Core
4 : * Purpose: Read metadata from Pleiades imagery.
5 : * Author: Alexander Lisovenko
6 : * Author: Dmitry Baryshnikov, polimax@mail.ru
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 2014-2018 NextGIS <info@nextgis.ru>
10 : *
11 : * SPDX-License-Identifier: MIT
12 : ****************************************************************************/
13 :
14 : #include "cpl_port.h"
15 : #include "reader_pleiades.h"
16 :
17 : #include <cstddef>
18 : #include <cstdio>
19 : #include <cstring>
20 : #include <ctime>
21 :
22 : #include <string>
23 :
24 : #include "cpl_conv.h"
25 : #include "cpl_error.h"
26 : #include "cpl_minixml.h"
27 : #include "cpl_string.h"
28 : #include "cpl_time.h"
29 :
30 : /**
31 : * GDALMDReaderPleiades()
32 : */
33 10710 : GDALMDReaderPleiades::GDALMDReaderPleiades(const char *pszPath,
34 10710 : char **papszSiblingFiles)
35 : : GDALMDReaderBase(pszPath, papszSiblingFiles), m_osBaseFilename(pszPath),
36 10710 : m_osIMDSourceFilename(CPLString()), m_osRPBSourceFilename(CPLString())
37 : {
38 10710 : const CPLString osBaseName = CPLGetBasenameSafe(pszPath);
39 10710 : const size_t nBaseNameLen = osBaseName.size();
40 10710 : if (nBaseNameLen < 4 || nBaseNameLen > 511)
41 1130 : return;
42 :
43 9580 : const CPLString osDirName = CPLGetDirnameSafe(pszPath);
44 :
45 : std::string osIMDSourceFilename = CPLFormFilenameSafe(
46 19160 : osDirName, (std::string("DIM_") + (osBaseName.c_str() + 4)).c_str(),
47 9580 : "XML");
48 : std::string osRPBSourceFilename = CPLFormFilenameSafe(
49 19160 : osDirName, (std::string("RPC_") + (osBaseName.c_str() + 4)).c_str(),
50 9580 : "XML");
51 :
52 : // find last underline
53 : char sBaseName[512];
54 9580 : size_t nLastUnderline = 0;
55 76843 : for (size_t i = 4; i < nBaseNameLen; i++)
56 : {
57 67263 : sBaseName[i - 4] = osBaseName[i];
58 67263 : if (osBaseName[i] == '_')
59 10681 : nLastUnderline = i - 4U;
60 : }
61 :
62 9580 : sBaseName[nLastUnderline] = 0;
63 :
64 : // Check if last 4 characters are fit in mask RjCj
65 : unsigned int iRow, iCol;
66 9580 : bool bHasRowColPart = nBaseNameLen > nLastUnderline + 5U;
67 9580 : if (!bHasRowColPart || sscanf(osBaseName.c_str() + nLastUnderline + 5U,
68 : "R%uC%u", &iRow, &iCol) != 2)
69 : {
70 9566 : return;
71 : }
72 :
73 : // Strip of suffix from PNEO products
74 14 : char *pszLastUnderScore = strrchr(sBaseName, '_');
75 14 : if (pszLastUnderScore &&
76 3 : (EQUAL(pszLastUnderScore, "_P") || EQUAL(pszLastUnderScore, "_RGB") ||
77 3 : EQUAL(pszLastUnderScore, "_NED")))
78 : {
79 0 : *pszLastUnderScore = 0;
80 : }
81 :
82 14 : if (CPLCheckForFile(&osIMDSourceFilename[0], papszSiblingFiles))
83 : {
84 0 : m_osIMDSourceFilename = std::move(osIMDSourceFilename);
85 : }
86 : else
87 : {
88 28 : osIMDSourceFilename = CPLFormFilenameSafe(
89 42 : osDirName, ("DIM_" + std::string(sBaseName)).c_str(), "XML");
90 14 : if (CPLCheckForFile(&osIMDSourceFilename[0], papszSiblingFiles))
91 : {
92 14 : m_osIMDSourceFilename = std::move(osIMDSourceFilename);
93 : }
94 : }
95 :
96 14 : if (CPLCheckForFile(&osRPBSourceFilename[0], papszSiblingFiles))
97 : {
98 0 : m_osRPBSourceFilename = std::move(osRPBSourceFilename);
99 : }
100 : else
101 : {
102 28 : osRPBSourceFilename = CPLFormFilenameSafe(
103 42 : osDirName, ("RPC_" + std::string(sBaseName)).c_str(), "XML");
104 14 : if (CPLCheckForFile(&osRPBSourceFilename[0], papszSiblingFiles))
105 : {
106 4 : m_osRPBSourceFilename = std::move(osRPBSourceFilename);
107 : }
108 : }
109 :
110 14 : if (!m_osIMDSourceFilename.empty())
111 14 : CPLDebug("MDReaderPleiades", "IMD Filename: %s",
112 : m_osIMDSourceFilename.c_str());
113 14 : if (!m_osRPBSourceFilename.empty())
114 4 : CPLDebug("MDReaderPleiades", "RPB Filename: %s",
115 : m_osRPBSourceFilename.c_str());
116 : }
117 :
118 6 : GDALMDReaderPleiades::GDALMDReaderPleiades()
119 6 : : GDALMDReaderBase(nullptr, nullptr)
120 : {
121 6 : }
122 :
123 : GDALMDReaderPleiades *
124 6 : GDALMDReaderPleiades::CreateReaderForRPC(const char *pszRPCSourceFilename)
125 : {
126 6 : GDALMDReaderPleiades *poReader = new GDALMDReaderPleiades();
127 6 : poReader->m_osRPBSourceFilename = pszRPCSourceFilename;
128 6 : return poReader;
129 : }
130 :
131 : /**
132 : * ~GDALMDReaderPleiades()
133 : */
134 16084 : GDALMDReaderPleiades::~GDALMDReaderPleiades()
135 : {
136 16084 : }
137 :
138 : /**
139 : * HasRequiredFiles()
140 : */
141 10710 : bool GDALMDReaderPleiades::HasRequiredFiles() const
142 : {
143 10710 : if (!m_osIMDSourceFilename.empty())
144 15 : return true;
145 10695 : if (!m_osRPBSourceFilename.empty())
146 0 : return true;
147 :
148 10695 : return false;
149 : }
150 :
151 : /**
152 : * GetMetadataFiles()
153 : */
154 15 : char **GDALMDReaderPleiades::GetMetadataFiles() const
155 : {
156 15 : char **papszFileList = nullptr;
157 15 : if (!m_osIMDSourceFilename.empty())
158 15 : papszFileList = CSLAddString(papszFileList, m_osIMDSourceFilename);
159 15 : if (!m_osRPBSourceFilename.empty())
160 4 : papszFileList = CSLAddString(papszFileList, m_osRPBSourceFilename);
161 :
162 15 : return papszFileList;
163 : }
164 :
165 : /**
166 : * LoadMetadata()
167 : */
168 25 : void GDALMDReaderPleiades::LoadMetadata()
169 : {
170 25 : if (m_bIsMetadataLoad)
171 11 : return;
172 :
173 14 : if (!m_osIMDSourceFilename.empty())
174 : {
175 14 : CPLXMLNode *psNode = CPLParseXMLFile(m_osIMDSourceFilename);
176 :
177 14 : if (psNode != nullptr)
178 : {
179 14 : CPLXMLNode *psisdNode = CPLSearchXMLNode(psNode, "=Dimap_Document");
180 :
181 14 : if (psisdNode != nullptr)
182 : {
183 14 : m_papszIMDMD = ReadXMLToList(psisdNode->psChild, m_papszIMDMD);
184 : }
185 14 : CPLDestroyXMLNode(psNode);
186 : }
187 : }
188 :
189 14 : if (!m_osRPBSourceFilename.empty())
190 : {
191 4 : m_papszRPCMD = LoadRPCXmlFile();
192 : }
193 :
194 14 : m_papszDEFAULTMD =
195 14 : CSLAddNameValue(m_papszDEFAULTMD, MD_NAME_MDTYPE, "DIMAP");
196 :
197 14 : m_bIsMetadataLoad = true;
198 :
199 14 : if (nullptr == m_papszIMDMD)
200 : {
201 0 : return;
202 : }
203 :
204 : // extract imagery metadata
205 14 : int nCounter = -1;
206 28 : const char *pszSatId1 = CSLFetchNameValue(
207 14 : m_papszIMDMD,
208 : "Dataset_Sources.Source_Identification.Strip_Source.MISSION");
209 14 : if (nullptr == pszSatId1)
210 : {
211 0 : nCounter = 1;
212 0 : for (int i = 0; i < 5; i++)
213 : {
214 0 : pszSatId1 = CSLFetchNameValue(
215 0 : m_papszIMDMD,
216 : CPLSPrintf("Dataset_Sources.Source_Identification_%d.Strip_"
217 : "Source.MISSION",
218 : nCounter));
219 0 : if (nullptr != pszSatId1)
220 0 : break;
221 0 : nCounter++;
222 : }
223 : }
224 :
225 : const char *pszSatId2;
226 14 : if (nCounter == -1)
227 14 : pszSatId2 = CSLFetchNameValue(
228 14 : m_papszIMDMD,
229 : "Dataset_Sources.Source_Identification.Strip_Source.MISSION_INDEX");
230 : else
231 0 : pszSatId2 = CSLFetchNameValue(
232 0 : m_papszIMDMD, CPLSPrintf("Dataset_Sources.Source_Identification_%d."
233 : "Strip_Source.MISSION_INDEX",
234 : nCounter));
235 :
236 14 : if (nullptr != pszSatId1 && nullptr != pszSatId2)
237 : {
238 28 : m_papszIMAGERYMD = CSLAddNameValue(
239 : m_papszIMAGERYMD, MD_NAME_SATELLITE,
240 28 : CPLSPrintf("%s %s", CPLStripQuotes(pszSatId1).c_str(),
241 28 : CPLStripQuotes(pszSatId2).c_str()));
242 : }
243 0 : else if (nullptr != pszSatId1 && nullptr == pszSatId2)
244 : {
245 0 : m_papszIMAGERYMD = CSLAddNameValue(m_papszIMAGERYMD, MD_NAME_SATELLITE,
246 0 : CPLStripQuotes(pszSatId1));
247 : }
248 0 : else if (nullptr == pszSatId1 && nullptr != pszSatId2)
249 : {
250 0 : m_papszIMAGERYMD = CSLAddNameValue(m_papszIMAGERYMD, MD_NAME_SATELLITE,
251 0 : CPLStripQuotes(pszSatId2));
252 : }
253 :
254 : const char *pszDate;
255 14 : if (nCounter == -1)
256 14 : pszDate = CSLFetchNameValue(
257 14 : m_papszIMDMD,
258 : "Dataset_Sources.Source_Identification.Strip_Source.IMAGING_DATE");
259 : else
260 0 : pszDate = CSLFetchNameValue(
261 0 : m_papszIMDMD, CPLSPrintf("Dataset_Sources.Source_Identification_%d."
262 : "Strip_Source.IMAGING_DATE",
263 : nCounter));
264 :
265 14 : if (nullptr != pszDate)
266 : {
267 : const char *pszTime;
268 14 : if (nCounter == -1)
269 14 : pszTime = CSLFetchNameValue(m_papszIMDMD,
270 : "Dataset_Sources.Source_Identification."
271 : "Strip_Source.IMAGING_TIME");
272 : else
273 0 : pszTime = CSLFetchNameValue(
274 0 : m_papszIMDMD,
275 : CPLSPrintf("Dataset_Sources.Source_Identification_%d.Strip_"
276 : "Source.IMAGING_TIME",
277 : nCounter));
278 :
279 14 : if (nullptr == pszTime)
280 0 : pszTime = "00:00:00.0Z";
281 :
282 : char buffer[80];
283 : GIntBig timeMid =
284 14 : GetAcquisitionTimeFromString(CPLSPrintf("%sT%s", pszDate, pszTime));
285 : struct tm tmBuf;
286 14 : strftime(buffer, 80, MD_DATETIMEFORMAT,
287 14 : CPLUnixTimeToYMDHMS(timeMid, &tmBuf));
288 14 : m_papszIMAGERYMD =
289 14 : CSLAddNameValue(m_papszIMAGERYMD, MD_NAME_ACQDATETIME, buffer);
290 : }
291 :
292 14 : m_papszIMAGERYMD =
293 14 : CSLAddNameValue(m_papszIMAGERYMD, MD_NAME_CLOUDCOVER, MD_CLOUDCOVER_NA);
294 : }
295 :
296 : /**
297 : * LoadRPCXmlFile()
298 : */
299 :
300 : static const char *const apszRPBMapPleiades[] = {
301 : RPC_LINE_OFF, "RFM_Validity.LINE_OFF", // do not change order !
302 : RPC_SAMP_OFF, "RFM_Validity.SAMP_OFF", // do not change order !
303 : RPC_LAT_OFF, "RFM_Validity.LAT_OFF",
304 : RPC_LONG_OFF, "RFM_Validity.LONG_OFF",
305 : RPC_HEIGHT_OFF, "RFM_Validity.HEIGHT_OFF",
306 : RPC_LINE_SCALE, "RFM_Validity.LINE_SCALE",
307 : RPC_SAMP_SCALE, "RFM_Validity.SAMP_SCALE",
308 : RPC_LAT_SCALE, "RFM_Validity.LAT_SCALE",
309 : RPC_LONG_SCALE, "RFM_Validity.LONG_SCALE",
310 : RPC_HEIGHT_SCALE, "RFM_Validity.HEIGHT_SCALE",
311 : nullptr, nullptr};
312 :
313 : static const char *const apszRPCTXT20ValItemsPleiades[] = {
314 : RPC_LINE_NUM_COEFF, RPC_LINE_DEN_COEFF, RPC_SAMP_NUM_COEFF,
315 : RPC_SAMP_DEN_COEFF, nullptr};
316 :
317 10 : char **GDALMDReaderPleiades::LoadRPCXmlFile()
318 : {
319 10 : CPLXMLNode *pNode = CPLParseXMLFile(m_osRPBSourceFilename);
320 :
321 10 : if (nullptr == pNode)
322 0 : return nullptr;
323 :
324 : // search Global_RFM
325 10 : char **papszRawRPCList = nullptr;
326 10 : CPLXMLNode *pGRFMNode = CPLSearchXMLNode(pNode, "=Global_RFM");
327 :
328 10 : if (pGRFMNode != nullptr)
329 : {
330 10 : papszRawRPCList = ReadXMLToList(pGRFMNode->psChild, papszRawRPCList);
331 : }
332 : else
333 : {
334 0 : pGRFMNode = CPLSearchXMLNode(pNode, "=Rational_Function_Model");
335 :
336 0 : if (pGRFMNode != nullptr)
337 : {
338 : papszRawRPCList =
339 0 : ReadXMLToList(pGRFMNode->psChild, papszRawRPCList);
340 : }
341 : }
342 :
343 10 : if (nullptr == papszRawRPCList)
344 : {
345 0 : CPLDestroyXMLNode(pNode);
346 0 : return nullptr;
347 : }
348 :
349 : // If we are not the top-left tile, then we must shift LINE_OFF and SAMP_OFF
350 10 : int nLineOffShift = 0;
351 10 : int nPixelOffShift = 0;
352 10 : for (int i = 1; TRUE; i++)
353 : {
354 11 : CPLString osKey;
355 : osKey.Printf("Raster_Data.Data_Access.Data_Files.Data_File_%d.DATA_"
356 : "FILE_PATH.href",
357 11 : i);
358 11 : const char *pszHref = CSLFetchNameValue(m_papszIMDMD, osKey);
359 11 : if (pszHref == nullptr)
360 9 : break;
361 2 : if (strcmp(CPLGetFilename(pszHref), CPLGetFilename(m_osBaseFilename)) ==
362 : 0)
363 : {
364 : osKey.Printf(
365 1 : "Raster_Data.Data_Access.Data_Files.Data_File_%d.tile_C", i);
366 1 : const char *pszC = CSLFetchNameValue(m_papszIMDMD, osKey);
367 : osKey.Printf(
368 1 : "Raster_Data.Data_Access.Data_Files.Data_File_%d.tile_R", i);
369 1 : const char *pszR = CSLFetchNameValue(m_papszIMDMD, osKey);
370 2 : const char *pszTileWidth = CSLFetchNameValue(
371 1 : m_papszIMDMD, "Raster_Data.Raster_Dimensions.Tile_Set.Regular_"
372 : "Tiling.NTILES_SIZE.ncols");
373 2 : const char *pszTileHeight = CSLFetchNameValue(
374 1 : m_papszIMDMD, "Raster_Data.Raster_Dimensions.Tile_Set.Regular_"
375 : "Tiling.NTILES_SIZE.nrows");
376 : const char *pszOVERLAP_COL =
377 1 : CSLFetchNameValueDef(m_papszIMDMD,
378 : "Raster_Data.Raster_Dimensions.Tile_Set."
379 : "Regular_Tiling.OVERLAP_COL",
380 : "0");
381 : const char *pszOVERLAP_ROW =
382 1 : CSLFetchNameValueDef(m_papszIMDMD,
383 : "Raster_Data.Raster_Dimensions.Tile_Set."
384 : "Regular_Tiling.OVERLAP_ROW",
385 : "0");
386 :
387 1 : if (pszC && pszR && pszTileWidth && pszTileHeight &&
388 1 : atoi(pszOVERLAP_COL) == 0 && atoi(pszOVERLAP_ROW) == 0)
389 : {
390 1 : nLineOffShift = -(atoi(pszR) - 1) * atoi(pszTileHeight);
391 1 : nPixelOffShift = -(atoi(pszC) - 1) * atoi(pszTileWidth);
392 : }
393 1 : break;
394 : }
395 1 : }
396 :
397 : // SPOT and PHR sensors use 1,1 as their upper left corner pixel convention
398 : // for RPCs which is non standard. This was fixed with PNEO which correctly
399 : // assumes 0,0.
400 : // Precompute the offset that will be applied to LINE_OFF and SAMP_OFF
401 : // in order to use the RPCs with the standard 0,0 convention
402 : double topleftOffset;
403 10 : CPLXMLNode *psDoc = CPLGetXMLNode(pNode, "=Dimap_Document");
404 10 : if (!psDoc)
405 0 : psDoc = CPLGetXMLNode(pNode, "=PHR_DIMAP_Document");
406 10 : const char *pszMetadataProfile = CPLGetXMLValue(
407 : psDoc, "Metadata_Identification.METADATA_PROFILE", "PHR_SENSOR");
408 10 : if (EQUAL(pszMetadataProfile, "PHR_SENSOR") ||
409 2 : EQUAL(pszMetadataProfile, "S7_SENSOR") ||
410 2 : EQUAL(pszMetadataProfile, "S6_SENSOR"))
411 : {
412 8 : topleftOffset = 1;
413 : }
414 2 : else if (EQUAL(pszMetadataProfile, "PNEO_SENSOR"))
415 : {
416 2 : topleftOffset = 0;
417 : }
418 : else
419 : {
420 : //CPLError(CE_Warning, CPLE_AppDefined,
421 : // "Unknown RPC Metadata Profile: %s. Assuming PHR_SENSOR",
422 : // pszMetadataProfile);
423 0 : topleftOffset = 1;
424 : }
425 :
426 : // format list
427 10 : char **papszRPB = nullptr;
428 110 : for (int i = 0; apszRPBMapPleiades[i] != nullptr; i += 2)
429 : {
430 : const char *pszValue =
431 100 : CSLFetchNameValue(papszRawRPCList, apszRPBMapPleiades[i + 1]);
432 100 : if ((i == 0 || i == 2) && pszValue) //i.e. LINE_OFF or SAMP_OFF
433 : {
434 20 : CPLString osField;
435 20 : double dfVal = CPLAtofM(pszValue) - topleftOffset;
436 20 : if (i == 0)
437 10 : dfVal += nLineOffShift;
438 : else
439 10 : dfVal += nPixelOffShift;
440 20 : osField.Printf("%.15g", dfVal);
441 : papszRPB =
442 20 : CSLAddNameValue(papszRPB, apszRPBMapPleiades[i], osField);
443 : }
444 : else
445 : {
446 : papszRPB =
447 80 : CSLAddNameValue(papszRPB, apszRPBMapPleiades[i], pszValue);
448 : }
449 : }
450 :
451 : // merge coefficients
452 50 : for (int i = 0; apszRPCTXT20ValItemsPleiades[i] != nullptr; i++)
453 : {
454 40 : CPLString value;
455 840 : for (int j = 1; j < 21; j++)
456 : {
457 : // We want to use the Inverse_Model
458 : // Quoting PleiadesUserGuideV2-1012.pdf:
459 : // """When using the inverse model (ground --> image), the user
460 : // supplies geographic coordinates (lon, lat) and an altitude
461 : // (alt)"""
462 800 : const char *pszValue = CSLFetchNameValue(
463 : papszRawRPCList,
464 : CPLSPrintf("Inverse_Model.%s_%d",
465 800 : apszRPCTXT20ValItemsPleiades[i], j));
466 800 : if (nullptr != pszValue)
467 : {
468 640 : value = value + " " + CPLString(pszValue);
469 : }
470 : else
471 : {
472 160 : pszValue = CSLFetchNameValue(
473 : papszRawRPCList,
474 : CPLSPrintf("GroundtoImage_Values.%s_%d",
475 160 : apszRPCTXT20ValItemsPleiades[i], j));
476 160 : if (nullptr != pszValue)
477 : {
478 160 : value = value + " " + CPLString(pszValue);
479 : }
480 : }
481 : }
482 : papszRPB =
483 40 : CSLAddNameValue(papszRPB, apszRPCTXT20ValItemsPleiades[i], value);
484 : }
485 :
486 10 : CSLDestroy(papszRawRPCList);
487 10 : CPLDestroyXMLNode(pNode);
488 10 : return papszRPB;
489 : }
|