lsgrmController.txx 12.2 KB
Newer Older
remicres's avatar
remicres committed
1
#include "lsgrmController.h"
2 3 4 5

namespace lsgrm
{

remicres's avatar
remicres committed
6 7 8
template<class TSegmenter>
Controller<TSegmenter>::Controller()
{
9
  m_TilingMode = LSGRM_TILING_AUTO;
remicres's avatar
remicres committed
10 11 12 13 14 15 16
  m_Margin = 0;
  m_NumberOfIterations = 0;
  m_NumberOfFirstIterations = 0;
  m_TileHeight = 0;
  m_TileWidth = 0;
  m_NbTilesX = 0;
  m_NbTilesY = 0;
17
  m_CleanTemporaryFiles = true;
remicres's avatar
remicres committed
18 19
  m_Threshold = 75;
  m_Memory = 0;
remicres's avatar
remicres committed
20 21 22 23 24 25 26 27 28 29 30
}

template<class TSegmenter>
Controller<TSegmenter>::~Controller()
{
}

template<class TSegmenter>
void Controller<TSegmenter>::RunSegmentation()
{

31
  if (m_TilingMode == LSGRM_TILING_AUTO || m_TilingMode == LSGRM_TILING_USER)
remicres's avatar
remicres committed
32
    {
remicres's avatar
remicres committed
33 34 35

    const unsigned int numberOfIterationsForPartialSegmentations = 3; // TODO: find a smart value
    unsigned int numberOfIterationsRemaining = m_NumberOfIterations;
36 37 38 39 40 41 42 43 44 45 46

    if(m_TilingMode == LSGRM_TILING_AUTO)
      {
      this->GetAutomaticConfiguration();
      }
    else if (m_TilingMode == LSGRM_TILING_USER)
      {
      m_Margin = static_cast<unsigned int>(pow(2, m_NumberOfFirstIterations + 1) - 2);
      }

    std::cout <<
remicres's avatar
remicres committed
47 48 49 50 51 52 53
        "--- Configuration: " <<
        "\n\tAvailable RAM: " << m_Memory <<
        "\n\tInput image dimensions: " << m_InputImage->GetLargestPossibleRegion().GetSize() <<
        "\n\tNumber of first iterations: " << m_NumberOfFirstIterations <<
        "\n\tStability margin: " << m_Margin <<
        "\n\tRegular tile size: " << m_TileWidth << " x " << m_TileHeight <<
        "\n\tTiling layout: " << m_NbTilesX << " x " << m_NbTilesY << std::endl;
remicres's avatar
remicres committed
54 55

    // Compute the splitting scheme
56 57
    m_Tiles = SplitOTBImage<ImageType>(m_InputImage, m_TileWidth, m_TileHeight, m_Margin,
        m_NbTilesX, m_NbTilesY, m_TemporaryFilesPrefix);
remicres's avatar
remicres committed
58

remicres's avatar
remicres committed
59 60 61 62 63 64 65 66 67 68
    // Boolean indicating if there are remaining fusions
    bool isFusion = false;

    // Run first partial segmentation
    boost::timer t; t.restart();
    auto accumulatedMemory = RunFirstPartialSegmentation<TSegmenter>(
        m_InputImage,
        m_SpecificParameters,
        m_Threshold,
        m_NumberOfFirstIterations,
69
        numberOfIterationsForPartialSegmentations,
remicres's avatar
remicres committed
70 71 72 73 74 75 76
        m_Tiles,
        m_NbTilesX,
        m_NbTilesY,
        m_TileWidth,
        m_TileHeight,
        isFusion);

77
    // Update the given number of iterations
remicres's avatar
remicres committed
78
    numberOfIterationsRemaining -= m_NumberOfFirstIterations;
79

remicres's avatar
remicres committed
80 81 82 83 84 85
    // Gathering useful variables
    GatherUsefulVariables(accumulatedMemory, isFusion);

    // Time monitoring
    ShowTime(t);

remicres's avatar
remicres committed
86 87
    std::cout << "accumulatedMemory=" << accumulatedMemory << std::endl;

remicres's avatar
remicres committed
88
    while(accumulatedMemory > m_Memory && isFusion)
remicres's avatar
remicres committed
89
      {
remicres's avatar
remicres committed
90 91 92

      isFusion = false;
      accumulatedMemory = RunPartialSegmentation<TSegmenter>(
remicres's avatar
remicres committed
93 94
          m_SpecificParameters,
          m_Threshold,
95
          numberOfIterationsForPartialSegmentations,
remicres's avatar
remicres committed
96 97 98 99 100 101
          m_Tiles,
          m_NbTilesX,
          m_NbTilesY,
          m_InputImage->GetLargestPossibleRegion().GetSize()[0],
          m_InputImage->GetLargestPossibleRegion().GetSize()[1],
          m_InputImage->GetNumberOfComponentsPerPixel(),
remicres's avatar
remicres committed
102 103
          isFusion);

104
      // Update the given number of iterations
remicres's avatar
remicres committed
105
      numberOfIterationsRemaining -= numberOfIterationsForPartialSegmentations;
106

remicres's avatar
remicres committed
107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130
      // Gathering useful variables
      GatherUsefulVariables(accumulatedMemory, isFusion);

      // Time monitoring
      ShowTime(t);
      }


#ifdef OTB_USE_MPI
    if (otb::MPIConfig::Instance()->GetMyRank() != 0)
      return;
#endif

    if(accumulatedMemory <= m_Memory)
      {
      // Merge all the graphs
      m_LabelImage = MergeAllGraphsAndAchieveSegmentation<TSegmenter>(
          m_SpecificParameters,
          m_Threshold,
          m_Tiles,
          m_NbTilesX,
          m_NbTilesY,
          m_InputImage->GetLargestPossibleRegion().GetSize()[0],
          m_InputImage->GetLargestPossibleRegion().GetSize()[1],
131
          m_InputImage->GetNumberOfComponentsPerPixel(),
remicres's avatar
remicres committed
132
          numberOfIterationsRemaining);
remicres's avatar
remicres committed
133

134
      ShowTime(t);
remicres's avatar
remicres committed
135 136

      }
137
    else // accumulatedMemory > m_Memory
remicres's avatar
remicres committed
138
      {
remicres's avatar
remicres committed
139
      // That means there are no more possible fusions but we can not store the output graph
remicres's avatar
remicres committed
140 141
      // Todo do not clean up temporary directory before copying resulting graph to the output directory
      // In the output directory add an info file to give the number of tiles.
142
      itkExceptionMacro(<< "No more possible fusions, but can not store the output graph");
remicres's avatar
remicres committed
143 144
      }
    }
145
  else // tiling_mode is none
remicres's avatar
remicres committed
146
    {
147
    // todo use classic grm
remicres's avatar
remicres committed
148 149 150 151 152
    }


}

153 154 155 156 157 158 159
/*
 * Compute the memory occupied by one node
 */
template<class TSegmenter>
unsigned int Controller<TSegmenter>::GetNodeMemory()
{
  // Create a unique node
160
  const unsigned int n = 100;
161 162 163 164
  typename ImageType::Pointer onePixelImage = ImageType::New();
  typename ImageType::IndexType start;
  start.Fill(0);
  typename ImageType::SizeType size;
165
  size.Fill(n);
166 167 168 169 170 171 172 173
  typename ImageType::RegionType region(start, size);
  onePixelImage->SetRegions(region);
  onePixelImage->SetNumberOfComponentsPerPixel(m_InputImage->GetNumberOfComponentsPerPixel());
  onePixelImage->Allocate();
  TSegmenter segmenter;
  segmenter.SetInput(onePixelImage);
  lsrm::GraphOperations<TSegmenter>::InitNodes(onePixelImage,segmenter,FOUR);

174
  unsigned int memory = segmenter.GetGraphMemory() / (n*n);
175
  itkDebugMacro(<<"Size of a node is " << memory);
176
  return memory;
177
}
remicres's avatar
remicres committed
178

179 180 181
/*
 * Compute the maximum number of nodes which can fit in the system memory
 */
remicres's avatar
remicres committed
182
template<class TSegmenter>
remicres's avatar
remicres committed
183
long unsigned int Controller<TSegmenter>::GetMaximumNumberOfNodesInMemory()
remicres's avatar
remicres committed
184
{
remicres's avatar
remicres committed
185 186
  itkDebugMacro(<< "Computing maximum number of nodes in memory");

remicres's avatar
remicres committed
187 188 189 190 191
  m_Memory = getMemorySize();
  assert(m_Memory > 0);

  m_Memory /= 2; // For safety and can prevent out of memory troubles

192
  return std::ceil(((float) m_Memory) / ((float) GetNodeMemory()));
remicres's avatar
remicres committed
193 194 195 196

}

template<class TSegmenter>
197 198
void
Controller<TSegmenter>::ComputeMaximumStabilityMargin(unsigned int width,
remicres's avatar
remicres committed
199
    unsigned int height, unsigned int &niter, unsigned int &margin)
200
    {
remicres's avatar
remicres committed
201 202 203 204
  itkDebugMacro(<< "Computing maximum stability margin");

  // Compute the stability margin. The naive strategy consider a margin value and a stable size equal.
  niter = 1;
205 206 207
  unsigned int maxMargin = std::min(width, height)/2;
  unsigned int currMargin = static_cast<unsigned int>(pow(2, niter + 1) - 2);
  margin = currMargin;
remicres's avatar
remicres committed
208

209
  while(currMargin < maxMargin)
remicres's avatar
remicres committed
210
    {
211
    margin = currMargin;
remicres's avatar
remicres committed
212
    niter++;
213
    currMargin = static_cast<unsigned int>(pow(2, niter + 1) - 2);
remicres's avatar
remicres committed
214 215 216 217
    }
  niter--;

  itkDebugMacro(<< "Number of iterations=" << niter << " margin=" << margin);
remicres's avatar
remicres committed
218

219
    }
remicres's avatar
remicres committed
220

221
/*
222
 * Compute a tiling layout which minimizes a criterion based on tile compactness
remicres's avatar
remicres committed
223
 * and memory usage
224 225
 *
 * TODO: use the lsgrmSplitter to truly compute the largest tile of a given layout
226
 */
remicres's avatar
remicres committed
227 228 229 230
template<class TSegmenter>
void Controller<TSegmenter>::GetAutomaticConfiguration()
{

remicres's avatar
remicres committed
231 232
  itkDebugMacro(<<"Get automatic configuration");

remicres's avatar
remicres committed
233 234
  // Compute the maximum number of nodes that can fit the memory
  unsigned long int maximumNumberOfNodesInMemory = GetMaximumNumberOfNodesInMemory();
remicres's avatar
remicres committed
235
  itkDebugMacro(<<"Maximum number of nodes in memory is " << maximumNumberOfNodesInMemory);
remicres's avatar
remicres committed
236

remicres's avatar
remicres committed
237
  // Number of nodes in the entire image
remicres's avatar
remicres committed
238 239
  const unsigned int imageWidth = m_InputImage->GetLargestPossibleRegion().GetSize()[0];
  const unsigned int imageHeight = m_InputImage->GetLargestPossibleRegion().GetSize()[1];
remicres's avatar
remicres committed
240 241
  const unsigned long int nbOfNodesInImage = imageWidth*imageHeight;

242
  // Default layout: 1x1
remicres's avatar
remicres committed
243 244
  m_NbTilesX = 1;
  m_NbTilesY = 1;
245

remicres's avatar
remicres committed
246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264
  // Without margins, the number of tiles maximizing memory use
  // is equal to: nbOfNodesInImage / maximumNumberOfNodesInMemory.
  // Actually, there is tile margins. And the best scenario is to have
  // square tiles with margin = width/2, that is tiles 4x larger.
  // Hence the number of tiles maximizing memory use is 4x larger.
  unsigned int minimumNumberOfTiles = std::ceil(4 * nbOfNodesInImage / ((float) maximumNumberOfNodesInMemory));
  itkDebugMacro(<<"Minimum number of tiles is " << minimumNumberOfTiles);

  // In the following steps, we will optimize tiling layout, starting from a number
  // of tiles equal to "minimumNumberOfTiles", up to a number of tiles equal to
  // twice the number of tiles (that is memory usage about 50%)
  unsigned int maximumNumberOfTiles = minimumNumberOfTiles * 4;

  // Search for layout which minimizes the criterion
  // The criterion is the ratio between compactness and memory usage
  // (i.e. tileWidth * tileHeight / maximumNumberOfNodesInMemory)
  itkDebugMacro(<<"Computing layouts properties:");
  float lowestCriterionValue = itk::NumericTraits<float>::max();
  for (unsigned int nbOfTiles = minimumNumberOfTiles ; nbOfTiles <= maximumNumberOfTiles ; nbOfTiles++)
265
    {
remicres's avatar
remicres committed
266 267
    // Get the multiples of k. For each one, compute the criterion of the tiling
    for (unsigned int layoutNCol = 1; layoutNCol<=nbOfTiles; layoutNCol++)
remicres's avatar
remicres committed
268
      {
remicres's avatar
remicres committed
269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317
#ifdef OTB_USE_MPI
      // We want number of tiles which is a multiple of the number of MPI processes
      if (nbOfTiles % layoutNCol == 0 && // Is it a multiple of the nb of ...Tiles?
          nbOfTiles % otb::MPIConfig::Instance()->GetNbProcs() == 0) // ...nProcs?
#else
        if (nbOfTiles % layoutNCol == 0) // Is it a multiple of the nb of Tiles?
#endif
          {
          // Tiling layout
          unsigned int layoutNRow = nbOfTiles / layoutNCol;
          unsigned int tileWidth = imageWidth / layoutNCol;
          unsigned int tileHeight = imageHeight / layoutNRow;

          // Compute margin for regular tiles of this layout
          unsigned int maxMargin, maxIter;
          ComputeMaximumStabilityMargin(tileWidth, tileHeight, maxIter, maxMargin);
          tileWidth += 2*maxMargin;
          tileHeight += 2*maxMargin;

          // Memory use efficiency
          float percentMemory = tileWidth * tileHeight / (float) maximumNumberOfNodesInMemory; // ]0, 1]

          // Compactness
          float perimeter = tileWidth + tileHeight;
          float surface = tileWidth * tileHeight;
          float compactness = perimeter / surface * (float) vcl_max(tileWidth,tileHeight); // [1,+inf]

          // Update minimum criterion
          float criterion = compactness / percentMemory; // ]0, +inf]

          itkDebugMacro(<< std::setprecision (2) << std::fixed
              << "Nb. tiles=" << nbOfTiles
              << " Layout: " << layoutNRow << "x" << layoutNCol
              << " Mem. use=" << percentMemory
              << " Compactness=" << compactness
              << " Criterion=" << criterion
              << " Size (no margin): " << (tileWidth-2*maxMargin)<< "x"<< (tileHeight-2*maxMargin)
              << " Size (with margin): " << tileWidth << "x" << tileHeight
              << " (margin=" << maxMargin << "/nb. iter=" << maxIter << ")" );

          if (criterion < lowestCriterionValue)
            {
            lowestCriterionValue = criterion;
            m_NbTilesX = layoutNCol;
            m_NbTilesY = layoutNRow;
            }
          }
      } // for each multiple of k
    }
remicres's avatar
remicres committed
318 319 320 321

  // Compute the tile size
  m_TileWidth = static_cast<unsigned int>(imageWidth/m_NbTilesX);
  m_TileHeight = static_cast<unsigned int>(imageHeight/m_NbTilesY);
remicres's avatar
remicres committed
322 323 324 325 326 327
  itkDebugMacro(<<"Selected layout: " << m_NbTilesX << "x" << m_NbTilesY
      << " (criterion=" << lowestCriterionValue << ")");

  // Compute the stability margin
  ComputeMaximumStabilityMargin(m_TileWidth, m_TileHeight,m_NumberOfFirstIterations, m_Margin);

328 329 330
  long long unsigned int memoryUsed = GetNodeMemory();
  memoryUsed *= static_cast<long long unsigned int>(m_TileHeight + 2*m_Margin);
  memoryUsed *= static_cast<long long unsigned int>(m_TileWidth + 2*m_Margin);
331
  itkDebugMacro(<< "An amount of " << memoryUsed/(1024.0*1024.0) << " Mbytes of RAM will be used for regular tiles of size "
remicres's avatar
remicres committed
332
      << (m_TileWidth + 2*m_Margin) << "x" << (m_TileHeight + 2*m_Margin) );
remicres's avatar
remicres committed
333

remicres's avatar
remicres committed
334 335
}

remicres's avatar
remicres committed
336 337 338 339 340 341 342 343
template <class TSegmenter>
void Controller<TSegmenter>::SetInternalMemoryAvailable(long long unsigned int v) // expecting a value in Mbytes.
{
  assert(v > 0);
  m_Memory = v * 1024ul * 1024ul;
}

template<class TSegmenter>
remicres's avatar
remicres committed
344
void Controller<TSegmenter>::SetInputImage(ImageType * inputImage)
remicres's avatar
remicres committed
345
{
remicres's avatar
remicres committed
346
  m_InputImage = inputImage;
remicres's avatar
remicres committed
347 348 349 350 351 352 353 354
}

template<class TSegmenter>
void Controller<TSegmenter>::SetSpecificParameters(const SegmentationParameterType& params)
{
  m_SpecificParameters = params;
}

remicres's avatar
remicres committed
355 356 357 358 359 360
template<class TSegmenter>
typename Controller<TSegmenter>::LabelImageType::Pointer
Controller<TSegmenter>::GetLabeledClusteredOutput()
{
  return m_LabelImage;
}
remicres's avatar
remicres committed
361
} // end of namespace lsgrm