ENH: Full MPI implementation.

app/otbLSGRM.cxx

@@ -132,6 +132,14 @@ private:
         tmpdir.append("/");
         }
       }
+
+    // Check that the directory exists
+    if (!itksys::SystemTools::FileExists(tmpdir.c_str(),false))
+      {
+      otbAppLogFATAL("The directory " << tmpdir << " does not exist.");
+      }
+
+    // Return the prefix
     std::string prefix = tmpdir + outbfname;
     return prefix;
   }
@@ -203,6 +211,7 @@ private:
     /*
         To add:
         the output directory in case the global graph cannot fit in memory
+        TODO: It seems to be a cast-related bug when changing the output pixel type (e.g. uint16) !
      */
 
     // Input image
@@ -245,6 +254,13 @@ private:
 
   void AfterExecuteAndWriteOutputs()
   {
+
+#ifdef OTB_USE_MPI
+    // When MPI is activated, only the master thread proceed
+    if (otb::MPIConfig::Instance()->GetMyRank() != 0)
+      return;
+#endif
+
     // Delete temporary files
     for (unsigned int i = 0 ; i < m_TemporaryFilesList.size() ; i++)
       {

include/lsgrmController.txx

@@ -37,13 +37,14 @@ void Controller<TSegmenter>::RunSegmentation()
 
   if (m_TilingMode == LSGRM_TILING_AUTO || m_TilingMode == LSGRM_TILING_USER)
     {
-
     if(m_TilingMode == LSGRM_TILING_AUTO)
       {
       this->GetAutomaticConfiguration();
       }
     else if (m_TilingMode == LSGRM_TILING_USER)
       {
+      m_NbTilesX = std::floor(m_InputImage->GetLargestPossibleRegion().GetSize()[0] / m_TileWidth);
+      m_NbTilesY = std::floor(m_InputImage->GetLargestPossibleRegion().GetSize()[1] / m_TileHeight);
       m_Margin = static_cast<unsigned int>(pow(2, m_NumberOfFirstIterations + 1) - 2);
       }
 
@@ -66,12 +67,10 @@ void Controller<TSegmenter>::RunSegmentation()
       std::cout << "Only one tile is needed. Fallback to tiling=none." << std::endl;
       SetTilingModeNone();
       }
-
     }
 
   if (m_TilingMode == LSGRM_TILING_AUTO || m_TilingMode == LSGRM_TILING_USER)
     {
-
     const unsigned int numberOfIterationsForPartialSegmentations = 3; // TODO: find a smart value
     unsigned int numberOfIterationsRemaining = m_NumberOfIterations;
 
@@ -108,7 +107,6 @@ void Controller<TSegmenter>::RunSegmentation()
 
     while(accumulatedMemory > m_Memory && isFusion)
       {
-
       isFusion = false;
       accumulatedMemory = RunPartialSegmentation<TSegmenter>(
           m_SpecificParameters,
@@ -135,13 +133,8 @@ void Controller<TSegmenter>::RunSegmentation()
         {
         numberOfIterationsRemaining -= numberOfIterationsForPartialSegmentations;
         }
-      else
-        {
-        break;
-        }
       }
 
-
 #ifdef OTB_USE_MPI
     // Only the master process is doing the next part
     // TODO: Use the MPI process wich has the largest amount of memory
@@ -176,6 +169,15 @@ void Controller<TSegmenter>::RunSegmentation()
     }
   else if (m_TilingMode == LSGRM_TILING_NONE)// tiling_mode is none
     {
+#ifdef OTB_USE_MPI
+    // Only the master process is doing the next part
+    if (otb::MPIConfig::Instance()->GetMyRank() > 0)
+      return;
+    else
+      // Warn that there is some unused MPI processes
+      if (otb::MPIConfig::Instance()->GetNbProcs() > 1)
+        itkWarningMacro(<< "Only 1 MPI process will be used");
+#endif
     // Update input image
     m_InputImage->Update();
 
@@ -410,6 +412,10 @@ template<class TSegmenter>
 typename Controller<TSegmenter>::LabelImageType::Pointer
 Controller<TSegmenter>::GetLabeledClusteredOutput()
 {
+#ifdef OTB_USE_MPI
+  // Get the label image from the master process (the one which achieves segmentation)
+  BroadcastImage<typename TSegmenter::LabelImageType>(m_LabelImage);
+#endif
   return m_LabelImage;
 }
 

include/lsgrmHeader.h

@@ -16,6 +16,9 @@
 #ifdef OTB_USE_MPI
 #include "otbMPIConfig.h"
 #include "mpi.h" // TODO: implement needed methods inside otbMPIConfig.h
+#include "otbExtractROI.h"
+#include "itkImageRegionIterator.h"
+#include "otbImageFileWriter.h"
 #endif
 
 /*
@@ -64,8 +67,112 @@ void GatherMe(T& x, MPI_Datatype dataType)
     MPI_Recv(&x, 1, dataType, 0, MPI_ANY_TAG, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
     }
 }
-#endif
 
+/*
+ * Function used to broadcast the label image to every MPI process
+ */
+template<class TImageType>
+void BroadcastImage(typename TImageType::Pointer & inPtr)
+{
+  otb::MPIConfig::Instance()->barrier();
+
+  unsigned int width;
+  unsigned int height;
+  unsigned int block_height;
+  unsigned int current_start_y;
+  if (otb::MPIConfig::Instance()->GetMyRank() == 0)
+    {
+    // Master process
+    width = inPtr->GetLargestPossibleRegion().GetSize()[0];
+    height = inPtr->GetLargestPossibleRegion().GetSize()[1];
+    }
+
+  // Broadcast width and height
+  MPI_Bcast(&width, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD);
+  MPI_Bcast(&height, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD);
+
+  // Slave processes do allocate image
+  typename TImageType::IndexType index;
+  index.Fill(0);
+  typename TImageType::SizeType size;
+  size[0] = width;
+  size[1] = height;
+  typename TImageType::RegionType region(index,size);
+  if (otb::MPIConfig::Instance()->GetMyRank() > 0)
+    {
+    inPtr = TImageType::New();
+    inPtr->SetRegions(region);
+    inPtr->SetNumberOfComponentsPerPixel(1);
+    inPtr->Allocate();
+    }
+
+  // Maximum data count that mpi can handle
+  unsigned int maximum_count = std::numeric_limits<int>::max();
+  block_height = std::floor((float) maximum_count / width);
+
+  // Broadcast array block by block (lines)
+  current_start_y = 0;
+  while (current_start_y < height)
+    {
+    if ( current_start_y + block_height > height )
+      block_height = height - current_start_y;
+
+    // Subregion of image
+    typename TImageType::Pointer tmpPtr = TImageType::New();
+
+    typename TImageType::IndexType subregion_index;
+    subregion_index[0] = 0;
+    subregion_index[1] = current_start_y;
+    typename TImageType::SizeType subregion_size;
+    subregion_size[0] = width;
+    subregion_size[1] = block_height;
+    typename TImageType::RegionType subregion(subregion_index, subregion_size);
+
+    // Slave processes do allocate subregion image
+    if (otb::MPIConfig::Instance()->GetMyRank() > 0)
+      {
+      tmpPtr->SetRegions(subregion);
+      tmpPtr->Allocate();
+      }
+    else
+      {
+      typedef typename otb::ExtractROI<typename TImageType::InternalPixelType,
+          typename TImageType::InternalPixelType> ExtractROIFilterType;
+      typename ExtractROIFilterType::Pointer filter = ExtractROIFilterType::New();
+      filter->SetInput(inPtr);
+      filter->SetStartX(0);
+      filter->SetStartY(current_start_y);
+      filter->SetSizeX(width);
+      filter->SetSizeY(block_height);
+      filter->SetReleaseDataFlag(false);
+      filter->Update();
+      tmpPtr = filter->GetOutput();
+      }
+
+    current_start_y += block_height;
+
+    // Broadcast buffer
+    MPI_Bcast(tmpPtr->GetBufferPointer(), width*block_height, MPI_UNSIGNED, 0, MPI_COMM_WORLD);
+
+    // Slave process must recopy the image
+    if (otb::MPIConfig::Instance()->GetMyRank() > 0)
+      {
+      typedef itk::ImageRegionIterator<TImageType> IteratorType;
+      IteratorType it1(inPtr, subregion);
+      IteratorType it2(tmpPtr, subregion);
+
+      for (it1.GoToBegin(), it2.GoToBegin(); !it1.IsAtEnd(); ++it1, ++it2)
+        {
+        it1.Set(it2.Get());
+        }
+      } // recopy image
+
+    } // while data to transmit
+
+}
+
+
+#endif
 
 /*
  * Gather accumulatedMemory and isFusion variables