From 339f79e96125111ade3a316ce2023d7a50908407 Mon Sep 17 00:00:00 2001
From: Julien Osman <julien.osman@csgroup.eu>
Date: Mon, 25 May 2020 15:06:27 +0200
Subject: [PATCH] BUG: #2046 Add an epsilon margin to compute the
 m_ReachableOutputRegion

The application Superimpose crashed with a segmentation fault during
resampling. Investigation revealed that the interpolated index was out
of image bounds. The application uses the PHR mode and resample the XS
image with a otb::GridResampleImageFilter. This filter doesn't check
explicitly that input indexes are in the input buffered region. It
uses a m_ReachableOutputRegion to crop the output region
processed. The problem appears because an output pixel is right on the
border of the input image extent, and when someone uses a nearest neighbor
interpolator, this one converts the continuous index to an
out-of-bound index.
---
 Data/Baseline/OTB/Images/apTvPrSuperimposePHR_nn.tif   |  3 +++
 .../Applications/AppProjection/app/otbSuperimpose.cxx  |  1 +
 Modules/Applications/AppProjection/test/CMakeLists.txt | 10 ++++++++++
 .../include/otbGridResampleImageFilter.h               |  5 +++++
 .../include/otbGridResampleImageFilter.hxx             | 10 ++++++----
 5 files changed, 25 insertions(+), 4 deletions(-)
 create mode 100644 Data/Baseline/OTB/Images/apTvPrSuperimposePHR_nn.tif

diff --git a/Data/Baseline/OTB/Images/apTvPrSuperimposePHR_nn.tif b/Data/Baseline/OTB/Images/apTvPrSuperimposePHR_nn.tif
new file mode 100644
index 0000000000..20bd7b1d3d
--- /dev/null
+++ b/Data/Baseline/OTB/Images/apTvPrSuperimposePHR_nn.tif
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:d6b8f0cc600d2b17c28cc7a3e7a2348bd986ed600f9ee86a3b96d7a83c75a88f
+size 8010952
diff --git a/Modules/Applications/AppProjection/app/otbSuperimpose.cxx b/Modules/Applications/AppProjection/app/otbSuperimpose.cxx
index 2edfa2aeeb..469fd1a1f6 100644
--- a/Modules/Applications/AppProjection/app/otbSuperimpose.cxx
+++ b/Modules/Applications/AppProjection/app/otbSuperimpose.cxx
@@ -192,6 +192,7 @@ private:
       NNInterpolatorType::Pointer interpolator = NNInterpolatorType::New();
       m_Resampler->SetInterpolator(interpolator);
       m_BasicResampler->SetInterpolator(interpolator);
+      m_BasicResampler->SetInterpolationMargin(1e-9);
     }
     break;
     case Interpolator_BCO:
diff --git a/Modules/Applications/AppProjection/test/CMakeLists.txt b/Modules/Applications/AppProjection/test/CMakeLists.txt
index 4637dacd5f..02265a85e5 100644
--- a/Modules/Applications/AppProjection/test/CMakeLists.txt
+++ b/Modules/Applications/AppProjection/test/CMakeLists.txt
@@ -367,3 +367,13 @@ otb_test_application(NAME apTvPrSuperimpose_phr
                      VALID  --compare-image ${EPSILON_7}
                         ${BASELINE}/apTvPrSuperimposePHR.tif
                         ${TEMP}/apTvPrSuperimposePHR.tif)
+
+otb_test_application(NAME apTvPrSuperimpose_phr_nn
+                     APP Superimpose
+                     OPTIONS -inr  ${INPUTDATA}/phr_pan.tif
+                             -inm ${INPUTDATA}/phr_xs.tif
+                             -interpolator nn
+                             -out ${TEMP}/apTvPrSuperimposePHR_nn.tif int16
+                     VALID  --compare-image ${EPSILON_7}
+                        ${BASELINE}/apTvPrSuperimposePHR_nn.tif
+                        ${TEMP}/apTvPrSuperimposePHR_nn.tif)
diff --git a/Modules/Filtering/ImageManipulation/include/otbGridResampleImageFilter.h b/Modules/Filtering/ImageManipulation/include/otbGridResampleImageFilter.h
index a45a37c15d..be41f30243 100644
--- a/Modules/Filtering/ImageManipulation/include/otbGridResampleImageFilter.h
+++ b/Modules/Filtering/ImageManipulation/include/otbGridResampleImageFilter.h
@@ -124,6 +124,9 @@ public:
   itkGetMacro(CheckOutputBounds, bool);
   itkBooleanMacro(CheckOutputBounds);
 
+  itkSetMacro(InterpolationMargin, double);
+  itkGetMacro(InterpolationMargin, double);
+
   itkSetObjectMacro(Interpolator, InterpolatorType);
   itkGetObjectMacro(Interpolator, InterpolatorType);
 
@@ -190,6 +193,8 @@ private:
 
   OutputPixelType m_EdgePaddingValue; // Default pixel value
 
+  double m_InterpolationMargin;
+
   bool m_CheckOutputBounds; // Shall we check
                             // output bounds when
                             // casting?
diff --git a/Modules/Filtering/ImageManipulation/include/otbGridResampleImageFilter.hxx b/Modules/Filtering/ImageManipulation/include/otbGridResampleImageFilter.hxx
index 34ecb1be90..7097361de8 100644
--- a/Modules/Filtering/ImageManipulation/include/otbGridResampleImageFilter.hxx
+++ b/Modules/Filtering/ImageManipulation/include/otbGridResampleImageFilter.hxx
@@ -42,6 +42,7 @@ GridResampleImageFilter<TInputImage, TOutputImage, TInterpolatorPrecision>::Grid
     m_OutputSpacing(),
     m_EdgePaddingValue(),
     m_CheckOutputBounds(true),
+    m_InterpolationMargin(0.0),
     m_Interpolator(),
     m_ReachableOutputRegion()
 {
@@ -227,7 +228,6 @@ void GridResampleImageFilter<TInputImage, TOutputImage, TInterpolatorPrecision>:
 
   // Compute the padding due to the interpolator
 
-
   IndexType inUL = this->GetInput()->GetBufferedRegion().GetIndex();
   IndexType inLR = this->GetInput()->GetBufferedRegion().GetIndex() + this->GetInput()->GetBufferedRegion().GetSize();
   inLR[0] -= 1;
@@ -246,9 +246,9 @@ void GridResampleImageFilter<TInputImage, TOutputImage, TInterpolatorPrecision>:
   this->GetInput()->TransformIndexToPhysicalPoint(inUL, inULp);
   this->GetInput()->TransformIndexToPhysicalPoint(inLR, inLRp);
 
-  inULp -= 0.5 * this->GetInput()->GetSignedSpacing();
-  inLRp += 0.5 * this->GetInput()->GetSignedSpacing();
-
+  inULp -= (0.5 - m_InterpolationMargin) * this->GetInput()->GetSignedSpacing();
+  inLRp += (0.5 - m_InterpolationMargin) * this->GetInput()->GetSignedSpacing();
+ 
   ContinuousInputIndexType outUL;
   ContinuousInputIndexType outLR;
   this->GetOutput()->TransformPhysicalPointToContinuousIndex(inULp, outUL);
@@ -265,6 +265,8 @@ void GridResampleImageFilter<TInputImage, TOutputImage, TInterpolatorPrecision>:
 
   m_ReachableOutputRegion.SetIndex(outputIndex);
   m_ReachableOutputRegion.SetSize(outputSize);
+
+  otbMsgDevMacro(<< "ReachableOutputRegion: " << m_ReachableOutputRegion);
 }
 
 template <typename TInputImage, typename TOutputImage, typename TInterpolatorPrecision>
-- 
GitLab