Commit 630f4fe7 authored by Cédric Traizet's avatar Cédric Traizet

ENH: parse geom metadata for SPOT5 and Formosat

parent 48e7d53e
Pipeline #6121 passed with stages
in 21 minutes and 29 seconds
SensorID FORMOSAT 2
Mission Formasat 2
GeometricLevel 1A
TileHintX 3000
TileHintY 500
DataType 1
......@@ -8,7 +7,7 @@ SunElevation 37.4423
SunAzimuth 160.749
SatElevation 57.9678
SatAzimuth 241.823
AcquisitionDate 2007-10-13T10:40:33.644225Z
AcquisitionDate 2007-10-13T10:40:31.793824Z
ProductionDate 2007-10-15T17:58:07Z
GCP <GCPParam>
BandName B1
......
......@@ -7,7 +7,7 @@ SunElevation 40.6549
SunAzimuth 156.82
SatElevation 82.7227
SatAzimuth -74.1736
AcquisitionDate 2009-03-16T11:03:49Z
AcquisitionDate 2009-03-16T11:03:49.972847Z
ProductionDate 2009-03-26T11:56:45Z
BandName XS1
PhysicalGain 2.24793
......
......@@ -7,7 +7,7 @@ SunElevation 67.9004
SunAzimuth 128.882
SatElevation 95.7372
SatAzimuth 103.464
AcquisitionDate 2003-07-23T07:27:43Z
AcquisitionDate 2003-07-23T07:27:43.889716Z
ProductionDate 2003-12-05T08:21:19Z
GCP <GCPParam>
BandName XS1
......
......@@ -24,7 +24,6 @@
#include "OTBMetadataExport.h"
#include "otbMetadataSupplierInterface.h"
namespace otb
{
......@@ -93,6 +92,71 @@ public:
return m_Data;
}
void ParseGeom(const MetadataSupplierInterface & mds)
{
// Read a vector value from
auto readVector = [&mds](const std::string & path)
{
auto strMetadata = mds.GetAs<std::string>(path);
strMetadata.erase(std::remove(strMetadata.begin(), strMetadata.end(), '\"'), strMetadata.end());
std::istringstream is( strMetadata );
return std::vector<double>(( std::istream_iterator<double>( is ) ),
( std::istream_iterator<double>() ) );
};
m_Data.mission = mds.GetAs<std::string>("sensor");
m_Data.ImageID = mds.GetAs<std::string>("support_data.image_id");
m_Data.ProductionDate = mds.GetAs<std::string>("support_data.production_date");
// Convert YYYY-MM-DD HH:mm:ss to YYYY-MM-DDTHH:mm:ss
auto pos = m_Data.ProductionDate.find(" ",10);
if (pos != std::string::npos)
{
m_Data.ProductionDate.replace(pos,1,"T");
}
m_Data.AcquisitionDate = mds.GetAs<std::string>("support_data.image_date");
pos = m_Data.AcquisitionDate.find(" ",10);
if (pos != std::string::npos)
{
m_Data.AcquisitionDate.replace(pos,1,"T");
}
m_Data.SunAzimuth = {mds.GetAs<double>("support_data.azimuth_angle")};
m_Data.SunElevation = {mds.GetAs<double>("support_data.elevation_angle")};
m_Data.IncidenceAngle = {mds.GetAs<double>("support_data.incident_angle") };
m_Data.SceneOrientation = {mds.GetAs<double>("support_data.scene_orientation")};
m_Data.StepCount = mds.GetAs<int>(0, "support_data.step_count");
m_Data.PhysicalGain = readVector("support_data.physical_gain");
m_Data.PhysicalBias = readVector("support_data.physical_bias");
m_Data.SolarIrradiance = readVector("support_data.solar_irradiance");
m_Data.softwareVersion = mds.GetAs<std::string>("", "support_data.software_version");
m_Data.SatAzimuth = mds.GetAs<double>(0, "support_data.sat_azimuth_angle");
auto acrossTrackViewingAngle = mds.GetAs<std::string>("", "viewing_angle_across_track");
auto alongTrackViewingAngle = mds.GetAs<std::string>("", "viewing_angle_along_track");
if (!acrossTrackViewingAngle.empty())
{
m_Data.AcrossTrackViewingAngle.push_back(std::stod("acrossTrackViewingAngle"));
}
if (!alongTrackViewingAngle.empty())
{
m_Data.AlongTrackViewingAngle.push_back(std::stod("alongTrackViewingAngle"));
}
}
void ParseDimapV1(const MetadataSupplierInterface & mds, const std::string prefix)
{
std::vector<double> defaultValue = {};
......@@ -171,6 +235,7 @@ public:
m_Data.ProductionDate = mds.GetAs<std::string>(path);
// Convert YYYY-MM-DD HH:mm:ss to YYYY-MM-DDTHH:mm:ss
auto pos = m_Data.ProductionDate.find(" ",10);
if (pos != std::string::npos)
{
......@@ -179,9 +244,19 @@ public:
std::vector<std::string> imagingDateVec;
auto imagingDate = GetSingleValueFromList<std::string>(mds, prefix + "Dataset_Sources.Source_Information", "Scene_Source.IMAGING_DATE" );
auto imagingTime = GetSingleValueFromList<std::string>(mds, prefix + "Dataset_Sources.Source_Information", "Scene_Source.IMAGING_TIME" );
m_Data.AcquisitionDate = imagingDate + "T" + imagingTime;
m_Data.AcquisitionDate = mds.GetAs<std::string>("", prefix + "Data_Strip.Sensor_Configuration.Time_Stamp.SCENE_CENTER_TIME");
// If SCENE_CENTER_TIME is not present, try /Dimap_Document/Data_Strip/Time_Stamp/REFERENCE_TIME instead (e.g. Formosat)
if (m_Data.AcquisitionDate.empty())
{
m_Data.AcquisitionDate = mds.GetAs<std::string>("", prefix + "Data_Strip.Time_Stamp.REFERENCE_TIME");
pos = m_Data.AcquisitionDate.find(" ",10);
if (pos != std::string::npos)
{
m_Data.AcquisitionDate.replace(pos,1,"T");
}
}
m_Data.Instrument = GetSingleValueFromList<std::string>(mds, prefix + "Dataset_Sources.Source_Information", "Scene_Source.INSTRUMENT" );
m_Data.InstrumentIndex = GetSingleValueFromList<std::string>(mds, prefix + "Dataset_Sources.Source_Information", "Scene_Source.INSTRUMENT_INDEX" );
......
......@@ -41,7 +41,7 @@ class OTBMetadata_EXPORT GeomMetadataSupplier
: public MetadataSupplierInterface
{
public:
using DictType = std::map<std::string, std::string>;
using DictType = std::unordered_map<std::string, std::string>;
GeomMetadataSupplier(const std::string &);
GeomMetadataSupplier(const GeomMetadataSupplier &) = delete;
......
......@@ -1373,44 +1373,49 @@ void FormosatImageMetadataInterface::FetchSpectralSensitivity()
void FormosatImageMetadataInterface::Parse(const MetadataSupplierInterface & mds)
{
assert(mds);
DimapMetadataHelper helper;
Fetch(MDStr::SensorID, mds, "IMAGERY/SATELLITEID");
if (m_Imd[MDStr::SensorID] == "FORMOSAT 2")
{
m_Imd.Add(MDStr::Mission, "Formasat 2");
}
else
// GDAL DIMAP metadata case
if (mds.GetAs<std::string>("", "IMAGERY/SATELLITEID") == "FORMOSAT 2")
{
otbGenericExceptionMacro(MissingMetadataException,<<"Not a Formosat product")
}
//Find the METADATA.DIM file (Dimap V1 metadata file)
auto resourceFiles = mds.GetResourceFiles();
std::string metadataFile;
for (const auto & file: resourceFiles)
{
if (file.find("METADATA.DIM")!=std::string::npos)
{
metadataFile = file;
break;
}
}
//Find the METADATA.DIM file (Dimap V1 metadata file)
auto resourceFiles = mds.GetResourceFiles();
std::string metadataFile;
for (const auto & file: resourceFiles)
{
if (file.find("METADATA.DIM")!=std::string::npos)
if (metadataFile.empty())
{
metadataFile = file;
break;
otbGenericExceptionMacro(MissingMetadataException,<<"Cannot find the METADATA.DIM file")
}
XMLMetadataSupplier xmlMds(metadataFile);
helper.ParseDimapV1(xmlMds, "Dimap_Document.");
}
// geom metadata case
else if (mds.GetAs<std::string>("", "support_data.sensorID" ) == "Formosat 2")
{
helper.ParseGeom(mds);
}
if (metadataFile.empty())
else
{
otbGenericExceptionMacro(MissingMetadataException,<<"Cannot find the METADATA.DIM file")
otbGenericExceptionMacro(MissingMetadataException,<<"Not a Formosat product")
}
XMLMetadataSupplier xmlMds(metadataFile);
m_Imd.Add(MDStr::Mission, "Formasat 2");
m_Imd.Add(MDStr::SensorID, "FORMOSAT 2");
DimapMetadataHelper helper;
helper.ParseDimapV1(xmlMds, "Dimap_Document.");
const auto & dimapData = helper.GetDimapData();
auto nbBands = m_Imd.Bands.size();
// Band names are not in the metadata. Two cases should be considered
......@@ -1481,15 +1486,8 @@ void FormosatImageMetadataInterface::Parse(const MetadataSupplierInterface & mds
<< "The number of bands in image metadatas is incoherent with the DIMAP product")
}
m_Imd.Add(MDStr::GeometricLevel, dimapData.ProcessingLevel);
FetchSpectralSensitivity();
// fill RPC model
if (m_Imd[MDStr::GeometricLevel] == "SENSOR")
{
FetchRPC(mds);
}
}
......
......@@ -138,10 +138,15 @@ void GeomMetadataSupplier::ReadGeomFile()
std::vector< std::string > keyVal;
while (std::getline(inputFile, line))
{
boost::split(keyVal, line, [](char c){return c == ':';});
boost::trim(keyVal[0]);
boost::trim(keyVal[1]);
this->m_MetadataDic[keyVal[0]] = keyVal[1];
auto pos = line.find(":");
if (pos != std::string::npos)
{
auto key = line.substr(0,pos);
auto value = line.substr(pos+1, line.size());
boost::trim(key);
boost::trim(value);
m_MetadataDic[key] = value;
}
}
}
......
......@@ -1745,41 +1745,49 @@ void SpotImageMetadataInterface::FetchSpectralSensitivity()
void SpotImageMetadataInterface::Parse(const MetadataSupplierInterface & mds)
{
Fetch(MDStr::SensorID, mds, "IMAGERY/SATELLITEID");
if (strncmp(m_Imd[MDStr::SensorID].c_str(), "SPOT 5", 6) == 0)
{
m_Imd.Add(MDStr::Mission, "SPOT 5");
}
else
{
otbGenericExceptionMacro(MissingMetadataException,<<"Not a Spot 5 product")
}
DimapMetadataHelper helper;
//Find the METADATA.DIM file (Dimap V1 metadata file)
auto resourceFiles = mds.GetResourceFiles();
std::string metadataFile;
for (const auto & file: resourceFiles)
// GDAL DIMAP metadata case
if (mds.GetAs<std::string>("", "IMAGERY/SATELLITEID") == "SPOT 5")
{
if (file.find("METADATA.DIM")!=std::string::npos)
//Find the METADATA.DIM file (Dimap V1 metadata file)
auto resourceFiles = mds.GetResourceFiles();
std::string metadataFile;
for (const auto & file: resourceFiles)
{
metadataFile = file;
break;
if (file.find("METADATA.DIM")!=std::string::npos)
{
metadataFile = file;
break;
}
}
}
if (metadataFile.empty())
{
otbGenericExceptionMacro(MissingMetadataException,<<"Cannot find the METADATA.DIM file")
}
if (metadataFile.empty())
{
otbGenericExceptionMacro(MissingMetadataException,<<"Cannot find the METADATA.DIM file")
}
XMLMetadataSupplier xmlMds(metadataFile);
XMLMetadataSupplier xmlMds(metadataFile);
DimapMetadataHelper helper;
helper.ParseDimapV1(xmlMds, "Dimap_Document.");
}
// geom metadata case
else if (mds.GetAs<std::string>("", "support_data.sensorID" ) == "Spot 5")
{
helper.ParseGeom(mds);
}
else
{
otbGenericExceptionMacro(MissingMetadataException,<<"Not a Spot 5 product")
}
helper.ParseDimapV1(xmlMds, "Dimap_Document.");
const auto & dimapData = helper.GetDimapData();
m_Imd.Add(MDStr::SensorID, "SPOT 5");
m_Imd.Add(MDStr::Mission, "SPOT 5");
m_Imd.Add(MDTime::ProductionDate,
boost::lexical_cast<MetaData::Time>(dimapData.ProductionDate));
m_Imd.Add(MDTime::AcquisitionDate,
......@@ -1801,24 +1809,20 @@ void SpotImageMetadataInterface::Parse(const MetadataSupplierInterface & mds)
auto nbBands = m_Imd.Bands.size();
if (dimapData.PhysicalBias.size() == nbBands
&& dimapData.PhysicalGain.size() == nbBands
&& dimapData.SolarIrradiance.size() == nbBands
&& dimapData.BandIDs.size() == nbBands)
&& dimapData.SolarIrradiance.size() == nbBands)
{
auto bias = dimapData.PhysicalBias.begin();
auto gain = dimapData.PhysicalGain.begin();
auto bandId = dimapData.BandIDs.begin();
auto solarIrradiance = dimapData.SolarIrradiance.begin();
for (auto & band: m_Imd.Bands)
{
band.Add(MDNum::PhysicalGain, *gain);
band.Add(MDNum::PhysicalBias, *bias);
band.Add(MDStr::BandName, *bandId);
band.Add(MDNum::SolarIrradiance, *solarIrradiance);
bias++;
gain++;
bandId++;
solarIrradiance++;
}
}
......@@ -1828,6 +1832,36 @@ void SpotImageMetadataInterface::Parse(const MetadataSupplierInterface & mds)
<< "The number of bands in image metadatas is incoherent with the DIMAP product")
}
// band IDs is not parsed in geom files. There are two possible cases: panchromatic and multispectral
if (dimapData.BandIDs.empty())
{
if (m_Imd.Bands.size() == 1)
{
m_Imd.Bands[0].Add(MDStr::BandName, "P");
}
else if (m_Imd.Bands.size() == 4)
{
m_Imd.Bands[0].Add(MDStr::BandName, "XS1");
m_Imd.Bands[1].Add(MDStr::BandName, "XS2");
m_Imd.Bands[2].Add(MDStr::BandName, "XS3");
m_Imd.Bands[3].Add(MDStr::BandName, "SWIR");
}
}
else if (dimapData.BandIDs.size() == nbBands )
{
auto bandId = dimapData.BandIDs.begin();
for (auto & band: m_Imd.Bands)
{
band.Add(MDStr::BandName, *bandId);
bandId++;
}
}
else
{
otbGenericExceptionMacro(MissingMetadataException,<<"Invalid band number in the SPOT 5 product")
}
if (nbBands == 4)
{
FetchSpectralSensitivity();
......
......@@ -431,6 +431,15 @@ otb_add_test(NAME ioTuotbXMLMetadataSupplierTest COMMAND otbMetadataTestDriver
${TEMP}/ioTuotbXMLMetadataSupplierTest.txt
)
# Image metadata interface tests.
# The image metadata are read from a product and dumped into a txt file, which is then compared
# to a baseline ${BASELINE_FILES}/ioTvImageMetadataInterfaceTest_${sensor}.txt.
# To create a test, add the name of the sensor in "sensor list" and define the ${sensor}_filename variable.
# Additionally, it is possible to compare the metadata of the product with the metadata of a geom file
# with the --compare-metadata utility by setting ${sensor}_geom_file
set(sensor_list pleiades spot6 spot5 spot5_2 formosat ikonos worldview2 quickbird Radarsat2 TerraSar CosmoSkyMed Sentinel1)
set(pleiades_filename LARGEINPUT{PLEIADES/MAIDO_JP2_DIMAPv2_PRIMARY_MS_lossy_12bits/IMG_PHR1A_MS_002/IMG_PHR1A_MS_201206050630064_SEN_559102101-002_R1C1.JP2})
......@@ -446,12 +455,33 @@ set(TerraSar_filename LARGEINPUT{TERRASARX/2008-03-10_GrandCanyon_SSC/TSX1_SAR__
set(CosmoSkyMed_filename LARGEINPUT{COSMOSKYMED/Toulouse_spotlight/CSKS3_SCS_B_S2_08_HH_RD_SF_20110418180325_20110418180332.h5})
set(Sentinel1_filename LARGEINPUT{SENTINEL1/S1A_S6_SLC__1SSV_20150619T195043})
set(spot5_2_geom_file ${INPUTDATA}/spot5/spot5-1.geom)
set(formosat_geom_file ${INPUTDATA}/FORMOSAT_ROI_1000_100.geom)
foreach(sensor ${sensor_list})
if (${sensor}_geom_file)
otb_add_test(NAME ioTvImageMetadataInterfaceTest_${sensor} COMMAND otbMetadataTestDriver
--compare-ascii 0
${BASELINE_FILES}/ioTvImageMetadataInterfaceTest_${sensor}.txt
${TEMP}/ioTvImageMetadataInterfaceTest_${sensor}.txt
--compare-metadata 0
${${sensor}_filename}
${${sensor}_filename}?&geom=${${sensor}_geom_file}
otbImageMetadataInterfaceTest
${${sensor}_filename}
${TEMP}/ioTvImageMetadataInterfaceTest_${sensor}.txt
)
else()
otb_add_test(NAME ioTvImageMetadataInterfaceTest_${sensor} COMMAND otbMetadataTestDriver
--compare-ascii 0 ${BASELINE_FILES}/ioTvImageMetadataInterfaceTest_${sensor}.txt
--compare-ascii 0
${BASELINE_FILES}/ioTvImageMetadataInterfaceTest_${sensor}.txt
${TEMP}/ioTvImageMetadataInterfaceTest_${sensor}.txt
otbImageMetadataInterfaceTest
${${sensor}_filename}
${TEMP}/ioTvImageMetadataInterfaceTest_${sensor}.txt
)
endif()
endforeach()
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment