Line data    Source code

       1             : /******************************************************************************
       2             :  *
       3             :  * Project:     TerraSAR-X XML Product Support
       4             :  * Purpose:     Support for TerraSAR-X XML Metadata files
       5             :  * Author:      Philippe Vachon <philippe@cowpig.ca>
       6             :  * Description: This driver adds support for reading metadata and georef data
       7             :  *              associated with TerraSAR-X products.
       8             :  *
       9             :  ******************************************************************************
      10             :  * Copyright (c) 2007, Philippe Vachon <philippe@cowpig.ca>
      11             :  * Copyright (c) 2009-2012, Even Rouault <even dot rouault at spatialys.com>
      12             :  *
      13             :  * SPDX-License-Identifier: MIT
      14             :  ****************************************************************************/
      15             : 
      16             : #include "cpl_minixml.h"
      17             : #include "gdal_frmts.h"
      18             : #include "gdal_pam.h"
      19             : #include "ogr_spatialref.h"
      20             : 
      21             : #define MAX_GCPS 5000  // this should be more than enough ground control points
      22             : 
      23             : namespace gdal::TSX
      24             : {
      25             : enum ePolarization
      26             : {
      27             :     HH = 0,
      28             :     HV,
      29             :     VH,
      30             :     VV
      31             : };
      32             : 
      33             : enum eProductType
      34             : {
      35             :     eSSC = 0,
      36             :     eMGD,
      37             :     eEEC,
      38             :     eGEC,
      39             :     eUnknown
      40             : };
      41             : }  // namespace gdal::TSX
      42             : 
      43             : using namespace gdal::TSX;
      44             : 
      45             : /************************************************************************/
      46             : /* Helper Functions                                                     */
      47             : /************************************************************************/
      48             : 
      49             : /* GetFilePath: return a relative path to a file within an XML node.
      50             :  * Returns Null on failure
      51             :  */
      52           0 : static CPLString GetFilePath(CPLXMLNode *psXMLNode, const char **pszNodeType)
      53             : {
      54             :     const char *pszDirectory =
      55           0 :         CPLGetXMLValue(psXMLNode, "file.location.path", "");
      56             :     const char *pszFilename =
      57           0 :         CPLGetXMLValue(psXMLNode, "file.location.filename", "");
      58           0 :     *pszNodeType = CPLGetXMLValue(psXMLNode, "type", " ");
      59             : 
      60           0 :     if (pszDirectory == nullptr || pszFilename == nullptr)
      61             :     {
      62           0 :         return "";
      63             :     }
      64             : 
      65           0 :     return CPLString(pszDirectory) + '/' + pszFilename;
      66             : }
      67             : 
      68             : /************************************************************************/
      69             : /* ==================================================================== */
      70             : /*                                TSXDataset                                 */
      71             : /* ==================================================================== */
      72             : /************************************************************************/
      73             : 
      74             : class TSXDataset final : public GDALPamDataset
      75             : {
      76             :     int nGCPCount;
      77             :     GDAL_GCP *pasGCPList;
      78             : 
      79             :     OGRSpatialReference m_oGCPSRS{};
      80             : 
      81             :     OGRSpatialReference m_oSRS{};
      82             :     double adfGeoTransform[6];
      83             :     bool bHaveGeoTransform;
      84             : 
      85             :     eProductType nProduct;
      86             : 
      87             :   public:
      88             :     TSXDataset();
      89             :     virtual ~TSXDataset();
      90             : 
      91             :     virtual int GetGCPCount() override;
      92             :     const OGRSpatialReference *GetGCPSpatialRef() const override;
      93             :     virtual const GDAL_GCP *GetGCPs() override;
      94             : 
      95             :     CPLErr GetGeoTransform(double *padfTransform) override;
      96             :     const OGRSpatialReference *GetSpatialRef() const override;
      97             : 
      98             :     static GDALDataset *Open(GDALOpenInfo *poOpenInfo);
      99             :     static int Identify(GDALOpenInfo *poOpenInfo);
     100             : 
     101             :   private:
     102             :     bool getGCPsFromGEOREF_XML(char *pszGeorefFilename);
     103             : };
     104             : 
     105             : /************************************************************************/
     106             : /* ==================================================================== */
     107             : /*                                TSXRasterBand                           */
     108             : /* ==================================================================== */
     109             : /************************************************************************/
     110             : 
     111             : class TSXRasterBand final : public GDALPamRasterBand
     112             : {
     113             :     GDALDataset *poBand;
     114             :     ePolarization ePol;
     115             : 
     116             :   public:
     117             :     TSXRasterBand(TSXDataset *poDSIn, GDALDataType eDataType,
     118             :                   ePolarization ePol, GDALDataset *poBand);
     119             :     virtual ~TSXRasterBand();
     120             : 
     121             :     virtual CPLErr IReadBlock(int nBlockXOff, int nBlockYOff,
     122             :                               void *pImage) override;
     123             : 
     124             :     static GDALDataset *Open(GDALOpenInfo *poOpenInfo);
     125             : };
     126             : 
     127             : /************************************************************************/
     128             : /*                            TSXRasterBand                             */
     129             : /************************************************************************/
     130             : 
     131           0 : TSXRasterBand::TSXRasterBand(TSXDataset *poDSIn, GDALDataType eDataTypeIn,
     132           0 :                              ePolarization ePolIn, GDALDataset *poBandIn)
     133           0 :     : poBand(poBandIn), ePol(ePolIn)
     134             : {
     135           0 :     poDS = poDSIn;
     136           0 :     eDataType = eDataTypeIn;
     137             : 
     138           0 :     switch (ePol)
     139             :     {
     140           0 :         case HH:
     141           0 :             SetMetadataItem("POLARIMETRIC_INTERP", "HH");
     142           0 :             break;
     143           0 :         case HV:
     144           0 :             SetMetadataItem("POLARIMETRIC_INTERP", "HV");
     145           0 :             break;
     146           0 :         case VH:
     147           0 :             SetMetadataItem("POLARIMETRIC_INTERP", "VH");
     148           0 :             break;
     149           0 :         case VV:
     150           0 :             SetMetadataItem("POLARIMETRIC_INTERP", "VV");
     151           0 :             break;
     152             :     }
     153             : 
     154             :     /* now setup the actual raster reader */
     155           0 :     GDALRasterBand *poSrcBand = poBandIn->GetRasterBand(1);
     156           0 :     poSrcBand->GetBlockSize(&nBlockXSize, &nBlockYSize);
     157           0 : }
     158             : 
     159             : /************************************************************************/
     160             : /*                            TSXRasterBand()                           */
     161             : /************************************************************************/
     162             : 
     163           0 : TSXRasterBand::~TSXRasterBand()
     164             : {
     165           0 :     if (poBand != nullptr)
     166           0 :         GDALClose(reinterpret_cast<GDALRasterBandH>(poBand));
     167           0 : }
     168             : 
     169             : /************************************************************************/
     170             : /*                             IReadBlock()                             */
     171             : /************************************************************************/
     172             : 
     173           0 : CPLErr TSXRasterBand::IReadBlock(int nBlockXOff, int nBlockYOff, void *pImage)
     174             : {
     175             :     int nRequestYSize;
     176             : 
     177             :     /* Check if the last strip is partial so we can avoid over-requesting */
     178           0 :     if ((nBlockYOff + 1) * nBlockYSize > nRasterYSize)
     179             :     {
     180           0 :         nRequestYSize = nRasterYSize - nBlockYOff * nBlockYSize;
     181           0 :         memset(pImage, 0,
     182           0 :                static_cast<size_t>(GDALGetDataTypeSizeBytes(eDataType)) *
     183           0 :                    nBlockXSize * nBlockYSize);
     184             :     }
     185             :     else
     186             :     {
     187           0 :         nRequestYSize = nBlockYSize;
     188             :     }
     189             : 
     190             :     /* Read Complex Data */
     191           0 :     if (eDataType == GDT_CInt16)
     192             :     {
     193           0 :         return poBand->RasterIO(
     194           0 :             GF_Read, nBlockXOff * nBlockXSize, nBlockYOff * nBlockYSize,
     195             :             nBlockXSize, nRequestYSize, pImage, nBlockXSize, nRequestYSize,
     196           0 :             GDT_CInt16, 1, nullptr, 4, nBlockXSize * 4, 0, nullptr);
     197             :     }
     198             : 
     199             :     // Detected Product
     200           0 :     return poBand->RasterIO(
     201           0 :         GF_Read, nBlockXOff * nBlockXSize, nBlockYOff * nBlockYSize,
     202             :         nBlockXSize, nRequestYSize, pImage, nBlockXSize, nRequestYSize,
     203           0 :         GDT_UInt16, 1, nullptr, 2, nBlockXSize * 2, 0, nullptr);
     204             : }
     205             : 
     206             : /************************************************************************/
     207             : /* ==================================================================== */
     208             : /*                                TSXDataset                                */
     209             : /* ==================================================================== */
     210             : /************************************************************************/
     211             : 
     212             : /************************************************************************/
     213             : /*                             TSXDataset()                             */
     214             : /************************************************************************/
     215             : 
     216           0 : TSXDataset::TSXDataset()
     217             :     : nGCPCount(0), pasGCPList(nullptr), bHaveGeoTransform(false),
     218           0 :       nProduct(eUnknown)
     219             : {
     220           0 :     m_oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
     221           0 :     m_oGCPSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
     222           0 :     adfGeoTransform[0] = 0.0;
     223           0 :     adfGeoTransform[1] = 1.0;
     224           0 :     adfGeoTransform[2] = 0.0;
     225           0 :     adfGeoTransform[3] = 0.0;
     226           0 :     adfGeoTransform[4] = 0.0;
     227           0 :     adfGeoTransform[5] = 1.0;
     228           0 : }
     229             : 
     230             : /************************************************************************/
     231             : /*                            ~TSXDataset()                             */
     232             : /************************************************************************/
     233             : 
     234           0 : TSXDataset::~TSXDataset()
     235             : {
     236           0 :     FlushCache(true);
     237             : 
     238           0 :     if (nGCPCount > 0)
     239             :     {
     240           0 :         GDALDeinitGCPs(nGCPCount, pasGCPList);
     241           0 :         CPLFree(pasGCPList);
     242             :     }
     243           0 : }
     244             : 
     245             : /************************************************************************/
     246             : /*                              Identify()                              */
     247             : /************************************************************************/
     248             : 
     249       53352 : int TSXDataset::Identify(GDALOpenInfo *poOpenInfo)
     250             : {
     251       53352 :     if (poOpenInfo->fpL == nullptr || poOpenInfo->nHeaderBytes < 260)
     252             :     {
     253       49655 :         if (poOpenInfo->bIsDirectory)
     254             :         {
     255           0 :             const CPLString osFilename = CPLFormCIFilenameSafe(
     256         420 :                 poOpenInfo->pszFilename,
     257         840 :                 CPLGetFilename(poOpenInfo->pszFilename), "xml");
     258             : 
     259             :             /* Check if the filename contains TSX1_SAR (TerraSAR-X) or TDX1_SAR
     260             :              * (TanDEM-X) or PAZ1_SAR (PAZ) */
     261        1260 :             if (!(STARTS_WITH_CI(CPLGetBasenameSafe(osFilename).c_str(),
     262         840 :                                  "TSX1_SAR") ||
     263         840 :                   STARTS_WITH_CI(CPLGetBasenameSafe(osFilename).c_str(),
     264             :                                  "TDX1_SAR") ||
     265         840 :                   STARTS_WITH_CI(CPLGetBasenameSafe(osFilename).c_str(),
     266             :                                  "PAZ1_SAR")))
     267         420 :                 return 0;
     268             : 
     269             :             VSIStatBufL sStat;
     270           0 :             if (VSIStatL(osFilename, &sStat) == 0)
     271           0 :                 return 1;
     272             :         }
     273             : 
     274       49235 :         return 0;
     275             :     }
     276             : 
     277             :     /* Check if the filename contains TSX1_SAR (TerraSAR-X) or TDX1_SAR
     278             :      * (TanDEM-X) or PAZ1_SAR (PAZ) */
     279       11091 :     if (!(STARTS_WITH_CI(CPLGetBasenameSafe(poOpenInfo->pszFilename).c_str(),
     280        7394 :                          "TSX1_SAR") ||
     281        7394 :           STARTS_WITH_CI(CPLGetBasenameSafe(poOpenInfo->pszFilename).c_str(),
     282             :                          "TDX1_SAR") ||
     283        7394 :           STARTS_WITH_CI(CPLGetBasenameSafe(poOpenInfo->pszFilename).c_str(),
     284             :                          "PAZ1_SAR")))
     285        3697 :         return 0;
     286             : 
     287             :     /* finally look for the <level1Product tag */
     288           0 :     if (!STARTS_WITH_CI(reinterpret_cast<char *>(poOpenInfo->pabyHeader),
     289             :                         "<level1Product"))
     290           0 :         return 0;
     291             : 
     292           0 :     return 1;
     293             : }
     294             : 
     295             : /************************************************************************/
     296             : /*                                getGCPsFromGEOREF_XML()               */
     297             : /*Reads georeferencing information from the TerraSAR-X GEOREF.XML file    */
     298             : /*and writes the information to the dataset's gcp list and projection     */
     299             : /*string.                                                                */
     300             : /*Returns true on success.                                                */
     301             : /************************************************************************/
     302           0 : bool TSXDataset::getGCPsFromGEOREF_XML(char *pszGeorefFilename)
     303             : {
     304             :     // open GEOREF.xml
     305           0 :     CPLXMLNode *psGeorefData = CPLParseXMLFile(pszGeorefFilename);
     306           0 :     if (psGeorefData == nullptr)
     307           0 :         return false;
     308             : 
     309             :     // get the ellipsoid and semi-major, semi-minor axes
     310           0 :     OGRSpatialReference osr;
     311             :     CPLXMLNode *psSphere =
     312           0 :         CPLGetXMLNode(psGeorefData, "=geoReference.referenceFrames.sphere");
     313           0 :     if (psSphere != nullptr)
     314             :     {
     315             :         const char *pszEllipsoidName =
     316           0 :             CPLGetXMLValue(psSphere, "ellipsoidID", "");
     317             :         const double minor_axis =
     318           0 :             CPLAtof(CPLGetXMLValue(psSphere, "semiMinorAxis", "0.0"));
     319             :         const double major_axis =
     320           0 :             CPLAtof(CPLGetXMLValue(psSphere, "semiMajorAxis", "0.0"));
     321             :         // save datum parameters to the spatial reference
     322           0 :         if (EQUAL(pszEllipsoidName, "") || minor_axis == 0.0 ||
     323             :             major_axis == 0.0)
     324             :         {
     325           0 :             CPLError(CE_Warning, CPLE_AppDefined,
     326             :                      "Warning- incomplete"
     327             :                      " ellipsoid information.  Using wgs-84 parameters.\n");
     328           0 :             osr.SetWellKnownGeogCS("WGS84");
     329             :         }
     330           0 :         else if (EQUAL(pszEllipsoidName, "WGS84"))
     331             :         {
     332           0 :             osr.SetWellKnownGeogCS("WGS84");
     333             :         }
     334             :         else
     335             :         {
     336           0 :             const double inv_flattening =
     337           0 :                 major_axis / (major_axis - minor_axis);
     338           0 :             osr.SetGeogCS("", "", pszEllipsoidName, major_axis, inv_flattening);
     339             :         }
     340             :     }
     341             : 
     342             :     // get gcps
     343             :     CPLXMLNode *psGeolocationGrid =
     344           0 :         CPLGetXMLNode(psGeorefData, "=geoReference.geolocationGrid");
     345           0 :     if (psGeolocationGrid == nullptr)
     346             :     {
     347           0 :         CPLDestroyXMLNode(psGeorefData);
     348           0 :         return false;
     349             :     }
     350           0 :     nGCPCount = atoi(
     351             :         CPLGetXMLValue(psGeolocationGrid, "numberOfGridPoints.total", "0"));
     352             :     // count the gcps if the given count value is invalid
     353           0 :     CPLXMLNode *psNode = nullptr;
     354           0 :     if (nGCPCount <= 0)
     355             :     {
     356           0 :         for (psNode = psGeolocationGrid->psChild; psNode != nullptr;
     357           0 :              psNode = psNode->psNext)
     358           0 :             if (EQUAL(psNode->pszValue, "gridPoint"))
     359           0 :                 nGCPCount++;
     360             :     }
     361             :     // if there are no gcps, fail
     362           0 :     if (nGCPCount <= 0)
     363             :     {
     364           0 :         CPLDestroyXMLNode(psGeorefData);
     365           0 :         return false;
     366             :     }
     367             : 
     368             :     // put some reasonable limits of the number of gcps
     369           0 :     if (nGCPCount > MAX_GCPS)
     370           0 :         nGCPCount = MAX_GCPS;
     371             :     // allocate memory for the gcps
     372           0 :     pasGCPList =
     373           0 :         reinterpret_cast<GDAL_GCP *>(CPLCalloc(sizeof(GDAL_GCP), nGCPCount));
     374             : 
     375             :     // loop through all gcps and set info
     376             : 
     377             :     // save the number allocated to ensure it does not run off the end of the
     378             :     // array
     379           0 :     const int gcps_allocated = nGCPCount;
     380           0 :     nGCPCount = 0;  // reset to zero and count
     381             :     // do a check on the grid point to make sure it has lat,long row, and column
     382             :     // it seems that only SSC products contain row, col - how to map lat long
     383             :     // otherwise?? for now fail if row and col are not present - just check the
     384             :     // first and assume the rest are the same
     385           0 :     for (psNode = psGeolocationGrid->psChild; psNode != nullptr;
     386           0 :          psNode = psNode->psNext)
     387             :     {
     388           0 :         if (!EQUAL(psNode->pszValue, "gridPoint"))
     389           0 :             continue;
     390             : 
     391           0 :         if (!strcmp(CPLGetXMLValue(psNode, "col", "error"), "error") ||
     392           0 :             !strcmp(CPLGetXMLValue(psNode, "row", "error"), "error") ||
     393           0 :             !strcmp(CPLGetXMLValue(psNode, "lon", "error"), "error") ||
     394           0 :             !strcmp(CPLGetXMLValue(psNode, "lat", "error"), "error"))
     395             :         {
     396           0 :             CPLDestroyXMLNode(psGeorefData);
     397           0 :             return false;
     398             :         }
     399             :     }
     400           0 :     for (psNode = psGeolocationGrid->psChild; psNode != nullptr;
     401           0 :          psNode = psNode->psNext)
     402             :     {
     403             :         // break out if the end of the array has been reached
     404           0 :         if (nGCPCount >= gcps_allocated)
     405             :         {
     406           0 :             CPLError(CE_Warning, CPLE_AppDefined,
     407             :                      "GDAL TSX driver: Truncating the number of GCPs.");
     408           0 :             break;
     409             :         }
     410             : 
     411           0 :         GDAL_GCP *psGCP = pasGCPList + nGCPCount;
     412             : 
     413           0 :         if (!EQUAL(psNode->pszValue, "gridPoint"))
     414           0 :             continue;
     415             : 
     416           0 :         nGCPCount++;
     417             : 
     418             :         char szID[32];
     419           0 :         snprintf(szID, sizeof(szID), "%d", nGCPCount);
     420           0 :         psGCP->pszId = CPLStrdup(szID);
     421           0 :         psGCP->pszInfo = CPLStrdup("");
     422           0 :         psGCP->dfGCPPixel = CPLAtof(CPLGetXMLValue(psNode, "col", "0"));
     423           0 :         psGCP->dfGCPLine = CPLAtof(CPLGetXMLValue(psNode, "row", "0"));
     424           0 :         psGCP->dfGCPX = CPLAtof(CPLGetXMLValue(psNode, "lon", ""));
     425           0 :         psGCP->dfGCPY = CPLAtof(CPLGetXMLValue(psNode, "lat", ""));
     426             :         // looks like height is in meters - should it be converted so xyz are
     427             :         // all on the same scale??
     428           0 :         psGCP->dfGCPZ = 0;
     429             :         // CPLAtof(CPLGetXMLValue(psNode,"height",""));
     430             :     }
     431             : 
     432           0 :     m_oGCPSRS = std::move(osr);
     433             : 
     434           0 :     CPLDestroyXMLNode(psGeorefData);
     435             : 
     436           0 :     return true;
     437             : }
     438             : 
     439             : /************************************************************************/
     440             : /*                                Open()                                */
     441             : /************************************************************************/
     442             : 
     443           0 : GDALDataset *TSXDataset::Open(GDALOpenInfo *poOpenInfo)
     444             : {
     445             :     /* -------------------------------------------------------------------- */
     446             :     /*      Is this a TerraSAR-X product file?                              */
     447             :     /* -------------------------------------------------------------------- */
     448           0 :     if (!TSXDataset::Identify(poOpenInfo))
     449             :     {
     450           0 :         return nullptr; /* nope */
     451             :     }
     452             : 
     453             :     /* -------------------------------------------------------------------- */
     454             :     /*      Confirm the requested access is supported.                      */
     455             :     /* -------------------------------------------------------------------- */
     456           0 :     if (poOpenInfo->eAccess == GA_Update)
     457             :     {
     458           0 :         CPLError(CE_Failure, CPLE_NotSupported,
     459             :                  "The TSX driver does not support update access to existing"
     460             :                  " datasets.\n");
     461           0 :         return nullptr;
     462             :     }
     463             : 
     464           0 :     CPLString osFilename;
     465             : 
     466           0 :     if (poOpenInfo->bIsDirectory)
     467             :     {
     468           0 :         osFilename = CPLFormCIFilenameSafe(
     469           0 :             poOpenInfo->pszFilename, CPLGetFilename(poOpenInfo->pszFilename),
     470           0 :             "xml");
     471             :     }
     472             :     else
     473           0 :         osFilename = poOpenInfo->pszFilename;
     474             : 
     475             :     /* Ingest the XML */
     476           0 :     CPLXMLNode *psData = CPLParseXMLFile(osFilename);
     477           0 :     if (psData == nullptr)
     478           0 :         return nullptr;
     479             : 
     480             :     /* find the product components */
     481             :     CPLXMLNode *psComponents =
     482           0 :         CPLGetXMLNode(psData, "=level1Product.productComponents");
     483           0 :     if (psComponents == nullptr)
     484             :     {
     485           0 :         CPLError(CE_Failure, CPLE_OpenFailed,
     486             :                  "Unable to find <productComponents> tag in file.\n");
     487           0 :         CPLDestroyXMLNode(psData);
     488           0 :         return nullptr;
     489             :     }
     490             : 
     491             :     /* find the product info tag */
     492             :     CPLXMLNode *psProductInfo =
     493           0 :         CPLGetXMLNode(psData, "=level1Product.productInfo");
     494           0 :     if (psProductInfo == nullptr)
     495             :     {
     496           0 :         CPLError(CE_Failure, CPLE_OpenFailed,
     497             :                  "Unable to find <productInfo> tag in file.\n");
     498           0 :         CPLDestroyXMLNode(psData);
     499           0 :         return nullptr;
     500             :     }
     501             : 
     502             :     /* -------------------------------------------------------------------- */
     503             :     /*      Create the dataset.                                             */
     504             :     /* -------------------------------------------------------------------- */
     505             : 
     506           0 :     TSXDataset *poDS = new TSXDataset();
     507             : 
     508             :     /* -------------------------------------------------------------------- */
     509             :     /*      Read in product info.                                           */
     510             :     /* -------------------------------------------------------------------- */
     511             : 
     512           0 :     poDS->SetMetadataItem(
     513             :         "SCENE_CENTRE_TIME",
     514             :         CPLGetXMLValue(psProductInfo,
     515           0 :                        "sceneInfo.sceneCenterCoord.azimuthTimeUTC", "unknown"));
     516           0 :     poDS->SetMetadataItem("OPERATIONAL_MODE",
     517             :                           CPLGetXMLValue(psProductInfo,
     518             :                                          "generationInfo.groundOperationsType",
     519           0 :                                          "unknown"));
     520           0 :     poDS->SetMetadataItem(
     521             :         "ORBIT_CYCLE",
     522           0 :         CPLGetXMLValue(psProductInfo, "missionInfo.orbitCycle", "unknown"));
     523           0 :     poDS->SetMetadataItem(
     524             :         "ABSOLUTE_ORBIT",
     525           0 :         CPLGetXMLValue(psProductInfo, "missionInfo.absOrbit", "unknown"));
     526           0 :     poDS->SetMetadataItem(
     527             :         "ORBIT_DIRECTION",
     528           0 :         CPLGetXMLValue(psProductInfo, "missionInfo.orbitDirection", "unknown"));
     529           0 :     poDS->SetMetadataItem("IMAGING_MODE",
     530             :                           CPLGetXMLValue(psProductInfo,
     531             :                                          "acquisitionInfo.imagingMode",
     532           0 :                                          "unknown"));
     533           0 :     poDS->SetMetadataItem("PRODUCT_VARIANT",
     534             :                           CPLGetXMLValue(psProductInfo,
     535             :                                          "productVariantInfo.productVariant",
     536           0 :                                          "unknown"));
     537           0 :     char *pszDataType = CPLStrdup(CPLGetXMLValue(
     538             :         psProductInfo, "imageDataInfo.imageDataType", "unknown"));
     539           0 :     poDS->SetMetadataItem("IMAGE_TYPE", pszDataType);
     540             : 
     541             :     /* Get raster information */
     542           0 :     int nRows = atoi(CPLGetXMLValue(
     543             :         psProductInfo, "imageDataInfo.imageRaster.numberOfRows", ""));
     544           0 :     int nCols = atoi(CPLGetXMLValue(
     545             :         psProductInfo, "imageDataInfo.imageRaster.numberOfColumns", ""));
     546             : 
     547           0 :     poDS->nRasterXSize = nCols;
     548           0 :     poDS->nRasterYSize = nRows;
     549             : 
     550           0 :     poDS->SetMetadataItem("ROW_SPACING",
     551             :                           CPLGetXMLValue(psProductInfo,
     552             :                                          "imageDataInfo.imageRaster.rowSpacing",
     553           0 :                                          "unknown"));
     554           0 :     poDS->SetMetadataItem(
     555             :         "COL_SPACING",
     556             :         CPLGetXMLValue(psProductInfo, "imageDataInfo.imageRaster.columnSpacing",
     557           0 :                        "unknown"));
     558           0 :     poDS->SetMetadataItem(
     559             :         "COL_SPACING_UNITS",
     560             :         CPLGetXMLValue(psProductInfo,
     561             :                        "imageDataInfo.imageRaster.columnSpacing.units",
     562           0 :                        "unknown"));
     563             : 
     564             :     /* Get equivalent number of looks */
     565           0 :     poDS->SetMetadataItem(
     566             :         "AZIMUTH_LOOKS",
     567             :         CPLGetXMLValue(psProductInfo, "imageDataInfo.imageRaster.azimuthLooks",
     568           0 :                        "unknown"));
     569           0 :     poDS->SetMetadataItem("RANGE_LOOKS",
     570             :                           CPLGetXMLValue(psProductInfo,
     571             :                                          "imageDataInfo.imageRaster.rangeLooks",
     572           0 :                                          "unknown"));
     573             : 
     574           0 :     const char *pszProductVariant = CPLGetXMLValue(
     575             :         psProductInfo, "productVariantInfo.productVariant", "unknown");
     576             : 
     577           0 :     poDS->SetMetadataItem("PRODUCT_VARIANT", pszProductVariant);
     578             : 
     579             :     /* Determine what product variant this is */
     580           0 :     if (STARTS_WITH_CI(pszProductVariant, "SSC"))
     581           0 :         poDS->nProduct = eSSC;
     582           0 :     else if (STARTS_WITH_CI(pszProductVariant, "MGD"))
     583           0 :         poDS->nProduct = eMGD;
     584           0 :     else if (STARTS_WITH_CI(pszProductVariant, "EEC"))
     585           0 :         poDS->nProduct = eEEC;
     586           0 :     else if (STARTS_WITH_CI(pszProductVariant, "GEC"))
     587           0 :         poDS->nProduct = eGEC;
     588             :     else
     589           0 :         poDS->nProduct = eUnknown;
     590             : 
     591             :     /* Start reading in the product components */
     592           0 :     char *pszGeorefFile = nullptr;
     593           0 :     CPLErr geoTransformErr = CE_Failure;
     594           0 :     for (CPLXMLNode *psComponent = psComponents->psChild;
     595           0 :          psComponent != nullptr; psComponent = psComponent->psNext)
     596             :     {
     597           0 :         const char *pszType = nullptr;
     598             :         const std::string osPath =
     599           0 :             CPLFormFilenameSafe(CPLGetDirnameSafe(osFilename).c_str(),
     600           0 :                                 GetFilePath(psComponent, &pszType).c_str(), "");
     601           0 :         const char *pszPolLayer = CPLGetXMLValue(psComponent, "polLayer", " ");
     602             : 
     603           0 :         if (!STARTS_WITH_CI(pszType, " "))
     604             :         {
     605           0 :             if (STARTS_WITH_CI(pszType, "MAPPING_GRID"))
     606             :             {
     607             :                 /* the mapping grid... save as a metadata item this path */
     608           0 :                 poDS->SetMetadataItem("MAPPING_GRID", osPath.c_str());
     609             :             }
     610           0 :             else if (STARTS_WITH_CI(pszType, "GEOREF"))
     611             :             {
     612             :                 /* save the path to the georef data for later use */
     613           0 :                 CPLFree(pszGeorefFile);
     614           0 :                 pszGeorefFile = CPLStrdup(osPath.c_str());
     615             :             }
     616             :         }
     617           0 :         else if (!STARTS_WITH_CI(pszPolLayer, " ") &&
     618           0 :                  STARTS_WITH_CI(psComponent->pszValue, "imageData"))
     619             :         {
     620             :             /* determine the polarization of this band */
     621             :             ePolarization ePol;
     622           0 :             if (STARTS_WITH_CI(pszPolLayer, "HH"))
     623             :             {
     624           0 :                 ePol = HH;
     625             :             }
     626           0 :             else if (STARTS_WITH_CI(pszPolLayer, "HV"))
     627             :             {
     628           0 :                 ePol = HV;
     629             :             }
     630           0 :             else if (STARTS_WITH_CI(pszPolLayer, "VH"))
     631             :             {
     632           0 :                 ePol = VH;
     633             :             }
     634             :             else
     635             :             {
     636           0 :                 ePol = VV;
     637             :             }
     638             : 
     639           0 :             GDALDataType eDataType = STARTS_WITH_CI(pszDataType, "COMPLEX")
     640           0 :                                          ? GDT_CInt16
     641             :                                          : GDT_UInt16;
     642             : 
     643             :             /* try opening the file that represents that band */
     644             :             GDALDataset *poBandData =
     645           0 :                 GDALDataset::FromHandle(GDALOpen(osPath.c_str(), GA_ReadOnly));
     646           0 :             if (poBandData != nullptr)
     647             :             {
     648             :                 TSXRasterBand *poBand =
     649           0 :                     new TSXRasterBand(poDS, eDataType, ePol, poBandData);
     650           0 :                 poDS->SetBand(poDS->GetRasterCount() + 1, poBand);
     651             : 
     652             :                 // copy georeferencing info from the band
     653             :                 // need error checking??
     654             :                 // it will just save the info from the last band
     655           0 :                 const auto poSrcSRS = poBandData->GetSpatialRef();
     656           0 :                 if (poSrcSRS)
     657           0 :                     poDS->m_oSRS = *poSrcSRS;
     658             : 
     659             :                 geoTransformErr =
     660           0 :                     poBandData->GetGeoTransform(poDS->adfGeoTransform);
     661             :             }
     662             :         }
     663             :     }
     664             : 
     665             :     // now check if there is a geotransform
     666           0 :     if (!poDS->m_oSRS.IsEmpty() && geoTransformErr == CE_None)
     667             :     {
     668           0 :         poDS->bHaveGeoTransform = TRUE;
     669             :     }
     670             :     else
     671             :     {
     672           0 :         poDS->bHaveGeoTransform = FALSE;
     673           0 :         poDS->m_oSRS.Clear();
     674           0 :         poDS->adfGeoTransform[0] = 0.0;
     675           0 :         poDS->adfGeoTransform[1] = 1.0;
     676           0 :         poDS->adfGeoTransform[2] = 0.0;
     677           0 :         poDS->adfGeoTransform[3] = 0.0;
     678           0 :         poDS->adfGeoTransform[4] = 0.0;
     679           0 :         poDS->adfGeoTransform[5] = 1.0;
     680             :     }
     681             : 
     682           0 :     CPLFree(pszDataType);
     683             : 
     684             :     /* -------------------------------------------------------------------- */
     685             :     /*      Check and set matrix representation.                            */
     686             :     /* -------------------------------------------------------------------- */
     687             : 
     688           0 :     if (poDS->GetRasterCount() == 4)
     689             :     {
     690           0 :         poDS->SetMetadataItem("MATRIX_REPRESENTATION", "SCATTERING");
     691             :     }
     692             : 
     693             :     /* -------------------------------------------------------------------- */
     694             :     /*      Read the four corners and centre GCPs in                        */
     695             :     /* -------------------------------------------------------------------- */
     696             : 
     697             :     CPLXMLNode *psSceneInfo =
     698           0 :         CPLGetXMLNode(psData, "=level1Product.productInfo.sceneInfo");
     699           0 :     if (psSceneInfo != nullptr)
     700             :     {
     701             :         /* extract the GCPs from the provided file */
     702           0 :         bool success = false;
     703           0 :         if (pszGeorefFile != nullptr)
     704           0 :             success = poDS->getGCPsFromGEOREF_XML(pszGeorefFile);
     705             : 
     706             :         // if the gcp's cannot be extracted from the georef file, try to get the
     707             :         // corner coordinates for now just SSC because the others don't have
     708             :         // refColumn and refRow
     709           0 :         if (!success && poDS->nProduct == eSSC)
     710             :         {
     711           0 :             int nGCP = 0;
     712           0 :             double dfAvgHeight = CPLAtof(
     713             :                 CPLGetXMLValue(psSceneInfo, "sceneAverageHeight", "0.0"));
     714             : 
     715             :             // count and allocate gcps - there should be five - 4 corners and a
     716             :             // centre
     717           0 :             poDS->nGCPCount = 0;
     718           0 :             CPLXMLNode *psNode = psSceneInfo->psChild;
     719           0 :             for (; psNode != nullptr; psNode = psNode->psNext)
     720             :             {
     721           0 :                 if (!EQUAL(psNode->pszValue, "sceneCenterCoord") &&
     722           0 :                     !EQUAL(psNode->pszValue, "sceneCornerCoord"))
     723           0 :                     continue;
     724             : 
     725           0 :                 poDS->nGCPCount++;
     726             :             }
     727           0 :             if (poDS->nGCPCount > 0)
     728             :             {
     729           0 :                 poDS->pasGCPList =
     730           0 :                     (GDAL_GCP *)CPLCalloc(sizeof(GDAL_GCP), poDS->nGCPCount);
     731             : 
     732             :                 /* iterate over GCPs */
     733           0 :                 for (psNode = psSceneInfo->psChild; psNode != nullptr;
     734           0 :                      psNode = psNode->psNext)
     735             :                 {
     736           0 :                     GDAL_GCP *psGCP = poDS->pasGCPList + nGCP;
     737             : 
     738           0 :                     if (!EQUAL(psNode->pszValue, "sceneCenterCoord") &&
     739           0 :                         !EQUAL(psNode->pszValue, "sceneCornerCoord"))
     740           0 :                         continue;
     741             : 
     742           0 :                     psGCP->dfGCPPixel =
     743           0 :                         CPLAtof(CPLGetXMLValue(psNode, "refColumn", "0.0"));
     744           0 :                     psGCP->dfGCPLine =
     745           0 :                         CPLAtof(CPLGetXMLValue(psNode, "refRow", "0.0"));
     746           0 :                     psGCP->dfGCPX =
     747           0 :                         CPLAtof(CPLGetXMLValue(psNode, "lon", "0.0"));
     748           0 :                     psGCP->dfGCPY =
     749           0 :                         CPLAtof(CPLGetXMLValue(psNode, "lat", "0.0"));
     750           0 :                     psGCP->dfGCPZ = dfAvgHeight;
     751           0 :                     psGCP->pszId = CPLStrdup(CPLSPrintf("%d", nGCP));
     752           0 :                     psGCP->pszInfo = CPLStrdup("");
     753             : 
     754           0 :                     nGCP++;
     755             :                 }
     756             : 
     757             :                 // set the projection string - the fields are lat/long - seems
     758             :                 // to be WGS84 datum
     759           0 :                 poDS->m_oGCPSRS.SetWellKnownGeogCS("WGS84");
     760             :             }
     761             :         }
     762             : 
     763             :         // gcps override geotransform - does it make sense to have both??
     764           0 :         if (poDS->nGCPCount > 0)
     765             :         {
     766           0 :             poDS->bHaveGeoTransform = FALSE;
     767           0 :             poDS->m_oSRS.Clear();
     768           0 :             poDS->adfGeoTransform[0] = 0.0;
     769           0 :             poDS->adfGeoTransform[1] = 1.0;
     770           0 :             poDS->adfGeoTransform[2] = 0.0;
     771           0 :             poDS->adfGeoTransform[3] = 0.0;
     772           0 :             poDS->adfGeoTransform[4] = 0.0;
     773           0 :             poDS->adfGeoTransform[5] = 1.0;
     774             :         }
     775             :     }
     776             :     else
     777             :     {
     778           0 :         CPLError(CE_Warning, CPLE_AppDefined,
     779             :                  "Unable to find sceneInfo tag in XML document. "
     780             :                  "Proceeding with caution.");
     781             :     }
     782             : 
     783           0 :     CPLFree(pszGeorefFile);
     784             : 
     785             :     /* -------------------------------------------------------------------- */
     786             :     /*      Initialize any PAM information.                                 */
     787             :     /* -------------------------------------------------------------------- */
     788           0 :     poDS->SetDescription(poOpenInfo->pszFilename);
     789           0 :     poDS->TryLoadXML();
     790             : 
     791             :     /* -------------------------------------------------------------------- */
     792             :     /*      Check for overviews.                                            */
     793             :     /* -------------------------------------------------------------------- */
     794           0 :     poDS->oOvManager.Initialize(poDS, poOpenInfo->pszFilename);
     795             : 
     796           0 :     CPLDestroyXMLNode(psData);
     797             : 
     798           0 :     return poDS;
     799             : }
     800             : 
     801             : /************************************************************************/
     802             : /*                            GetGCPCount()                             */
     803             : /************************************************************************/
     804             : 
     805           0 : int TSXDataset::GetGCPCount()
     806             : {
     807           0 :     return nGCPCount;
     808             : }
     809             : 
     810             : /************************************************************************/
     811             : /*                          GetGCPSpatialRef()                          */
     812             : /************************************************************************/
     813             : 
     814           0 : const OGRSpatialReference *TSXDataset::GetGCPSpatialRef() const
     815             : {
     816           0 :     return m_oGCPSRS.IsEmpty() ? nullptr : &m_oGCPSRS;
     817             : }
     818             : 
     819             : /************************************************************************/
     820             : /*                               GetGCPs()                              */
     821             : /************************************************************************/
     822             : 
     823           0 : const GDAL_GCP *TSXDataset::GetGCPs()
     824             : {
     825           0 :     return pasGCPList;
     826             : }
     827             : 
     828             : /************************************************************************/
     829             : /*                          GetSpatialRef()                             */
     830             : /************************************************************************/
     831             : 
     832           0 : const OGRSpatialReference *TSXDataset::GetSpatialRef() const
     833             : 
     834             : {
     835           0 :     return m_oSRS.IsEmpty() ? nullptr : &m_oSRS;
     836             : }
     837             : 
     838             : /************************************************************************/
     839             : /*                               GetGeotransform()                      */
     840             : /************************************************************************/
     841           0 : CPLErr TSXDataset::GetGeoTransform(double *padfTransform)
     842             : {
     843           0 :     memcpy(padfTransform, adfGeoTransform, sizeof(double) * 6);
     844             : 
     845           0 :     if (bHaveGeoTransform)
     846           0 :         return CE_None;
     847             : 
     848           0 :     return CE_Failure;
     849             : }
     850             : 
     851             : /************************************************************************/
     852             : /*                         GDALRegister_TSX()                           */
     853             : /************************************************************************/
     854             : 
     855        1682 : void GDALRegister_TSX()
     856             : {
     857        1682 :     if (GDALGetDriverByName("TSX") != nullptr)
     858         301 :         return;
     859             : 
     860        1381 :     GDALDriver *poDriver = new GDALDriver();
     861             : 
     862        1381 :     poDriver->SetDescription("TSX");
     863        1381 :     poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
     864        1381 :     poDriver->SetMetadataItem(GDAL_DMD_LONGNAME, "TerraSAR-X Product");
     865        1381 :     poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/tsx.html");
     866        1381 :     poDriver->SetMetadataItem(GDAL_DCAP_VIRTUALIO, "YES");
     867             : 
     868        1381 :     poDriver->pfnOpen = TSXDataset::Open;
     869        1381 :     poDriver->pfnIdentify = TSXDataset::Identify;
     870             : 
     871        1381 :     GetGDALDriverManager()->RegisterDriver(poDriver);
     872             : }