From bdd042df373a4b77a381c393eb827a7da769b208 Mon Sep 17 00:00:00 2001
From: Julien Michel <julien.michel@orfeo-toolbox.org>
Date: Mon, 21 Mar 2011 18:39:14 +0100
Subject: [PATCH] TEST: Adding a test to assess performances of sensor models
---
Testing/Code/Projections/CMakeLists.txt | 14 ++
.../otbGenericRSTransformFromImage.cxx | 125 ++++++++++++++++++
.../Code/Projections/otbProjectionsTests3.cxx | 1 +
3 files changed, 140 insertions(+)
diff --git a/Testing/Code/Projections/CMakeLists.txt b/Testing/Code/Projections/CMakeLists.txt
index 7b94e96cfd..c2d6ebdaf1 100644
--- a/Testing/Code/Projections/CMakeLists.txt
+++ b/Testing/Code/Projections/CMakeLists.txt
@@ -842,6 +842,20 @@ ADD_TEST(prTvGenericRSTransformFromImage ${PROJECTIONS_TESTS3}
otbGenericRSTransformFromImage
${INPUTDATA}/WithoutProjRefWithKeywordlist.tif)
+IF(OTB_DATA_USE_LARGEINPUT)
+ADD_TEST(prTvGenericRSTransformQuickbirdToulouseGeodesicPointChecking ${PROJECTIONS_TESTS3}
+otbGenericRSTransformImageAndMNTToWGS84ConversionChecking
+${LARGEINPUT}/QUICKBIRD/TOULOUSE/000000128955_01_P001_PAN/02APR01105228-P1BS-000000128955_01_P001.TIF
+${INPUTDATA}/DEM/srtm_directory/
+16271 15647 # IGN reference point in front of st sernin church pointed in image
+1.48353 43.55968 # IGN reference point lon/lat
+192.205 # IGN reference point elevation (prec < 10 cm)
+5
+0.00001
+)
+
+ENDIF(OTB_DATA_USE_LARGEINPUT)
+
ADD_TEST(prTuVectorDataProjectionFilterNew ${PROJECTIONS_TESTS3} otbVectorDataProjectionFilterNew )
#test points
diff --git a/Testing/Code/Projections/otbGenericRSTransformFromImage.cxx b/Testing/Code/Projections/otbGenericRSTransformFromImage.cxx
index 415d40aefd..92707d539c 100644
--- a/Testing/Code/Projections/otbGenericRSTransformFromImage.cxx
+++ b/Testing/Code/Projections/otbGenericRSTransformFromImage.cxx
@@ -24,11 +24,14 @@
#include "otbGenericRSTransform.h"
#include <ogr_spatialref.h>
#include <fstream>
+#include "itkEuclideanDistance.h"
+
typedef otb::VectorImage<double, 2> ImageType;
typedef otb::ImageFileReader<ImageType> ReaderType;
typedef otb::GenericRSTransform<> TransformType;
typedef TransformType::InputPointType PointType;
+typedef itk::Statistics::EuclideanDistance<ImageType::PointType> DistanceType;
int otbGenericRSTransformFromImage(int argc, char* argv[])
{
@@ -64,3 +67,125 @@ int otbGenericRSTransformFromImage(int argc, char* argv[])
return EXIT_SUCCESS;
}
+
+int otbGenericRSTransformImageAndMNTToWGS84ConversionChecking(int argc, char* argv[])
+{
+ std::string infname = argv[1];
+ std::string demdir = argv[2];
+
+ ImageType::PointType refImgPt,refGeoPt, estimatedImgPt, estimatedGeoPt;
+ refImgPt[0] = atof(argv[3]);
+ refImgPt[1] = atof(argv[4]);
+ refGeoPt[0] = atof(argv[5]);
+ refGeoPt[1] = atof(argv[6]);
+ double refHeight = atof(argv[7]);
+ double imgTol = atof(argv[8]);
+ double geoTol = atof(argv[9]);
+
+ bool pass = true;
+
+// Reader
+ ReaderType::Pointer reader = ReaderType::New();
+ reader->SetFileName(argv[1]);
+ reader->UpdateOutputInformation();
+
+ // Build wgs ref
+ OGRSpatialReference oSRS;
+ oSRS.SetWellKnownGeogCS("WGS84");
+ char * wgsRef = NULL;
+ oSRS.exportToWkt(&wgsRef);
+
+ DistanceType::Pointer distance = DistanceType::New();
+
+ // Instanciate WGS->Image transform
+ TransformType::Pointer wgs2img = TransformType::New();
+ wgs2img->SetInputProjectionRef(wgsRef);
+ wgs2img->SetOutputProjectionRef(reader->GetOutput()->GetProjectionRef());
+ wgs2img->SetOutputKeywordList(reader->GetOutput()->GetImageKeywordlist());
+ wgs2img->SetDEMDirectory(argv[2]);
+ wgs2img->InstanciateTransform();
+
+ // Instanciate Image->WGS transform
+ TransformType::Pointer img2wgs = TransformType::New();
+ img2wgs->SetInputProjectionRef(reader->GetOutput()->GetProjectionRef());
+ img2wgs->SetInputKeywordList(reader->GetOutput()->GetImageKeywordlist());
+ img2wgs->SetOutputProjectionRef(wgsRef);
+ img2wgs->SetDEMDirectory(argv[2]);
+ img2wgs->InstanciateTransform();
+
+ std::cout<< std::setprecision(8) << std::endl;
+
+ std::cout<<"References: "<<std::endl;
+ std::cout<<"Geo (wgs84): "<<refGeoPt<<std::endl;
+ std::cout<<"Image: "<<refImgPt<<std::endl;
+ std::cout<<"Height: "<<refHeight<<std::endl;
+
+ std::cout<<"Using elevation from "<<demdir<<std::endl;
+
+ estimatedImgPt = wgs2img->TransformPoint(refGeoPt);
+ estimatedGeoPt = img2wgs->TransformPoint(refImgPt);
+
+ std::cout<<"Forward(refImg): "<<refImgPt<<" -> " << estimatedGeoPt<<std::endl;
+ std::cout<<"Inverse(refGeo): "<<refGeoPt<<" -> "<< estimatedImgPt <<std::endl;
+
+ double imgRes = distance->Evaluate(refImgPt,estimatedImgPt);
+ double geoRes = distance->Evaluate(refGeoPt,estimatedGeoPt);
+
+ if(imgRes > imgTol)
+ {
+ pass = false;
+ std::cerr<<"Image residual is too high: "<<imgRes<<std::endl;
+ }
+
+ if(geoRes > geoTol)
+ {
+ pass = false;
+ std::cerr<<"Geographic (WGS84) residual is too high: "<<geoRes<<std::endl;
+ }
+
+ std::cout<<"Using reference elevation: "<<std::endl;
+ wgs2img = TransformType::New();
+ wgs2img->SetInputProjectionRef(wgsRef);
+ wgs2img->SetOutputProjectionRef(reader->GetOutput()->GetProjectionRef());
+ wgs2img->SetOutputKeywordList(reader->GetOutput()->GetImageKeywordlist());
+ wgs2img->SetAverageElevation(refHeight);
+ wgs2img->InstanciateTransform();
+
+ img2wgs = TransformType::New();
+ img2wgs->SetInputProjectionRef(reader->GetOutput()->GetProjectionRef());
+ img2wgs->SetInputKeywordList(reader->GetOutput()->GetImageKeywordlist());
+ img2wgs->SetOutputProjectionRef(wgsRef);
+ img2wgs->SetAverageElevation(refHeight);
+ img2wgs->InstanciateTransform();
+
+ estimatedImgPt = wgs2img->TransformPoint(refGeoPt);
+ estimatedGeoPt = img2wgs->TransformPoint(refImgPt);
+
+ std::cout<<"Forward(refImg): "<<refImgPt<<" -> " << estimatedGeoPt<<std::endl;
+ std::cout<<"Inverse(refGeo): "<<refGeoPt<<" -> "<< estimatedImgPt <<std::endl;
+
+ imgRes = distance->Evaluate(refImgPt,estimatedImgPt);
+ geoRes = distance->Evaluate(refGeoPt,estimatedGeoPt);
+
+ if(imgRes > imgTol)
+ {
+ pass = false;
+ std::cerr<<"Image residual is too high: "<<imgRes<<std::endl;
+ }
+
+ if(geoRes > geoTol)
+ {
+ pass = false;
+ std::cerr<<"Geographic (WGS84) residual is too high: "<<geoRes<<std::endl;
+ }
+
+ if(pass)
+ {
+ return EXIT_SUCCESS;
+ }
+ else
+ {
+ std::cerr<<"There were imprecise results."<<std::endl;
+ return EXIT_FAILURE;
+ }
+}
diff --git a/Testing/Code/Projections/otbProjectionsTests3.cxx b/Testing/Code/Projections/otbProjectionsTests3.cxx
index cfac112ba3..418a62ad10 100644
--- a/Testing/Code/Projections/otbProjectionsTests3.cxx
+++ b/Testing/Code/Projections/otbProjectionsTests3.cxx
@@ -34,6 +34,7 @@ void RegisterTests()
REGISTER_TEST(otbGenericRSTransform);
REGISTER_TEST(otbGenericRSTransformWithSRID);
REGISTER_TEST(otbGenericRSTransformFromImage);
+ REGISTER_TEST(otbGenericRSTransformImageAndMNTToWGS84ConversionChecking);
REGISTER_TEST(otbVectorDataProjectionFilterNew);
REGISTER_TEST(otbVectorDataProjectionFilter);
REGISTER_TEST(otbVectorDataProjectionFilterFromMapToSensor);
--
GitLab