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