From 3dd1b0e530950990844ffce772de6ad0cb567c48 Mon Sep 17 00:00:00 2001
From: Cyrille Valladeau <cyrille.valladeau@c-s.fr>
Date: Mon, 10 Dec 2007 16:32:39 +0000
Subject: [PATCH] Test sur la radio.
---
...UnaryImageFunctorWithVectorImageFilter.txx | 95 +++-------------
.../otbAtmosphericCorrectionParameters.h | 12 ++
...arametersTo6SAtmosphericRadiativeTerms.cxx | 107 ++++++++----------
.../otbAtmosphericRadiativeTerms.cxx | 45 ++++++--
.../otbImageToReflectanceImageFilter.h | 4 +-
.../otbLuminanceToReflectanceImageFilter.h | 8 +-
Code/Radiometry/otbSIXSTraits.cxx | 6 +-
Testing/Code/Radiometry/CMakeLists.txt | 77 ++++++++-----
...arametersTo6SAtmosphericRadiativeTerms.cxx | 33 +-----
.../Code/Radiometry/otbRadiometryTests.cxx | 1 +
...ectanceToSurfaceReflectanceImageFilter.cxx | 22 +++-
11 files changed, 189 insertions(+), 221 deletions(-)
diff --git a/Code/BasicFilters/otbUnaryImageFunctorWithVectorImageFilter.txx b/Code/BasicFilters/otbUnaryImageFunctorWithVectorImageFilter.txx
index 8bd63038d7..3e992424d0 100644
--- a/Code/BasicFilters/otbUnaryImageFunctorWithVectorImageFilter.txx
+++ b/Code/BasicFilters/otbUnaryImageFunctorWithVectorImageFilter.txx
@@ -22,6 +22,7 @@
#include "itkImageRegionIterator.h"
#include "itkImageRegionConstIterator.h"
#include "itkProgressReporter.h"
+#include "itkNumericTraits.h"
namespace otb
{
@@ -63,86 +64,9 @@ UnaryImageFunctorWithVectorImageFilter<TInputImage,TOutputImage,TFunction>
{
return;
}
-
-// // Set the output image largest possible region. Use a RegionCopier
-// // so that the input and output images can be different dimensions.
-// OutputImageRegionType outputLargestPossibleRegion;
-// this->CallCopyInputRegionToOutputRegion(outputLargestPossibleRegion,
-// inputPtr->GetLargestPossibleRegion());
-// outputPtr->SetLargestPossibleRegion( outputLargestPossibleRegion );
-
-// // Set the output spacing and origin
-// const itk::ImageBase<Superclass::InputImageDimension> *phyData;
-
-// phyData
-// = dynamic_cast<const itk::ImageBase<Superclass::InputImageDimension>*>(this->GetInput());
-
-// if (phyData)
-// {
-// // Copy what we can from the image from spacing and origin of the input
-// // This logic needs to be augmented with logic that select which
-// // dimensions to copy
-// unsigned int i, j;
-// const typename InputImageType::SpacingType&
-// inputSpacing = inputPtr->GetSpacing();
-// const typename InputImageType::PointType&
-// inputOrigin = inputPtr->GetOrigin();
-// const typename InputImageType::DirectionType&
-// inputDirection = inputPtr->GetDirection();
-
-// typename OutputImageType::SpacingType outputSpacing;
-// typename OutputImageType::PointType outputOrigin;
-// typename OutputImageType::DirectionType outputDirection;
-
-// // copy the input to the output and fill the rest of the
-// // output with zeros.
-// for (i=0; i < Superclass::InputImageDimension; ++i)
-// {
-// outputSpacing[i] = inputSpacing[i];
-// outputOrigin[i] = inputOrigin[i];
-// for (j=0; j < Superclass::OutputImageDimension; j++)
-// {
-// if (j < Superclass::InputImageDimension)
-// {
-// outputDirection[j][i] = inputDirection[j][i];
-// }
-// else
-// {
-// outputDirection[j][i] = 0.0;
-// }
-// }
-// }
-// for (; i < Superclass::OutputImageDimension; ++i)
-// {
-// outputSpacing[i] = 1.0;
-// outputOrigin[i] = 0.0;
-// for (j=0; j < Superclass::OutputImageDimension; j++)
-// {
-// if (j == i)
-// {
-// outputDirection[j][i] = 1.0;
-// }
-// else
-// {
-// outputDirection[j][i] = 0.0;
-// }
-// }
-// }
-
-// // set the spacing and origin
-// outputPtr->SetSpacing( outputSpacing );
-// outputPtr->SetOrigin( outputOrigin );
-// outputPtr->SetDirection( outputDirection );
outputPtr->SetNumberOfComponentsPerPixel( // propagate vector length info
inputPtr->GetNumberOfComponentsPerPixel());
- // }
-// else
-// {
-// // pointer could not be cast back down
-// itkExceptionMacro(<< "otb::UnaryImageFunctorWithVectorImageFilter::GenerateOutputInformation "
-// << "cannot cast input to "
-// << typeid(itk::ImageBase<Superclass::InputImageDimension>*).name() );
-// }
+
}
/**
@@ -171,16 +95,25 @@ UnaryImageFunctorWithVectorImageFilter<TInputImage,TOutputImage,TFunction>
inputIt.GoToBegin();
outputIt.GoToBegin();
-
+
+ // Null piuxel construction
+ InputPixelType nullPixel;
+ nullPixel.SetSize( outputIt.Get().GetSize() );
+ nullPixel.Fill(itk::NumericTraits<OutputInternalPixelType>::Zero);
+
+
while( !inputIt.IsAtEnd() )
{
InputPixelType inPixel = inputIt.Get();
- OutputPixelType outPixel = outputIt.Get();
-
+ OutputPixelType outPixel = nullPixel;//outputIt.Get();
+ // if the input pixel in null, the output is considered as null ( no sensor informations )
+ if( inPixel!= nullPixel)
+ {
for (unsigned int j=0; j<inPixel.GetSize(); j++)
{
outPixel[j] = m_FunctorVector[j]( inPixel[j] );
}
+ }
outputIt.Set(outPixel);
++inputIt;
++outputIt;
diff --git a/Code/Radiometry/otbAtmosphericCorrectionParameters.h b/Code/Radiometry/otbAtmosphericCorrectionParameters.h
index e07c19f7ae..dfb068cfda 100644
--- a/Code/Radiometry/otbAtmosphericCorrectionParameters.h
+++ b/Code/Radiometry/otbAtmosphericCorrectionParameters.h
@@ -206,6 +206,18 @@ public:
* Set/Get the wavelenght spectral band.
*/
void SetWavelenghtSpectralBand( const WavelenghtSpectralBandVectorType & waveband){ m_WavelenghtSpectralBand = waveband; };
+ void SetWavelenghtSpectralBandWithIndex( unsigned int id, const FilterFunctionValues::Pointer & function)
+ {
+ if (m_WavelenghtSpectralBand.size() < id+1)
+ {
+ for(unsigned int j=0; j<(id+1-m_WavelenghtSpectralBand.size());j++)
+ {
+ FilterFunctionValues::Pointer temp;
+ m_WavelenghtSpectralBand.push_back(temp);
+ }
+ }
+ m_WavelenghtSpectralBand[id] = function;
+ };
WavelenghtSpectralBandVectorType GetWavelenghtSpectralBand(){ return m_WavelenghtSpectralBand; };
WavelenghtSpectralBandVectorType * GetWavelenghtSpectralBandRef(){ return &m_WavelenghtSpectralBand; };
diff --git a/Code/Radiometry/otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms.cxx b/Code/Radiometry/otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms.cxx
index 7dadbf40b0..df73348d65 100644
--- a/Code/Radiometry/otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms.cxx
+++ b/Code/Radiometry/otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms.cxx
@@ -138,65 +138,54 @@ void
AtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms
::GenerateData()
{
-
- AtmosphericCorrectionParametersPointer input = this->GetInput();
- AtmosphericRadiativeTermsPointer output = this->GetOutput();
-
- output->GetValues().clear();
- typedef AtmosphericCorrectionParameters::WavelenghtSpectralBandVectorType WavelenghtSpectralBandVectorType;
- WavelenghtSpectralBandVectorType WavelenghtSpectralBandVector = input->GetWavelenghtSpectralBand();
- unsigned int NbBand = WavelenghtSpectralBandVector.size();
-
- double atmosphericReflectance(0.);
- double atmosphericSphericalAlbedo(0.);
- double totalGaseousTransmission(0.);
- double downwardTransmittance(0.);
- double upwardTransmittance(0.);
-
- for(unsigned int i=0 ; i<NbBand ; i++)
- {
- atmosphericReflectance = 0.;
- atmosphericSphericalAlbedo = 0.;
- totalGaseousTransmission = 0.;
- downwardTransmittance = 0.;
- upwardTransmittance = 0.;
- SIXSTraits::ComputeAtmosphericParameters(
- input->GetSolarZenithalAngle(), /** The Solar zenithal angle */
- input->GetSolarAzimutalAngle(), /** The Solar azimutal angle */
- input->GetViewingZenithalAngle(), /** The Viewing zenithal angle */
- input->GetViewingAzimutalAngle(), /** The Viewing azimutal angle */
- input->GetMonth(), /** The Month */
- input->GetDay(), /** The Day (in the month) */
- input->GetAtmosphericPressure(), /** The Atmospheric pressure */
- input->GetWaterVaporAmount(), /** The Water vapor amount (Total water vapor content over vertical atmospheric column) */
- input->GetOzoneAmount(), /** The Ozone amount (Stratospheric ozone layer content) */
- input->GetAerosolModel(), /** The Aerosol model */
- input->GetAerosolOptical(), /** The Aerosol optical (radiative impact of aerosol for the reference wavelenght 550-nm) */
- input->GetWavelenghtSpectralBand()[i], /** Wavelenght for the spectral band definition */
- /** Note : The Max wavelenght spectral band value must be updated ! */
- atmosphericReflectance, /** Atmospheric reflectance */
- atmosphericSphericalAlbedo, /** atmospheric spherical albedo */
- totalGaseousTransmission, /** Total gaseous transmission */
- downwardTransmittance, /** downward transmittance */
- upwardTransmittance ); /** upward transmittance */
-
- /*
- atmRadiativeTerm = New();
- atmRadiativeTerm->SetAtmosphericReflectances(AtmosphericReflectance);
- atmRadiativeTerm->SetAtmosphericSphericalAlbedo(AtmosphericSphericalAlbedo);
- atmRadiativeTerm->SetTotalGaseousTransmission(TotalGaseousTransmission);
- atmRadiativeTerm->SetDownwardTransmittance(DownwardTransmittance);
- atmRadiativeTerm->SetUpwardTransmittance(UpwardTransmittance);
- //output->Get...().push_back(atmRadiativeTerm);
- output->SetValueByIndex(i, atmRadiativeTerm);
- */
- output->SetIntrinsicAtmosphericReflectances(i, atmosphericReflectance);
- output->SetSphericalAlbedos(i, atmosphericSphericalAlbedo);
- output->SetTotalGaseousTransmissions(i, totalGaseousTransmission);
- output->SetDownwardTransmittances(i, downwardTransmittance);
- output->SetUpwardTransmittances(i, upwardTransmittance);
-
- }
+ AtmosphericCorrectionParametersPointer input = this->GetInput();
+ AtmosphericRadiativeTermsPointer output = this->GetOutput();
+
+ output->GetValues().clear();
+ typedef AtmosphericCorrectionParameters::WavelenghtSpectralBandVectorType WavelenghtSpectralBandVectorType;
+ WavelenghtSpectralBandVectorType WavelenghtSpectralBandVector = input->GetWavelenghtSpectralBand();
+ unsigned int NbBand = WavelenghtSpectralBandVector.size();
+
+ double atmosphericReflectance(0.);
+ double atmosphericSphericalAlbedo(0.);
+ double totalGaseousTransmission(0.);
+ double downwardTransmittance(0.);
+ double upwardTransmittance(0.);
+
+ for(unsigned int i=0 ; i<NbBand ; i++)
+ {
+ atmosphericReflectance = 0.;
+ atmosphericSphericalAlbedo = 0.;
+ totalGaseousTransmission = 0.;
+ downwardTransmittance = 0.;
+ upwardTransmittance = 0.;
+ SIXSTraits::ComputeAtmosphericParameters(
+ input->GetSolarZenithalAngle(), /** The Solar zenithal angle */
+ input->GetSolarAzimutalAngle(), /** The Solar azimutal angle */
+ input->GetViewingZenithalAngle(), /** The Viewing zenithal angle */
+ input->GetViewingAzimutalAngle(), /** The Viewing azimutal angle */
+ input->GetMonth(), /** The Month */
+ input->GetDay(), /** The Day (in the month) */
+ input->GetAtmosphericPressure(), /** The Atmospheric pressure */
+ input->GetWaterVaporAmount(), /** The Water vapor amount (Total water vapor content over vertical atmospheric column) */
+ input->GetOzoneAmount(), /** The Ozone amount (Stratospheric ozone layer content) */
+ input->GetAerosolModel(), /** The Aerosol model */
+ input->GetAerosolOptical(), /** The Aerosol optical (radiative impact of aerosol for the reference wavelenght 550-nm) */
+ input->GetWavelenghtSpectralBand()[i], /** Wavelenght for the spectral band definition */
+ /** Note : The Max wavelenght spectral band value must be updated ! */
+ atmosphericReflectance, /** Atmospheric reflectance */
+ atmosphericSphericalAlbedo, /** atmospheric spherical albedo */
+ totalGaseousTransmission, /** Total gaseous transmission */
+ downwardTransmittance, /** downward transmittance */
+ upwardTransmittance ); /** upward transmittance */
+
+ output->SetIntrinsicAtmosphericReflectances(i, atmosphericReflectance);
+ output->SetSphericalAlbedos(i, atmosphericSphericalAlbedo);
+ output->SetTotalGaseousTransmissions(i, totalGaseousTransmission);
+ output->SetDownwardTransmittances(i, downwardTransmittance);
+ output->SetUpwardTransmittances(i, upwardTransmittance);
+
+ }
}
diff --git a/Code/Radiometry/otbAtmosphericRadiativeTerms.cxx b/Code/Radiometry/otbAtmosphericRadiativeTerms.cxx
index 739945fde5..3cf15bed6d 100644
--- a/Code/Radiometry/otbAtmosphericRadiativeTerms.cxx
+++ b/Code/Radiometry/otbAtmosphericRadiativeTerms.cxx
@@ -134,6 +134,14 @@ AtmosphericRadiativeTerms
{
if ( m_IsInitialized )
{
+ if ( m_Values.size()<id+1 )
+ {
+ for(unsigned int j=0; j<(id+1-m_Values.size());j++)
+ {
+ ValueType temp = AtmosphericRadiativeTermsSingleChannel::New();
+ m_Values.push_back(temp);
+ }
+ }
m_Values[id] = val;
}
else
@@ -147,8 +155,11 @@ AtmosphericRadiativeTerms
{
if ( m_Values.size()<id+1 )
{
- ValueType temp = AtmosphericRadiativeTermsSingleChannel::New();
- m_Values.push_back(temp);
+ for(unsigned int j=0; j<(id+1-m_Values.size());j++)
+ {
+ ValueType temp = AtmosphericRadiativeTermsSingleChannel::New();
+ m_Values.push_back(temp);
+ }
}
m_Values[id]->SetIntrinsicAtmosphericReflectance(val);
}
@@ -157,9 +168,12 @@ AtmosphericRadiativeTerms
::SetSphericalAlbedos(unsigned int id, const double & val)
{
if ( m_Values.size()<id+1 )
- {
- ValueType temp = AtmosphericRadiativeTermsSingleChannel::New();
- m_Values.push_back(temp);
+ {
+ for(unsigned int j=0; j<(id+1-m_Values.size());j++)
+ {
+ ValueType temp = AtmosphericRadiativeTermsSingleChannel::New();
+ m_Values.push_back(temp);
+ }
}
m_Values[id]->SetSphericalAlbedo(val);
}
@@ -169,8 +183,11 @@ AtmosphericRadiativeTerms
{
if ( m_Values.size()<id+1 )
{
- ValueType temp = AtmosphericRadiativeTermsSingleChannel::New();
- m_Values.push_back(temp);
+ for(unsigned int j=0; j<(id+1-m_Values.size());j++)
+ {
+ ValueType temp = AtmosphericRadiativeTermsSingleChannel::New();
+ m_Values.push_back(temp);
+ }
}
m_Values[id]->SetTotalGaseousTransmission(val);
}
@@ -180,8 +197,11 @@ AtmosphericRadiativeTerms
{
if ( m_Values.size()<id+1 )
{
- ValueType temp = AtmosphericRadiativeTermsSingleChannel::New();
- m_Values.push_back(temp);
+ for(unsigned int j=0; j<(id+1-m_Values.size());j++)
+ {
+ ValueType temp = AtmosphericRadiativeTermsSingleChannel::New();
+ m_Values.push_back(temp);
+ }
}
m_Values[id]->SetDownwardTransmittance(val);
}
@@ -191,8 +211,11 @@ AtmosphericRadiativeTerms
{
if ( m_Values.size()<id+1 )
{
- ValueType temp = AtmosphericRadiativeTermsSingleChannel::New();
- m_Values.push_back(temp);
+ for(unsigned int j=0; j<(id+1-m_Values.size());j++)
+ {
+ ValueType temp = AtmosphericRadiativeTermsSingleChannel::New();
+ m_Values.push_back(temp);
+ }
}
m_Values[id]->SetUpwardTransmittance(val);
}
diff --git a/Code/Radiometry/otbImageToReflectanceImageFilter.h b/Code/Radiometry/otbImageToReflectanceImageFilter.h
index 9128e45504..88d9b6d79d 100755
--- a/Code/Radiometry/otbImageToReflectanceImageFilter.h
+++ b/Code/Radiometry/otbImageToReflectanceImageFilter.h
@@ -207,7 +207,7 @@ public:
otb_6s_integer mounth = static_cast<otb_6s_integer>(m_Month);
int cr(0);
cr = otb_6s_varsol_(&day, &mounth, &dsol);
- coefTemp = vcl_cos(m_ZenithalSolarRadius)*static_cast<double>(dsol);
+ coefTemp = vcl_cos(m_ZenithalSolarRadius*M_PI/180)*static_cast<double>(dsol);
}
else
{
@@ -216,7 +216,7 @@ public:
}
else
{
- coefTemp = vcl_cos(m_ZenithalSolarRadius)*m_FluxNormalizationCoefficient*m_FluxNormalizationCoefficient;
+ coefTemp = vcl_cos(m_ZenithalSolarRadius*M_PI/180)*m_FluxNormalizationCoefficient*m_FluxNormalizationCoefficient;
}
functor.SetIlluminationCorrectionCoefficient(1. / coefTemp);
functor.SetAlpha(m_Alpha[i]);
diff --git a/Code/Radiometry/otbLuminanceToReflectanceImageFilter.h b/Code/Radiometry/otbLuminanceToReflectanceImageFilter.h
index 61745ee395..0107133b42 100644
--- a/Code/Radiometry/otbLuminanceToReflectanceImageFilter.h
+++ b/Code/Radiometry/otbLuminanceToReflectanceImageFilter.h
@@ -177,8 +177,9 @@ public:
/** Update the functor list */
virtual void BeforeThreadedGenerateData(void)
- {
+ {std::cout<<"BeforeThreadedGenerateData LumToRef"<<std::endl;
this->GetFunctorVector().clear();
+ std::cout<<"1. / coefTemp"<<this->GetInput()->GetNumberOfComponentsPerPixel()<<std::endl;
for(unsigned int i = 0;i<this->GetInput()->GetNumberOfComponentsPerPixel();++i)
{
FunctorType functor;
@@ -192,7 +193,7 @@ public:
otb_6s_integer month = static_cast<otb_6s_integer>(m_Month);
int cr(0);
cr = otb_6s_varsol_(&day, &month, &dsol);
- coefTemp = vcl_cos(m_ZenithalSolarRadius)*static_cast<double>(dsol);
+ coefTemp = vcl_cos(m_ZenithalSolarRadius*M_PI/180)*static_cast<double>(dsol);
}
else
{
@@ -201,10 +202,11 @@ public:
}
else
{
- coefTemp = vcl_cos(m_ZenithalSolarRadius)*m_FluxNormalizationCoefficient*m_FluxNormalizationCoefficient;
+ coefTemp = vcl_cos(m_ZenithalSolarRadius*M_PI/180)*m_FluxNormalizationCoefficient*m_FluxNormalizationCoefficient;
}
functor.SetIlluminationCorrectionCoefficient(1. / coefTemp);
functor.SetSolarIllumination(static_cast<double>(m_SolarIllumination[i]));
+
this->GetFunctorVector().push_back(functor);
}
}
diff --git a/Code/Radiometry/otbSIXSTraits.cxx b/Code/Radiometry/otbSIXSTraits.cxx
index d0cec0b1da..5984af6c37 100644
--- a/Code/Radiometry/otbSIXSTraits.cxx
+++ b/Code/Radiometry/otbSIXSTraits.cxx
@@ -60,8 +60,8 @@ SIXSTraits::ComputeAtmosphericParameters(
// atmospheric model
otb_6s_integer iaer(static_cast<otb_6s_integer>(AerosolModel));
otb_6s_real taer55(static_cast<otb_6s_real>(AerosolOptical));
-
-
+
+
// Init output parameters
AtmosphericReflectance = 0.;
AtmosphericSphericalAlbedo = 0.;
@@ -131,7 +131,7 @@ SIXSTraits::ComputeAtmosphericParameters(
AtmosphericSphericalAlbedo = static_cast<double>(sast);
TotalGaseousTransmission = static_cast<double>(tgasm);
DownwardTransmittance = static_cast<double>(sdtott);
- UpwardTransmittance = static_cast<double>(sutott);
+ UpwardTransmittance = static_cast<double>(sutott);
}
diff --git a/Testing/Code/Radiometry/CMakeLists.txt b/Testing/Code/Radiometry/CMakeLists.txt
index 1e5f95af91..4ffe661ca7 100755
--- a/Testing/Code/Radiometry/CMakeLists.txt
+++ b/Testing/Code/Radiometry/CMakeLists.txt
@@ -262,11 +262,13 @@ ADD_TEST(raTvImageToReflectanceImageFilter ${RADIOMETRY_TESTS}
)
ADD_TEST(raTvRomaniaImageToReflectance ${RADIOMETRY_TESTS}
- --compare-image ${EPSILON} ${BASELINE}/raTvRomaniaReflectanceImage.tif
+ --compare-image ${EPSILON} #${BASELINE}/raTvS4_20020518.img.hdr
+ ${BASELINE}/raTvRomania_S4_20020518_Reflectance.tif
${TEMP}/raTvRomaniaReflectanceImage.tif
otbImageToReflectanceImageFilter
${IMAGEDATA}/ROMANIA/S4_20020518.tif
${TEMP}/raTvRomaniaReflectanceImage.tif
+#
27.3 # = 90-62.70 : elevation et azimuth solaire
1.981247824 #channel 1 alpha = 0.881338*2.24800
4.332207085 #channel 2 alpha = 0.858713*5.04500
@@ -384,14 +386,34 @@ ADD_TEST(raTvReflectanceToSurfaceReflectanceImageFilter ${RADIOMETRY_TESTS}
${INPUTDATA}/poupees_sub.png
${TEMP}/raTvReflectanceToSurfaceReflectanceImageFilter.tif
# SAME VALUE FOR EACH CHANNEL
- 1 # intrinsic atmospheric reflectance
- 2 # spherical albedo of the atmosphere
- 3 # total transmission
- 2 # downward transmittance
- 3 # upward transmittance
+ 1 1 1 1 # intrinsic atmospheric reflectance
+ 2 2 2 2 # spherical albedo of the atmosphere
+ 3 3 3 3 # total transmission
+ 2 2 2 2 # downward transmittance
+ 3 3 3 3 # upward transmittance
)
-
+ADD_TEST(raTvRomaniaReflectanceToRomaniaSurfaceReflectanceImageFilter ${RADIOMETRY_TESTS}
+ --compare-image ${EPSILON}
+ ${IMAGEDATA}/ROMANIA/raTvS4_20020518_Corrected.img.hdr
+ ${TEMP}/raTvRomaniaReflectanceToRomaniaSurfaceReflectanceImageFilter.tif
+ otbRomaniaReflectanceToRomaniaSurfaceReflectanceImageFilter
+ #${IMAGEDATA}/ROMANIA/raTvS4_20020518_Reflectance.img.hdr #
+ ${BASELINE}/raTvRomania_S4_20020518_Reflectance.tif
+ ${TEMP}/raTvRomaniaReflectanceToRomaniaSurfaceReflectanceImageFilter.tif
+ ${INPUTDATA}/testRomania.txt # atmo param;
+ ${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_B3.txt # wavelenghts, channel 3
+ ${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_B2.txt # wavelenghts, channel 2
+ ${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_B1.txt # wavelenghts, channel 1
+ ${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_MIR.txt # wavelenghts, channel 4
+
+ # SAME VALUE FOR EACH CHANNEL
+ # 0.02571013197 0.02571013197 0.02571013197 0.02571013197 # intrinsic atmospheric reflectance
+ # 0.06615033001 0.06615033001 0.06615033001 0.06615033001 # spherical albedo of the atmosphere
+ # 0.9276281595 0.9276281595 0.9276281595 0.9276281595 # total transmission
+ # 0.9506403804 0.9506403804 0.9506403804 0.9506403804 # downward transmittance
+ # 0.9373859167 0.9373859167 0.9373859167 0.9373859167 # upward transmittance
+ )
# ------- otb::SurfaceAdjencyEffect6SCorrectionSchemeFilter ------------------------------
ADD_TEST(raTuSurfaceAdjencyEffect6SCorrectionSchemeFilterNew ${RADIOMETRY_TESTS}
@@ -414,7 +436,7 @@ ADD_TEST(raTuSurfaceAdjencyEffect6SCorrectionSchemeFilter ${RADIOMETRY_TESTS}
0 0 0 0 0 1 0 1 0 # Channel 4 ponderation values;
)
-# ------- otb::AtmosphericCorrectionParametersTo6SAtmosphericRadiativeTermsNew ------------------------------
+# ------- otb::AtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms ------------------------------
ADD_TEST(raTuAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTermsNew ${RADIOMETRY_TESTS}
otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTermsNew
)
@@ -425,27 +447,27 @@ ADD_TEST(raTvAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms ${RADI
otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms
${INPUTDATA}/in6S_otb
${TEMP}/raTvCorrectionTo6SRadiative.txt
-
- #${IMAGEDATA}/Romania/S4_2002518.tif # input image
- #${TEMP}/raTvAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms.tif # output image
- #${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_B1.txt # wavelet length file for channel 1
- #${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_B2.txt # wavelet length file for channel 2
- #${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_B3.txt # wavelet length file for channel 3
- #${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_MIR.txt # wavelet length file for channel 4
- #27.3 # = 90-62.7, solar zenithal angle
- #152.7 # solar azimutal angle
- #2 # = 90-88, viewing zenithal angle
- #-77.0 # viewing azimutal angle
- #12 # acquisition month
- #17 # acquisition day
- #1013.00 # atmospheric pressure
- #2.48134 # water vapor amount
- #0.344000 # ozone amount
- #1 # aerosol model type (0 = NO_AEROSOL=0, 1 = CONTINENTAL, 2 = MARITIME, 3 =)NO_AEROSOL=0,CONTINENTAL=1,MARITIME=2,URBAN=3,DESERTIC=5
- #1.99854 # aerosol opical thichness
+)
+
+# ${IMAGEDATA}/Romania/S4_2002518.tif # input image
+# ${TEMP}/raTvAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms.tif # output image
+# ${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_B1.txt # wavelet length file for channel 1
+# ${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_B2.txt # wavelet length file for channel 2
+# ${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_B3.txt # wavelet length file for channel 3
+# ${INPUTDATA}/RADIO_WAVELENGHT_SPECTRAL_BAND_SPOT4_1_MIR.txt # wavelet length file for channel 4
+# 27.3 # = 90-62.7, solar zenithal angle
+# 152.7 # solar azimutal angle
+# 2 # = 90-88, viewing zenithal angle
+# -77.0 # viewing azimutal angle
+# 12 # acquisition month
+# 17 # acquisition day
+# 1013.00 # atmospheric pressure
+# 2.48134 # water vapor amount
+# 0.344000 # ozone amount
+# 1 # aerosol model type (0 = NO_AEROSOL=0, 1 = CONTINENTAL, 2 = MARITIME, 3 =)NO_AEROSOL=0,CONTINENTAL=1,MARITIME=2,URBAN=3,DESERTIC=5
+# 1.99854 # aerosol opical thichness
-)
@@ -479,6 +501,7 @@ otbSurfaceAdjencyEffect6SCorrectionSchemeFilterNew.cxx
otbSurfaceAdjencyEffect6SCorrectionSchemeFilter.cxx
otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTermsNew.cxx
otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms.cxx
+otbRomaniaReflectanceToRomaniaSurfaceReflectanceImageFilter.cxx
)
diff --git a/Testing/Code/Radiometry/otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms.cxx b/Testing/Code/Radiometry/otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms.cxx
index 2e0c4e6342..5a53f80a00 100644
--- a/Testing/Code/Radiometry/otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms.cxx
+++ b/Testing/Code/Radiometry/otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms.cxx
@@ -26,13 +26,6 @@
int otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms(int argc, char * argv[])
{
- /*
- std::vector<const char *> wavelenghFiles;
- wavelenghFiles.push_back( argv[3] );
- wavelenghFiles.push_back( argv[4] );
- wavelenghFiles.push_back( argv[5] );
- wavelenghFiles.push_back( argv[6] );
- */
const char * wavelenghFile = argv[1];
const char * outputFile = argv[2];
@@ -102,18 +95,16 @@ int otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms(int argc, ch
fin >> value;
vect.push_back(value);
}
- // Remove the last vector element which is added by fin, and not contains in the original file.
- //vect.pop_back();
+
fin.close();
functionValues->SetFilterFunctionValues(vect);
functionValues->SetMinSpectralValue(minSpectralValue);
functionValues->SetMaxSpectralValue(maxSpectralValue);
functionValues->SetUserStep( val );
- param->GetWavelenghtSpectralBandRef()->push_back(functionValues);
-
+ param->SetWavelenghtSpectralBandWithIndex(0, functionValues);
//}
- aerosolModel = static_cast<AerosolModelType>(::atoi(argv[16]));
+ //aerosolModel = static_cast<AerosolModelType>(::atoi(argv[16]));
// Set parameters
param->SetSolarZenithalAngle(static_cast<double>(solarZenithalAngle));
@@ -128,24 +119,8 @@ int otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms(int argc, ch
param->SetAerosolModel(aerosolModel);
param->SetAerosolOptical(static_cast<double>(aerosolOptical));
- /*
- aerosolModel = static_cast<AerosolModelType>(::atoi(argv[16]));
-
- // Set parameters
- param->SetSolarZenithalAngle(static_cast<double>(::atof(argv[7])));
- param->SetSolarAzimutalAngle(static_cast<double>(::atof(argv[8])));
- param->SetViewingZenithalAngle(static_cast<double>(::atof(argv[9])));
- param->SetViewingAzimutalAngle(static_cast<double>(::atof(argv[10])));
- param->SetMonth(::atoi(argv[11]));
- param->SetDay(::atoi(argv[12]));
- param->SetAtmosphericPressure(static_cast<double>(::atof(argv[13])));
- param->SetWaterVaporAmount(static_cast<double>(::atof(argv[14])));
- param->SetOzoneAmount(static_cast<double>(::atof(argv[15])));
- param->SetAerosolModel(aerosolModel);
- param->SetAerosolOptical(static_cast<double>(::atof(argv[17])));
- */
object->SetInput( param );
- object->GenerateData();
+ object->Update();
radiative = object->GetOutput();
fout.open(outputFile);
diff --git a/Testing/Code/Radiometry/otbRadiometryTests.cxx b/Testing/Code/Radiometry/otbRadiometryTests.cxx
index 018ee38197..7cab6af0f0 100755
--- a/Testing/Code/Radiometry/otbRadiometryTests.cxx
+++ b/Testing/Code/Radiometry/otbRadiometryTests.cxx
@@ -55,5 +55,6 @@ REGISTER_TEST(otbSurfaceAdjencyEffect6SCorrectionSchemeFilterNew);
REGISTER_TEST(otbSurfaceAdjencyEffect6SCorrectionSchemeFilter);
REGISTER_TEST(otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTermsNew);
REGISTER_TEST(otbAtmosphericCorrectionParametersTo6SAtmosphericRadiativeTerms);
+REGISTER_TEST(otbRomaniaReflectanceToRomaniaSurfaceReflectanceImageFilter);
}
diff --git a/Testing/Code/Radiometry/otbReflectanceToSurfaceReflectanceImageFilter.cxx b/Testing/Code/Radiometry/otbReflectanceToSurfaceReflectanceImageFilter.cxx
index 81ae200f63..db248efb0f 100644
--- a/Testing/Code/Radiometry/otbReflectanceToSurfaceReflectanceImageFilter.cxx
+++ b/Testing/Code/Radiometry/otbReflectanceToSurfaceReflectanceImageFilter.cxx
@@ -49,12 +49,22 @@ int otbReflectanceToSurfaceReflectanceImageFilter(int argc, char * argv[])
reader->GenerateOutputInformation();
unsigned int nbChannel = reader->GetOutput()->GetNumberOfComponentsPerPixel();
-
- DataVectorType intrinsic(nbChannel, static_cast<double>(atof(argv[3])));
- DataVectorType albedo(nbChannel, static_cast<double>(atof(argv[4])));
- DataVectorType gaseous(nbChannel, static_cast<double>(atof(argv[5])));
- DataVectorType downTrans(nbChannel, static_cast<double>(atof(argv[6])));
- DataVectorType upTrans(nbChannel, static_cast<double>(atof(argv[7])));
+
+ DataVectorType intrinsic;
+ DataVectorType albedo;
+ DataVectorType gaseous;
+ DataVectorType downTrans;
+ DataVectorType upTrans;
+
+ std::cout<<nbChannel<<std::endl;
+ for( unsigned int j=0; j<nbChannel; j++)
+ {
+ intrinsic.push_back(static_cast<double>(atof(argv[3+j])));
+ albedo.push_back(static_cast<double>(atof(argv[3+j+nbChannel])));
+ gaseous.push_back(static_cast<double>(atof(argv[3+j+2*nbChannel])));
+ downTrans.push_back(static_cast<double>(atof(argv[3+j+3*nbChannel])));
+ upTrans.push_back(static_cast<double>(atof(argv[3+j+4*nbChannel])));
+ }
atmo->SetIntrinsicAtmosphericReflectances(intrinsic);
atmo->SetSphericalAlbedos(albedo);
--
GitLab