OTB  9.0.0
Orfeo Toolbox
otbNRIBandImagesToOneNComplexBandsImage.hxx
Go to the documentation of this file.
1 /*
2  * Copyright (C) 2005-2022 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 otbNRIBandImagesToOneNComplexBandsImage_hxx
22 #define otbNRIBandImagesToOneNComplexBandsImage_hxx
23 
25 
26 #include "itkImageRegionIterator.h"
27 #include "itkImageRegionConstIterator.h"
28 #include "itkProgressReporter.h"
29 #include "itkVariableLengthVector.h"
30 
31 
32 namespace otb
33 {
34 
38 template <class TInputImage, class TOutputImage>
40 {
41  // this->SetNumberOfThreads(1);
42 }
43 
47 template <class TInputImage, class TOutputImage>
49 {
50  Superclass::GenerateOutputInformation();
51 
52  unsigned int nbInputs = this->GetNumberOfInputs();
53  this->GetOutput()->SetNumberOfComponentsPerPixel(nbInputs);
54 }
55 
59 template <class TInputImage, class TOutputImage>
61 {
62  unsigned int nbInputs = this->GetNumberOfInputs();
63 
64  for (unsigned int i = 0; i < nbInputs; i++)
65  if (this->GetInput(i)->GetNumberOfComponentsPerPixel() != 2)
66  itkExceptionMacro("Input images must be made of two bands and only two (see input #" << i << ").");
67 }
68 
72 template <class TInputImage, class TOutputImage>
74  itk::ThreadIdType threadId)
75 {
76 
77  // support progress methods/callbacks
78  itk::ProgressReporter progress(this, threadId, outputRegionForThread.GetNumberOfPixels());
79 
80  unsigned int nbInputs = this->GetNumberOfInputs();
81 
82  itk::VariableLengthVector<std::complex<typename InputPixelType::ValueType>> vlv(nbInputs);
83 
84  std::vector<typename itk::ImageRegionConstIterator<TInputImage>> vInIt;
85  for (unsigned int i = 0; i < nbInputs; i++)
86  vInIt.push_back(itk::ImageRegionConstIterator<TInputImage>(this->GetInput(i), outputRegionForThread));
87 
88 
89  itk::ImageRegionIterator<OutputImageType> outIt;
90  outIt = itk::ImageRegionIterator<OutputImageType>(this->GetOutput(), outputRegionForThread);
91 
92  for (unsigned int i = 0; i < nbInputs; i++)
93  vInIt[i].GoToBegin();
94  outIt.GoToBegin();
95 
96  while (!outIt.IsAtEnd())
97  {
98 
99 
100  for (unsigned int i = 0; i < nbInputs; i++)
101  {
102  vlv[i] = std::complex<typename InputPixelType::ValueType>(vInIt[i].Get()[0], vInIt[i].Get()[1]);
103  // std::cout << "i = " << i << " " << vInIt[i].Get()[0] << " " << vInIt[i].Get()[1] << std::endl;
104  }
105 
106 
107  outIt.Set(vlv);
108 
109  // std::cout << "outIt.Get() = " << outIt.Get() << std::endl;
110 
111  for (unsigned int i = 0; i < nbInputs; i++)
112  ++vInIt[i];
113  ++outIt;
114 
115  progress.CompletedPixel();
116  }
117 }
118 
122 template <class TInputImage, class TOutput>
123 void NRIBandImagesToOneNComplexBandsImage<TInputImage, TOutput>::PrintSelf(std::ostream& os, itk::Indent indent) const
124 {
125  Superclass::PrintSelf(os, indent);
126 }
127 
128 } // end namespace otb
129 
130 #endif
otb::NRIBandImagesToOneNComplexBandsImage::ThreadedGenerateData
void ThreadedGenerateData(const OutputImageRegionType &outputRegionForThread, itk::ThreadIdType threadId) override
Definition: otbNRIBandImagesToOneNComplexBandsImage.hxx:73
otb
The "otb" namespace contains all Orfeo Toolbox (OTB) classes.
Definition: otbJoinContainer.h:32
otb::NRIBandImagesToOneNComplexBandsImage::NRIBandImagesToOneNComplexBandsImage
NRIBandImagesToOneNComplexBandsImage()
Definition: otbNRIBandImagesToOneNComplexBandsImage.hxx:39
otb::NRIBandImagesToOneNComplexBandsImage::BeforeThreadedGenerateData
void BeforeThreadedGenerateData(void) override
Definition: otbNRIBandImagesToOneNComplexBandsImage.hxx:60
otb::NRIBandImagesToOneNComplexBandsImage::OutputImageRegionType
OutputImageType::RegionType OutputImageRegionType
Definition: otbNRIBandImagesToOneNComplexBandsImage.h:67
otb::NRIBandImagesToOneNComplexBandsImage::GenerateOutputInformation
void GenerateOutputInformation(void) override
Definition: otbNRIBandImagesToOneNComplexBandsImage.hxx:48
otbNRIBandImagesToOneNComplexBandsImage.h
otb::NRIBandImagesToOneNComplexBandsImage::PrintSelf
void PrintSelf(std::ostream &os, itk::Indent indent) const override
Definition: otbNRIBandImagesToOneNComplexBandsImage.hxx:123