OTB  9.0.0
Orfeo Toolbox
otbMultiChannelExtractROI.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 otbMultiChannelExtractROI_hxx
22 #define otbMultiChannelExtractROI_hxx
23 
25 
26 #include "itkImageRegionIterator.h"
27 #include "itkObjectFactory.h"
28 #include "itkExtractImageFilterRegionCopier.h"
29 #include "itkProgressReporter.h"
30 
31 namespace otb
32 {
33 
37 template <class TInputPixelType, class TOutputPixelType>
39  : ExtractROIBase<VectorImage<TInputPixelType, 2>, VectorImage<TOutputPixelType, 2>>(), m_FirstChannel(0), m_LastChannel(0), m_ChannelsKind(0)
40 {
41  ClearChannels();
42 }
43 
47 template <class TInputPixelType, class TOutputPixelType>
49 {
50  if (m_ChannelsKind == 1)
51  {
52  itkExceptionMacro(<< "m_Channels already set using channels interval.");
53  }
54  m_Channels.push_back(channel);
55  if (m_ChannelsKind == 0)
56  {
57  m_ChannelsKind = 2;
58  }
59  this->Modified();
60 }
62 
63 template <class TInputPixelType, class TOutputPixelType>
65 {
66  if (m_ChannelsKind == 2)
67  {
68  itkExceptionMacro(<< "m_Channels already set using SetChannels method.");
69  }
70  m_FirstChannel = id;
71  if (m_ChannelsKind == 0)
72  {
73  m_ChannelsKind = 1;
74  }
75  this->Modified();
76 }
77 
78 template <class TInputPixelType, class TOutputPixelType>
80 {
81  if (m_ChannelsKind == 2)
82  {
83  itkExceptionMacro(<< "m_Channels already set using SetChannels method.");
84  }
85  m_LastChannel = id;
86  if (m_ChannelsKind == 0)
87  {
88  m_ChannelsKind = 1;
89  }
90  this->Modified();
91 }
92 
96 template <class TInputPixelType, class TOutputPixelType>
98 {
99  m_FirstChannel = 0;
100  m_LastChannel = 0;
101  m_Channels.clear();
102  m_ChannelsKind = 0;
103  m_ChannelsWorks.clear();
104 }
106 
110 template <class TInputPixelType, class TOutputPixelType>
111 void MultiChannelExtractROI<TInputPixelType, TOutputPixelType>::PrintSelf(std::ostream& os, itk::Indent indent) const
112 {
113  Superclass::PrintSelf(os, indent);
114 }
115 
116 template <class TInputPixelType, class TOutputPixelType>
118 {
119  // The following conditions can be gathered but we'd loose in comprehension
120  m_ChannelsWorks.clear();
121  // First passage in the method:
122  if (m_Channels.empty() == true)
123  {
124  // - User SetFirst/LastChannel()
125  if (m_ChannelsKind == 1)
126  {
127  this->SetChannelsWorkWithLimits();
128  }
129  else
130  {
131  // - User called SetChannels()
132  if (m_Channels.empty() == true && m_ChannelsKind == 2)
133  {
134  m_ChannelsWorks = m_Channels;
135  }
136  }
137  }
138  // Second passage in the method: Update already donne
139  else
140  {
141  // - User SetFirst/LastChannel()
142  if (m_ChannelsKind == 1)
143  {
144  m_Channels.clear();
145  this->SetChannelsWorkWithLimits();
146  }
147  else
148  {
149  // - User called SetChannels()
150  if (m_ChannelsKind == 2)
151  {
152  m_ChannelsWorks = m_Channels;
153  }
154  }
155  }
156 }
157 
158 template <class TInputPixelType, class TOutputPixelType>
160 {
161  if ((m_FirstChannel == 0) || (m_LastChannel == 0))
162  {
163  itkExceptionMacro(<< "otb::ExtractImageFilter::GenerateOutputInformation "
164  << "Channels must reside into [1...] " << typeid(itk::ImageBase<InputImageDimension>*).name());
165  }
166  if (m_FirstChannel > m_LastChannel)
167  {
168  itkExceptionMacro(<< "otb::ExtractImageFilter::GenerateOutputInformation "
169  << "FirstChannel is greater than LastChannel" << typeid(itk::ImageBase<InputImageDimension>*).name());
170  }
171 
172  for (unsigned int channel = m_FirstChannel; channel <= m_LastChannel; channel++)
173  {
174  m_ChannelsWorks.push_back(channel);
175  }
176 
177  m_Channels = m_ChannelsWorks;
178 }
179 
189 template <class TInputPixelType, class TOutputPixelType>
191 {
192  // Call to the superclass implementation
193  Superclass::GenerateOutputInformation();
194  this->ChannelsReInitialization();
196 
197  typename Superclass::InputImageConstPointer inputPtr = this->GetInput();
198  typename Superclass::OutputImagePointer outputPtr = this->GetOutput();
199 
200  unsigned int nbComponentsPerPixel = inputPtr->GetNumberOfComponentsPerPixel();
201  if (m_ChannelsKind != 0)
202  {
203  // Test if the asked channels index exists in the input image
204  ChannelsType m_BadChannels;
205  m_BadChannels.clear();
206  for (unsigned int i = 0; i < m_ChannelsWorks.size(); ++i)
207  {
208  if ((m_ChannelsWorks[i] < 1) || (m_ChannelsWorks[i] > nbComponentsPerPixel))
209  {
210  bool isInsideBadChannels = false;
211  for (unsigned int j = 0; j < m_BadChannels.size(); ++j)
212  {
213  if (m_BadChannels[j] == m_ChannelsWorks[i])
214  isInsideBadChannels = true;
215  }
216  if (!isInsideBadChannels)
217  m_BadChannels.push_back(m_ChannelsWorks[i]);
218  }
219  }
220  if (m_BadChannels.empty() == false)
221  {
222  std::ostringstream oss;
223  oss << "otb::ExtractImageFilter::GenerateOutputInformation : ";
224  oss << "Channel(s) [ ";
225  for (unsigned int i = 0; i < m_BadChannels.size(); ++i)
226  {
227  oss << m_BadChannels[i] << " ";
228  }
229  oss << "] not authorized.";
230  oss << " Each channel index has to be in [1," << nbComponentsPerPixel << "].";
231  itkExceptionMacro(<< oss.str().c_str());
232  }
233  nbComponentsPerPixel = m_ChannelsWorks.size();
234  }
235 
236  // initialize the number of channels of the output image
237  outputPtr->SetNumberOfComponentsPerPixel(nbComponentsPerPixel);
238 }
239 
240 template <class TInputPixelType, class TOutputPixelType>
242  itk::ThreadIdType threadId)
243 {
244  itkDebugMacro(<< "Actually executing");
245  // Get the input and output pointers
246  typename Superclass::InputImageConstPointer inputPtr = this->GetInput();
247  typename Superclass::OutputImagePointer outputPtr = this->GetOutput();
248 
249  // support progress methods/callbacks
250  itk::ProgressReporter progress(this, threadId, outputRegionForThread.GetNumberOfPixels());
251 
252  // Define the portion of the input to walk for this thread
253  InputImageRegionType inputRegionForThread;
254  this->CallCopyOutputRegionToInputRegion(inputRegionForThread, outputRegionForThread);
255 
256  // Define the iterators.
257  typedef itk::ImageRegionIterator<OutputImageType> OutputIterator;
258  typedef itk::ImageRegionConstIterator<InputImageType> InputIterator;
259 
260  OutputIterator outIt(outputPtr, outputRegionForThread);
261  InputIterator inIt(inputPtr, inputRegionForThread);
262 
263  outIt.GoToBegin();
264  inIt.GoToBegin();
265 
266  // if default behaviour
267  if (m_ChannelsKind == 0)
268  {
269  // walk the output region, and sample the input image
270  while (!outIt.IsAtEnd())
271  {
272  outIt.Set(inIt.Get());
273  ++outIt;
274  ++inIt;
275  progress.CompletedPixel();
276  }
277  }
278  // Specific behaviour
279  else
280  {
281  // for each channel to process
282  unsigned int channelIn(0);
283  unsigned int channelOut(0);
284  unsigned int nbChannels(0);
285 
286  InputImagePixelType pixelInput;
287  while (!outIt.IsAtEnd())
288  {
289  OutputImagePixelType pixelOutput;
290  pixelOutput.Reserve(outputPtr->GetVectorLength());
291  pixelInput = inIt.Get();
292  channelOut = 0;
293  for (nbChannels = 0; nbChannels < m_ChannelsWorks.size(); ++nbChannels)
294  {
295  channelIn = m_ChannelsWorks[nbChannels] - 1;
296  pixelOutput[channelOut] = static_cast<OutputValueType>(pixelInput[channelIn]);
297  channelOut++;
298  }
299  outIt.Set(pixelOutput);
300  ++outIt;
301  ++inIt;
302  progress.CompletedPixel();
303  }
304  }
305 }
306 
307 } // end namespace otb
308 
309 #endif
otb::MultiChannelExtractROI::SetChannel
void SetChannel(unsigned int channel)
Definition: otbMultiChannelExtractROI.hxx:48
otb::MultiChannelExtractROI::SetLastChannel
void SetLastChannel(unsigned int id)
Definition: otbMultiChannelExtractROI.hxx:79
otb::MultiChannelExtractROI::InputImagePixelType
InputImageType::PixelType InputImagePixelType
Definition: otbMultiChannelExtractROI.h:75
otbMultiChannelExtractROI.h
otb::MultiChannelExtractROI::ClearChannels
void ClearChannels(void)
Definition: otbMultiChannelExtractROI.hxx:97
otb::MultiChannelExtractROI::PrintSelf
void PrintSelf(std::ostream &os, itk::Indent indent) const override
Definition: otbMultiChannelExtractROI.hxx:111
otb
The "otb" namespace contains all Orfeo Toolbox (OTB) classes.
Definition: otbJoinContainer.h:32
otb::MultiChannelExtractROI::ThreadedGenerateData
void ThreadedGenerateData(const OutputImageRegionType &outputRegionForThread, itk::ThreadIdType threadId) override
Definition: otbMultiChannelExtractROI.hxx:241
otb::MultiChannelExtractROI::OutputValueType
OutputImageType::InternalPixelType OutputValueType
Definition: otbMultiChannelExtractROI.h:67
otb::MultiChannelExtractROI::SetFirstChannel
void SetFirstChannel(unsigned int id)
Definition: otbMultiChannelExtractROI.hxx:64
otb::MultiChannelExtractROI::ChannelsType
std::vector< unsigned int > ChannelsType
Definition: otbMultiChannelExtractROI.h:86
otb::MultiChannelExtractROI::GenerateOutputInformation
void GenerateOutputInformation() override
Definition: otbMultiChannelExtractROI.hxx:190
otb::MultiChannelExtractROI::MultiChannelExtractROI
MultiChannelExtractROI()
Definition: otbMultiChannelExtractROI.hxx:38
otb::MultiChannelExtractROI::SetChannelsWorkWithLimits
void SetChannelsWorkWithLimits()
Definition: otbMultiChannelExtractROI.hxx:159
otb::ExtractROIBase
Base class to extract area of images.
Definition: otbExtractROIBase.h:48
otb::MultiChannelExtractROI::OutputImageRegionType
OutputImageType::RegionType OutputImageRegionType
Definition: otbMultiChannelExtractROI.h:70
otb::MultiChannelExtractROI::ChannelsReInitialization
void ChannelsReInitialization()
Definition: otbMultiChannelExtractROI.hxx:117
otb::MultiChannelExtractROI::InputImageRegionType
InputImageType::RegionType InputImageRegionType
Definition: otbMultiChannelExtractROI.h:71
otb::MultiChannelExtractROI::OutputImagePixelType
OutputImageType::PixelType OutputImagePixelType
Definition: otbMultiChannelExtractROI.h:74
otb::VectorImage
Creation of an "otb" vector image which contains metadata.
Definition: otbVectorImage.h:45