Commit 753dd44b authored by Cédric Traizet's avatar Cédric Traizet

Merge branch 'app_mosaic' into 'develop'

Move Mosaic in OTB

See merge request !596
parents 1a84c84e c490d43f
Pipeline #2684 failed with stages
in 50 minutes and 55 seconds
......@@ -70,7 +70,7 @@ if(XDK_PATH)
set(cmake_configure_option
"${cmake_configure_option}
CMAKE_PREFIX_PATH=${XDK_PATH}")
foreach(remote_module OTBTemporalGapFilling Mosaic SertitObject DiapOTBModule)#otbGRM
foreach(remote_module OTBTemporalGapFilling SertitObject DiapOTBModule)#otbGRM
set(cmake_configure_option
"${cmake_configure_option}
Module_${remote_module}:BOOL=ON")
......
......@@ -185,6 +185,7 @@ r'''(.*
)+(# Copyright \(C\) 1999-2011 Insight Software Consortium
|# Copyright \(C\) 20\d\d(-20\d\d)? Centre National d'Etudes Spatiales \(CNES\)
|# Copyright \(C\) 20\d\d(-20\d\d)? CS Systemes d'Information \(CS SI\)
|# Copyright \(C\) 20\d\d(-20\d\d)? IRSTEA
|# Copyright \(C\) 2007-2012 Institut Mines Telecom / Telecom Bretagne
)+#
# This file is part of Orfeo Toolbox
......
......@@ -87,3 +87,8 @@ otb_create_application(
NAME DynamicConvert
SOURCES otbDynamicConvert.cxx
LINK_LIBRARIES ${${otb-module}_LIBRARIES})
OTB_CREATE_APPLICATION(
NAME Mosaic
SOURCES otbMosaic.cxx
LINK_LIBRARIES ${${otb-module}_LIBRARIES})
This diff is collapsed.
......@@ -39,6 +39,7 @@ otb_module(OTBAppImageUtils
OTBStreaming
OTBTransform
OTBFunctor
OTBMosaic
TEST_DEPENDS
OTBCommandLine
......
......@@ -455,3 +455,63 @@ otb_test_application(NAME apTvUtSplitImage
${INPUTDATA}/poupees_sub_c3.png
${TEMP}/apTvUtSplitImageOutput_2.tif)
#----------- Mosaic TESTS ----------------
otb_test_application(NAME MosaicTestLargeFeathering
APP Mosaic
OPTIONS -il ${INPUTDATA}/SP67_FR_subset_1.tif ${INPUTDATA}/SP67_FR_subset_2.tif
-out ${TEMP}/apTvMosaicTestLargeFeathering.tif uint8
-comp.feather large
VALID --compare-image ${EPSILON_8}
${BASELINE}/apTvMosaicTestLargeFeathering.tif
${TEMP}/apTvMosaicTestLargeFeathering.tif)
otb_test_application(NAME MosaicTestSlimFeathering
APP Mosaic
OPTIONS -il ${INPUTDATA}/SP67_FR_subset_1.tif ${INPUTDATA}/SP67_FR_subset_2.tif
-out ${TEMP}/apTvMosaicTestSlimFeathering.tif uint8
-comp.feather slim
-comp.feather.slim.lenght 100
VALID --compare-image ${EPSILON_8}
${BASELINE}/apTvMosaicTestSlimFeathering.tif
${TEMP}/apTvMosaicTestSlimFeathering.tif)
otb_test_application(NAME MosaicTestSimpleWithHarmoBandRmse
APP Mosaic
OPTIONS -il ${INPUTDATA}/SP67_FR_subset_1.tif ${INPUTDATA}/SP67_FR_subset_2.tif
-out ${TEMP}/apTvMosaicTestSimpleWithHarmoBandRmse.tif uint8
-harmo.method band
-harmo.cost rmse
VALID --compare-image ${EPSILON_8}
${BASELINE}/apTvMosaicTestSimpleWithHarmoBandRmse.tif
${TEMP}/apTvMosaicTestSimpleWithHarmoBandRmse.tif)
otb_test_application(NAME MosaicTestSimpleWithHarmoRgbRmse
APP Mosaic
OPTIONS -il ${INPUTDATA}/SP67_FR_subset_1.tif ${INPUTDATA}/SP67_FR_subset_2.tif
-out ${TEMP}/apTvMosaicTestSimpleWithHarmoRgbRmse.tif uint8
-harmo.method rgb
-harmo.cost rmse
VALID --compare-image ${EPSILON_8}
${BASELINE}/apTvMosaicTestSimpleWithHarmoRgbRmse.tif
${TEMP}/apTvMosaicTestSimpleWithHarmoRgbRmse.tif)
otb_test_application(NAME MosaicTestSimpleWithCutline
APP Mosaic
OPTIONS -il ${INPUTDATA}/SP67_FR_subset_1.tif ${INPUTDATA}/SP67_FR_subset_2.tif
-out ${TEMP}/apTvMosaicTestSimpleWithCutline.tif uint8
-vdcut ${INPUTDATA}/SP67_FR_subset_1_cutline.shp ${INPUTDATA}/SP67_FR_subset_2_cutline.shp
VALID --compare-image ${EPSILON_8}
${BASELINE}/apTvMosaicTestSimpleWithCutline.tif
${TEMP}/apTvMosaicTestSimpleWithCutline.tif)
otb_test_application(NAME MosaicTestSimpleWithVdstats
APP Mosaic
OPTIONS -il ${INPUTDATA}/SP67_FR_subset_1.tif ${INPUTDATA}/SP67_FR_subset_2.tif
-out ${TEMP}/apTvMosaicTestSimpleWithVdstats.tif uint8
-vdstats ${INPUTDATA}/SP67_FR_subset_1_cutline.shp ${INPUTDATA}/SP67_FR_subset_2_cutline.shp
VALID --compare-image ${EPSILON_8}
${BASELINE}/apTvMosaicTestSimpleWithVdstats.tif
${TEMP}/apTvMosaicTestSimpleWithVdstats.tif)
#
# Copyright (C) 2005-2019 Centre National d'Etudes Spatiales (CNES)
# Copyright (C) 2019 IRSTEA
#
# This file is part of Orfeo Toolbox
#
# https://www.orfeo-toolbox.org/
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
# http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
#
project(OTBMosaic)
otb_module_impl()
/*
* Copyright (C) 1999-2011 Insight Software Consortium
* Copyright (C) 2005-2019 Centre National d'Etudes Spatiales (CNES)
* Copyright (C) 2016-2019 IRSTEA
*
* This file is part of Orfeo Toolbox
*
* https://www.orfeo-toolbox.org/
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
#ifndef MODULES_REMOTE_MOSAIC_INCLUDE_OTBMOSAICFUNCTORS_H_
#define MODULES_REMOTE_MOSAIC_INCLUDE_OTBMOSAICFUNCTORS_H_
#include "vnl/vnl_matrix.h"
#include "vcl_compiler.h"
namespace otb
{
namespace Functor
{
/**
* \class RGB2LAB Functor
* \brief Base class for converting RGB into LAB color space (Ruderman et al.)
*
* \ingroup OTBMosaic
*/
template <class TInput, class TOutput>
class RGB2LAB
{
public:
RGB2LAB()
{
M.set_size(3, 3);
M[0][0] = 0.3811;
M[0][1] = 0.5783;
M[0][2] = 0.0406;
M[1][0] = 0.1967;
M[1][1] = 0.7244;
M[1][2] = 0.0790;
M[2][0] = 0.0241;
M[2][1] = 0.1288;
M[2][2] = 0.8531;
D1.set_size(3, 3);
D1.fill(0.0);
D1[0][0] = 1.0 / vcl_sqrt(3.0);
D1[1][1] = 1.0 / vcl_sqrt(6.0);
D1[2][2] = 1.0 / vcl_sqrt(2.0);
D2.set_size(3, 3);
D2.fill(1.0);
D2[1][2] = -2.0;
D2[2][1] = -1.0;
D2[2][2] = 0.0;
}
~RGB2LAB()
{
}
bool operator!=(const RGB2LAB&) const
{
return false;
}
bool operator==(const RGB2LAB& other) const
{
return !(*this != other);
}
inline TOutput operator()(const TInput& A) const
{
TOutput output;
output.SetSize(3);
if (A[0] == 0 && A[1] == 0 && A[2] == 0)
{
output.Fill(0);
return output;
}
// RGB
vnl_matrix<double> rgb(3, 1);
rgb[0][0] = A[0];
rgb[1][0] = A[1];
rgb[2][0] = A[2];
// LMS
vnl_matrix<double> lms(3, 1);
lms = M * rgb;
// LMS (log10)
const double log10 = vcl_log(10);
lms[0][0] = vcl_log(lms[0][0]) / log10;
lms[1][0] = vcl_log(lms[1][0]) / log10;
lms[2][0] = vcl_log(lms[2][0]) / log10;
// LAB
vnl_matrix<double> lab(3, 1);
lab = D1 * (D2 * lms);
output[0] = lab[0][0];
output[1] = lab[1][0];
output[2] = lab[2][0];
return output;
}
size_t OutputSize(const std::array<size_t, 1>&) const
{
return 3;
}
private:
vnl_matrix<double> M;
vnl_matrix<double> D1;
vnl_matrix<double> D2;
};
/**
* \class LAB2RGB Functor
* \brief Base class for converting LAB into RGB color space (Ruderman et al.)
*
* TODO: invert the function RGB2LAB than using the hardcoded one
*
* \ingroup OTBMosaic
*/
template <class TInput, class TOutput>
class LAB2RGB
{
public:
LAB2RGB()
{
M.set_size(3, 3);
M[0][0] = 4.4687;
M[0][1] = -3.5887;
M[0][2] = 0.1197;
M[1][0] = -1.2197;
M[1][1] = 2.3831;
M[1][2] = -0.1626;
M[2][0] = 0.0579;
M[2][1] = -0.2584;
M[2][2] = 1.1934;
D1.set_size(3, 3);
D1.fill(0.0);
D1[0][0] = 1.0 / vcl_sqrt(3.0);
D1[1][1] = 1.0 / vcl_sqrt(6.0);
D1[2][2] = 1.0 / vcl_sqrt(2.0);
D2.set_size(3, 3);
D2.fill(1.0);
D2[1][2] = -1.0;
D2[2][1] = -2.0;
D2[2][2] = 0.0;
}
~LAB2RGB()
{
}
bool operator!=(const LAB2RGB&) const
{
return false;
}
bool operator==(const LAB2RGB& other) const
{
return !(*this != other);
}
inline TOutput operator()(const TInput& A) const
{
TOutput output;
output.SetSize(3);
if (A[0] == 0 && A[1] == 0 && A[2] == 0)
{
output.Fill(0);
return output;
}
// LAB
vnl_matrix<double> lab(3, 1);
lab[0][0] = A[0];
lab[1][0] = A[1];
lab[2][0] = A[2];
// LMS
vnl_matrix<double> lms(3, 1);
lms = D2 * (D1 * lab);
lms[0][0] = vcl_pow(10.0, lms[0][0]);
lms[1][0] = vcl_pow(10.0, lms[1][0]);
lms[2][0] = vcl_pow(10.0, lms[2][0]);
// RGB
vnl_matrix<double> rgb(3, 1);
rgb = M * lms;
output[0] = rgb[0][0];
output[1] = rgb[1][0];
output[2] = rgb[2][0];
return output;
}
inline size_t OutputSize(const std::array<size_t, 1>&) const
{
return 3;
}
private:
vnl_matrix<double> M;
vnl_matrix<double> D1;
vnl_matrix<double> D2;
};
} // namespace functor
} // namespace otb
#endif /* MODULES_REMOTE_MOSAIC_INCLUDE_OTBMOSAICFUNCTORS_H_ */
/*
* Copyright (C) 1999-2011 Insight Software Consortium
* Copyright (C) 2005-2019 Centre National d'Etudes Spatiales (CNES)
* Copyright (C) 2016-2019 IRSTEA
*
* This file is part of Orfeo Toolbox
*
* https://www.orfeo-toolbox.org/
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
#ifndef MODULES_REMOTE_MOSAIC_INCLUDE_OTBPERSISTENTMOSAICFILTER_H_
#define MODULES_REMOTE_MOSAIC_INCLUDE_OTBPERSISTENTMOSAICFILTER_H_
#include "otbStreamingMosaicFilterBase.h"
namespace otb
{
/** \class PersistentMosaicFilter
* \brief This filter is the base class for all mosaic filter persisting data through multiple
* update.
* For instance, a filter computing global statistics on an mosaic with streaming
* capabilities will have to keep the temporary results for each streamed piece of the
* image in order to synthesize the global statistics at the end. This filter is an
* itk::ImageToImageFilter, providing two additional methods. The first one, Synthetize(),
* allows the user to synthesize temporary data produced by the multiple updates on different
* pieces of the image to the global result. The second one, Reset(), allows the user to
* reset the temporary data for a new input image for instance.
*
* \note This class contains pure virtual method, and can not be instantiated.
*
* \sa PersistentMosaicFilter
* \sa StreamingStatisticsMosaicFilter
*
* \ingroup OTBMosaic
*/
template <class TInputImage, class TOutputImage, class TPrecisionType>
class ITK_EXPORT PersistentMosaicFilter : public otb::StreamingMosaicFilterBase<TInputImage, TOutputImage, TPrecisionType>
{
public:
/** Standard typedefs */
typedef PersistentMosaicFilter Self;
typedef otb::StreamingMosaicFilterBase<TInputImage, TOutputImage, TPrecisionType> Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
/** Creation through object factory macro */
itkTypeMacro(PersistentMosaicFilter, StreamingMosaicFilterBase);
/** Input image typedefs. */
typedef typename Superclass::InputImageType InputImageType;
typedef typename Superclass::InputImagePointer InputImagePointer;
typedef typename Superclass::InputImagePointType InputImagePointType;
typedef typename Superclass::InputImagePixelType InputImagePixelType;
typedef typename Superclass::InputImageIndexType InputImageIndexType;
typedef typename Superclass::InputImageSizeType InputImageSizeType;
typedef typename Superclass::InputImageSpacingType InputImageSpacingType;
typedef typename Superclass::InputImageInternalPixelType InputImageInternalPixelType;
typedef typename Superclass::InputImageRegionType InputImageRegionType;
/** Output image typedefs. */
typedef typename Superclass::OutputImageType OutputImageType;
typedef typename Superclass::OutputImagePointer OutputImagePointer;
typedef typename Superclass::OutputImagePointType OutputImagePointType;
typedef typename Superclass::OutputImagePixelType OutputImagePixelType;
typedef typename Superclass::OutputImageIndexType OutputImageIndexType;
typedef typename Superclass::OutputImageSizeType OutputImageSizeType;
typedef typename Superclass::OutputImageSpacingType OutputImageSpacingType;
typedef typename Superclass::OutputImageInternalPixelType OutputImageInternalPixelType;
typedef typename Superclass::OutputImageRegionType OutputImageRegionType;
/** Internal computing typedef support. */
typedef typename Superclass::InternalValueType InternalValueType;
typedef typename Superclass::ContinuousIndexType ContinuousIndexType;
typedef typename Superclass::InterpolatorType InterpolatorType;
typedef typename Superclass::InterpolatorPointerType InterpolatorPointerType;
typedef typename Superclass::DefaultInterpolatorType DefaultInterpolatorType;
typedef typename Superclass::InternalImageType InternalImageType;
typedef typename Superclass::InternalPixelType InternalPixelType;
typedef typename Superclass::IteratorType IteratorType;
typedef typename Superclass::ConstIteratorType ConstIteratorType;
typedef typename Superclass::StreamingTraitsType StreamingTraitsType;
/**
* Reset the persistent data of the filter.
*/
virtual void Reset(void) = 0;
/**
* Synthesize the persistent data of the filter.
*/
virtual void Synthetize(void) = 0;
protected:
/** Constructor */
PersistentMosaicFilter()
{
}
/** Destructor */
~PersistentMosaicFilter() override
{
}
/**PrintSelf method */
void PrintSelf(std::ostream& os, itk::Indent indent) const override
{
Superclass::PrintSelf(os, indent);
}
private:
PersistentMosaicFilter(const Self&); // purposely not implemented
void operator=(const Self&); // purposely not implemented
};
} // End namespace otb
#endif /* MODULES_REMOTE_MOSAIC_INCLUDE_OTBPERSISTENTMOSAICFILTER_H_ */
/*
* Copyright (C) 1999-2011 Insight Software Consortium
* Copyright (C) 2005-2019 Centre National d'Etudes Spatiales (CNES)
* Copyright (C) 2016-2019 IRSTEA
*
* This file is part of Orfeo Toolbox
*
* https://www.orfeo-toolbox.org/
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
#ifndef QuadraticallyConstrainedSimpleSolver_H_
#define QuadraticallyConstrainedSimpleSolver_H_
#include "itkObjectFactory.h"
#include "itkLightObject.h"
#include "itkNumericTraits.h"
#include <vnl/vnl_matrix.h>
#include "vnl/algo/vnl_solve_qp.h"
namespace otb
{
/**
* \class QuadraticallyConstrainedSimpleSolver
* \brief Solves the optimisation problem for radiometric harmonisation of multiple
* overlapping images.
*
* This solver inputs statistics of the overlapping images, and produces a
* zero-y intercept linear correction model. Various cost functions can be
* employed: RMSE based, Mean based, Mean+Standard deviation based, and Mean
* + weighted Standard deviation bases.
*
* Inputs:
* -N x N Matrix of mean of overlaps ij
* -N x N Matrix of standard deviation of overlaps ij
* -N x N Matrix of area of overlaps ij
* -N x N Matrix of mean of pixels products of overlaps ij
*
* For all i and j, m_{ij} = stats of image i in overlap ij
*
* Output:
* N x 1 Vector of scales to apply to images
*
* For more details, see Cresson, Remi, and Nathalie Saint-Geours.
* "Natural color satellite image mosaicking using quadratic programming in decorrelated color space."
* IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing 8.8 (2015): 4151-4162.
*
* https://doi.org/10.1109/JSTARS.2015.2449233
* https://hal.archives-ouvertes.fr/hal-01373314/file/cresson2015.pdf
*
* \ingroup OTBMosaic
*/
template <class ValueType>
class ITK_EXPORT QuadraticallyConstrainedSimpleSolver : public itk::LightObject
{
public:
/** Standard class typedef */
typedef QuadraticallyConstrainedSimpleSolver Self;
typedef itk::LightObject Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
/** Runtime information support. */
itkTypeMacro(QuadraticallyConstrainedSimpleSolver, LightObject);
/** Method for creation through the object factory. */
itkNewMacro(Self);
/** Typedefs */
typedef vnl_matrix<ValueType> RealMatrixType;
typedef vnl_vector<ValueType> RealVectorType;
typedef vnl_matrix<double> DoubleMatrixType;
typedef vnl_vector<double> DoubleVectorType;
typedef std::vector<unsigned int> ListIndexType;
/** Enum for objective function type */
enum ObjectiveFunctionType
{
Cost_Function_rmse, // Root mean square error based
Cost_Function_musig, // Mean and standard deviation based
Cost_Function_mu, // Mean based
Cost_Function_weighted_musig // Mean and weighted standard deviation based
};
/** Mean-in-overlaps matrix */
void SetMeanInOverlaps(const RealMatrixType& matrix)
{
m_MeanInOverlaps = RealMatrixType(matrix);
}
const RealMatrixType GetMeanInOverlaps()
{
return m_MeanInOverlaps;
}
/** Standard-deviation-in-overlaps matrix */
void SetStandardDeviationInOverlaps(const RealMatrixType& matrix)
{
m_StandardDeviationInOverlaps = RealMatrixType(matrix);
}
const RealMatrixType GetStandardDeviationInOverlaps()
{
return m_StandardDeviationInOverlaps;
}
/** Area-in-overlaps matrix */
void SetAreaInOverlaps(const RealMatrixType& matrix)
{
m_AreaInOverlaps = RealMatrixType(matrix);
}
const RealMatrixType GetAreaInOverlaps()
{
return m_AreaInOverlaps;
}
/** Mean-of-pixels-products-in-overlaps matrix */
void SetMeanOfProductsInOverlaps(const RealMatrixType& matrix)
{
m_MeanOfProductsInOverlaps = RealMatrixType(matrix);
}
const RealMatrixType GetMeanOfProductsInOverlaps()
{
return m_MeanOfProductsInOverlaps;
}
/** Output correction model */
const RealVectorType GetOutputCorrectionModel()
{
return m_OutputCorrectionModel;
}
/**
* STD weight in harmonization
* if value is near 0, importance is accorded to MEAN
* if value is 1, importance is the same than MEAN
* if value is higher than 1, importance is accorder to STD
*/
void SetWeightOfStandardDeviationTerm(ValueType weight)
{
m_WeightOfStandardDeviationTerm = weight;
}
ValueType GetWeightOfStandardDeviationTerm()
{
return m_WeightOfStandardDeviationTerm;
}
/** Solving routine */
void Solve();
/** Set the cost function type */
void SetMeanBased()
{
oft = Cost_Function_mu;
}
void SetMeanAndStandardDeviationBased()
{
oft = Cost_Function_musig;
}
void SetRMSEBased()
{
oft = Cost_Function_rmse;
}
void SetWeightedMeanAndStandardDeviationBased()
{
oft = Cost_Function_weighted_musig;
}
protected:
QuadraticallyConstrainedSimpleSolver();
virtual ~QuadraticallyConstrainedSimpleSolver();
private:
// Check inputs
void CheckInputs(void) const;
// Deep First Search
void DFS(std::vector<bool>& marked, unsigned int s) const;
// Compute the objective matrix
const DoubleMatrixType GetQuadraticObjectiveMatrix(const DoubleMatrixType& areas, const DoubleMatrixType& means, const DoubleMatrixType& stds,
const DoubleMatrixType& mops);
// Extract a sub matrix from indices list
const DoubleMatrixType ExtractMatrix(const RealMatrixType& mat, const ListIndexType& idx);
// Input
RealMatrixType m_MeanInOverlaps;
RealMatrixType m_StandardDeviationInOverlaps;
RealMatrixType m_AreaInOverlaps;
RealMatrixType m_MeanOfProductsInOverlaps;
// Params
</