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