OTB  7.4.0
Orfeo Toolbox
otbSARStreamingGridHistogramFilter.h
Go to the documentation of this file.
1 /*
2  * Copyright (C) 2005-2018 Centre National d'Etudes Spatiales (CNES)
3  *
4  * This file is part of Orfeo Toolbox
5  *
6  * https://www.orfeo-toolbox.org/
7  *
8  * Licensed under the Apache License, Version 2.0 (the "License");
9  * you may not use this file except in compliance with the License.
10  * You may obtain a copy of the License at
11  *
12  * http://www.apache.org/licenses/LICENSE-2.0
13  *
14  * Unless required by applicable law or agreed to in writing, software
15  * distributed under the License is distributed on an "AS IS" BASIS,
16  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
17  * See the License for the specific language governing permissions and
18  * limitations under the License.
19  */
20 
21 #ifndef otbSARStreamingGridHistogramFilter_h
22 #define otbSARStreamingGridHistogramFilter_h
23 
25 #include "itkNumericTraits.h"
26 #include "itkArray.h"
27 #include "itkSimpleDataObjectDecorator.h"
28 #include "itkPoint.h"
29 
31 
33 #include "otbImageKeywordlist.h"
34 
35 #include <complex>
36 #include <cmath>
37 
38 namespace otb
39 {
40 
55 template<class TInputImage>
57  public PersistentImageFilter<TInputImage, TInputImage>
58 {
59 public:
63  typedef itk::SmartPointer<Self> Pointer;
64  typedef itk::SmartPointer<const Self> ConstPointer;
65 
67  itkNewMacro(Self);
68 
71 
73  typedef TInputImage ImageType;
74  typedef typename TInputImage::Pointer InputImagePointer;
75 
76  typedef typename TInputImage::RegionType RegionType;
77  typedef typename TInputImage::SizeType SizeType;
79  typedef typename TInputImage::PixelType PixelType;
80  typedef typename ImageType::IndexValueType IndexValueType;
81  typedef typename ImageType::SizeValueType SizeValueType;
82 
83  itkStaticConstMacro(InputImageDimension, unsigned int,
84  TInputImage::ImageDimension);
85 
87  itkStaticConstMacro(ImageDimension, unsigned int,
88  TInputImage::ImageDimension);
89 
91  // typedef typename itk::NumericTraits<PixelType>::RealType RealType;
93 
95  typedef typename itk::DataObject::Pointer DataObjectPointer;
96  typedef itk::ProcessObject::DataObjectPointerArraySizeType DataObjectPointerArraySizeType;
97 
99  typedef itk::SimpleDataObjectDecorator<RealType> RealObjectType;
100  typedef itk::SimpleDataObjectDecorator<std::vector<RealType> *> RealVectorObjectType;
101  typedef itk::SimpleDataObjectDecorator<long> LongObjectType;
102  typedef itk::SimpleDataObjectDecorator<PixelType> PixelObjectType;
103 
105  std::vector<RealType> * GetVectorShifts() const
106  {
107  return this->GetVectorShifts_Output()->Get();
108  }
109  RealVectorObjectType* GetVectorShifts_Output();
110  const RealVectorObjectType* GetVectorShifts_Output() const;
112 
113 
116  DataObjectPointer MakeOutput(DataObjectPointerArraySizeType idx) ITK_OVERRIDE;
117  using Superclass::MakeOutput;
118 
122  void Synthetize(void) ITK_OVERRIDE;
123  void Reset(void) ITK_OVERRIDE;
125 
126  // Getter/Setter
127  itkSetMacro(StepHist, RealType);
128  itkGetMacro(StepHist, RealType);
129  itkSetMacro(MeanShift, RealType);
130  itkGetMacro(MeanShift, RealType);
131  itkSetMacro(Dimension, int);
132  itkGetMacro(Dimension, int);
133  itkSetMacro(Threshold, float);
134  itkGetMacro(Threshold, float);
135  itkSetMacro(NumberHist, int);
136  itkGetMacro(NumberHist, int);
137 
138 protected:
141  void PrintSelf(std::ostream& os, itk::Indent indent) const ITK_OVERRIDE;
142 
143 
147  virtual void GenerateInputRequestedRegion() ITK_OVERRIDE;
148 
150  void ThreadedGenerateData(const RegionType&
151  outputRegionForThread,
152  itk::ThreadIdType threadId) ITK_OVERRIDE;
153 
157  void AllocateOutputs() ITK_OVERRIDE;
158  void GenerateOutputInformation() ITK_OVERRIDE;
160 
161 private:
162  PersistentGridHistogramFilter(const Self &); //purposely not implemented
163  void operator =(const Self&); //purposely not implemented
164 
165  std::vector<std::vector<RealType> * > m_Thread_Shifts;
166 
167  // Step for the histogram
168  RealType m_StepHist;
169  RealType m_MeanShift;
170  // Dimension for the selection : 0 = range, 1 = azimut
172  // Threshold on correlation rate
173  float m_Threshold;
174  // Number of values range
176 
177 
178 }; // end of class PersistentGridHistogramFilter
179 
180 /*===========================================================================*/
181 
218 template<class TInputImage>
220  public PersistentFilterStreamingDecorator<PersistentGridHistogramFilter<TInputImage> >
221 {
222 public:
227  typedef itk::SmartPointer<Self> Pointer;
228  typedef itk::SmartPointer<const Self> ConstPointer;
229 
231  itkNewMacro(Self);
232 
234  itkTypeMacro(StreamingGridHistogramFilter, PersistentFilterStreamingDecorator);
235 
236  typedef typename Superclass::FilterType StatFilterType;
239  typedef TInputImage InputImageType;
240 
242  typedef itk::SimpleDataObjectDecorator<RealType> RealObjectType;
243  typedef itk::SimpleDataObjectDecorator<long> LongObjectType;
244  typedef itk::SimpleDataObjectDecorator<PixelType> PixelObjectType;
245  typedef itk::SimpleDataObjectDecorator<std::vector<RealType> *> RealVectorObjectType;
246 
248  typedef TInputImage ImageType;
249  typedef typename TInputImage::Pointer InputImagePointer;
250 
251  typedef typename TInputImage::RegionType RegionType;
252  typedef typename TInputImage::SizeType SizeType;
254 
255  // Define Point2DType and Point3DType
256  using Point2DType = itk::Point<double,2>;
257  using Point3DType = itk::Point<double,3>;
258 
259  using Superclass::SetInput;
260  void SetInput(InputImageType * input)
261  {
262  this->GetFilter()->SetInput(input);
263  }
264  const InputImageType * GetInput()
265  {
266  return this->GetFilter()->GetInput();
267  }
268 
269  // Setter for Streaming and Persistent Filter
270  void SetStepHist(RealType stepHist)
271  {
272  this->GetFilter()->SetStepHist(stepHist);
273  }
274  void SetMeanShift(RealType meanShift)
275  {
276  this->GetFilter()->SetMeanShift(meanShift);
277  }
278  void SetThreshold(float threshold)
279  {
280  this->GetFilter()->SetThreshold(threshold);
281  }
282  void SetNumberHist(int numberHist)
283  {
284  this->GetFilter()->SetNumberHist(numberHist);
285  m_NumberHist = numberHist;
286  }
287  void SetDimension(int dimension)
288  {
289  if (dimension == 0 || dimension == 1)
290  {
291  this->GetFilter()->SetDimension(dimension);
292  }
293  else
294  {
295  std::cout << "Wrong dimension : Only 0 for range and 1 for azimut" << std::endl;
296  this->GetFilter()->SetDimension(0);
297  }
298  }
299 
300 
302  std::vector<RealType> * GetVectorShifts() const
303  {
304  return this->GetFilter()->GetVectorShifts_Output()->Get();
305  }
306  RealVectorObjectType* GetVectorShifts_Output()
307  {
308  return this->GetFilter()->GetVectorShifts_Output();
309  }
310  const RealVectorObjectType* GetVectorShifts_Output() const
311  {
312  return this->GetFilter()->GetVectorShifts_Output();
313  }
315 
316 
317  // Get step and max occurence for the current histogram
318  void GetGridHistogramInformation(double & step, double & meanShift)
319  {
320  // Retrieve selected shifts
321  std::vector<RealType> * shifts = this->GetVectorShifts();
322 
323  // Estimate min, max and mean for each dimension
324  double minShift = 1E20;
325  double maxShift = -1E20;
326 
327  double sumShift = 0.;
328 
329  for (long unsigned int i = 0; i < shifts->size(); i++)
330  {
331  sumShift += shifts->at(i);
332  minShift = std::min(minShift, shifts->at(i));
333  maxShift = std::max(maxShift, shifts->at(i));
334  }
335 
336  meanShift = sumShift/shifts->size();
337 
338  // Assignate the step
339  step = (maxShift - minShift) / (float)(m_NumberHist - 1.);
340 
341  // Initialize the histogram
342  float *histo = new float[m_NumberHist];
343  float *plages = new float[m_NumberHist];
344  for (int index = 0; index < m_NumberHist; index++)
345  {
346  histo[index] = 0;
347  plages[index] = minShift + index * step;
348  }
349 
350  // Estimate the histogram
351  for (long unsigned int i = 0; i < shifts->size(); i++)
352  {
353  int index = int((shifts->at(i) - minShift) / step);
354  if ((index >= 0) && (index < m_NumberHist))
355  {
356  ++histo[index];
357  }
358  }
359 
360  // Determinate the maximum
361  int index_top = 0;
362 
363  for (int index = 1; index < m_NumberHist; index++)
364  {
365  if (histo[index] > histo[index_top])
366  {
367  index_top = index;
368  }
369  }
370 
371  // Assignate meanShift
372  meanShift = plages[index_top];
373 
374  delete [] histo;
375  delete [] plages;
376  }
377 
378 
379 protected:
382 
385 
386 private:
387  SARStreamingGridHistogramFilter(const Self &); //purposely not implemented
388  void operator =(const Self&); //purposely not implemented
389 
391 
392 };
393 
394 } // end namespace otb
395 
396 #ifndef OTB_MANUAL_INSTANTIATION
397 #include "otbSARStreamingGridHistogramFilter.txx"
398 #endif
399 
400 #endif
itk::Size< Monteverdi_DIMENSION > SizeType
Definition: mvdTypes.h:137
itk::Index< Monteverdi_DIMENSION > IndexType
Definition: mvdTypes.h:133
void GetGridHistogramInformation(double &step, double &meanShift)
itk::SimpleDataObjectDecorator< std::vector< RealType > * > RealVectorObjectType
PersistentImageFilter< TInputImage, TInputImage > Superclass
Retrive some shifts from the input grid. The grid is a vector image of three components : shift_range...
itk::SimpleDataObjectDecorator< PixelType > PixelObjectType
itk::SimpleDataObjectDecorator< long > LongObjectType
PersistentFilterStreamingDecorator< PersistentGridHistogramFilter< TInputImage > > Superclass
This class streams the whole input image through the PersistentGridHistogramFilter.
std::vector< RealType > * GetVectorShifts() const
itk::SimpleDataObjectDecorator< std::vector< RealType > * > RealVectorObjectType
Monteverdi_FLOATING_TYPE RealType
Definition: mvdTypes.h:84
itk::SimpleDataObjectDecorator< RealType > RealObjectType
The "otb" namespace contains all Orfeo Toolbox (OTB) classes.
itk::ProcessObject::DataObjectPointerArraySizeType DataObjectPointerArraySizeType
const RealVectorObjectType * GetVectorShifts_Output() const
itk::SimpleDataObjectDecorator< RealType > RealObjectType
itk::SimpleDataObjectDecorator< PixelType > PixelObjectType
itk::NumericTraits< double >::RealType RealType
std::vector< std::vector< RealType > *> m_Thread_Shifts
itk::SimpleDataObjectDecorator< long > LongObjectType
This filter is the base class for all filter persisting data through multiple update. For instance, a filter computing global statistics on an image with streaming capabilities will have to keep the temporary results for each streamed piece of the image in order to synthesize the global statistics at the end. This filter is an itk::ImageToImageFilter, providing two additional methods. The first one, Synthetize(), allows the user to synthesize temporary data produced by the multiple updates on different pieces of the image to the global result. The second one, Reset(), allows the user to reset the temporary data for a new input image for instance.
This filter link a persistent filter with a StreamingImageVirtualWriter.