otbPipelineMemoryPrintCalculator.cxx 8.64 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
/*
 * Copyright (C) 2005-2017 Centre National d'Etudes Spatiales (CNES)
 *
 * This file is part of Orfeo Toolbox
 *
 *     https://www.orfeo-toolbox.org/
 *
 * Licensed under the Apache License, Version 2.0 (the "License");
 * you may not use this file except in compliance with the License.
 * You may obtain a copy of the License at
 *
 *     http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the License is distributed on an "AS IS" BASIS,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the License for the specific language governing permissions and
 * limitations under the License.
 */
20

21

22 23 24
#include <complex>
#include <ostream>

25 26
#include "otbPipelineMemoryPrintCalculator.h"

27
#include "otbMacro.h"
28
#include "otbMath.h"
29 30
#include "otbImage.h"
#include "otbVectorImage.h"
31
#include "itkFixedArray.h"
32
#include "otbImageList.h"
33

34 35
namespace otb
{
Emmanuel Christophe's avatar
Emmanuel Christophe committed
36
const double PipelineMemoryPrintCalculator::ByteToMegabyte = 1./vcl_pow(2.0, 20);
37
const double PipelineMemoryPrintCalculator::MegabyteToByte = vcl_pow(2.0, 20);
38 39 40 41

PipelineMemoryPrintCalculator
::PipelineMemoryPrintCalculator()
  : m_MemoryPrint(0),
42
    m_DataToWrite(nullptr),
43 44 45 46 47 48
    m_BiasCorrectionFactor(1.),
    m_VisitedProcessObjects()
{}

PipelineMemoryPrintCalculator
::~PipelineMemoryPrintCalculator()
49 50
{}

51 52 53 54 55 56 57 58 59 60 61
// [static]
unsigned long
PipelineMemoryPrintCalculator
::EstimateOptimalNumberOfStreamDivisions(MemoryPrintType memoryPrint, MemoryPrintType availableMemory)
{
  unsigned long divisions;
  divisions = vcl_ceil(static_cast<double>(memoryPrint)
                       / availableMemory);
  return divisions;
}

62
void
63 64 65 66
PipelineMemoryPrintCalculator
::PrintSelf(std::ostream& os, itk::Indent indent) const
{
  // Call superclass implementation
OTB Bot's avatar
STYLE  
OTB Bot committed
67
  Superclass::PrintSelf(os, indent);
68 69 70

  // Display parameters
  os<<indent<<"Data to write:                      "<<m_DataToWrite<<std::endl;
71
  os<<indent<<"Memory print of whole pipeline:     "<<m_MemoryPrint * ByteToMegabyte <<" Mb"<<std::endl;
72 73 74
  os<<indent<<"Bias correction factor applied:     "<<m_BiasCorrectionFactor<<std::endl;
}

75
void
76
PipelineMemoryPrintCalculator
77
::Compute(bool propagate)
78
{
79
  // Clear the visited process objects set
80 81
  m_VisitedProcessObjects.clear();

82
  // Dry run of pipeline synchronisation
83 84 85 86 87 88
  if (propagate)
    {
    m_DataToWrite->UpdateOutputInformation();
    m_DataToWrite->SetRequestedRegionToLargestPossibleRegion();
    m_DataToWrite->PropagateRequestedRegion();
    }
89 90 91 92 93 94 95 96

  // Get the source process object
  ProcessObjectType * source = m_DataToWrite->GetSource();

  // Check if source exists
  if(source)
    {
    // Call the recursive memory print evaluation
97
    m_MemoryPrint = EvaluateProcessObjectPrintRecursive(source);
98 99 100 101 102 103 104 105
    }
  else
    {
    // Get memory print for this data only
    m_MemoryPrint = EvaluateDataObjectPrint(m_DataToWrite);
    }

  // Apply bias correction factor
Emmanuel Christophe's avatar
Emmanuel Christophe committed
106
  m_MemoryPrint *= m_BiasCorrectionFactor;
107 108 109

}

110
PipelineMemoryPrintCalculator::MemoryPrintType
111
PipelineMemoryPrintCalculator
112
::EvaluateProcessObjectPrintRecursive(ProcessObjectType * process)
113
{
114
  otbLogMacro(Debug,<<"Recursive evaluation of memory print for ProcessObject" << process->GetNameOfClass() << " (" << process << ")");
115
  // This variable will store the final print
116
  MemoryPrintType print = 0;
117

118 119 120 121 122 123 124 125 126 127 128
  // Check if this process object has already been visited
  if(m_VisitedProcessObjects.count(process))
    {
    return print;
    }
  // Else register it as a visited process object
  else
    {
    m_VisitedProcessObjects.insert(process);
    }

129 130 131
  // Retrieve the array of inputs
  ProcessObjectType::DataObjectPointerArray inputs = process->GetInputs();
  // First, recurse on each input source
Emmanuel Christophe's avatar
Emmanuel Christophe committed
132
  for(unsigned int i = 0; i < process->GetNumberOfInputs(); ++i)
133 134 135
    {
    // Retrieve the data object
    DataObjectType * input = inputs[i];
136

137
    if( input )
138
      {
139 140
        // Retrieve possible source
        ProcessObjectType * source = input->GetSource();
141

142 143 144
        // If data object has a source
        if(source)
          {
145
            print += this->EvaluateProcessObjectPrintRecursive(source);
146 147 148 149 150 151
          }
        else
          {
            MemoryPrintType localPrint = this->EvaluateDataObjectPrint(input);
            print += localPrint;
          }
152 153
      }
    }
154

155 156 157 158
  // Retrieve the output array
  ProcessObjectType::DataObjectPointerArray outputs = process->GetOutputs();

  // Now, evaluate the current object print
Emmanuel Christophe's avatar
Emmanuel Christophe committed
159
  for(unsigned int i = 0; i < process->GetNumberOfOutputs(); ++i)
160
    {
161
      MemoryPrintType localPrint = this->EvaluateDataObjectPrint(outputs[i]);
162
      print += localPrint;
163 164 165 166 167 168
    }

  // Finally, return the total print
  return print;
}

169
PipelineMemoryPrintCalculator::MemoryPrintType
170
PipelineMemoryPrintCalculator
171
::EvaluateDataObjectPrint(DataObjectType * data)
172
{
173 174 175
    
  otbLogMacro(Debug,<<"Evaluation of memory print for DataObject " << data->GetNameOfClass() << " (" << data << ")");
    
176
#define OTB_IMAGE_SIZE_BLOCK(type)                                      \
OTB Bot's avatar
STYLE  
OTB Bot committed
177
  if(dynamic_cast<itk::Image<type, 2> *>(data) != NULL)                  \
178
    {                                                                   \
OTB Bot's avatar
STYLE  
OTB Bot committed
179
    itk::Image<type, 2> * image = dynamic_cast<itk::Image<type, 2> *>(data); \
180
    return image->GetRequestedRegion().GetNumberOfPixels()              \
181
      * image->GetNumberOfComponentsPerPixel() * sizeof(type); \
182
    }                                                                   \
OTB Bot's avatar
STYLE  
OTB Bot committed
183
  if(dynamic_cast<itk::VectorImage<type, 2> * >(data) != NULL)           \
184
    {                                                                   \
OTB Bot's avatar
STYLE  
OTB Bot committed
185
    itk::VectorImage<type, 2> * image = dynamic_cast<itk::VectorImage<type, 2> *>(data); \
186
    return image->GetRequestedRegion().GetNumberOfPixels()              \
187
      * image->GetNumberOfComponentsPerPixel() * sizeof(type); \
188
    }                                                                   \
OTB Bot's avatar
STYLE  
OTB Bot committed
189
  if(dynamic_cast<ImageList<Image<type, 2> > *>(data) != NULL)   \
190
    {                                                                   \
OTB Bot's avatar
STYLE  
OTB Bot committed
191
    ImageList<Image<type, 2> > * imageList = dynamic_cast<otb::ImageList<otb::Image<type, 2> > *>(data); \
192
    MemoryPrintType print(0);                                         \
193
    for(ImageList<Image<type, 2> >::Iterator it = imageList->Begin(); \
OTB Bot's avatar
STYLE  
OTB Bot committed
194
       it != imageList->End(); ++it)                                    \
195
       {                                                             \
196 197 198 199
       if(it.Get()->GetSource())                                        \
         print += this->EvaluateProcessObjectPrintRecursive(it.Get()->GetSource());\
       else                                                             \
         print += this->EvaluateDataObjectPrint(it.Get());              \
200 201 202
       }                                                           \
    return print;                                                  \
    }                                                              \
OTB Bot's avatar
STYLE  
OTB Bot committed
203
  if(dynamic_cast<ImageList<VectorImage<type, 2> > *>(data) != NULL)   \
204
    {                                                                   \
OTB Bot's avatar
STYLE  
OTB Bot committed
205
    ImageList<VectorImage<type, 2> > * imageList = dynamic_cast<otb::ImageList<otb::VectorImage<type, 2> > *>(data); \
206
    MemoryPrintType print(0);                                         \
OTB Bot's avatar
STYLE  
OTB Bot committed
207 208
    for(ImageList<VectorImage<type, 2> >::ConstIterator it = imageList->Begin(); \
       it != imageList->End(); ++it)                                    \
209
       {                                                             \
210 211 212 213
       if(it.Get()->GetSource())                                        \
         print += this->EvaluateProcessObjectPrintRecursive(it.Get()->GetSource());\
       else                                                             \
         print += this->EvaluateDataObjectPrint(it.Get());              \
214 215 216
       }                                                           \
    return print;                                                  \
    }                                                              \
217

218

219 220
  // Call the macro for each pixel type
  OTB_IMAGE_SIZE_BLOCK(unsigned char)
Emmanuel Christophe's avatar
Emmanuel Christophe committed
221 222 223 224 225 226 227 228 229 230 231
  OTB_IMAGE_SIZE_BLOCK(char)
  OTB_IMAGE_SIZE_BLOCK(unsigned short)
  OTB_IMAGE_SIZE_BLOCK(short)
  OTB_IMAGE_SIZE_BLOCK(unsigned int)
  OTB_IMAGE_SIZE_BLOCK(int)
  OTB_IMAGE_SIZE_BLOCK(unsigned long)
  OTB_IMAGE_SIZE_BLOCK(long)
  OTB_IMAGE_SIZE_BLOCK(float)
  OTB_IMAGE_SIZE_BLOCK(double)
  OTB_IMAGE_SIZE_BLOCK(std::complex<float>)
  OTB_IMAGE_SIZE_BLOCK(std::complex<double>)
OTB Bot's avatar
STYLE  
OTB Bot committed
232 233
  typedef itk::FixedArray<float, 2> FloatFixedArray2Type;
  typedef itk::FixedArray<float, 2> DoubleFixedArray2Type;
Emmanuel Christophe's avatar
Emmanuel Christophe committed
234 235 236 237 238
  OTB_IMAGE_SIZE_BLOCK(FloatFixedArray2Type)
  OTB_IMAGE_SIZE_BLOCK(DoubleFixedArray2Type)

  // If we are still here, none of the macro call succeed
  return 0;
239 240 241
}

} // End namespace otb