Commit 2f586d7a authored by Cédric Traizet's avatar Cédric Traizet

ENH: add parsing for WorldView-2 metadata from geom

parent f1158f5a
Pipeline #6177 failed with stages
in 121 minutes and 9 seconds
SensorID WV02
Mission WorldView-2
ProductType Multi
GeometricLevel LV2A
TileHintX 1150
TileHintY 455
DataType 3
SunElevation 24.9
SunAzimuth 171.7
SatElevation 299.2
SatAzimuth 60.5
SatElevation 60.5
SatAzimuth 299.2
AcquisitionDate 2009-12-10T10:30:18.142149Z
ProductionDate 2010-01-11T19:51:30Z
ProjectionWKT PROJCS["WGS 84 / UTM zone 33N",
......
adjustment_0.adj_param_0.center: 0
adjustment_0.adj_param_0.description: intrack_offset
adjustment_0.adj_param_0.lock_flag: 0
adjustment_0.adj_param_0.parameter: 0
adjustment_0.adj_param_0.sigma: 50
adjustment_0.adj_param_0.units: pixel
adjustment_0.adj_param_1.center: 0
adjustment_0.adj_param_1.description: crtrack_offset
adjustment_0.adj_param_1.lock_flag: 0
adjustment_0.adj_param_1.parameter: 0
adjustment_0.adj_param_1.sigma: 50
adjustment_0.adj_param_1.units: pixel
adjustment_0.adj_param_2.center: 0
adjustment_0.adj_param_2.description: intrack_scale
adjustment_0.adj_param_2.lock_flag: 0
adjustment_0.adj_param_2.parameter: 0
adjustment_0.adj_param_2.sigma: 50
adjustment_0.adj_param_2.units: unknown
adjustment_0.adj_param_3.center: 0
adjustment_0.adj_param_3.description: crtrack_scale
adjustment_0.adj_param_3.lock_flag: 0
adjustment_0.adj_param_3.parameter: 0
adjustment_0.adj_param_3.sigma: 50
adjustment_0.adj_param_3.units: unknown
adjustment_0.adj_param_4.center: 0
adjustment_0.adj_param_4.description: map_rotation
adjustment_0.adj_param_4.lock_flag: 0
adjustment_0.adj_param_4.parameter: 0
adjustment_0.adj_param_4.sigma: 0.1
adjustment_0.adj_param_4.units: degrees
adjustment_0.description: Initial adjustment
adjustment_0.dirty_flag: 0
adjustment_0.number_of_params: 5
bias_error: 0
ce90_absolute: 0
ce90_relative: 0
current_adjustment: 0
height_off: 94
height_scale: 500
image_id: 09DEC10103019-M2AS-052298844010_01_P001
lat_off: 41.8901
lat_scale: 0.0111
line_den_coeff_00: 1
line_den_coeff_01: -7.80273e-05
line_den_coeff_02: 0.003246576
line_den_coeff_03: 0
line_den_coeff_04: -5.185663e-07
line_den_coeff_05: 0
line_den_coeff_06: 0
line_den_coeff_07: 3.258874e-07
line_den_coeff_08: 1.052847e-05
line_den_coeff_09: -2.694567e-07
line_den_coeff_10: 0
line_den_coeff_11: 0
line_den_coeff_12: 0
line_den_coeff_13: 0
line_den_coeff_14: 0
line_den_coeff_15: 0
line_den_coeff_16: 0
line_den_coeff_17: 0
line_den_coeff_18: 0
line_den_coeff_19: 0
line_num_coeff_00: 0.006527383
line_num_coeff_01: 0.06434966
line_num_coeff_02: -1.071741
line_num_coeff_03: 0
line_num_coeff_04: 0.0002952761
line_num_coeff_05: 0
line_num_coeff_06: 0
line_num_coeff_07: -0.0003969196
line_num_coeff_08: -0.00348033
line_num_coeff_09: 0
line_num_coeff_10: 0
line_num_coeff_11: 5.867179e-08
line_num_coeff_12: 1.237231e-06
line_num_coeff_13: -1.734225e-08
line_num_coeff_14: -1.671568e-06
line_num_coeff_15: -1.128591e-05
line_num_coeff_16: 2.888045e-07
line_num_coeff_17: 0
line_num_coeff_18: 0
line_num_coeff_19: 0
line_off: 575
line_scale: 575
ll_lat: 41.8795189359076
ll_lon: 12.4788418224102
long_off: 12.4922
long_scale: 0.0305
lr_lat: 41.8801236343259
lr_lon: 12.5065153441922
meters_per_pixel_x: 1.99973246942754
meters_per_pixel_y: 1.99973454897226
number_lines: 1151
number_of_adjustments: 1
number_samples: 1150
polynomial_format: B
rand_error: 0
rect: 0 0 1149 1150
ref_point_hgt: 94
ref_point_lat: 41.8901692107072
ref_point_line: 575
ref_point_lon: 12.4922735919611
ref_point_samp: 574.5
samp_den_coeff_00: 1
samp_den_coeff_01: 0.001067082
samp_den_coeff_02: 9.526734e-05
samp_den_coeff_03: 0
samp_den_coeff_04: 1.607511e-08
samp_den_coeff_05: 0
samp_den_coeff_06: 0
samp_den_coeff_07: 1.133002e-06
samp_den_coeff_08: 1.612854e-08
samp_den_coeff_09: 0
samp_den_coeff_10: 0
samp_den_coeff_11: 6.6949e-08
samp_den_coeff_12: 0
samp_den_coeff_13: 0
samp_den_coeff_14: 0
samp_den_coeff_15: 0
samp_den_coeff_16: 0
samp_den_coeff_17: 0
samp_den_coeff_18: 0
samp_den_coeff_19: 0
samp_num_coeff_00: -0.002132228
samp_num_coeff_01: 1.012216
samp_num_coeff_02: 0.01441729
samp_num_coeff_03: 0
samp_num_coeff_04: -6.379023e-05
samp_num_coeff_05: 0
samp_num_coeff_06: 0
samp_num_coeff_07: 0.001078912
samp_num_coeff_08: 2.941825e-06
samp_num_coeff_09: 0
samp_num_coeff_10: 0
samp_num_coeff_11: 1.140182e-06
samp_num_coeff_12: -1.754211e-08
samp_num_coeff_13: 0
samp_num_coeff_14: -1.505595e-07
samp_num_coeff_15: 0
samp_num_coeff_16: 0
samp_num_coeff_17: 0
samp_num_coeff_18: 0
samp_num_coeff_19: 0
samp_off: 574
samp_scale: 1250
sensor: WV02
support_data.B_band_absCalFactor: 0.007291212
support_data.C_band_absCalFactor: 0.009295654
support_data.G_band_absCalFactor: 0.005654403
support_data.N2_band_absCalFactor: 0.009042234
support_data.N_band_absCalFactor: 0.005076947
support_data.RE_band_absCalFactor: 0.004539619
support_data.R_band_absCalFactor: 0.004569992
support_data.TDI_level: 24
support_data.Y_band_absCalFactor: 0.005101088
support_data.azimuth_angle: 171.7
support_data.band_id: Multi
support_data.band_name_list: C B G Y R RE N N2
support_data.bits_per_pixel: 16
support_data.elevation_angle: 24.9
support_data.generation_date: 2010-01-11T19:51:30.000000Z
support_data.sat_azimuth_angle: 299.2
support_data.sat_elevation_angle: 60.5
support_data.sat_id: WV02
support_data.tlc_date: 2009-12-10T10:30:18.142149Z;
support_data.type: ossimQuickbirdMetaData
type: ossimQuickbirdRpcModel
ul_lat: 41.9002147870885
ul_lon: 12.47803183973
ur_lat: 41.9008194855069
ur_lon: 12.505705361512
......@@ -150,9 +150,7 @@ private:
void FetchPhysicalBias();
void FetchSolarIrradiance();
void FetchPhysicalGain(const MetadataSupplierInterface &);
void FetchDates(const MetadataSupplierInterface &);
void FetchWavelengths();
void FetchWavelengths();
void FetchSpectralSensitivity();
};
......
......@@ -116,7 +116,6 @@ ImageMetadataInterfaceFactory
{
// silent catch of MissingMetadataException
// just means that this IMI can't parse the file
std::cout << e.what() << std::endl;
}
}
else
......
......@@ -87,39 +87,6 @@ void WorldView2ImageMetadataInterface::FetchSolarIrradiance()
}
}
void WorldView2ImageMetadataInterface::FetchDates(const MetadataSupplierInterface & mds)
{
if(m_Imd[MDStr::SensorID] != "WV02")
{
itkExceptionMacro(<< "Invalid Metadata, no WorldView2 Image");
}
bool hasValue = 0;
// Use TLC time in priority, only available for level 1B products.
// If it is not available use firstLineTime instead.
auto tlcTimeStr = mds.GetMetadataValue("IMAGE_1.TLCTime", hasValue);
if (hasValue)
{
try
{
m_Imd.Add(MDTime::AcquisitionDate,
boost::lexical_cast<MetaData::Time>(tlcTimeStr));
}
catch (boost::bad_lexical_cast&)
{
otbGenericExceptionMacro(MissingMetadataException,
<<"Bad metadata value for 'IMAGE_1.TLCTime', got: "
<<tlcTimeStr);
}
}
else
{
Fetch(MDTime::AcquisitionDate, mds, "IMD/IMAGE_1.firstLineTime" );
}
Fetch(MDTime::ProductionDate, mds, "IMD/generationTime" );
}
void WorldView2ImageMetadataInterface::FetchWavelengths()
{
if(m_Imd[MDStr::SensorID] != "WV02")
......@@ -695,53 +662,6 @@ WorldView2ImageMetadataInterface::VariableLengthVectorType WorldView2ImageMetada
return outputValuesVariableLengthVector;
}
void WorldView2ImageMetadataInterface::FetchPhysicalGain(const MetadataSupplierInterface &mds)
{
if(m_Imd[MDStr::SensorID] != "WV02")
{
itkExceptionMacro(<< "Invalid Metadata, no WorldView2 Image");
}
auto bitsPerPixel = mds.GetAs<int>("IMD/bitsPerPixel");
if (bitsPerPixel != 16)
{
itkExceptionMacro(<< "Invalid bitsPerPixel " << bitsPerPixel);
}
auto productType = m_Imd[MDStr::ProductType];
std::string panchro("P");
std::string multi("Multi");
std::string ms1("MS1");
if (productType != panchro && productType != multi && productType != ms1)
{
itkExceptionMacro(<< "Invalid bandID " << productType);
}
for (auto & band: m_Imd.Bands)
{
auto prefix = "IMD/BAND_" + band[MDStr::BandName] + ".";
auto TDILevel = mds.GetAs<int>(prefix + "TDILevel");
// Verify if TDI level is correct, according to the document
// Radiometric_Use_of_WorldView-2_Imagery, Technical Note, 2010 (Digital Globe)
if ( (productType == panchro && (TDILevel <8 || TDILevel > 64))
|| ( (productType == multi || productType == ms1 ) && (TDILevel <3 || TDILevel > 24)) )
{
otbGenericExceptionMacro(MissingMetadataException,
<< "Invalid TDI settings for band " << band[MDStr::BandName]);
}
// In previous version of OTB the effective bandwith was tabulated
// because Ossim did not read this specific metadata.
auto effectiveBandwidth = mds.GetAs<double>(prefix + "effectiveBandwidth");
auto absCalFactor = mds.GetAs<double>(prefix + "absCalFactor");
band.Add(MDNum::PhysicalGain, effectiveBandwidth / absCalFactor);
}
}
double WorldView2ImageMetadataInterface::GetSatElevation() const
{
const MetaDataDictionaryType& dict = this->GetMetaDataDictionary();
......@@ -1770,112 +1690,255 @@ namespace
s.erase( s.size() - 1 ); // erase the last character
}
}
}
void WorldView2ImageMetadataInterface::Parse(const MetadataSupplierInterface & mds)
{
// Check if there is DG metadatas
bool hasValue = false;
auto metadatatype = mds.GetMetadataValue("METADATATYPE", hasValue);
if (!hasValue || metadatatype != "DG")
struct WorldView2Metadata
{
otbGenericExceptionMacro(MissingMetadataException,
<< "No Digital Globe metadata has been found")
}
std::string sensorId;
std::string bandId;
std::string acquisitionDateStr;
std::string productionDateStr;
// Check if the sensor is WorldView 2
auto sensorID = mds.GetMetadataValue("IMD/IMAGE_1.satId", hasValue);
double sunElevation;
double sunAzimuth;
double satElevation;
double satAzimuth;
int bitsPerPixel;
std::vector<std::string> bandNameList;
std::unordered_map<std::string, double> absCalFactor;
std::unordered_map<std::string, int> TDILevel;
std::unordered_map<std::string, double> effectiveBandwidth;
};
if (sensorID.find("WV02") != std::string::npos)
void ParseGeomMetadata(const MetadataSupplierInterface & mds, WorldView2Metadata & outMetadata)
{
m_Imd.Add(MDStr::SensorID, "WV02");
m_Imd.Add(MDStr::Mission, "WorldView-2");
}
else
{
otbGenericExceptionMacro(MissingMetadataException, << "Not a Worldview 2 sensor")
}
outMetadata.sensorId = mds.GetAs<std::string>("sensor");
outMetadata.bandId = mds.GetAs<std::string>("support_data.band_id");
// the date might have a semicolon at its end in the geom, e.g.
// support_data.tlc_date: 2009-12-10T10:30:18.142149Z;
outMetadata.acquisitionDateStr = mds.GetAs<std::string>("support_data.tlc_date");
if (outMetadata.acquisitionDateStr.back() == ';')
{
outMetadata.acquisitionDateStr.pop_back();
}
auto ressourceFiles = mds.GetResourceFiles();
outMetadata.productionDateStr = mds.GetAs<std::string>("support_data.generation_date");
outMetadata.sunElevation = mds.GetAs<double>("support_data.elevation_angle");
outMetadata.sunAzimuth = mds.GetAs<double>("support_data.azimuth_angle");
outMetadata.satElevation = mds.GetAs<double>("support_data.sat_elevation_angle");
outMetadata.satAzimuth = mds.GetAs<double>("support_data.sat_azimuth_angle");
// return true if the file extension is .IMD
auto lambda = [](const std::string & filename){return itksys::SystemTools::GetFilenameLastExtension(filename) == ".IMD";};
auto filename = std::find_if(ressourceFiles.begin(), ressourceFiles.end(), lambda);
if (filename == ressourceFiles.end())
{
otbGenericExceptionMacro(MissingMetadataException,
<< "Cannot find the .IMD file associated with the Digital Globe dataset");
outMetadata.bitsPerPixel = mds.GetAs<int>("support_data.bits_per_pixel");
outMetadata.bandNameList = mds.GetAsVector<std::string>("support_data.band_name_list");
// In the case of WorldView-2, we need to divide the absolute calibration
// factor by the effective bandwidth see for details:
// http://www.digitalglobe.com/downloads/Radiometric_Use_of_WorldView-2_Imagery.pdf
// These values are not written in the geom file by OSSIM as
// there are specific to WV2 We did not retrieve those values in the ossim
// class and consider them as constant values
outMetadata.effectiveBandwidth = {
{"P", 2.846000e-01},
{"C", 4.730000e-02},
{"B", 5.430000e-02},
{"G", 6.300000e-02},
{"Y", 3.740000e-02},
{"R", 5.740000e-02},
{"RE", 3.930000e-02},
{"N", 9.890000e-02},
{"N2", 9.960000e-02}
};
auto TDILevel = mds.GetAs<int>("support_data.TDI_level");
for (const auto & bandName : outMetadata.bandNameList)
{
outMetadata.TDILevel[bandName] = TDILevel;
outMetadata.absCalFactor[bandName] = mds.GetAs<double>("support_data." + bandName + "_band_absCalFactor");
}
}
std::ifstream file(*filename);
auto bandIt = m_Imd.Bands.begin();
if (file.is_open())
void ParseProductMetadata(const MetadataSupplierInterface & mds, WorldView2Metadata & outMetadata)
{
// Find band names at group definitions :
// BEGIN_GROUP = BAND_N
const auto key = "BEGIN_GROUP";
const auto key_size = strlen(key);
// Open the metadata file to find the band names, as it is not parsed by GDAL
auto ressourceFiles = mds.GetResourceFiles();
const auto prefix = "BAND_";
const auto prefix_size = strlen(prefix);
// return true if the file extension is .IMD
auto lambda = [](const std::string & filename){return itksys::SystemTools::GetFilenameLastExtension(filename) == ".IMD";};
auto filename = std::find_if(ressourceFiles.begin(), ressourceFiles.end(), lambda);
if (filename == ressourceFiles.end())
{
otbGenericExceptionMacro(MissingMetadataException,
<< "Cannot find the .IMD file associated with the Digital Globe dataset");
}
std::ifstream file(*filename);
std::string line;
while (std::getline(file, line))
if (file.is_open())
{
if (starts_with(line, "BEGIN_GROUP"))
// Find band names at group definitions :
// BEGIN_GROUP = BAND_N
const auto key = "BEGIN_GROUP";
const auto key_size = strlen(key);
const auto prefix = "BAND_";
const auto prefix_size = strlen(prefix);
std::string line;
while (std::getline(file, line))
{
auto pos = line.find(prefix, key_size);
if (pos != std::string::npos)
if (starts_with(line, "BEGIN_GROUP"))
{
if (bandIt == m_Imd.Bands.end())
auto pos = line.find(prefix, key_size);
if (pos != std::string::npos)
{
otbGenericExceptionMacro(MissingMetadataException,
<< "Number of bands in "
<< *filename
<< " inconsistent with the number of bands of the image");
outMetadata.bandNameList.push_back(line.substr(pos+prefix_size, line.size()-pos-prefix_size));
}
bandIt->Add(MDStr::BandName,
line.substr(pos+prefix_size, line.size()-pos-prefix_size));
bandIt++;
}
}
file.close();
}
outMetadata.sensorId = mds.GetAs<std::string>("IMD/IMAGE_1.satId");
outMetadata.bandId = mds.GetAs<std::string>("IMD/bandId");
outMetadata.sunElevation = mds.GetAs<double>("IMD/IMAGE_1.meanSunEl");
outMetadata.sunAzimuth = mds.GetAs<double>("IMD/IMAGE_1.meanSunAz");
outMetadata.satElevation = mds.GetAs<double>("IMD/IMAGE_1.meanSatEl");
outMetadata.satAzimuth = mds.GetAs<double>("IMD/IMAGE_1.meanSatAz");
outMetadata.productionDateStr = mds.GetAs<std::string>("IMD/generationTime");
// Use TLC time in priority, only available for level 1B products.
// If it is not available use firstLineTime instead.
outMetadata.acquisitionDateStr = mds.GetAs<std::string>("", "IMAGE_1.TLCTime");
if (outMetadata.acquisitionDateStr == "")
{
outMetadata.acquisitionDateStr = mds.GetAs<std::string>("IMD/IMAGE_1.firstLineTime");
}
outMetadata.bitsPerPixel = mds.GetAs<int>("IMD/bitsPerPixel");
if (outMetadata.bitsPerPixel != 16)
{
otbGenericExceptionMacro(MissingMetadataException,<< "Invalid bitsPerPixel " << outMetadata.bitsPerPixel);
}
for (auto & bandName: outMetadata.bandNameList)
{
auto prefix = "IMD/BAND_" + bandName + ".";
auto TDILevel = mds.GetAs<int>(prefix + "TDILevel");
outMetadata.TDILevel[bandName] = TDILevel;
// Verify if TDI level is correct, according to the document
// Radiometric_Use_of_WorldView-2_Imagery, Technical Note, 2010 (Digital Globe)
if ( (outMetadata.bandId == "P" && (TDILevel <8 || TDILevel > 64))
|| ( (outMetadata.bandId == "Multi" || outMetadata.bandId == "MS1" ) && (TDILevel <3 || TDILevel > 24)) )
{
otbGenericExceptionMacro(MissingMetadataException,
<< "Invalid TDI settings for band " << bandName);
}
// In previous version of OTB the effective bandwith was tabulated
// because Ossim did not read this specific metadata.
outMetadata.effectiveBandwidth[bandName] = mds.GetAs<double>(prefix + "effectiveBandwidth");
outMetadata.absCalFactor[bandName] = mds.GetAs<double>(prefix + "absCalFactor");
}
file.close();
}
auto geometricLevel = mds.GetMetadataValue("IMD/productLevel", hasValue);
// throw missing
Unquote(geometricLevel);
m_Imd.Add(MDStr::GeometricLevel, geometricLevel);
auto productType = mds.GetMetadataValue("IMD/bandId", hasValue);
// throw missing
Unquote(productType);
m_Imd.Add(MDStr::ProductType, productType);
}
void WorldView2ImageMetadataInterface::Parse(const MetadataSupplierInterface & mds)
{
WorldView2Metadata metadata;
// Check if the sensor is WorldView 2
if (mds.GetAs<std::string>("", "IMD/IMAGE_1.satId").find("WV02") != std::string::npos)
{
ParseProductMetadata(mds, metadata);
FetchRPC(mds);
}
else if (mds.GetAs<std::string>("", "support_data.sat_id").find("WV02") != std::string::npos)
{
ParseGeomMetadata(mds, metadata);
}
else
{
otbGenericExceptionMacro(MissingMetadataException, << "Not a Worldview 2 sensor")
}
m_Imd.Add(MDStr::SensorID, "WV02");
m_Imd.Add(MDStr::Mission, "WorldView-2");
// Product type (Multi, MS1 or P)
Unquote(metadata.bandId);
m_Imd.Add(MDStr::ProductType, metadata.bandId);
if (m_Imd.Bands.size() == metadata.bandNameList.size())
{
auto bandIt = m_Imd.Bands.begin();
for (const auto & bandName : metadata.bandNameList)
{
bandIt->Add(MDStr::BandName, bandName);
bandIt++;
}
}
else
{
otbGenericExceptionMacro(MissingMetadataException,
<< "The Number of bands in the metadata is inconsistent with the number of bands of the image");
}
// Acquisition and production dates
FetchDates(mds);
try
{
m_Imd.Add(MDTime::AcquisitionDate,
boost::lexical_cast<MetaData::Time>(metadata.acquisitionDateStr));
}
catch (boost::bad_lexical_cast&)
{
otbGenericExceptionMacro(MissingMetadataException,
<< "Invalid value for the acquisition date, got: "
<< metadata.acquisitionDateStr);
}
try
{
m_Imd.Add(MDTime::ProductionDate,
boost::lexical_cast<MetaData::Time>(metadata.productionDateStr));
}
catch (boost::bad_lexical_cast&)
{
otbGenericExceptionMacro(MissingMetadataException,
<< "Invalid value for the production date, got: "
<< metadata.productionDateStr);
}
//Radiometry
Fetch(MDNum::SunElevation, mds, "IMD/IMAGE_1.meanSunEl");
Fetch(MDNum::SunAzimuth, mds, "IMD/IMAGE_1.meanSunAz");
Fetch(MDNum::SatElevation, mds, "IMD/IMAGE_1.meanSatAz");
Fetch(MDNum::SatAzimuth , mds, "IMD/IMAGE_1.meanSatEl");