LCOV - code coverage report
Current view: top level - gcore/mdreader - reader_pleiades.cpp (source / functions) Hit Total Coverage
Test: gdal_filtered.info Lines: 170 203 83.7 %
Date: 2025-01-18 12:42:00 Functions: 9 9 100.0 %

          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             : }

Generated by: LCOV version 1.14