OTB  5.0.0
Orfeo Toolbox
otbWaveletFilterBank.txx
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: ORFEO Toolbox
4  Language: C++
5  Date: $Date$
6  Version: $Revision$
7 
8 
9  Copyright (c) Centre National d'Etudes Spatiales. All rights reserved.
10  See OTBCopyright.txt for details.
11 
12  Copyright (c) Institut Telecom / Telecom Bretagne. All rights reserved.
13  See ITCopyright.txt for details.
14 
15  This software is distributed WITHOUT ANY WARRANTY; without even
16  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
17  PURPOSE. See the above copyright notices for more information.
18 
19 =========================================================================*/
20 
21 #ifndef __otbWaveletFilterBank_txx
22 #define __otbWaveletFilterBank_txx
23 #include "otbWaveletFilterBank.h"
24 
27 
29 
30 // FIXME
31 #define __myDebug__ 0
32 
33 namespace otb {
34 
39 template <class TInputImage, class TOutputImage, class TWaveletOperator>
42 {
43  this->SetNumberOfRequiredInputs(1);
44  this->SetNumberOfRequiredInputs(1);
45 
46  unsigned int numOfOutputs = 1 << InputImageDimension;
47 
48  this->SetNumberOfRequiredOutputs(numOfOutputs);
49  for (unsigned int i = 0; i < numOfOutputs; ++i)
50  {
51  this->SetNthOutput(i, OutputImageType::New());
52  }
53 
54  m_UpSampleFilterFactor = 0;
55  m_SubsampleImageFactor = 1;
56 
57  //this->SetNumberOfThreads(1);
58 }
59 
60 template <class TInputImage, class TOutputImage, class TWaveletOperator>
61 void
64 {
65  Superclass::GenerateOutputInformation();
66 
67  if (GetSubsampleImageFactor() == 1) return;
68 
69 #if __myDebug__
70  otbGenericMsgDebugMacro(<< " down sampling output regions by a factor of " << GetSubsampleImageFactor());
71  otbGenericMsgDebugMacro(<< "initial region " << this->GetInput()->GetLargestPossibleRegion().GetSize()[0]
72  << "," << this->GetInput()->GetLargestPossibleRegion().GetSize()[1]);
73 #endif
74 
75  OutputImageRegionType newRegion;
76  this->CallCopyInputRegionToOutputRegion(newRegion, this->GetInput()->GetLargestPossibleRegion());
77 
78  for (unsigned int i = 0; i < this->GetNumberOfOutputs(); ++i)
79  {
80  this->GetOutput(i)->SetRegions(newRegion);
81  }
82 
83 #if __myDebug__
84  otbGenericMsgDebugMacro(<< "new region output " << newRegion.GetSize()[0] << "," << newRegion.GetSize()[1]);
85 #endif
86 }
87 
88 template <class TInputImage, class TOutputImage, class TWaveletOperator>
89 void
92 throw (itk::InvalidRequestedRegionError)
93  {
94  Superclass::GenerateInputRequestedRegion();
95 
96  // Filter length calculation
97  LowPassOperatorType lowPassOperator;
98  lowPassOperator.SetDirection(0);
99  lowPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
100  lowPassOperator.CreateDirectional();
101 
102  unsigned int radius = lowPassOperator.GetRadius()[0];
103 
104  HighPassOperatorType highPassOperator;
105  highPassOperator.SetDirection(0);
106  highPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
107  highPassOperator.CreateDirectional();
108 
109  if (radius < highPassOperator.GetRadius()[0]) radius = highPassOperator.GetRadius()[0];
110 
111  // Get the requested region and pad it
112  InputImagePointerType input = const_cast<InputImageType*>(this->GetInput());
113  InputImageRegionType inputRegion = input->GetRequestedRegion();
114  inputRegion.PadByRadius(radius);
115 
116  if (inputRegion.Crop(input->GetLargestPossibleRegion()))
117  {
118  input->SetRequestedRegion(inputRegion);
119  }
120  else
121  {
122  input->SetRequestedRegion(inputRegion);
123  itk::InvalidRequestedRegionError err(__FILE__, __LINE__);
124  err.SetLocation(ITK_LOCATION);
125  err.SetDescription("Requested region is (at least partially) outside the largest possible region.");
126  err.SetDataObject(input);
127  throw err;
128  }
129  }
130 
131 template <class TInputImage, class TOutputImage, class TWaveletOperator>
132 void
135 {
136  if (m_SubsampleImageFactor > 1)
137  {
138  // Check the dimension
139  for (unsigned int i = 0; i < InputImageDimension; ++i)
140  {
141  if ((m_SubsampleImageFactor
142  * (this->GetInput()->GetRequestedRegion().GetSize()[i] / m_SubsampleImageFactor))
143  != this->GetInput()->GetRequestedRegion().GetSize()[i])
144  {
145  itk::InvalidRequestedRegionError err(__FILE__, __LINE__);
146  err.SetLocation(ITK_LOCATION);
147  err.SetDescription("Requested region dimension cannot be used in multiresolution analysis (crop it).");
148  err.SetDataObject(const_cast<InputImageType*>(this->GetInput()));
149  throw err;
150  }
151  }
152 
153  if (InputImageDimension > 1)
154  {
155  // Internal images will be used only if m_SubsampledInputImages != 1
156  m_InternalImages.resize(InputImageDimension - 1);
157  for (unsigned int i = 0; i < m_InternalImages.size(); ++i)
158  {
159  // the size is linked to the SubsampleImageFactor that is assume to be 2!!!
160  m_InternalImages[InputImageDimension - 2 - i].resize(1 << (i + 1));
161  }
162 
163  OutputImageRegionType intermediateRegion;
164  this->Superclass::CallCopyInputRegionToOutputRegion(intermediateRegion,
165  this->GetInput()->GetLargestPossibleRegion());
166 
167  AllocateInternalData(intermediateRegion);
168  }
169  }
170 }
171 
172 template <class TInputImage, class TOutputImage, class TWaveletOperator>
173 void
176  (const OutputImageRegionType& outputRegion)
177 {
178  OutputImageRegionType smallerRegion;
179  OutputImageRegionType largerRegion = outputRegion;
180 
181  for (unsigned int direction = 0; direction < InputImageDimension - 1; direction++)
182  {
183  this->CallCopyInputRegionToOutputRegion(InputImageDimension - 1 - direction,
184  smallerRegion, largerRegion);
185 
186  for (unsigned int i = 0; i < m_InternalImages[direction].size(); ++i)
187  {
188  m_InternalImages[InputImageDimension - 2 - direction][i] = OutputImageType::New();
189  m_InternalImages[InputImageDimension - 2 - direction][i]->SetRegions(smallerRegion);
190  m_InternalImages[InputImageDimension - 2 - direction][i]->Allocate();
191  m_InternalImages[InputImageDimension - 2 - direction][i]->FillBuffer(0);
192  }
193 
194  largerRegion = smallerRegion;
195  }
196 }
197 
198 template <class TInputImage, class TOutputImage, class TWaveletOperator>
199 void
202 {
203  if (m_SubsampleImageFactor > 1 && InputImageDimension > 1)
204  {
205  m_InternalImages.clear();
206  }
207 }
208 
209 template <class TInputImage, class TOutputImage, class TWaveletOperator>
210 void
213  (InputImageRegionType& destRegion, const OutputImageRegionType& srcRegion)
214 {
215  Superclass::CallCopyOutputRegionToInputRegion(destRegion, srcRegion);
216 
217  if (GetSubsampleImageFactor() > 1)
218  {
219  OutputIndexType srcIndex = srcRegion.GetIndex();
220  OutputSizeType srcSize = srcRegion.GetSize();
221 
222  InputIndexType destIndex;
223  InputSizeType destSize;
224 
225  for (unsigned int i = 0; i < InputImageDimension; ++i)
226  {
227  destIndex[i] = srcIndex[i] * GetSubsampleImageFactor();
228  destSize[i] = srcSize[i] * GetSubsampleImageFactor();
229  }
230 
231  destRegion.SetIndex(destIndex);
232  destRegion.SetSize(destSize);
233 
234 #if 0
235  // Contrairement a Wavelet::INVERSE, ici ca ne sera a rien apparemment...
236 
237  // Region Padding
238  LowPassOperatorType lowPassOperator;
239  lowPassOperator.SetDirection(0);
240  lowPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
241  lowPassOperator.CreateDirectional();
242 
243  unsigned long radius[InputImageDimension];
244  radius[0] = lowPassOperator.GetRadius()[0];
245 
246  HighPassOperatorType highPassOperator;
247  highPassOperator.SetDirection(0);
248  highPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
249  highPassOperator.CreateDirectional();
250 
251  if (radius[0] < highPassOperator.GetRadius()[0]) radius[0] = highPassOperator.GetRadius()[0];
252 
253  for (unsigned int i = 1; i < InputImageDimension; ++i)
254  radius[i] = 0;
255 
256  InputImageRegionType paddedRegion = destRegion;
257  paddedRegion.PadByRadius(radius);
258 
259  if (paddedRegion.Crop(this->GetInput()->GetLargestPossibleRegion()))
260  {
261  destRegion = paddedRegion;
262  }
263 #endif
264 
265  }
266 }
267 
268 template <class TInputImage, class TOutputImage, class TWaveletOperator>
269 void
272  (unsigned int direction,
273  InputImageRegionType& destRegion, const OutputImageRegionType& srcRegion)
274 {
275  Superclass::CallCopyOutputRegionToInputRegion(destRegion, srcRegion);
276 
277  if (GetSubsampleImageFactor() > 1)
278  {
279  OutputIndexType srcIndex = srcRegion.GetIndex();
280  OutputSizeType srcSize = srcRegion.GetSize();
281 
282  InputIndexType destIndex;
283  InputSizeType destSize;
284 
285  for (unsigned int i = 0; i < InputImageDimension; ++i)
286  {
287  if (i == direction)
288  {
289  destIndex[i] = srcIndex[i] * GetSubsampleImageFactor();
290  destSize[i] = srcSize[i] * GetSubsampleImageFactor();
291  }
292  else
293  {
294  destIndex[i] = srcIndex[i];
295  destSize[i] = srcSize[i];
296  }
297  }
298 
299  destRegion.SetIndex(destIndex);
300  destRegion.SetSize(destSize);
301  }
302 }
303 
304 template <class TInputImage, class TOutputImage, class TWaveletOperator>
305 void
308  (OutputImageRegionType& destRegion, const InputImageRegionType& srcRegion)
309 {
310  Superclass::CallCopyInputRegionToOutputRegion(destRegion, srcRegion);
311 
312  if (GetSubsampleImageFactor() > 1)
313  {
314  typename InputImageRegionType::IndexType srcIndex = srcRegion.GetIndex();
315  typename InputImageRegionType::SizeType srcSize = srcRegion.GetSize();
316 
317  typename OutputImageRegionType::IndexType destIndex;
318  typename OutputImageRegionType::SizeType destSize;
319 
320  for (unsigned int i = 0; i < InputImageDimension; ++i)
321  {
322  // TODO: This seems not right in odd index cases
323  destIndex[i] = srcIndex[i] / GetSubsampleImageFactor();
324  destSize[i] = srcSize[i] / GetSubsampleImageFactor();
325  }
326 
327  destRegion.SetIndex(destIndex);
328  destRegion.SetSize(destSize);
329  }
330 }
331 
332 template <class TInputImage, class TOutputImage, class TWaveletOperator>
333 void
336  (unsigned int direction,
337  OutputImageRegionType& destRegion, const InputImageRegionType& srcRegion)
338 {
339  Superclass::CallCopyInputRegionToOutputRegion(destRegion, srcRegion);
340 
341  if (GetSubsampleImageFactor() > 1)
342  {
343  typename InputImageRegionType::IndexType srcIndex = srcRegion.GetIndex();
344  typename InputImageRegionType::SizeType srcSize = srcRegion.GetSize();
345 
346  typename OutputImageRegionType::IndexType destIndex;
347  typename OutputImageRegionType::SizeType destSize;
348 
349  for (unsigned int i = 0; i < InputImageDimension; ++i)
350  {
351  if (i == direction)
352  {
353  // TODO: This seems not right in odd index cases
354  destIndex[i] = srcIndex[i] / GetSubsampleImageFactor();
355  destSize[i] = srcSize[i] / GetSubsampleImageFactor();
356  }
357  else
358  {
359  destIndex[i] = srcIndex[i];
360  destSize[i] = srcSize[i];
361  }
362  }
363 
364  destRegion.SetIndex(destIndex);
365  destRegion.SetSize(destSize);
366  }
367 }
368 
369 template <class TInputImage, class TOutputImage, class TWaveletOperator>
370 void
373  (const OutputImageRegionType& outputRegionForThread, itk::ThreadIdType threadId)
374 {
375  unsigned int dir = InputImageDimension - 1;
376 
377  if ((1 << dir) >= static_cast<int>(this->GetNumberOfOutputs()))
378  {
379  std::ostringstream msg;
380  msg << "Output number 1<<" << dir << " = " << (1 << dir) << " not allocated\n";
381  msg << "Number of excpected outputs " << this->GetNumberOfOutputs() << "\n";
382  throw itk::ExceptionObject(__FILE__, __LINE__, msg.str().c_str(), ITK_LOCATION);
383  }
384 
385  itk::ProgressReporter reporter(this, threadId,
386  outputRegionForThread.GetNumberOfPixels() * this->GetNumberOfOutputs() * 2);
387 
388  const InputImageType * input = this->GetInput();
389  InputImageRegionType inputRegionForThread;
390  this->CallCopyOutputRegionToInputRegion(inputRegionForThread, outputRegionForThread);
391 
392  OutputImagePointerType outputLowPass = this->GetOutput(0);
393  OutputImagePointerType outputHighPass = this->GetOutput(1 << dir);
394 
395  // On multidimensional case, if m_SubsampleImageFactor != 1, we need internal images of different size
396  if (dir != 0 && m_SubsampleImageFactor > 1)
397  {
398  outputLowPass = m_InternalImages[dir - 1][0];
399  outputHighPass = m_InternalImages[dir - 1][1];
400  }
401 
402  // typedef for the iterations over the input image
403  typedef itk::ConstNeighborhoodIterator<InputImageType> NeighborhoodIteratorType;
404  typedef itk::NeighborhoodInnerProduct<InputImageType> InnerProductType;
405  typedef itk::ImageRegionIterator<OutputImageType> IteratorType;
406 
407  // Prepare the subsampling image factor, if any.
408  typedef SubsampledImageRegionConstIterator<InputImageType> SubsampleIteratorType;
409  typename SubsampleIteratorType::IndexType subsampling;
410  subsampling.Fill(1);
411  subsampling[dir] = GetSubsampleImageFactor();
412 
413  // Inner product
414  InnerProductType innerProduct;
415 
416  // High pass part calculation
417  HighPassOperatorType highPassOperator;
418  highPassOperator.SetDirection(dir);
419  highPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
420  highPassOperator.CreateDirectional();
421 
422  SubsampleIteratorType subItHighPass(input, inputRegionForThread);
423  subItHighPass.SetSubsampleFactor(subsampling);
424  subItHighPass.GoToBegin();
425 
426  NeighborhoodIteratorType itHighPass(highPassOperator.GetRadius(), input, inputRegionForThread);
428  itHighPass.OverrideBoundaryCondition(&boundaryCondition);
429 
430  IteratorType outHighPass(outputHighPass, subItHighPass.GenerateOutputInformation());
431  outHighPass.GoToBegin();
432 
433  while (!subItHighPass.IsAtEnd() && !outHighPass.IsAtEnd())
434  {
435  itHighPass.SetLocation(subItHighPass.GetIndex());
436  outHighPass.Set(innerProduct(itHighPass, highPassOperator));
437 
438  ++subItHighPass;
439  ++outHighPass;
440 
441  reporter.CompletedPixel();
442  }
443 
444  // Low pass part calculation
445  LowPassOperatorType lowPassOperator;
446  lowPassOperator.SetDirection(dir);
447  lowPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
448  lowPassOperator.CreateDirectional();
449 
450  SubsampleIteratorType subItLowPass(input, inputRegionForThread);
451  subItLowPass.SetSubsampleFactor(subsampling);
452  subItLowPass.GoToBegin();
453 
454  NeighborhoodIteratorType itLowPass(lowPassOperator.GetRadius(), input, inputRegionForThread);
455  itLowPass.OverrideBoundaryCondition(&boundaryCondition);
456 
457  IteratorType outLowPass(outputLowPass, subItLowPass.GenerateOutputInformation());
458  outLowPass.GoToBegin();
459 
460  while (!subItLowPass.IsAtEnd() && !outLowPass.IsAtEnd())
461  {
462  itLowPass.SetLocation(subItLowPass.GetIndex());
463  outLowPass.Set(innerProduct(itLowPass, lowPassOperator));
464 
465  ++subItLowPass;
466  ++outLowPass;
467 
468  reporter.CompletedPixel();
469  }
470 
471  if (dir > 0)
472  {
473  // Note that outputImageRegion correspond to the actual region of (local) input !
474  OutputImageRegionType outputImageRegion;
475  this->CallCopyInputRegionToOutputRegion(dir, outputImageRegion, inputRegionForThread);
476 
477  ThreadedGenerateDataAtDimensionN(0, dir - 1, reporter, outputImageRegion, threadId);
478  ThreadedGenerateDataAtDimensionN(1 << dir, dir - 1, reporter, outputImageRegion, threadId);
479  }
480 }
481 
482 template <class TInputImage, class TOutputImage, class TWaveletOperator>
483 void
486  (unsigned int idx, unsigned int direction, itk::ProgressReporter& reporter,
487  const OutputImageRegionType& outputRegionForThread, itk::ThreadIdType threadId)
488 {
489  // Note that outputRegionForThread correspond to the actual region of input !
490  OutputImagePointerType input = this->GetOutput(idx);
491  OutputImagePointerType outputHighPass = this->GetOutput(idx + (1 << direction));
492  OutputImagePointerType outputLowPass = OutputImageType::New();
493 
494  OutputImageRegionType outputImageRegion;
495  this->CallCopyInputRegionToOutputRegion(direction, outputImageRegion, outputRegionForThread);
496 
497  if (m_SubsampleImageFactor > 1)
498  {
499  input = m_InternalImages[direction][idx / (1 << (direction + 1))];
500 
501  if (direction != 0)
502  {
503  outputLowPass = m_InternalImages[direction - 1][idx / (1 << direction)];
504  outputHighPass = m_InternalImages[direction - 1][idx / (1 << direction) + 1];
505  }
506  }
507 
508  if (direction == 0)
509  {
510  // The output image has to be allocated
511  // May be not valid on multithreaded process ???
512  outputLowPass->SetRegions(outputImageRegion);
513  outputLowPass->Allocate(); // FIXME Check this line...
514  outputLowPass->FillBuffer(0);
515  }
516 
518  typedef itk::ConstNeighborhoodIterator<OutputImageType> NeighborhoodIteratorType;
519  typedef itk::NeighborhoodInnerProduct<OutputImageType> InnerProductType;
520  typedef itk::ImageRegionIterator<OutputImageType> IteratorType;
521 
522  // Prepare the subsampling image factor, if any.
523  typedef SubsampledImageRegionConstIterator<InputImageType> SubsampleIteratorType;
524  typename SubsampleIteratorType::IndexType subsampling;
525  subsampling.Fill(1);
526  subsampling[direction] = GetSubsampleImageFactor();
527 
528  // Inner products
529  InnerProductType innerProduct;
530 
531  // High pass part calculation
532  HighPassOperatorType highPassOperator;
533  highPassOperator.SetDirection(direction);
534  highPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
535  highPassOperator.CreateDirectional();
536 
537  SubsampleIteratorType subItHighPass(input, outputRegionForThread);
538  subItHighPass.SetSubsampleFactor(subsampling);
539  subItHighPass.GoToBegin();
540 
541  NeighborhoodIteratorType itHighPass(highPassOperator.GetRadius(), input, outputRegionForThread);
543  itHighPass.OverrideBoundaryCondition(&boundaryCondition);
544 
545  IteratorType outHighPass(outputHighPass, subItHighPass.GenerateOutputInformation());
546  outHighPass.GoToBegin();
547 
548  while (!subItHighPass.IsAtEnd() && !outHighPass.IsAtEnd())
549  {
550  itHighPass.SetLocation(subItHighPass.GetIndex());
551  outHighPass.Set(innerProduct(itHighPass, highPassOperator));
552 
553  ++subItHighPass;
554  ++outHighPass;
555 
556  reporter.CompletedPixel();
557  }
558 
559  // Low pass part calculation
560  LowPassOperatorType lowPassOperator;
561  lowPassOperator.SetDirection(direction);
562  lowPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
563  lowPassOperator.CreateDirectional();
564 
565  SubsampleIteratorType subItLowPass(input, outputRegionForThread);
566  subItLowPass.SetSubsampleFactor(subsampling);
567  subItLowPass.GoToBegin();
568 
569  NeighborhoodIteratorType itLowPass(lowPassOperator.GetRadius(), input, outputRegionForThread);
570  itLowPass.OverrideBoundaryCondition(&boundaryCondition);
571 
572  IteratorType outLowPass(outputLowPass, subItLowPass.GenerateOutputInformation());
573  outLowPass.GoToBegin();
574 
575  while (!subItLowPass.IsAtEnd() && !outLowPass.IsAtEnd())
576  {
577  itLowPass.SetLocation(subItLowPass.GetIndex());
578  outLowPass.Set(innerProduct(itLowPass, lowPassOperator));
579 
580  ++subItLowPass;
581  ++outLowPass;
582 
583  reporter.CompletedPixel();
584  }
585 
586  // Swap input and lowPassOutput
587  itk::ImageRegionConstIterator<OutputImageType> inIt(outputLowPass, outputImageRegion);
588  IteratorType outIt((direction != 0 && m_SubsampleImageFactor > 1) ?
589  static_cast<OutputImageType*>(m_InternalImages[direction - 2][idx / (1 << (direction - 1))])
590  : this->GetOutput(idx),
591  outputImageRegion);
592 
593  for (inIt.GoToBegin(), outIt.GoToBegin(); !inIt.IsAtEnd(); ++inIt, ++outIt)
594  {
595  outIt.Set(inIt.Get());
596  }
597 
598  if (direction != 0)
599  {
600  ThreadedGenerateDataAtDimensionN(idx, direction - 1, reporter, outputImageRegion, threadId);
601  ThreadedGenerateDataAtDimensionN(idx + (1 << direction), direction - 1, reporter, outputImageRegion, threadId);
602  }
603 }
604 
609 template <class TInputImage, class TOutputImage, class TWaveletOperator>
612 {
613  this->SetNumberOfRequiredInputs(1 << InputImageDimension);
614 
615  m_UpSampleFilterFactor = 0;
616  m_SubsampleImageFactor = 1;
617 
618  // TODO: For now, we force the number threads to 1 because there is a bug with multithreading in INVERSE transform
619  // Resulting in discontinuities in the reconstructed images
620  this->SetNumberOfThreads(1);
621 }
622 
623 template <class TInputImage, class TOutputImage, class TWaveletOperator>
624 void
627 {
628  Superclass::GenerateOutputInformation();
629 
630  for (unsigned int i = 1; i < this->GetNumberOfInputs(); ++i)
631  {
632  for (unsigned int dim = 0; dim < InputImageDimension; dim++)
633  {
634  if (this->GetInput(0)->GetLargestPossibleRegion().GetSize()[dim]
635  != this->GetInput(i)->GetLargestPossibleRegion().GetSize()[dim])
636  {
637  throw itk::ExceptionObject(__FILE__, __LINE__,
638  "Input images are not of the same dimension", ITK_LOCATION);
639  }
640  }
641  }
642 
643 #if __myDebug__
644  otbGenericMsgDebugMacro(<< " up sampling output regions by a factor of " << GetSubsampleImageFactor());
645 
646  otbGenericMsgDebugMacro(<< "initial region "
647  << this->GetInput(0)->GetLargestPossibleRegion().GetSize()[0]
648  << "," << this->GetInput(0)->GetLargestPossibleRegion().GetSize()[1]);
649 #endif
650 
651  OutputImageRegionType newRegion;
652  this->CallCopyInputRegionToOutputRegion(newRegion, this->GetInput(0)->GetLargestPossibleRegion());
653  this->GetOutput()->SetRegions(newRegion);
654 
655 #if __myDebug__
656  otbGenericMsgDebugMacro(<< "new region output " << newRegion.GetSize()[0] << "," << newRegion.GetSize()[1]);
657 #endif
658 }
659 
660 template <class TInputImage, class TOutputImage, class TWaveletOperator>
661 void
664 throw (itk::InvalidRequestedRegionError)
665  {
666  Superclass::GenerateInputRequestedRegion();
667 
668  // Filter length calculation
669  LowPassOperatorType lowPassOperator;
670  lowPassOperator.SetDirection(0);
671  lowPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
672  lowPassOperator.CreateDirectional();
673 
674  unsigned int radius = lowPassOperator.GetRadius()[0];
675 
676  HighPassOperatorType highPassOperator;
677  highPassOperator.SetDirection(0);
678  highPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
679  highPassOperator.CreateDirectional();
680 
681  if (radius < highPassOperator.GetRadius()[0]) radius = highPassOperator.GetRadius()[0];
682 
683  // Get the requested regionand pad it
684  for (unsigned int idx = 0; idx < this->GetNumberOfInputs(); ++idx)
685  {
686  InputImagePointerType input = const_cast<InputImageType*>(this->GetInput(idx));
687  InputImageRegionType inputRegion = input->GetRequestedRegion();
688  inputRegion.PadByRadius(radius);
689 
690  if (inputRegion.Crop(input->GetLargestPossibleRegion()))
691  {
692  input->SetRequestedRegion(inputRegion);
693  }
694  else
695  {
696  input->SetRequestedRegion(inputRegion);
697  itk::InvalidRequestedRegionError err(__FILE__, __LINE__);
698  err.SetLocation(ITK_LOCATION);
699  err.SetDescription("Requested region is (at least partially) outside the largest possible region.");
700  err.SetDataObject(input);
701  throw err;
702  }
703  }
704  }
705 
706 template <class TInputImage, class TOutputImage, class TWaveletOperator>
707 void
710 {
711  if (InputImageDimension > 1)
712  {
713  // Internal images will be used only if m_SubsampleImageFactor != 1
714  m_InternalImages.resize(InputImageDimension - 1);
715  for (unsigned int i = 0; i < m_InternalImages.size(); ++i)
716  {
717  // the size is linked to the SubsampleImageFactor that is assume to be 2!!!
718  m_InternalImages[i].resize(1 << (i + 1));
719  }
720 
721  OutputImageRegionType intermediateRegion;
722  Superclass::CallCopyInputRegionToOutputRegion(intermediateRegion,
723  this->GetInput(0)->GetLargestPossibleRegion());
724 
725  AllocateInternalData(intermediateRegion);
726  }
727 }
728 
729 template <class TInputImage, class TOutputImage, class TWaveletOperator>
730 void
733  (const OutputImageRegionType& outputRegion)
734 {
735  OutputImageRegionType largerRegion;
736  OutputImageRegionType smallerRegion = outputRegion;
737 
738  for (unsigned int direction = 0; direction < InputImageDimension - 1; direction++)
739  {
740  this->CallCopyInputRegionToOutputRegion(direction,
741  largerRegion, smallerRegion);
742 
743  for (unsigned int i = 0; i < m_InternalImages[direction].size(); ++i)
744  {
745  m_InternalImages[direction][i] = OutputImageType::New();
746  m_InternalImages[direction][i]->SetRegions(largerRegion);
747  m_InternalImages[direction][i]->Allocate();
748  m_InternalImages[direction][i]->FillBuffer(0);
749  }
750 
751  smallerRegion = largerRegion;
752  }
753 }
754 
755 template <class TInputImage, class TOutputImage, class TWaveletOperator>
756 void
759 {
760  if (m_SubsampleImageFactor > 1 && InputImageDimension > 1)
761  {
762  m_InternalImages.clear();
763  }
764 }
765 
766 template <class TInputImage, class TOutputImage, class TWaveletOperator>
767 void
770  (InputImageRegionType& destRegion, const OutputImageRegionType& srcRegion)
771 {
772  Superclass::CallCopyOutputRegionToInputRegion(destRegion, srcRegion);
773 
774  if (GetSubsampleImageFactor() > 1)
775  {
776  OutputIndexType srcIndex = srcRegion.GetIndex();
777  OutputSizeType srcSize = srcRegion.GetSize();
778 
779  InputIndexType destIndex;
780  InputSizeType destSize;
781 
782  for (unsigned int i = 0; i < InputImageDimension; ++i)
783  {
784  // TODO: This seems not right in odd index cases
785  destIndex[i] = srcIndex[i] / GetSubsampleImageFactor();
786  destSize[i] = srcSize[i] / GetSubsampleImageFactor();
787  }
788 
789  destRegion.SetIndex(destIndex);
790  destRegion.SetSize(destSize);
791 
792 #if 1
793  // Region Padding
794  LowPassOperatorType lowPassOperator;
795  lowPassOperator.SetDirection(0);
796  lowPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
797  lowPassOperator.CreateDirectional();
798 
799  typename InputImageRegionType::SizeType radius;
800  radius[0] = lowPassOperator.GetRadius()[0];
801 
802  HighPassOperatorType highPassOperator;
803  highPassOperator.SetDirection(0);
804  highPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
805  highPassOperator.CreateDirectional();
806 
807  if (radius[0] < highPassOperator.GetRadius()[0]) radius[0] = highPassOperator.GetRadius()[0];
808 
809  for (unsigned int i = 1; i < InputImageDimension; ++i)
810  radius[i] = 0;
811 
812 // for ( unsigned int i = 0; i < InputImageDimension; ++i )
813 // {
814 // radius[i] = lowPassOperator.GetRadius()[i];
815 // if ( radius[i] < highPassOperator.GetRadius()[i] )
816 // radius[i] = highPassOperator.GetRadius()[i];
817 // }
818 
819  InputImageRegionType paddedRegion = destRegion;
820  paddedRegion.PadByRadius(radius);
821 
822  if (paddedRegion.Crop(this->GetInput(0)->GetLargestPossibleRegion()))
823  {
824  destRegion = paddedRegion;
825  }
826 #endif
827  }
828 }
829 
830 template <class TInputImage, class TOutputImage, class TWaveletOperator>
831 void
834  (OutputImageRegionType& destRegion, const InputImageRegionType& srcRegion)
835 {
836  Superclass::CallCopyInputRegionToOutputRegion(destRegion, srcRegion);
837 
838  if (GetSubsampleImageFactor() > 1)
839  {
840  OutputIndexType srcIndex = srcRegion.GetIndex();
841  OutputSizeType srcSize = srcRegion.GetSize();
842 
843  InputIndexType destIndex;
844  InputSizeType destSize;
845 
846  for (unsigned int i = 0; i < InputImageDimension; ++i)
847  {
848  destIndex[i] = srcIndex[i] * GetSubsampleImageFactor();
849  destSize[i] = srcSize[i] * GetSubsampleImageFactor();
850  }
851 
852  destRegion.SetIndex(destIndex);
853  destRegion.SetSize(destSize);
854  }
855 }
856 
857 template <class TInputImage, class TOutputImage, class TWaveletOperator>
858 void
861  (unsigned int direction,
862  InputImageRegionType& destRegion, const OutputImageRegionType& srcRegion)
863 {
864  Superclass::CallCopyOutputRegionToInputRegion(destRegion, srcRegion);
865 
866  if (GetSubsampleImageFactor() > 1)
867  {
868  OutputIndexType srcIndex = srcRegion.GetIndex();
869  OutputSizeType srcSize = srcRegion.GetSize();
870 
871  InputIndexType destIndex;
872  InputSizeType destSize;
873 
874  for (unsigned int i = 0; i < InputImageDimension; ++i)
875  {
876  if (i == direction)
877  {
878  // TODO: This seems not right in odd index cases
879  destIndex[i] = srcIndex[i] / GetSubsampleImageFactor();
880  destSize[i] = srcSize[i] / GetSubsampleImageFactor();
881  }
882  else
883  {
884  destIndex[i] = srcIndex[i];
885  destSize[i] = srcSize[i];
886  }
887  }
888 
889  destRegion.SetIndex(destIndex);
890  destRegion.SetSize(destSize);
891  }
892 }
893 
894 template <class TInputImage, class TOutputImage, class TWaveletOperator>
895 void
898  (unsigned int direction,
899  OutputImageRegionType& destRegion, const InputImageRegionType& srcRegion)
900 {
901  Superclass::CallCopyInputRegionToOutputRegion(destRegion, srcRegion);
902 
903  if (GetSubsampleImageFactor() > 1)
904  {
905  typename InputImageRegionType::IndexType srcIndex = srcRegion.GetIndex();
906  typename InputImageRegionType::SizeType srcSize = srcRegion.GetSize();
907 
908  typename OutputImageRegionType::IndexType destIndex;
909  typename OutputImageRegionType::SizeType destSize;
910 
911  for (unsigned int i = 0; i < InputImageDimension; ++i)
912  {
913  if (i == direction)
914  {
915  destIndex[i] = srcIndex[i] * GetSubsampleImageFactor();
916  destSize[i] = srcSize[i] * GetSubsampleImageFactor();
917  }
918  else
919  {
920  destIndex[i] = srcIndex[i];
921  destSize[i] = srcSize[i];
922  }
923  }
924 
925  destRegion.SetIndex(destIndex);
926  destRegion.SetSize(destSize);
927  }
928 }
929 
930 template <class TInputImage, class TOutputImage, class TWaveletOperator>
931 void
934  (const OutputImageRegionType& outputRegionForThread, itk::ThreadIdType threadId)
935 {
936  itk::ProgressReporter reporter(this, threadId,
937  outputRegionForThread.GetNumberOfPixels() * this->GetNumberOfInputs());
938 
939  InputImageRegionType inputRegionForThread;
940  this->CallCopyOutputRegionToInputRegion(inputRegionForThread, outputRegionForThread);
941 
942  unsigned int dir = 0;
943 
944  // Low pass part calculation
945  LowPassOperatorType lowPassOperator;
946  lowPassOperator.SetDirection(dir);
947  lowPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
948  lowPassOperator.CreateDirectional();
949 
950  // High pass part calculation
951  HighPassOperatorType highPassOperator;
952  highPassOperator.SetDirection(dir);
953  highPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
954  highPassOperator.CreateDirectional();
955 
956  // typedef for the iterations over the input image
957  typedef itk::ConstNeighborhoodIterator<OutputImageType> NeighborhoodIteratorType;
958  typedef itk::NeighborhoodInnerProduct<OutputImageType> InnerProductType;
959 
960  // Faces iterations
961  typename NeighborhoodIteratorType::RadiusType radiusMax;
962  for (unsigned int idx = 0; idx < OutputImageDimension; ++idx)
963  {
964  radiusMax[idx] = lowPassOperator.GetRadius(idx);
965  if (radiusMax[idx] < highPassOperator.GetRadius(idx)) radiusMax[idx] = highPassOperator.GetRadius(idx);
966  }
967 
968  // The multiresolution case requires a SubsampleImageFilter step
969  if (m_SubsampleImageFactor > 1)
970  {
971  for (unsigned int i = 0; i < this->GetNumberOfInputs(); i += 2)
972  {
973  InputImagePointerType imgLowPass = const_cast<InputImageType*>(this->GetInput(i));
974  InputImagePointerType imgHighPass = const_cast<InputImageType*>(this->GetInput(i + 1));
975 
976  OutputImagePointerType outputImage = this->GetOutput();
977  if (dir != InputImageDimension - 1)
978  {
979  outputImage = m_InternalImages[0][i / 2];
980  }
981 
983  typename FilterType::InputImageIndexType delta;
984  delta.Fill(1);
985  delta[dir] = this->GetSubsampleImageFactor();
986 
987  InputImagePointerType cropedLowPass = InputImageType::New();
988  cropedLowPass->SetRegions(inputRegionForThread);
989  cropedLowPass->Allocate();
990  cropedLowPass->FillBuffer(0.);
991  itk::ImageRegionIterator<InputImageType> cropedLowPassIt(cropedLowPass, inputRegionForThread);
992  itk::ImageRegionIterator<InputImageType> imgLowPassIt(imgLowPass, inputRegionForThread);
993  for (cropedLowPassIt.GoToBegin(), imgLowPassIt.GoToBegin();
994  !cropedLowPassIt.IsAtEnd() && !imgLowPassIt.IsAtEnd();
995  ++cropedLowPassIt, ++imgLowPassIt)
996  {
997  cropedLowPassIt.Set(imgLowPassIt.Get());
998  }
999 
1000  typename FilterType::Pointer overSampledLowPass = FilterType::New();
1001  overSampledLowPass->SetInput(cropedLowPass);
1002  overSampledLowPass->SetSubsampleFactor(delta);
1003  overSampledLowPass->Update();
1004 
1005  InputImagePointerType cropedHighPass = InputImageType::New();
1006  cropedHighPass->SetRegions(inputRegionForThread);
1007  cropedHighPass->Allocate();
1008  cropedHighPass->FillBuffer(0.);
1009  itk::ImageRegionIterator<InputImageType> cropedHighPassIt(cropedHighPass, inputRegionForThread);
1010  itk::ImageRegionIterator<InputImageType> imgHighPassIt(imgHighPass, inputRegionForThread);
1011  for (cropedHighPassIt.GoToBegin(), imgHighPassIt.GoToBegin();
1012  !cropedHighPassIt.IsAtEnd() && !imgHighPassIt.IsAtEnd();
1013  ++cropedHighPassIt, ++imgHighPassIt)
1014  {
1015  cropedHighPassIt.Set(imgHighPassIt.Get());
1016  }
1017 
1018  typename FilterType::Pointer overSampledHighPass = FilterType::New();
1019  overSampledHighPass->SetInput(cropedHighPass);
1020  overSampledHighPass->SetSubsampleFactor(delta);
1021  overSampledHighPass->Update();
1022 
1023  InnerProductType innerProduct;
1024 
1026  (outputImage, overSampledLowPass->GetOutput()->GetRequestedRegion());
1027 
1028  NeighborhoodIteratorType lowIter(lowPassOperator.GetRadius(),
1029  overSampledLowPass->GetOutput(),
1030  overSampledLowPass->GetOutput()->GetRequestedRegion());
1032  lowIter.OverrideBoundaryCondition(&boundaryCondition);
1033 
1034  NeighborhoodIteratorType highIter(highPassOperator.GetRadius(),
1035  overSampledHighPass->GetOutput(),
1036  overSampledHighPass->GetOutput()->GetRequestedRegion());
1037  highIter.OverrideBoundaryCondition(&boundaryCondition);
1038 
1039  lowIter.GoToBegin();
1040  highIter.GoToBegin();
1041  out.GoToBegin();
1042 
1043  while (!out.IsAtEnd())
1044  {
1045  out.Set(innerProduct(lowIter, lowPassOperator)
1046  + innerProduct(highIter, highPassOperator));
1047 
1048  ++lowIter;
1049  ++highIter;
1050  ++out;
1051 
1052  reporter.CompletedPixel();
1053  }
1054  } // end for each overSampledLowPass/overSampledhighPass pair of entry
1055  }
1056  else // multiscale case
1057  {
1058  for (unsigned int i = 0; i < this->GetNumberOfInputs(); i += 2)
1059  {
1060  InputImagePointerType imgLowPass = const_cast<InputImageType*>(this->GetInput(i));
1061  InputImagePointerType imgHighPass = const_cast<InputImageType*>(this->GetInput(i + 1));
1062 
1063  OutputImagePointerType outputImage = this->GetOutput();
1064  if (dir != InputImageDimension - 1)
1065  {
1066  outputImage = m_InternalImages[0][i / 2];
1067  }
1068 
1069  InnerProductType innerProduct;
1070 
1071  itk::ImageRegionIterator<OutputImageType> out(outputImage, imgLowPass->GetRequestedRegion());
1072 
1073  NeighborhoodIteratorType lowIter(lowPassOperator.GetRadius(), imgLowPass, imgLowPass->GetRequestedRegion());
1075  lowIter.OverrideBoundaryCondition(&boundaryCondition);
1076 
1077  NeighborhoodIteratorType highIter(highPassOperator.GetRadius(), imgHighPass, imgLowPass->GetRequestedRegion());
1078  highIter.OverrideBoundaryCondition(&boundaryCondition);
1079 
1080  lowIter.GoToBegin();
1081  highIter.GoToBegin();
1082  out.GoToBegin();
1083 
1084  while (!out.IsAtEnd())
1085  {
1086  out.Set((innerProduct(lowIter, lowPassOperator)
1087  + innerProduct(highIter, highPassOperator)) / 2.);
1088 
1089  ++lowIter;
1090  ++highIter;
1091  ++out;
1092 
1093  reporter.CompletedPixel();
1094  }
1095  } // end for each imgLowPass/imghighPass pair of entry
1096  } // end multiscale case
1097 
1098  if (dir != InputImageDimension - 1)
1099  {
1100  // Note that outputImageRegion correspond to the actual region of (local) input !
1101  OutputImageRegionType outputImageRegion;
1102  this->CallCopyInputRegionToOutputRegion(dir, outputImageRegion, inputRegionForThread);
1103 
1104  ThreadedGenerateDataAtDimensionN(dir + 1, reporter, outputImageRegion, threadId);
1105  }
1106 
1107 }
1108 
1109 template <class TInputImage, class TOutputImage, class TWaveletOperator>
1110 void
1113  (unsigned int direction,
1114  itk::ProgressReporter& reporter,
1115  const OutputImageRegionType& outputRegionForThread, itk::ThreadIdType threadId)
1116 {
1117  OutputImageRegionType outputImageRegion;
1118  this->CallCopyInputRegionToOutputRegion(direction, outputImageRegion, outputRegionForThread);
1119 
1120  // Low pass part calculation
1121  LowPassOperatorType lowPassOperator;
1122  lowPassOperator.SetDirection(direction);
1123  lowPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
1124  lowPassOperator.CreateDirectional();
1125 
1126  // High pass part calculation
1127  HighPassOperatorType highPassOperator;
1128  highPassOperator.SetDirection(direction);
1129  highPassOperator.SetUpSampleFactor(this->GetUpSampleFilterFactor());
1130  highPassOperator.CreateDirectional();
1131 
1132  // typedef for the iterations over the input image
1133  typedef itk::ConstNeighborhoodIterator<OutputImageType> NeighborhoodIteratorType;
1134  typedef itk::NeighborhoodInnerProduct<OutputImageType> InnerProductType;
1135  // Faces iterations
1136  typename NeighborhoodIteratorType::RadiusType radiusMax;
1137  for (unsigned int i = 0; i < InputImageDimension; ++i)
1138  {
1139  radiusMax[i] = lowPassOperator.GetRadius(i);
1140  if (radiusMax[i] < highPassOperator.GetRadius(i)) radiusMax[i] = highPassOperator.GetRadius(i);
1141  }
1142 
1143  // The multiresolution case requires a SubsampleImageFilter step
1144  if (m_SubsampleImageFactor > 1)
1145  {
1146  for (unsigned int i = 0; i < m_InternalImages[direction - 1].size(); i += 2)
1147  {
1148  InputImagePointerType imgLowPass = m_InternalImages[direction - 1][i];
1149  InputImagePointerType imgHighPass = m_InternalImages[direction - 1][i + 1];
1150 
1151  OutputImagePointerType outputImage = this->GetOutput();
1152  if (direction < InputImageDimension - 1)
1153  {
1154  outputImage = m_InternalImages[direction][i / 2];
1155  }
1156 
1158  typename FilterType::InputImageIndexType delta;
1159  delta.Fill(1);
1160  delta[direction] = this->GetSubsampleImageFactor();
1161 
1162  InputImagePointerType cropedLowPass = InputImageType::New();
1163  cropedLowPass->SetRegions(outputRegionForThread);
1164  cropedLowPass->Allocate();
1165  cropedLowPass->FillBuffer(0.);
1166  itk::ImageRegionIterator<InputImageType> cropedLowPassIt(cropedLowPass, outputRegionForThread);
1167  itk::ImageRegionIterator<InputImageType> imgLowPassIt(imgLowPass, outputRegionForThread);
1168  for (cropedLowPassIt.GoToBegin(), imgLowPassIt.GoToBegin();
1169  !cropedLowPassIt.IsAtEnd() && !imgLowPassIt.IsAtEnd();
1170  ++cropedLowPassIt, ++imgLowPassIt)
1171  {
1172  cropedLowPassIt.Set(imgLowPassIt.Get());
1173  }
1174 
1175  typename FilterType::Pointer overSampledLowPass = FilterType::New();
1176  overSampledLowPass->SetInput(cropedLowPass);
1177  overSampledLowPass->SetSubsampleFactor(delta);
1178  overSampledLowPass->Update();
1179 
1180  InputImagePointerType cropedHighPass = InputImageType::New();
1181  cropedHighPass->SetRegions(outputRegionForThread);
1182  cropedHighPass->Allocate();
1183  cropedHighPass->FillBuffer(0.);
1184  itk::ImageRegionIterator<InputImageType> cropedHighPassIt(cropedHighPass, outputRegionForThread);
1185  itk::ImageRegionIterator<InputImageType> imgHighPassIt(imgHighPass, outputRegionForThread);
1186  for (cropedHighPassIt.GoToBegin(), imgHighPassIt.GoToBegin();
1187  !cropedHighPassIt.IsAtEnd() && !imgHighPassIt.IsAtEnd();
1188  ++cropedHighPassIt, ++imgHighPassIt)
1189  {
1190  cropedHighPassIt.Set(imgHighPassIt.Get());
1191  }
1192 
1193  typename FilterType::Pointer overSampledHighPass = FilterType::New();
1194  overSampledHighPass->SetInput(cropedHighPass);
1195  overSampledHighPass->SetSubsampleFactor(delta);
1196  overSampledHighPass->Update();
1197 
1198  InnerProductType innerProduct;
1199 
1201  overSampledLowPass->GetOutput()->GetRequestedRegion());
1202 
1203  // TODO: This might be the cause of the multithreading bug : we use a neighborhood iterator on cropped data
1204  // Are we sure that we have cropped enough data to access the neighborhood ?
1205  NeighborhoodIteratorType lowIter(lowPassOperator.GetRadius(),
1206  overSampledLowPass->GetOutput(),
1207  overSampledLowPass->GetOutput()->GetRequestedRegion());
1209  lowIter.OverrideBoundaryCondition(&boundaryCondition);
1210 
1211  NeighborhoodIteratorType highIter(highPassOperator.GetRadius(),
1212  overSampledHighPass->GetOutput(),
1213  overSampledHighPass->GetOutput()->GetRequestedRegion());
1214  highIter.OverrideBoundaryCondition(&boundaryCondition);
1215 
1216  lowIter.GoToBegin();
1217  highIter.GoToBegin();
1218  out.GoToBegin();
1219 
1220  while (!out.IsAtEnd())
1221  {
1222  out.Set(innerProduct(lowIter, lowPassOperator)
1223  + innerProduct(highIter, highPassOperator));
1224 
1225  ++lowIter;
1226  ++highIter;
1227  ++out;
1228 
1229  reporter.CompletedPixel();
1230  }
1231  } // end for each overSampledLowPass/overSampledhighPass pair of entry
1232  }
1233  else // multiscale case
1234  {
1235  for (unsigned int i = 0; i < m_InternalImages[direction - 1].size(); i += 2)
1236  {
1237  InputImagePointerType imgLowPass = m_InternalImages[direction - 1][i];
1238  InputImagePointerType imgHighPass = m_InternalImages[direction - 1][i + 1];
1239 
1240  OutputImagePointerType outputImage = this->GetOutput();
1241  if (direction < InputImageDimension - 1)
1242  {
1243  outputImage = m_InternalImages[direction][i / 2];
1244  }
1245 
1246  InnerProductType innerProduct;
1247 
1248  itk::ImageRegionIterator<OutputImageType> out(outputImage, imgLowPass->GetRequestedRegion());
1249 
1250  NeighborhoodIteratorType lowIter(lowPassOperator.GetRadius(), imgLowPass, imgLowPass->GetRequestedRegion());
1252  lowIter.OverrideBoundaryCondition(&boundaryCondition);
1253 
1254  NeighborhoodIteratorType highIter(highPassOperator.GetRadius(), imgHighPass, imgLowPass->GetRequestedRegion());
1255  highIter.OverrideBoundaryCondition(&boundaryCondition);
1256 
1257  lowIter.GoToBegin();
1258  highIter.GoToBegin();
1259  out.GoToBegin();
1260 
1261  while (!out.IsAtEnd())
1262  {
1263  out.Set((innerProduct(lowIter, lowPassOperator)
1264  + innerProduct(highIter, highPassOperator)) / 2.);
1265 
1266  ++lowIter;
1267  ++highIter;
1268  ++out;
1269 
1270  reporter.CompletedPixel();
1271  }
1272  } // end for each imgLowPass/imghighPass pair of entry
1273  }
1274 
1275  if (direction < InputImageDimension - 1)
1276  {
1277  ThreadedGenerateDataAtDimensionN(direction + 1, reporter, outputImageRegion, threadId);
1278  }
1279 }
1280 
1281 } // end of namespace
1282 
1283 #endif
void SetDataObject(DataObject *dobj)
virtual void SetDescription(const std::string &s)
virtual void SetLocation(const std::string &s)
void Set(const PixelType &value) const
virtual void ThreadedGenerateData(const OutputImageRegionType &outputRegionForThread, ThreadIdType threadId)
Performs a down sampling of an image.
virtual void CallCopyOutputRegionToInputRegion(InputImageRegionType &destRegion, const OutputImageRegionType &srcRegion)
virtual void BeforeThreadedGenerateData()
#define otbGenericMsgDebugMacro(x)
Definition: otbMacro.h:67
Regular subsample iterator over an image.
InputImageType::RegionType InputImageRegionType
TInputImage InputImageType
virtual void GenerateOutputInformation()
virtual void AfterThreadedGenerateData()
virtual void GenerateInputRequestedRegion()
OutputImageType::RegionType OutputImageRegionType
virtual void CallCopyInputRegionToOutputRegion(OutputImageRegionType &destRegion, const InputImageRegionType &srcRegion)
bool IsAtEnd(void) const
PixelType Get(void) const
One level stationary wavelet transform.
unsigned int ThreadIdType