lsgrmController.txx 14 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
17
18
  m_Margin = 0;
  m_NumberOfIterations = 0;
  m_NumberOfFirstIterations = 0;
  m_TileHeight = 0;
  m_TileWidth = 0;
  m_NbTilesX = 0;
  m_NbTilesY = 0;
  m_Threshold = 75;
  m_Memory = 0;
remicres's avatar
remicres committed
19
20
21
22
23
24
25
}

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

26
27
28
29
30
31
32
template<class TSegmenter>
void Controller<TSegmenter>::Modified()
{
  Superclass::Modified();
  m_Tiles.clear();
}

remicres's avatar
remicres committed
33
34
35
template<class TSegmenter>
void Controller<TSegmenter>::RunSegmentation()
{
remicres's avatar
remicres committed
36
  itkDebugMacro(<< "Entering RunSegmentation()");
remicres's avatar
remicres committed
37

38
  if (m_TilingMode == LSGRM_TILING_AUTO || m_TilingMode == LSGRM_TILING_USER)
remicres's avatar
remicres committed
39
    {
remicres's avatar
remicres committed
40

41
42
43
44
45
46
47
48
49
50
    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
51
52
53
54
55
56
57
        "--- 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
58
59

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

63
    // If there is only one tile, then fallback to LSGRM_TILING_NONE case
64
65
66
67
68
69
    if (m_Tiles.size() == 1)
      {
      std::cout << "Only one tile is needed. Fallback to tiling=none." << std::endl;
      SetTilingModeNone();
      }

70
71
72
73
74
75
76
77
    }

  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;

remicres's avatar
remicres committed
78
79
80
81
82
83
84
85
86
87
    // 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,
88
        numberOfIterationsForPartialSegmentations,
remicres's avatar
remicres committed
89
90
91
92
93
94
95
        m_Tiles,
        m_NbTilesX,
        m_NbTilesY,
        m_TileWidth,
        m_TileHeight,
        isFusion);

96
    // Update the given number of iterations
remicres's avatar
remicres committed
97
    numberOfIterationsRemaining -= m_NumberOfFirstIterations;
98

99
#ifdef OTB_USE_MPI
remicres's avatar
remicres committed
100
101
    // Gathering useful variables
    GatherUsefulVariables(accumulatedMemory, isFusion);
102
#endif
remicres's avatar
remicres committed
103
104
105
106

    // Time monitoring
    ShowTime(t);

remicres's avatar
remicres committed
107
108
    std::cout << "accumulatedMemory=" << accumulatedMemory << std::endl;

remicres's avatar
remicres committed
109
    while(accumulatedMemory > m_Memory && isFusion)
remicres's avatar
remicres committed
110
      {
remicres's avatar
remicres committed
111
112
113

      isFusion = false;
      accumulatedMemory = RunPartialSegmentation<TSegmenter>(
remicres's avatar
remicres committed
114
115
          m_SpecificParameters,
          m_Threshold,
116
          numberOfIterationsForPartialSegmentations,
remicres's avatar
remicres committed
117
118
119
120
121
122
          m_Tiles,
          m_NbTilesX,
          m_NbTilesY,
          m_InputImage->GetLargestPossibleRegion().GetSize()[0],
          m_InputImage->GetLargestPossibleRegion().GetSize()[1],
          m_InputImage->GetNumberOfComponentsPerPixel(),
remicres's avatar
remicres committed
123
124
          isFusion);

125
#ifdef OTB_USE_MPI
remicres's avatar
remicres committed
126
127
      // Gathering useful variables
      GatherUsefulVariables(accumulatedMemory, isFusion);
128
#endif
129

remicres's avatar
remicres committed
130
131
      // Time monitoring
      ShowTime(t);
132
133
134
135
136
137
138
139
140
141

      // Update the number of iterations remaining
      if (numberOfIterationsRemaining > numberOfIterationsForPartialSegmentations)
        {
        numberOfIterationsRemaining -= numberOfIterationsForPartialSegmentations;
        }
      else
        {
        break;
        }
remicres's avatar
remicres committed
142
143
144
145
      }


#ifdef OTB_USE_MPI
146
    // Only the master process is doing the next part
remicres's avatar
remicres committed
147
    // TODO: Use the MPI process wich has the largest amount of memory
remicres's avatar
remicres committed
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
    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],
163
          m_InputImage->GetNumberOfComponentsPerPixel(),
remicres's avatar
remicres committed
164
          numberOfIterationsRemaining);
remicres's avatar
remicres committed
165

166
      ShowTime(t);
remicres's avatar
remicres committed
167
168

      }
169
    else // accumulatedMemory > m_Memory
remicres's avatar
remicres committed
170
      {
remicres's avatar
remicres committed
171
      // That means there are no more possible fusions but we can not store the output graph
remicres's avatar
remicres committed
172
173
      // 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.
174
      itkExceptionMacro(<< "No more possible fusions, but can not store the output graph");
remicres's avatar
remicres committed
175
176
      }
    }
remicres's avatar
remicres committed
177
  else if (m_TilingMode == LSGRM_TILING_NONE)// tiling_mode is none
remicres's avatar
remicres committed
178
    {
179
180
181
    // Update input image
    m_InputImage->Update();

remicres's avatar
remicres committed
182
183
184
185
186
187
188
189
    // Use classic grm
    TSegmenter segmenter;
    segmenter.SetParam(m_SpecificParameters);
    segmenter.SetThreshold(m_Threshold);
    segmenter.SetDoFastSegmentation(true);
    segmenter.SetNumberOfIterations(m_NumberOfIterations);
    segmenter.SetInput(m_InputImage);
    segmenter.Update();
190
191
192

    // Get label image
    m_LabelImage = segmenter.GetLabeledClusteredOutput();
remicres's avatar
remicres committed
193
194
195
196
    }
  else
    {
    itkExceptionMacro(<<"Unknow tiling mode!");
remicres's avatar
remicres committed
197
198
199
200
    }

}

201
202
/*
 * Compute the memory occupied by one node
203
 * TODO: compute the exact value
204
205
206
207
 */
template<class TSegmenter>
unsigned int Controller<TSegmenter>::GetNodeMemory()
{
208
  // Create a n*n image
209
  const unsigned int n = 100;
210
211
212
213
  typename ImageType::Pointer onePixelImage = ImageType::New();
  typename ImageType::IndexType start;
  start.Fill(0);
  typename ImageType::SizeType size;
214
  size.Fill(n);
215
216
217
218
  typename ImageType::RegionType region(start, size);
  onePixelImage->SetRegions(region);
  onePixelImage->SetNumberOfComponentsPerPixel(m_InputImage->GetNumberOfComponentsPerPixel());
  onePixelImage->Allocate();
219
220

  // Instanciate and initialize a segmenter
221
222
  TSegmenter segmenter;
  segmenter.SetInput(onePixelImage);
223
  grm::GraphOperations<TSegmenter>::InitNodes(onePixelImage,segmenter,FOUR);
224

225
  // Get the memory occupied by the graph, normalize it by n*n
226
  unsigned int memory = segmenter.GetGraphMemory() / (n*n);
227

228
  itkDebugMacro(<<"Size of a node is " << memory);
229

230
  return memory;
231
}
remicres's avatar
remicres committed
232

233
234
235
/*
 * Compute the maximum number of nodes which can fit in the system memory
 */
remicres's avatar
remicres committed
236
template<class TSegmenter>
remicres's avatar
remicres committed
237
long unsigned int Controller<TSegmenter>::GetMaximumNumberOfNodesInMemory()
remicres's avatar
remicres committed
238
{
remicres's avatar
remicres committed
239
240
  itkDebugMacro(<< "Computing maximum number of nodes in memory");

remicres's avatar
remicres committed
241
242
243
244
245
  m_Memory = getMemorySize();
  assert(m_Memory > 0);

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

246
  return std::ceil(((float) m_Memory) / ((float) GetNodeMemory()));
remicres's avatar
remicres committed
247
248
249
250

}

template<class TSegmenter>
251
void Controller<TSegmenter>::ComputeMaximumStabilityMargin(unsigned int width,
remicres's avatar
remicres committed
252
    unsigned int height, unsigned int &niter, unsigned int &margin)
253
    {
remicres's avatar
remicres committed
254
255
256
257
  itkDebugMacro(<< "Computing maximum stability margin");

  // Compute the stability margin. The naive strategy consider a margin value and a stable size equal.
  niter = 1;
258
259
260
  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
261

262
  while(currMargin < maxMargin)
remicres's avatar
remicres committed
263
    {
264
    margin = currMargin;
remicres's avatar
remicres committed
265
    niter++;
266
    currMargin = static_cast<unsigned int>(pow(2, niter + 1) - 2);
remicres's avatar
remicres committed
267
268
269
270
    }
  niter--;

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

272
    }
remicres's avatar
remicres committed
273

274
/*
275
 * Compute a tiling layout which minimizes a criterion based on tile compactness
remicres's avatar
remicres committed
276
 * and memory usage
277
278
 *
 * TODO: use the lsgrmSplitter to truly compute the largest tile of a given layout
279
 */
remicres's avatar
remicres committed
280
281
282
283
template<class TSegmenter>
void Controller<TSegmenter>::GetAutomaticConfiguration()
{

remicres's avatar
remicres committed
284
285
  itkDebugMacro(<<"Get automatic configuration");

remicres's avatar
remicres committed
286
  // Compute the maximum number of nodes that can fit the memory
remicres's avatar
remicres committed
287
  // TODO: Use the smallest number amongst MPI processes
remicres's avatar
remicres committed
288
  unsigned long int maximumNumberOfNodesInMemory = GetMaximumNumberOfNodesInMemory();
remicres's avatar
remicres committed
289
  itkDebugMacro(<<"Maximum number of nodes in memory is " << maximumNumberOfNodesInMemory);
remicres's avatar
remicres committed
290

remicres's avatar
remicres committed
291
  // Number of nodes in the entire image
remicres's avatar
remicres committed
292
293
  const unsigned int imageWidth = m_InputImage->GetLargestPossibleRegion().GetSize()[0];
  const unsigned int imageHeight = m_InputImage->GetLargestPossibleRegion().GetSize()[1];
remicres's avatar
remicres committed
294
295
  const unsigned long int nbOfNodesInImage = imageWidth*imageHeight;

296
  // Default layout: 1x1
remicres's avatar
remicres committed
297
298
  m_NbTilesX = 1;
  m_NbTilesY = 1;
299

remicres's avatar
remicres committed
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
  // 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++)
319
    {
remicres's avatar
remicres committed
320
321
    // 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
322
      {
remicres's avatar
remicres committed
323
324
#ifdef OTB_USE_MPI
      // We want number of tiles which is a multiple of the number of MPI processes
325
326
      if (nbOfTiles % layoutNCol == 0 && // Is it a multiple of the nb of Tiles and nProcs?
          nbOfTiles % otb::MPIConfig::Instance()->GetNbProcs() == 0)
remicres's avatar
remicres committed
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
#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;
348
          float compactness = perimeter / surface * (float) std::max(tileWidth,tileHeight); // [1,+inf]
remicres's avatar
remicres committed
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371

          // 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
372
373
374
375

  // 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
376
377
378
379
380
381
  itkDebugMacro(<<"Selected layout: " << m_NbTilesX << "x" << m_NbTilesY
      << " (criterion=" << lowestCriterionValue << ")");

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

382
383
384
  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);
385
  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
386
      << (m_TileWidth + 2*m_Margin) << "x" << (m_TileHeight + 2*m_Margin) );
remicres's avatar
remicres committed
387

remicres's avatar
remicres committed
388
389
}

remicres's avatar
remicres committed
390
391
392
393
394
395
396
397
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
398
void Controller<TSegmenter>::SetInputImage(ImageType * inputImage)
remicres's avatar
remicres committed
399
{
remicres's avatar
remicres committed
400
  m_InputImage = inputImage;
remicres's avatar
remicres committed
401
402
403
404
405
406
407
408
}

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

remicres's avatar
remicres committed
409
410
411
412
413
414
template<class TSegmenter>
typename Controller<TSegmenter>::LabelImageType::Pointer
Controller<TSegmenter>::GetLabeledClusteredOutput()
{
  return m_LabelImage;
}
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429

template <class TSegmenter>
std::vector<std::string> Controller<TSegmenter>::GetTemporaryFilesList()
{
  std::vector<std::string> list;
  for (unsigned int i = 0; i < m_Tiles.size(); i++)
    {
    list.push_back(m_Tiles[i].edgeFileName);
    list.push_back(m_Tiles[i].edgeMarginFileName);
    list.push_back(m_Tiles[i].nodeFileName);
    list.push_back(m_Tiles[i].nodeMarginFileName);
    }
  return list;
}

remicres's avatar
remicres committed
430
} // end of namespace lsgrm