Commit 541e730d authored by Luc Hermitte's avatar Luc Hermitte

Move ClampROI and Synthetize to OTB

parent 7a9cb051
......@@ -12,9 +12,6 @@ This module can be built as a remote module or as a standalone module
The module provides the following applications:
- ClampROI: that sets margins to 0 outside the ROI,
- Synthetize: that iterates overs a list of overlapping images to keep the
first non null pixel
- MultitempFilteringFilter: that implements Quegan speckle filter for
SAR Images.
- MultitempFilteringOutcore: that computes the outcore of the filter.
......
......@@ -18,16 +18,6 @@
# limitations under the License.
#
otb_create_application(
NAME ClampROI
SOURCES otbClampROI.cxx
LINK_LIBRARIES ${${otb-module}_LIBRARIES})
otb_create_application(
NAME Synthetize
SOURCES otbSynthetize.cxx
LINK_LIBRARIES ${${otb-module}_LIBRARIES})
otb_create_application(
NAME MultitempFilteringOutcore
SOURCES otbMultitempFilteringOutcore.cxx
......
/*
* Copyright(C) 2005-2020 Centre National d'Etudes Spatiales(CNES)
*
* This file is part of S1Tiling remote module for 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.
*/
#include "otbClampROIFilter.h"
#include "otbWrapperApplication.h"
#include "otbWrapperApplicationFactory.h"
namespace otb
{
namespace Wrapper
{
/**
* Application that fills margins to 0.
*
* This application is similar to ExtractROI with the difference the margin is
* kept, and filled with 0.
*
* This application is used to implement the _cut_ processing in S1Tiling
* chain.
*
* \author Luc Hermitte (CS Group)
* \copyright CNES
*/
class ClampROI : public Application
{
public:
using Self = ClampROI;
using Pointer = itk::SmartPointer<Self>;
itkNewMacro(Self);
itkTypeMacro(ClampROI, otb::Wrapper::Application);
private:
void DoInit() override
{
SetName("ClampROI");
SetDescription("This is the ClampROI application, version X.X.X");
SetDocLongDescription(
"This application is similar to ExtractROI in the sense it extracts a Region of Interrest.\n"
"However, the region outside of the ROI isn't trimmed, but set to 0.\n"
"\n"
"The filter set lines of index < threshold.y, and of index >= threshold.y to 0\n"
"The filter set columns of index < threshold.x, and of index >= threshold.x to 0");
SetDocLimitations("None");
SetDocAuthors("Luc Hermitte (CS Group)");
SetDocSeeAlso("");
AddDocTag("S1Tiling"); // any better idea?
AddParameter(ParameterType_InputImage, "in", "Input image");
SetParameterDescription("in", "Input image");
AddParameter(ParameterType_OutputImage, "out", "Output Image");
SetParameterDescription("out","Output image");
AddParameter(ParameterType_Group, "threshold", "threshold group");
AddParameter(ParameterType_Group, "threshold.y", "threshold group");
MandatoryOff("threshold");
MandatoryOff("threshold.y");
AddParameter(ParameterType_Int, "threshold.x", "Column index threshold");
SetParameterDescription("threshold.x", "Column index threshold");
SetDefaultParameterInt("threshold.x", 0);
AddParameter(ParameterType_Int, "threshold.y.start", "Top line index threshold");
SetParameterDescription("threshold.y.start", "Top line index threshold");
SetDefaultParameterInt("threshold.y.start", 0);
AddParameter(ParameterType_Int, "threshold.y.end", "Bottom line index threshold");
SetParameterDescription("threshold.y.end", "Bottom line index threshold");
SetDefaultParameterInt("threshold.y.end", 0);
AddRAMParameter();
}
void DoUpdateParameters() override
{}
void DoExecute() override
{
auto const thrX = GetParameterInt("threshold.x");
auto const thrYtop = GetParameterInt("threshold.y.start");
auto const thrYbot = GetParameterInt("threshold.y.end");
if (thrX < 0)
itkExceptionMacro("The column threshold is expected to be positive");
if (thrYtop < 0)
itkExceptionMacro("The top line threshold is expected to be positive");
if (thrYbot < 0)
itkExceptionMacro("The bottom line threshold is expected to be positive");
if (thrX == 0 && thrYtop == 0 && thrYbot == 0)
itkExceptionMacro("Don't use ClampROI to clamp nothing!");
auto filter = ClampROIFilter<FloatImageType>::New();
assert(thrX >= 0);
assert(thrYtop >= 0);
assert(thrYbot >= 0);
filter->SetThresholdX(thrX);
filter->SetThresholdYtop(thrYtop);
filter->SetThresholdYbot(thrYbot);
filter->SetInput(GetParameterFloatImage("in"));
SetParameterOutputImage("out", filter->GetOutput());
RegisterPipeline();
}
};
} // otb::Wrapper namespace
} // otb namespace
OTB_APPLICATION_EXPORT(otb::Wrapper::ClampROI)
/*
* Copyright(C) 2005-2020 Centre National d'Etudes Spatiales(CNES)
*
* This file is part of S1Tiling remote module for 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.
*/
#include "otbSynthetizeFilter.h"
#include "otbWrapperApplication.h"
#include "otbWrapperApplicationFactory.h"
#include "otbImageFileReader.h"
#include <set>
namespace otb
{
namespace Wrapper
{
/**
* This application synthetizes/reduces multiple inputs into a single one.
* In that particular case, for each output pixel, this application will
* consider the corresponding pixels from all the input images, and keep the
* first one that isn't equal to 0.
*
* This application is used to implement the _concatenate_ processing in
* S1Tiling chain.
*
* \author Luc Hermitte (CS Group)
* \copyright CNES
* \todo find a better name for the application. Alas `otbConcatenate` is
* already used...
*/
class Synthetize : public Application
{
public:
using Self = Synthetize;
using Pointer = itk::SmartPointer<Self>;
itkNewMacro(Self);
itkTypeMacro(Synthetize, otb::Wrapper::Application);
private:
using ReaderType = otb::ImageFileReader<FloatImageType>;
void DoInit() override
{
SetName("Synthetize");
SetDescription("This is the Synthetize application, version X.X.X");
SetDocLongDescription("Concatenate a list of images of the same size into a single single-channel image.\n\
It keeps the first non-null pixel value found in the input list.");
SetDocLimitations("None");
SetDocAuthors("Luc Hermitte (CS Group)");
SetDocSeeAlso("");
AddDocTag("S1Tiling"); // any better idea?
AddParameter(ParameterType_StringList, "il", "Input images list");
SetParameterDescription("il", "Input image list");
AddParameter(ParameterType_OutputImage, "out", "Output Image");
SetParameterDescription("out","Output image.");
AddRAMParameter();
}
void DoUpdateParameters() override
{}
void DoExecute() override
{
// Get the input image list
auto inNameList = GetParameterStringList("il");
// checking the input images list validity
auto const nbImages = inNameList.size();
if (nbImages == 0)
{
itkExceptionMacro("No input Image set...; please set at least one input image");
}
auto functor = [](auto input) {
assert(!input.empty());
auto const wh = std::find_if(
input.begin(), input.end()-1,
[](auto v){ return v != 0;});
return *wh;
};
auto filter = MakeSynthetizeFilter<FloatImageType, FloatImageType>(functor);
for (unsigned int i = 0; i < nbImages; i++)
{
// Given the explicit use of a Reader, this application cannot be used in
// a in-memory pipeline
auto reader = ReaderType::New();
// currentImage->SetExtendedFileName(inNameList[i]);
reader->SetFileName(inNameList[i]);
auto currentImage = reader->GetOutput();
currentImage->UpdateOutputInformation();
otbAppLogINFO(<< "Image #" << i + 1 << " has " << currentImage->GetNumberOfComponentsPerPixel() << " components");
filter->SetInput(i, currentImage);
m_Cache.insert(reader);
}
SetParameterOutputImage("out", filter->GetOutput());
RegisterPipeline(); // TODO: check!!
}
// Needed to register the inputs handled manually
// and not with a VectorImageList through GetParameterImageList
std::set<ReaderType::Pointer> m_Cache;
};
} // otb::Wrapper namespace
} // otb namespace
OTB_APPLICATION_EXPORT(otb::Wrapper::Synthetize)
/*
* Copyright(C) 2005-2020 Centre National d'Etudes Spatiales(CNES)
*
* This file is part of S1Tiling remote module for 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 otbClampROIFilter_h
#define otbClampROIFilter_h
#include "itkImageToImageFilter.h"
namespace otb
{
/**
* Region clamping filter.
* This filter is a kind of ROI pass filter. Data within the ROI is kept with
* its original value. Data outside ROI is forced to 0.
*
* Also, this filter propagate the exact ROI upstream in the pipeline. This
* way, if it's piped after another filter, the upstream filter isn't executed
* on the data outside the ROI.
*
* \tparam TImage Image type.
* \sa `otb::ExtractROI<>`
* \author Luc Hermitte (CS Group)
* \copyright CNES
*/
template <typename TImage>
class ClampROIFilter : public itk::ImageToImageFilter<TImage, TImage>
{
public:
/**\name Convenient typedefs for simplifying declarations */
//@{
using InputImageType = TImage;
using OutputImageType = TImage;
//@}
/**\name Extract dimension from input and output images */
//@{
itkStaticConstMacro(InputImageDimension, unsigned int, InputImageType::ImageDimension);
itkStaticConstMacro(OutputImageDimension, unsigned int, OutputImageType::ImageDimension);
//@}
/**\name Standard class typedefs */
//@{
using Self = ClampROIFilter;
using Superclass = itk::ImageToImageFilter<InputImageType, OutputImageType>;
using Pointer = itk::SmartPointer<Self>;
using ConstPointer = itk::SmartPointer<const Self>;
//@}
/** Method for creation through the object factory. */
itkNewMacro(Self);
/** Run-time type information (and related methods). */
itkTypeMacro(ClampROIFilter, unused);
/**\name Image typedef support */
//@{
using InputPixelType = typename InputImageType::PixelType;
using OutputPixelType = typename OutputImageType::PixelType;
using InputRealType = typename itk::NumericTraits<InputPixelType>::RealType;
using InputImageRegionType = typename InputImageType::RegionType;
using OutputImageRegionType = typename OutputImageType::RegionType;
using InputIndexType = typename InputImageType::IndexType;
using InputSizeType = typename InputImageType::SizeType;
using OutputIndexType = typename OutputImageType::IndexType;
using OutputSizeType = typename OutputImageType::SizeType;
static_assert(InputImageDimension == OutputImageDimension, "Images have the same number of components");
//@}
/** Column threshold setter. */
void SetThresholdX(long threshold) noexcept
{ m_thresholdX = threshold; }
/** Column threshold getter. */
long GetThresholdX() const noexcept
{ return m_thresholdX;}
/** Top line threshold setter. */
void SetThresholdYtop(long threshold) noexcept
{ m_thresholdYtop = threshold; }
/** Top line threshold getter. */
long GetThresholdYtop() const noexcept
{ return m_thresholdYtop;}
/** Bottom line threshold setter. */
void SetThresholdYbot(long threshold) noexcept
{ m_thresholdYbot = threshold; }
/** Bottom line threshold getter. */
long GetThresholdYbot() const noexcept
{ return m_thresholdYbot;}
protected:
/// Hidden constructor
ClampROIFilter() = default;
InputImageType * GetInputImage() { return const_cast<InputImageType*>(this->GetInput()); }
InputImageType const* GetInputImage() const { return this->GetInput(); }
/** otbClampROIFilter doesn't need an input requested region as large as the
* output requested region.
* \sa ImageToImageFilter::GenerateInputRequestedRegion()
*/
void CallCopyOutputRegionToInputRegion(
InputImageRegionType & destRegion,
OutputImageRegionType const& srcRegion) override
{
destRegion = OutputRegionToInputRegion(srcRegion);
}
/**
* Functional implementation of `CallCopyOutputRegionToInputRegion()`.
*/
InputImageRegionType OutputRegionToInputRegion(
OutputImageRegionType const& srcRegion);
/**
* Main computation function called by each thread.
* \param[in] outputRegionForThread Specified output region to compute
* \param[in] threadId Id of the computing threads
*/
void ThreadedGenerateData(
OutputImageRegionType const& outputRegionForThread,
itk::ThreadIdType threadId) override;
private:
long m_thresholdX = 0;
long m_thresholdYtop = 0;
long m_thresholdYbot = 0;
};
} // otb namespace
#ifndef OTB_MANUAL_INSTANTIATION
#include "otbClampROIFilter.hxx"
#endif
#endif // otbClampROIFilter_h
/*
* Copyright(C) 2005-2020 Centre National d'Etudes Spatiales(CNES)
*
* This file is part of S1Tiling remote module for 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 otbClampROIFilter_hxx
#define otbClampROIFilter_hxx
#include "otbClampROIFilter.h"
#include "otbInterval.h"
#include "otbMacro.h"
#include "otbLogHelpers.h"
#include "itkImageScanlineConstIterator.h"
#include "itkImageScanlineIterator.h"
#include "itkProgressReporter.h"
#include <boost/numeric/interval.hpp>
#include <algorithm>
#include <cassert>
#include <ostream>
template <typename T, typename P>
inline
std::ostream & operator<<(std::ostream & os, boost::numeric::interval<T,P> const& v)
{
return os << '[' << v.lower() << ".." << v.upper() << '[';
}
namespace otb
{
template<typename TImage>
void
ClampROIFilter<TImage>
::ThreadedGenerateData(
OutputImageRegionType const& outputRegionForThread,
itk::ThreadIdType threadId)
{
// otbMsgDevMacro("ThreadedGenerateData begin("<<NeatRegionLogger{outputRegionForThread}<<")");
using InputIterator = itk::ImageScanlineConstIterator<InputImageType const>;
using OutputIterator = itk::ImageScanlineIterator<OutputImageType>;
auto const* input = this->GetInput();
auto * output = this->GetOutput();
assert(input);
assert(output);
InputIterator inputIterator (input, OutputRegionToInputRegion(outputRegionForThread));
OutputIterator outputIterator(output, outputRegionForThread);
auto const& imgRegion = output->GetLargestPossibleRegion();
auto const imgSizeX = imgRegion.GetSize()[0];
auto const imgSizeY = imgRegion.GetSize()[1];
itk::IndexValueType const imgEndX = imgRegion.GetIndex()[0] + imgSizeX;
itk::IndexValueType const imgEndY = imgRegion.GetIndex()[1] + imgSizeY;
auto const& size = outputRegionForThread.GetSize();
auto const& index = outputRegionForThread.GetIndex();
auto const sizeX = size[0];
auto const sizeY = size[1];
auto const startX = index[0];
auto const startY = index[1];
itk::IndexValueType const endX = startX + sizeX;
itk::IndexValueType const endY = startY + sizeY;
auto const thrX1 = std::min(endX, m_thresholdX);
auto const thrX2 = std::min(endX, imgEndX - m_thresholdX);
auto const thrY1 = std::min(endY, m_thresholdYtop);
auto const thrY2 = std::min(endY, imgEndY - m_thresholdYbot);
assert(thrX1 <= endX && "Iterations shall stay within requested region");
assert(thrX2 <= endX && "Iterations shall stay within requested region");
assert(thrY1 <= endY && "Iterations shall stay within requested region");
assert(thrY2 <= endY && "Iterations shall stay within requested region");
// using interval_t = boost::numeric::interval<long>;
using interval_t = Interval;
auto const region = interval_t{startX, endX};
auto const zero_left = intersect(interval_t{startX, thrX1}, region);
auto const copy_middle = intersect(interval_t{thrX1, thrX2}, region);
auto const zero_right = intersect(interval_t{thrX2, endX}, region);
otbMsgDevMacro("X in " << zero_left << " <<-- 0");
otbMsgDevMacro("X in " << copy_middle << " <<-- copy input");
otbMsgDevMacro("X in " << zero_right << " <<-- 0");
otbMsgDevMacro("Y in ["<<startY<<".."<<thrY1<<"[ <<--- 0");
auto const nb_z_l = zero_left.upper() - zero_left.lower();
auto const nb_c_m = copy_middle.upper() - copy_middle.lower();
auto const nb_z_r = zero_right.upper() - zero_right.lower();
assert(nb_z_l + nb_c_m + nb_z_r == sizeX);
itk::ProgressReporter progress( this, threadId, outputRegionForThread.GetNumberOfPixels() / sizeY );
outputIterator.GoToBegin();
// TODO: Can we consider that lines are contiguous in memory?
// If so, we could have only 2 fill_n and one copy_n!
auto y = startY;
for (
; y < thrY1
; ++y, outputIterator.NextLine())
{
// If there is any trimming of first lines, the inputIterator iterator will
// directly point to the right region. we shall not increment it!
// otbMsgDevMacro("o(" << y << ") <-- 0");
assert(! outputIterator.IsAtEnd());
outputIterator.GoToBeginOfLine();
std::fill_n(&outputIterator.Value(), sizeX, OutputPixelType{});
progress.CompletedPixel(); // Completed...Line()
}
assert(y == thrY1 || y == startY);
otbMsgDevMacro("Y in ["<<thrY1<<".."<<thrY2<<"[ <<--- Input");
inputIterator.GoToBegin();
for (
; y < thrY2
; ++y, inputIterator.NextLine(), outputIterator.NextLine())
{
// otbMsgDevMacro("o(" << y << ") <-- Input");
assert(! inputIterator.IsAtEnd());
assert(! outputIterator.IsAtEnd());
inputIterator.GoToBeginOfLine();
outputIterator.GoToBeginOfLine();
auto const t1 = std::fill_n(&outputIterator.Value(), nb_z_l, OutputPixelType{});
// If there is any trimming of first columns, the inputIterator iterator
// will directly point to the right region. we shall not apply an offset!
auto const t2 = std::copy_n(&inputIterator.Value(), nb_c_m, t1);
std::fill_n(t2, nb_z_r, OutputPixelType{});
progress.CompletedPixel(); // Completed...Line()
}
assert(y == thrY2 || y == startY);
otbMsgDevMacro("Y in ["<<thrY2<<".."<<endY<<"[ <<--- 0");
for (
; y < endY
; ++y, outputIterator.NextLine())
{
// If there is any trimming of last lines, the inputIterator iterator will
// directly point to the right region. we shall not increment it!
// otbMsgDevMacro("o(" << y << ") <-- 0");
assert(! outputIterator.IsAtEnd());
outputIterator.GoToBeginOfLine();
std::fill_n(&outputIterator.Value(), sizeX, OutputPixelType{});
progress.CompletedPixel(); // Completed...Line()
}
assert(y == endY);
otbMsgDevMacro("ThreadedGenerateData end");
}
template<typename TImage>
typename ClampROIFilter<TImage>::InputImageRegionType
ClampROIFilter<TImage>
::OutputRegionToInputRegion(OutputImageRegionType const& srcRegion)
{
auto const* output = this->GetOutput();
assert(output);
auto const& maxRegion = output->GetLargestPossibleRegion();
auto const& maxSize = maxRegion.GetSize();
auto const& maxStart = maxRegion.GetIndex();
auto const& reqRegion = srcRegion;
auto const& reqSize = reqRegion.GetSize();
auto const& reqStart = reqRegion.GetIndex();
// using interval_t = boost::numeric::interval<long>;
using interval_t = Interval;
auto const maxRegionX = interval_t{
maxStart[0]+m_thresholdX,
static_cast<itk::IndexValueType>(maxStart[0]+maxSize[0]-m_thresholdX)
};
auto const maxRegionY = interval_t{
maxStart[1]+m_thresholdYtop,
static_cast<itk::IndexValueType>(maxStart[1]+maxSize[1]-m_thresholdYbot)
};
auto const reqRegionX = interval_t::OfLength(reqStart[0], reqSize[0]);
auto const reqRegionY = interval_t::OfLength(reqStart[1], reqSize[1]);
#if 0
otbMsgDevMacro("OutputRegionToInputRegion: "
<< "out="<< NeatRegionLogger{reqRegion}
<< "; max: x="<<maxRegionX << " y="<<maxRegionY
<< "; req: x="<<reqRegionX << " y="<<reqRegionY
);
#endif