OTB  9.0.0
Orfeo Toolbox
otbStreamingSimpleMosaicFilter.hxx
Go to the documentation of this file.
1 /*
2  * Copyright (C) 1999-2011 Insight Software Consortium
3  * Copyright (C) 2005-2022 Centre National d'Etudes Spatiales (CNES)
4  * Copyright (C) 2016-2019 IRSTEA
5  *
6  * This file is part of Orfeo Toolbox
7  *
8  * https://www.orfeo-toolbox.org/
9  *
10  * Licensed under the Apache License, Version 2.0 (the "License");
11  * you may not use this file except in compliance with the License.
12  * You may obtain a copy of the License at
13  *
14  * http://www.apache.org/licenses/LICENSE-2.0
15  *
16  * Unless required by applicable law or agreed to in writing, software
17  * distributed under the License is distributed on an "AS IS" BASIS,
18  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
19  * See the License for the specific language governing permissions and
20  * limitations under the License.
21  */
22 #ifndef __StreamingSimpleMosaicFilter_hxx
23 #define __StreamingSimpleMosaicFilter_hxx
24 
26 
27 namespace otb
28 {
29 
33 template <class TInputImage, class TOutputImage, class TInternalValueType>
35  itk::ThreadIdType threadId)
36 {
37 
38  // Debug info
39  itkDebugMacro(<< "Actually executing thread " << threadId << " in region " << outputRegionForThread);
40 
41  // Support progress methods/callbacks
42  itk::ProgressReporter progress(this, threadId, outputRegionForThread.GetNumberOfPixels());
43 
44  // Get output pointer
45  OutputImageType* mosaicImage = this->GetOutput();
46 
47  // Get number of used inputs
48  const unsigned int nbOfUsedInputImages = Superclass::GetNumberOfUsedInputImages();
49 
50  // Get number of bands
51  const unsigned int nBands = Superclass::GetNumberOfBands();
52 
53  // Iterate through the thread region
54  IteratorType outputIt(mosaicImage, outputRegionForThread);
55 
56  // Prepare interpolated pixel
57  InternalPixelType interpolatedMathPixel;
58  interpolatedMathPixel.SetSize(nBands);
59 
60  // Prepare input pointers, interpolators, and valid regions (input images)
61  typename std::vector<InputImageType*> currentImage;
62  typename std::vector<InterpolatorPointerType> interp;
63  Superclass::PrepareImageAccessors(currentImage, interp);
64 
65  // Container for geo coordinates
66  OutputImagePointType geoPoint;
67 
68  for (outputIt.GoToBegin(); !outputIt.IsAtEnd(); ++outputIt)
69  {
70  // Prepare output pixel
71  OutputImagePixelType outputPixel(Superclass::GetNoDataOutputPixel());
72 
73  // Current pixel --> Geographical point
74  mosaicImage->TransformIndexToPhysicalPoint(outputIt.GetIndex(), geoPoint);
75 
76  // Loop on used input images
77  for (unsigned int i = 0; i < nbOfUsedInputImages; i++)
78  {
79  // Get the input image pointer
80  unsigned int imgIndex = Superclass::GetUsedInputImageIndice(i);
81 
82  // Check if the point is inside the transformed thread region
83  // (i.e. the region in the current input image which match the thread
84  // region)
85  if (interp[i]->IsInsideBuffer(geoPoint))
86  {
87 
88  // Compute the interpolated pixel value
89  InputImagePixelType interpolatedPixel = interp[i]->Evaluate(geoPoint);
90 
91  // Check that interpolated pixel is not empty
92  if (Superclass::IsPixelNotEmpty(interpolatedPixel))
93  {
94  // Update the output pixel
95  for (unsigned int band = 0; band < nBands; band++)
96  {
97  if (this->GetShiftScaleInputImages())
98  {
99  InternalValueType pixelValue = static_cast<InternalValueType>(interpolatedPixel[band]);
100  this->ShiftScaleValue(pixelValue, imgIndex, band);
101  outputPixel[band] = static_cast<OutputImageInternalPixelType>(pixelValue);
102  }
103  else
104  {
105  outputPixel[band] = static_cast<OutputImageInternalPixelType>(interpolatedPixel[band]);
106  }
107  }
108  } // Interpolated pixel is not empty
109  } // point inside buffer
110  } // next image
111 
112  // Update output pixel value
113  outputIt.Set(outputPixel);
114 
115  // Update progress
116  progress.CompletedPixel();
117 
118  } // next output pixel
119 }
120 
121 } // end namespace gtb
122 
123 #endif
otb::StreamingSimpleMosaicFilter::InternalValueType
Superclass::InternalValueType InternalValueType
Definition: otbStreamingSimpleMosaicFilter.h:74
otb::StreamingSimpleMosaicFilter::OutputImagePixelType
Superclass::OutputImagePixelType OutputImagePixelType
Definition: otbStreamingSimpleMosaicFilter.h:69
otb::StreamingSimpleMosaicFilter::OutputImageType
Superclass::OutputImageType OutputImageType
Definition: otbStreamingSimpleMosaicFilter.h:67
otb
The "otb" namespace contains all Orfeo Toolbox (OTB) classes.
Definition: otbJoinContainer.h:32
otb::StreamingSimpleMosaicFilter::OutputImagePointType
Superclass::OutputImagePointType OutputImagePointType
Definition: otbStreamingSimpleMosaicFilter.h:68
otbStreamingSimpleMosaicFilter.h
otb::StreamingSimpleMosaicFilter::IteratorType
Superclass::IteratorType IteratorType
Definition: otbStreamingSimpleMosaicFilter.h:62
otb::StreamingSimpleMosaicFilter::ThreadedGenerateData
void ThreadedGenerateData(const OutputImageRegionType &outputRegionForThread, itk::ThreadIdType threadId) override
Definition: otbStreamingSimpleMosaicFilter.hxx:34
otb::StreamingSimpleMosaicFilter::OutputImageInternalPixelType
Superclass::OutputImageInternalPixelType OutputImageInternalPixelType
Definition: otbStreamingSimpleMosaicFilter.h:70
otb::StreamingSimpleMosaicFilter::OutputImageRegionType
Superclass::OutputImageRegionType OutputImageRegionType
Definition: otbStreamingSimpleMosaicFilter.h:71
otb::StreamingSimpleMosaicFilter::InternalPixelType
Superclass::InternalPixelType InternalPixelType
Definition: otbStreamingSimpleMosaicFilter.h:75
otb::StreamingSimpleMosaicFilter::InputImagePixelType
Superclass::InputImagePixelType InputImagePixelType
Definition: otbStreamingSimpleMosaicFilter.h:61