WIP/BUG: Trying to fix end-bound @ EOL

The dichotomic search doesn't assign a correct value to "last" end-bound when
it's close the end of image line

include/otbSARDEMPolygonsAnalysisImageFilter2.hxx

@@ -28,6 +28,7 @@
 #include "otbStringUtilities.h"
 #include "otbVectorImage.h"
 #include "otbPositionHelpers.h"
+#include "otbLogHelpers.h"
 
 #include "itkImageScanlineConstIterator.h"
 #include "itkImageScanlineIterator.h"
@@ -501,6 +502,7 @@ SARDEMPolygonsAnalysisImageFilter< TImageIn, TImageOut, TImageDEM, TImageSAR, TF
   inputRequestedRegionSize[1] = static_cast<ImageInIndexValueType>(lastL - firstL)+1;
 
   ImageInRegionType inputRequestedRegion(inputRequestedRegionIndex, inputRequestedRegionSize);
+  // assert(this->GetInput()->GetBufferedRegion().IsInside(inputRequestedRegion));
 
   return inputRequestedRegion;
 }
@@ -887,15 +889,17 @@ SARDEMPolygonsAnalysisImageFilter< TImageIn, TImageOut, TImageDEM, TImageSAR, TF
   auto point = [=](unsigned dx, unsigned y) {
     return VLVPointIterator<ImageInSubPixelType const>(in_data + nBands*dx + (y-in_startY)*in_line_off, nBands);
   };  // ------------------------------[ point
-  auto line = [](VLVPointIterator<ImageInSubPixelType const> it)
-  {
+  auto line = [](VLVPointIterator<ImageInSubPixelType const> it) {
     return (*it)[1];
   };  // ------------------------------[ line
 
+  // For 1 minute of execution, this function takes less than 0.5 seconds.
+  // => there is no need to optimize it more (parallelization, or skipping of
+  // lines already analysed)
   m_all_lines_are_sorted = true;
+  bool last_was_ascending, last_was_descending;
   for (unsigned y = in_startY ; y < in_endY ; ++y)
-  { // TODO: parallelize!
-    // TODO: don't analyse lines already analysed...
+  {
     int nb_up = 0;
     int nb_dw = 0;
     auto prev = line(point(0, y));
@@ -904,11 +908,25 @@ SARDEMPolygonsAnalysisImageFilter< TImageIn, TImageOut, TImageDEM, TImageSAR, TF
       auto const UL_line = line(point(dx, y));
       if (UL_line > prev) ++nb_up;
       else                ++nb_dw;
+      if (nb_up > 0 && nb_dw > 0)
+      { // no need to continue to scan the line!
+        break;
+      }
       prev = UL_line;
     }
-    if (nb_up == 0 || nb_dw == 0) m_is_lineX_sorted[y] = true;
-    else                          m_all_lines_are_sorted = false;
+    if (nb_up == 0 || nb_dw == 0) {
+      m_is_lineX_sorted[y] = true;
+      last_was_ascending = nb_up > 0 && nb_dw == 0;
+      last_was_descending = nb_up == 0 && nb_dw > 0;
+    } else {
+      m_all_lines_are_sorted = false;
+    }
   }
+  otbMsgDevMacro("BeforeThreadedGenerateData: in current block: all lines are sorted: " << m_all_lines_are_sorted);
+  otbMsgDevMacro("BeforeThreadedGenerateData: order is " << (last_was_ascending ? "ASC" : last_was_descending ? "DSC" : "???"));
+
+  otbMsgDevMacro("BeforeThreadedGenerateData: buffered region: " << NeatRegionLogger{input->GetBufferedRegion()});
+  otbMsgDevMacro("BeforeThreadedGenerateData: input    region: " << NeatRegionLogger{inputRegion});
 }
 
 /**
@@ -936,6 +954,8 @@ SARDEMPolygonsAnalysisImageFilter< TImageIn, TImageOut, TImageDEM, TImageSAR, TF
   ImageInSizeType const & inputRegionForThreadSize = inputRegionForThread.GetSize();
   unsigned int    const   nbDEMLines   =  static_cast<unsigned int>(inputRegionForThreadSize[1]);
   unsigned int    const   nbDEMColumns =  static_cast<unsigned int>(inputRegionForThreadSize[0]);
+  otbMsgDevMacro("ThreadedGenerateData: input region #"<<threadId<<" : " << NeatRegionLogger{inputRegionForThread});
+  assert(this->GetInput()->GetBufferedRegion().IsInside(inputRegionForThread));
 
 
   // Allocate outValueTab (size = nbCol)
@@ -1102,7 +1122,7 @@ SARDEMPolygonsAnalysisImageFilter< TImageIn, TImageOut, TImageDEM, TImageSAR, TF
         {
           auto find_lo_hi_window = [lo_margin, hi_margin, line, is_in_window](auto line_start, auto line_end)
           {
-            if (line(line_start) < line(line_end-1))
+            if (line(line_start) < line(line_end))
             {
               return equal_range_interval(
                   line_start, line_end,
@@ -1111,7 +1131,7 @@ SARDEMPolygonsAnalysisImageFilter< TImageIn, TImageOut, TImageDEM, TImageSAR, TF
             else
             {
               auto const lo_hi = equal_range_interval(
-                  std::make_reverse_iterator(line_end),
+                  std::make_reverse_iterator(line_end-1),
                   std::make_reverse_iterator(line_start),
                   lo_margin,  hi_margin, details::CmpVLVComponent(1));
               auto const lo = lo_hi.first;
@@ -1138,40 +1158,87 @@ SARDEMPolygonsAnalysisImageFilter< TImageIn, TImageOut, TImageDEM, TImageSAR, TF
             ++first;
           }
           //
-          while (first < last && !is_in_window(line(last-1)))
+          if (last == line_end - 2 && is_in_window(line(last)))
+          { // make "last" past-end
+            ++last;
+            assert(last-line_start == in_sizeX - 1);
+          }
+          else
           {
-            --last;
+            while (first < last && !is_in_window(line(last-1)))
+            {
+              --last;
+            }
           }
+          assert(is_in_window(line(point(in_sizeX-2,dy)))
+              // == ((in_sizeX-1) == (last-line_start)));
+              == (last == line_end - 1));
 
           // std::cout << "#y="<<(dy+in_startY) << "\t; start[0]="<<(*point(0, dy))[1]
           //   << "\t; end[" << (in_sizeX-1) << "]= " << (*point(in_sizeX-1, dy))[1]
           //   << "\t VS ["<<lo_margin<<", " << hi_margin << "]"
           //   << "\t -> " << (first - line_start) << " -- " << (last - line_start)
           //   << "\n";
-
           has_match = first.data() != last.data();
           if (has_match)
           { // a window has been found!
-            if (first.data() > line_start.data())
+#if 0
+            // Check bounds are consistent with generic code
+            auto first2 = line_start;
+            auto last2 = line_start;
+            bool was_in = false;
+            for (unsigned dx = 0 ; dx < in_sizeX -1 ; ++dx)
+            {
+              auto const pt_UL = *point(dx, dy);
+              auto const UL_line = pt_UL[1];
+
+              // Protect Memory access and input access
+              if (is_in_window(UL_line)) // true 2% of the time...
+              {
+                if (! was_in)
+                {
+                  first2 = point(dx, dy);
+                  was_in = true;
+                }
+              } else if (was_in) {
+                last2 = point(dx, dy);
+                was_in = false;
+              }
+            }
+            if (was_in)
+            {
+              last2 += in_sizeX - 1;
+            }
+            ASSERT(first == first2, "")(first.data())(first2.data())(first-first2)(line(first))(line(first2))(dy)(lo_margin)(hi_margin)(line(line_start))(line(line_end-1));
+            ASSERT(last == last2, "")(last.data())(last2.data())(last-last2)(line(last))(line(last2))(dy)(lo_margin)(hi_margin)(line(line_start))(line(line_end-1))(line_start.data());
+#endif
+
+            if (first > line_start)
             { // assert before is outside
               ASSERT(!is_in_window(line(first-1)), "")(in_data)(line(first-1))(dy)(lo_margin)(hi_margin)(line(line_start))(line(line_end-1));
             }
-            if (last.data() < line_end.data())
+            if (last < (line_end-1))
             { // assert after is outside
               ASSERT(!is_in_window(line(last)), "")(in_data)(line(last))(dy)(lo_margin)(hi_margin)(line(line_start))(line(line_end-1));
             }
             // after within is inside
             ASSERT(is_in_window(line(first)), "")(in_data)(line(first))(dy)(lo_margin)(hi_margin)(line(line_start))(line(line_end-1));
-            if (last.data() > line_start.data())
+            if (last > line_start)
               ASSERT(is_in_window(line(last-1)), "")(in_data)(line(last-1))(line(last))(dy)(lo_margin)(hi_margin)(line(line_start))(line(line_end-1));
 
+            assert((last-line_start) <= (in_sizeX - 1));
             for (auto dx = (first-line_start); dx < (last-line_start) ; ++dx)
             {
+              assert(dx < (in_sizeX - 1));
               auto const pt_UL = *point(dx, dy);
               auto const UL_line = pt_UL[1];
               assert (is_in_window(UL_line));
               register_pixel(dx, dy, pt_UL);
             }
+            // comment asserter que l'on va au moins jusqu'à sX-1???
+            assert(is_in_window(line(point(in_sizeX-2,dy)))
+              // == ((in_sizeX-1) == (last-line_start)));
+              == (last == line_end - 1));
           } // has_match
           // cpt_out_window += (hi - lo);
         }