ENH: Set Sentinel1 IMI with new framwork

Data/Baseline/OTB/Files/ioTvImageMetadataInterfaceTest_Sentinel1.txt

@@ -11,12 +11,6 @@ TileHintX 17663
 TileHintY 1
 DataType 11
 OrbitNumber 6447
-NumberOfLines 31106
-NumberOfColumns 17663
-AverageSceneHeight 19
-CalScale 1
-PRF 1663.48
-RadarFrequency 5.405e+09
 LineSpacing 4.08568
 PixelSpacing 3.19483
 AcquisitionDate 2015-06-19T19:50:43.223221Z
@@ -27,6 +21,12 @@ GCP <GCPParam>
 Extra.FACILITY_IDENTIFIER UPA_
 Polarization VV
 Swath S6
+NumberOfLines 31106
+NumberOfColumns 17663
+AverageSceneHeight 19
+CalScale 1
+PRF 1663.48
+RadarFrequency 5.405e+09
 SAR <SARParam>
 {"Projection": "GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]", 
 [{"GCP_Id": "1", "GCP_Info": "", "GCP_Row": "0", "GCP_Col": "0", "GCP_X": "-24.4551", "GCP_Y": "14.4996", "GCP_Z": "0", }, 

Modules/Core/Metadata/include/otbCosmoImageMetadataInterface.h

@@ -92,6 +92,10 @@ public:
 
   void Parse(const MetadataSupplierInterface &) override;
 
+  void ParseGdal(const MetadataSupplierInterface &);
+  
+  void ParseGeom(const MetadataSupplierInterface &) override;
+
   std::vector<std::map<std::string, std::string> > saveMetadataBands(std::string file) ;
   std::vector<Orbit> getOrbits(const MetadataSupplierInterface &, const std::string & referenceTime);
 

Modules/Core/Metadata/include/otbRadarsat2ImageMetadataInterface.h

@@ -98,6 +98,11 @@ public:
   /*get lookup data for calculating backscatter */
   void CreateCalibrationLookupData(const short type) override;
 
+  virtual void ParseGdal(const MetadataSupplierInterface &) override
+  {
+    itkExceptionMacro("ParseGdal not implemented in Radarsat2ImageMetadataInterface.");
+  }
+
   void Parse(const MetadataSupplierInterface &) override;
 
 

Modules/Core/Metadata/include/otbSarDefaultImageMetadataInterface.h

@@ -205,6 +205,11 @@ public:
     itkExceptionMacro("GetCenterIncidenceAngle not implemented in SarDefaultImageMetadataInterface, no captor type found");
   }
 
+  virtual void ParseGdal(const MetadataSupplierInterface &) override
+  {
+    itkExceptionMacro("ParseGdal not implemented in SarDefaultImageMetadataInterface, no captor type found");
+  }
+
   bool CanRead() const override
   {
     // This class is the default one, it has to be able to call every metadata

Modules/Core/Metadata/include/otbSarImageMetadataInterface.h

@@ -112,7 +112,11 @@ public:
    * then translate it into ImageMetadata. Handles most SAR sensors.
    * Returns true if succeed. */
   bool ConvertImageKeywordlistToImageMetadata() override;
-  
+
+  virtual void ParseGdal(const MetadataSupplierInterface &) =0;
+
+  virtual void ParseGeom(const MetadataSupplierInterface &);
+
 protected:
   SarImageMetadataInterface();
   ~SarImageMetadataInterface() override

Modules/Core/Metadata/include/otbSentinel1ImageMetadataInterface.h

@@ -100,6 +100,10 @@ public:
   /*get lookup data for calculating backscatter */
   void CreateCalibrationLookupData(const short type) override;
 
+  void ParseGdal(const MetadataSupplierInterface &);
+
+  void ParseGeom(const MetadataSupplierInterface &);
+
   void Parse(const MetadataSupplierInterface &) override;
 
 protected:
@@ -114,15 +118,19 @@ protected:
 
   /* Fetch the DopplerCentroid metadata */
   std::vector<DopplerCentroid> GetDopplerCentroid(const XMLMetadataSupplier&) const;
+  std::vector<DopplerCentroid> GetDopplerCentroidGeom(const MetadataSupplierInterface&) const;
 
   /* Fetch the Orbits metadata */
   std::vector<Orbit> GetOrbits(const XMLMetadataSupplier&) const;
+  std::vector<Orbit> GetOrbitsGeom(const MetadataSupplierInterface&) const;
 
   /* Fetch the Calibration metadata */
   std::vector<CalibrationVector> GetCalibrationVector(const XMLMetadataSupplier&) const;
+  std::vector<CalibrationVector> GetCalibrationVectorGeom(const MetadataSupplierInterface&) const;
 
   /* Fetch the noise LUTs */
   std::vector<SARNoise> GetNoiseVector(const XMLMetadataSupplier&) const;
+  std::vector<SARNoise> GetNoiseVectorGeom(const MetadataSupplierInterface&) const;
 
   /* Compute the mean terrain elevation */
   double getBandTerrainHeight(const XMLMetadataSupplier&) const;

Modules/Core/Metadata/include/otbTerraSarImageMetadataInterface.h

@@ -174,6 +174,11 @@ public:
    *  in the order R, G, B */
   std::vector<unsigned int> GetDefaultDisplay() const override;
 
+  virtual void ParseGdal(const MetadataSupplierInterface &) override
+  {
+    itkExceptionMacro("ParseGdal not implemented in TerraSarImageMetadataInterface.");
+  }
+
   void Parse(const MetadataSupplierInterface &) override;
 
 protected:

Modules/Core/Metadata/src/otbCosmoImageMetadataInterface.cxx

@@ -255,8 +255,6 @@ double CosmoImageMetadataInterface::GetCenterIncidenceAngle() const
   return 0;
 }
 
-
-
 std::vector<std::map<std::string, std::string> > CosmoImageMetadataInterface::saveMetadataBands(std::string file)
 {
   // Create GDALImageIO to retrieve all metadata (from .h5 input file)
@@ -298,9 +296,6 @@ std::vector<std::map<std::string, std::string> > CosmoImageMetadataInterface::sa
 
 }
 
-
-
-
 std::vector<Orbit> CosmoImageMetadataInterface::getOrbits(const MetadataSupplierInterface & mds, const std::string & reference_UTC)
 {
     ////////////////// Add Orbit List ////////////////
@@ -353,19 +348,9 @@ std::vector<Orbit> CosmoImageMetadataInterface::getOrbits(const MetadataSupplier
   return orbitVector ;
 }
 
-
-
-void CosmoImageMetadataInterface::Parse(const MetadataSupplierInterface & mds)
+void CosmoImageMetadataInterface::ParseGdal(const MetadataSupplierInterface & mds)
 {
-  // Check Mission Id, acquisition mode and product type
-  Fetch(MDStr::Mission, mds, "MISSION_ID");
-
-  if (m_Imd[MDStr::Mission] != "CSK" )
-  {
-    otbGenericExceptionMacro(MissingMetadataException,
-      << "Not a CosmoSkyMed product");
-  }
-
+  // Check acquisition mode and product type
   Fetch(MDStr::Mode, mds, "Acquisition_Mode");
   if((m_Imd[MDStr::Mode] != "HIMAGE") &&
         (m_Imd[MDStr::Mode] != "SPOTLIGHT") &&
@@ -431,11 +416,23 @@ void CosmoImageMetadataInterface::Parse(const MetadataSupplierInterface & mds)
   m_Imd.Bands[0].Add(MDGeom::SAR, sarParam);
 }
 
+void CosmoImageMetadataInterface::ParseGeom(const MetadataSupplierInterface & mds)
+{
+  Superclass::ParseGeom(mds);
+}
 
-
-
-
-
-
+void CosmoImageMetadataInterface::Parse(const MetadataSupplierInterface & mds)
+{
+  // Try to fetch the metadata from GDAL Metadata Supplier
+  if (mds.GetAs<std::string>("", "MISSION_ID") == "CSK")
+    this->ParseGdal(mds);
+  // Try to fetch the metadata from GEOM file
+  else if (mds.GetAs<std::string>("", "support_data.sensor") == "CSK")
+    this->ParseGeom(mds);
+  // Failed to fetch the metadata
+  else
+    otbGenericExceptionMacro(MissingMetadataException,
+			     << "Not a CosmoSkyMed product");
+}
 
 } // end namespace otb

Modules/Core/Metadata/src/otbSarImageMetadataInterface.cxx

@@ -24,6 +24,7 @@
 #include "otbMath.h"
 #include "itkMetaDataObject.h"
 #include "otbImageKeywordlist.h"
+#include "otbSARMetadata.h"
 
 namespace otb
 {
@@ -189,6 +190,22 @@ SarImageMetadataInterface::IndexType SarImageMetadataInterface::GetRadiometricCa
 }
 
 
+void SarImageMetadataInterface::ParseGeom(const MetadataSupplierInterface & mds)
+{
+  Fetch(MDTime::AcquisitionStartTime, mds, "support_data.first_line_time");
+  Fetch(MDTime::AcquisitionStopTime, mds, "support_data.last_line_time");
+  Fetch(MDStr::BeamMode, mds, "manifest_data.acquisition_mode");
+  Fetch(MDStr::BeamSwath, mds, "manifest_data.swath");
+  Fetch(MDNum::LineSpacing, mds, "support_data.azimuth_spacing");
+  Fetch(MDStr::Mission, mds, "manifest_data.instrument");
+  Fetch(MDStr::OrbitDirection, mds, "manifest_data.orbit_pass");
+  Fetch(MDNum::OrbitNumber, mds, "manifest_data.abs_orbit");
+  Fetch(MDNum::PixelSpacing, mds, "support_data.range_spacing");
+  Fetch(MDStr::ProductType, mds, "manifest_data.product_type");
+  Fetch(MDTime::ProductionDate, mds, "manifest_data.date");
+  Fetch(MDTime::AcquisitionDate, mds, "manifest_data.image_date");
+}
+
 void SarImageMetadataInterface::PrintSelf(std::ostream& os, itk::Indent indent) const
 {
   Superclass::PrintSelf(os, indent);

Modules/Core/Metadata/src/otbSentinel1ImageMetadataInterface.cxx

@@ -392,10 +392,30 @@ std::vector<DopplerCentroid> Sentinel1ImageMetadataInterface::GetDopplerCentroid
     std::istringstream(xmlMS.GetAs<std::string>(path_root + ".azimuthTime")) >> dopplerCent.azimuthTime;
     dopplerCent.t0 = xmlMS.GetAs<double>(path_root + ".t0");
     dopplerCent.dopCoef = xmlMS.GetAsVector<double>(path_root + ".dataDcPolynomial",
-    		' ', xmlMS.GetAs<int>(path_root + ".dataDcPolynomial.count"));
-	dopplerCent.geoDopCoef = xmlMS.GetAsVector<double>(path_root + ".geometryDcPolynomial",
-    		' ', xmlMS.GetAs<int>(path_root + ".geometryDcPolynomial.count"));
-	dopplerCentroidVector.push_back(dopplerCent);
+						    ' ', xmlMS.GetAs<int>(path_root + ".dataDcPolynomial.count"));
+    dopplerCent.geoDopCoef = xmlMS.GetAsVector<double>(path_root + ".geometryDcPolynomial",
+						       ' ', xmlMS.GetAs<int>(path_root + ".geometryDcPolynomial.count"));
+    dopplerCentroidVector.push_back(dopplerCent);
+  }
+  return dopplerCentroidVector;
+}
+
+std::vector<DopplerCentroid> Sentinel1ImageMetadataInterface::GetDopplerCentroidGeom(const MetadataSupplierInterface & mds) const
+{
+  std::vector<DopplerCentroid> dopplerCentroidVector;
+  // Path: dopplerCentroid.dop_coef_list<listId>.{dop_coef_time,slant_range_time,{1,2,3}.dop_coef}
+  // This streams wild hold the iteration number
+  std::ostringstream oss;
+  for (int listId = 1 ; mds.GetAs<std::string>("", std::string("dopplerCentroid.dop_coef_list")+std::to_string(listId)+std::string(".slant_range_time")) != "" ; ++listId)
+  {
+    oss.str("");
+    oss << listId;
+    // Base path to the data, that depends on the iteration number
+    std::string path_root = "dopplerCentroid.dop_coef_list" + oss.str();
+    DopplerCentroid dopplerCent;
+    std::istringstream(mds.GetAs<std::string>(path_root + ".dop_coef_time")) >> dopplerCent.azimuthTime;
+    dopplerCent.t0 = mds.GetAs<double>(path_root + ".slant_range_time");
+    dopplerCentroidVector.push_back(dopplerCent);
   }
   return dopplerCentroidVector;
 }
@@ -426,6 +446,32 @@ std::vector<Orbit> Sentinel1ImageMetadataInterface::GetOrbits(const XMLMetadataS
   return orbitVector;
 }
 
+std::vector<Orbit> Sentinel1ImageMetadataInterface::GetOrbitsGeom(const MetadataSupplierInterface & mds) const
+{
+  std::vector<Orbit> orbitVector;
+  // Number of entries in the vector
+  int listCount = mds.GetAs<int>("orbitList.nb_orbits");
+  // This streams wild hold the iteration number
+  std::ostringstream oss;
+  for (int listId = 0 ; listId <= listCount - 1 ; ++listId)
+  {
+    oss.str("");
+    oss << listId;
+    // Base path to the data, that depends on the iteration number
+    std::string path_root = "orbitList.orbit[" + oss.str() + "]";
+    Orbit orbit;
+    std::istringstream(mds.GetAs<std::string>(path_root + ".time")) >> orbit.time;
+    orbit.posX = mds.GetAs<double>(path_root + ".x_pos");
+    orbit.posY = mds.GetAs<double>(path_root + ".y_pos");
+    orbit.posZ = mds.GetAs<double>(path_root + ".z_pos");
+    orbit.velX = mds.GetAs<double>(path_root + ".x_vel");
+    orbit.velY = mds.GetAs<double>(path_root + ".y_vel");
+    orbit.velZ = mds.GetAs<double>(path_root + ".z_vel");
+    orbitVector.push_back(orbit);
+  }
+  return orbitVector;
+}
+
 std::vector<CalibrationVector> Sentinel1ImageMetadataInterface::GetCalibrationVector(const XMLMetadataSupplier &xmlMS) const
 {
   std::vector<CalibrationVector> calibrationVector;
@@ -478,6 +524,54 @@ std::vector<CalibrationVector> Sentinel1ImageMetadataInterface::GetCalibrationVe
   return calibrationVector;
 }
 
+std::vector<CalibrationVector> Sentinel1ImageMetadataInterface::GetCalibrationVectorGeom(const MetadataSupplierInterface & mds) const
+{
+  std::vector<CalibrationVector> calibrationVector;
+  // Number of entries in the vector
+  int listCount = mds.GetAs<int>("calibration.count");
+  // This streams wild hold the iteration number
+  std::ostringstream oss;
+  for (int listId = 0 ; listId <= listCount - 1 ; ++listId)
+  {
+    oss.str("");
+    oss << listId;
+    // Base path to the data, that depends on the iteration number
+    std::string path_root = "calibration.calibrationVector[" + oss.str() + "]";
+
+    CalibrationVector calVect;
+    std::istringstream(mds.GetAs<std::string>(path_root + ".azimuthTime")) >> calVect.azimuthTime;
+    calVect.line = mds.GetAs<int>(path_root + ".line");
+
+    // Same axe for all LUTs
+    MetaData::LUTAxis ax1;
+    ax1.Size = mds.GetAs<int>(path_root + ".pixel_count");
+    ax1.Values = mds.GetAsVector<double>(path_root + ".pixel", ' ', ax1.Size);
+
+    MetaData::LUT1D sigmaNoughtLut;
+    sigmaNoughtLut.Axis[0] = ax1;
+    sigmaNoughtLut.Array = mds.GetAsVector<double>(path_root + ".sigmaNought", ' ', ax1.Size);
+    calVect.sigmaNought = sigmaNoughtLut;
+
+    MetaData::LUT1D betaNoughtLut;
+    betaNoughtLut.Axis[0] = ax1;
+    betaNoughtLut.Array = mds.GetAsVector<double>(path_root + ".betaNought", ' ', ax1.Size);
+    calVect.betaNought = betaNoughtLut;
+
+    MetaData::LUT1D gammaLut;
+    gammaLut.Axis[0] = ax1;
+    gammaLut.Array = mds.GetAsVector<double>(path_root + ".gamma", ' ', ax1.Size);
+    calVect.gamma = gammaLut;
+
+    MetaData::LUT1D dnLut;
+    dnLut.Axis[0] = ax1;
+    dnLut.Array = mds.GetAsVector<double>(path_root + ".dn", ' ', ax1.Size);
+    calVect.dn = dnLut;
+
+    calibrationVector.push_back(calVect);
+  }
+  return calibrationVector;
+}
+
 std::vector<SARNoise> Sentinel1ImageMetadataInterface::GetNoiseVector(const XMLMetadataSupplier &xmlMS) const
 {
   std::vector<SARNoise> noiseVector;
@@ -506,6 +600,32 @@ std::vector<SARNoise> Sentinel1ImageMetadataInterface::GetNoiseVector(const XMLM
   return noiseVector;
 }
 
+std::vector<SARNoise> Sentinel1ImageMetadataInterface::GetNoiseVectorGeom(const MetadataSupplierInterface & mds) const
+{
+  std::vector<SARNoise> noiseVector;
+  // This streams wild hold the iteration number
+  std::ostringstream oss;
+  // Path: noise.noiseVector[<listId>].{azimuthTime,line,noiseLut,pixel,pixel_count}
+  for (int listId = 0 ; mds.GetAs<std::string>("", std::string("noise.noiseVector[")+std::to_string(listId)+std::string("].pixel_count")) != "" ; ++listId)
+  {
+    oss.str("");
+    oss << listId;
+    // Base path to the data, that depends on the iteration number
+    std::string path_root = "noise.noiseVector[" + oss.str() + "]";
+    SARNoise noiseVect;
+    std::istringstream(mds.GetAs<std::string>(path_root + ".azimuthTime")) >> noiseVect.azimuthTime;
+    MetaData::LUT1D noiseLut;
+    MetaData::LUTAxis ax1;
+    ax1.Size = mds.GetAs<int>(path_root + ".pixel_count");
+    ax1.Values = mds.GetAsVector<double>(path_root + ".pixel", ' ', ax1.Size);
+    noiseLut.Axis[0] = ax1;
+    noiseLut.Array = mds.GetAsVector<double>(path_root + ".noiseLut", ' ', ax1.Size);
+    noiseVect.noiseLut = noiseLut;
+    noiseVector.push_back(noiseVect);
+  }
+  return noiseVector;
+}
+
 double Sentinel1ImageMetadataInterface::getBandTerrainHeight(const XMLMetadataSupplier &xmlMS) const
 {
   double heightSum = 0.0;
@@ -524,7 +644,7 @@ double Sentinel1ImageMetadataInterface::getBandTerrainHeight(const XMLMetadataSu
   return heightSum / (double)listCount;
 }
 
-void Sentinel1ImageMetadataInterface::Parse(const MetadataSupplierInterface & mds)
+void Sentinel1ImageMetadataInterface::ParseGdal(const MetadataSupplierInterface & mds)
 {
   // Metadata read by GDAL
   Fetch(MDTime::AcquisitionStartTime, mds, "ACQUISITION_START_TIME");
@@ -578,7 +698,6 @@ void Sentinel1ImageMetadataInterface::Parse(const MetadataSupplierInterface & md
     m_Imd.Bands[bandId].Add(MDNum::AverageSceneHeight, this->getBandTerrainHeight(AnnotationFilePath));
     m_Imd.Bands[bandId].Add(MDNum::RadarFrequency, AnnotationMS.GetAs<double>("product.generalAnnotation.productInformation.radarFrequency"));
     m_Imd.Bands[bandId].Add(MDNum::PRF, AnnotationMS.GetAs<double>("product.imageAnnotation.imageInformation.azimuthFrequency"));
-    m_Imd.Bands[bandId].Add(MDNum::CenterIncidenceAngle, AnnotationMS.GetAs<double>("product.imageAnnotation.imageInformation.incidenceAngleMidSwath"));
 
     // Calibration file
     std::string CalibrationFilePath = itksys::SystemTools::GetFilenamePath(AnnotationFilePath)
@@ -605,4 +724,46 @@ void Sentinel1ImageMetadataInterface::Parse(const MetadataSupplierInterface & md
   }
 }
 
+void Sentinel1ImageMetadataInterface::ParseGeom(const MetadataSupplierInterface & mds)
+{
+  Superclass::ParseGeom(mds);
+
+  SARParam sarParam;
+  Fetch("FACILITY_IDENTIFIER", mds, "manifest_data.Processing_system_identifier");
+  // TODO: sarParam.azimuthFmRates = ???
+  sarParam.dopplerCentroids = this->GetDopplerCentroidGeom(mds);
+  sarParam.orbits = this->GetOrbitsGeom(mds);
+  m_Imd.Add(MDStr::SensorID, "SAR");
+  auto instrument = mds.GetAs<std::string>("sensor");
+  instrument.pop_back();
+  m_Imd.Add(MDStr::Instrument, instrument);
+  Fetch(MDStr::Mode, mds, "header.swath");
+  Fetch(MDStr::Swath, mds, "manifest_data.swath");
+  Fetch(MDNum::NumberOfLines, mds, "number_lines", 0);
+  Fetch(MDNum::NumberOfColumns, mds, "number_samples", 0);
+  Fetch(MDNum::AverageSceneHeight, mds, "support_data.avg_scene_height", 0);
+  Fetch(MDNum::RadarFrequency, mds, "support_data.radar_frequency", 0);
+  Fetch(MDNum::PRF, mds, "support_data.pulse_repetition_frequency", 0);
+  Fetch(MDNum::CalScale, mds, "calibration.absoluteCalibrationConstant", 0);
+  sarParam.calibrationVectors = this->GetCalibrationVectorGeom(mds);
+  std::istringstream(mds.GetAs<std::string>("calibration.startTime")) >> sarParam.calibrationStartTime;
+  std::istringstream(mds.GetAs<std::string>("calibration.stopTime")) >> sarParam.calibrationStopTime;
+  sarParam.noiseVector = this->GetNoiseVectorGeom(mds);
+  m_Imd.Bands[0].Add(MDGeom::SAR, sarParam);
+}
+
+void Sentinel1ImageMetadataInterface::Parse(const MetadataSupplierInterface & mds)
+{
+  // Try to fetch the metadata from GDAL Metadata Supplier
+  if (mds.GetAs<std::string>("", "MISSION_ID") == "S1A" || mds.GetAs<std::string>("", "MISSION_ID") == "S1B")
+    this->ParseGdal(mds);
+  // Try to fetch the metadata from GEOM file
+  else if (mds.GetAs<std::string>("", "sensor") == "SENTINEL-1A" || mds.GetAs<std::string>("", "sensor") == "SENTINEL-1B")
+    this->ParseGeom(mds);
+  // Failed to fetch the metadata
+  else
+    otbGenericExceptionMacro(MissingMetadataException,
+			     << "Not a Sentinel1 product");
+}
+
 } // end namespace otb