OTB  9.0.0
Orfeo Toolbox
otbDisparityMapEstimationMethod.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 otbDisparityMapEstimationMethod_hxx
22 #define otbDisparityMapEstimationMethod_hxx
23 
25 #include "itkImageRegistrationMethod.h"
26 #include "otbExtractROI.h"
27 #include "otbMacro.h"
28 
29 // #include "otbImageFileWriter.h"
30 //
31 
32 namespace otb
33 {
34 /*
35  * Constructor.
36  */
37 template <typename TFixedImage, typename TMovingImage, class TPointSet>
39 {
40  this->SetNumberOfRequiredInputs(3);
41  // this->SetNumberOfRequiredOutputs(1);
42  // this->SetReleaseDataFlag(false);
43  this->SetReleaseDataBeforeUpdateFlag(false);
44  m_Transform = nullptr; // has to be provided by the user.
45  m_Interpolator = nullptr; // has to be provided by the user.
46  m_Metric = nullptr; // has to be provided by the user.
47  m_Optimizer = nullptr; // has to be provided by the user.
48  m_WinSize.Fill(15);
49  m_ExploSize.Fill(10);
50  m_InitialTransformParameters = ParametersType(1);
51  m_InitialTransformParameters.Fill(0.0f);
52 }
53 /*
54  * Destructor.
55  */
56 template <class TFixedImage, class TMovingImage, class TPointSet>
58 {
59 }
64 template <class TFixedImage, class TMovingImage, class TPointSet>
66 {
67  this->itk::ProcessObject::SetNthInput(0, const_cast<PointSetType*>(pointset));
68 }
69 
74 template <class TFixedImage, class TMovingImage, class TPointSet>
76 {
77  return static_cast<const PointSetType*>(this->itk::ProcessObject::GetInput(0));
78 }
79 
84 template <class TFixedImage, class TMovingImage, class TPointSet>
86 {
87  this->itk::ProcessObject::SetNthInput(1, const_cast<FixedImageType*>(image));
88 }
89 
94 template <class TFixedImage, class TMovingImage, class TPointSet>
96 {
97  return static_cast<const FixedImageType*>(this->itk::ProcessObject::GetInput(1));
98 }
99 
104 template <class TFixedImage, class TMovingImage, class TPointSet>
106 {
107  this->itk::ProcessObject::SetNthInput(2, const_cast<MovingImageType*>(image));
108 }
109 
114 template <class TFixedImage, class TMovingImage, class TPointSet>
116 {
117  return static_cast<const MovingImageType*>(this->itk::ProcessObject::GetInput(2));
118 }
119 
123 template <class TFixedImage, class TMovingImage, class TPointSet>
125 {
126  // inputs pointers
127  const FixedImageType* fixed = this->GetFixedImage();
128  const MovingImageType* moving = this->GetMovingImage();
129  const PointSetType* pointSet = this->GetPointSet();
130  PointSetType* output = this->GetOutput();
132 
133  // Typedefs
134  typedef typename PointSetType::PointsContainer PointsContainer;
135  typedef typename PointsContainer::ConstIterator PointsIterator;
136  typedef itk::ImageRegistrationMethod<FixedImageType, MovingImageType> RegistrationType;
137  typedef otb::ExtractROI<FixedPixelType, FixedPixelType> FixedExtractType;
138  typedef otb::ExtractROI<MovingPixelType, MovingPixelType> MovingExtractType;
139 
140  // points retrieving
141  typename PointsContainer::ConstPointer points = pointSet->GetPoints();
142 
143  // Iterator set up
144  PointsIterator pointIterator = points->Begin();
145  PointsIterator end = points->End();
146  unsigned int dataId = 0;
147 
148  otbMsgDevMacro(<< "Starting registration");
149 
151  while (pointIterator != end)
152  {
153  typename PointSetType::PointType p = pointIterator.Value(); // access the point
154 
155  // Extract the needed sub-images
156  typename FixedExtractType::Pointer fixedExtractor = FixedExtractType::New();
157  typename MovingExtractType::Pointer movingExtractor = MovingExtractType::New();
158  fixedExtractor->SetInput(fixed);
159  movingExtractor->SetInput(moving);
160 
161  // Fixed extractor setup
162  fixedExtractor->SetStartX(static_cast<unsigned long>(p[0] - m_ExploSize[0]));
163  fixedExtractor->SetStartY(static_cast<unsigned long>(p[1] - m_ExploSize[1]));
164  fixedExtractor->SetSizeX(static_cast<unsigned long>(2 * m_ExploSize[0] + 1));
165  fixedExtractor->SetSizeY(static_cast<unsigned long>(2 * m_ExploSize[1] + 1));
166  otbMsgDevMacro(<< "Point id: " << dataId);
167  otbMsgDevMacro(<< "Fixed region: origin(" << p[0] - m_ExploSize[0] << ", " << p[1] - m_ExploSize[1] << ") size(" << 2 * m_ExploSize[0] + 1 << ", "
168  << 2 * m_ExploSize[1] + 1 << ")");
169  // Moving extractor setup
170  movingExtractor->SetStartX(static_cast<unsigned long>(p[0] - m_WinSize[0]));
171  movingExtractor->SetStartY(static_cast<unsigned long>(p[1] - m_WinSize[1]));
172  movingExtractor->SetSizeX(static_cast<unsigned long>(2 * m_WinSize[0] + 1));
173  movingExtractor->SetSizeY(static_cast<unsigned long>(2 * m_WinSize[1] + 1));
174  otbMsgDevMacro(<< "Moving region: origin(" << p[0] - m_WinSize[0] << ", " << p[1] - m_WinSize[1] << ") size(" << 2 * m_WinSize[0] + 1 << ", "
175  << 2 * m_WinSize[1] + 1 << ")");
176  // update the extractors
177  fixedExtractor->Update();
178  movingExtractor->Update();
179 
180  // std::cout<<"Fixed extract origin: "<<fixedExtractor->GetOutput()->GetOrigin()<<std::endl;
181  // std::cout<<"Fixed extract spacing: "<<fixedExtractor->GetOutput()->GetSignedSpacing()<<std::endl;
182  // std::cout<<"Moving extract origin: "<<movingExtractor->GetOutput()->GetOrigin()<<std::endl;
183  // std::cout<<"Moving extract spacing: "<<movingExtractor->GetOutput()->GetSignedSpacing()<<std::endl;
184 
185  // typedef otb::ImageFileWriter<FixedImageType> FixedWriterType;
186  // typedef otb::ImageFileWriter<MovingImageType> MovingWriterType;
187 
188  // typename FixedWriterType::Pointer fwriter = FixedWriterType::New();
189  // typename MovingWriterType::Pointer mwriter = MovingWriterType::New();
190 
191  // std::ostringstream oss;
192  // oss.str("");
193  // oss<<"Fixed"<<dataId<<".tif";
194  // fwriter->SetInput(fixedExtractor->GetOutput());
195  // fwriter->SetFileName(oss.str());
196  // fwriter->Update();
197 
198  // oss.str("");
199  // oss<<"Moving"<<dataId<<".tif";
200  // mwriter->SetInput(movingExtractor->GetOutput());
201  // mwriter->SetFileName(oss.str());
202  // mwriter->Update();
203 
204  // Registration filter definition
205  typename RegistrationType::Pointer registration = RegistrationType::New();
206 
207  // Registration filter setup
208  registration->SetOptimizer(m_Optimizer);
209  registration->SetTransform(m_Transform);
210  registration->SetInterpolator(m_Interpolator);
211  registration->SetMetric(m_Metric);
212  registration->SetFixedImage(fixedExtractor->GetOutput());
213  registration->SetMovingImage(movingExtractor->GetOutput());
214 
215  // initial transform parameters setup
216  registration->SetInitialTransformParameters(m_InitialTransformParameters);
217  m_Interpolator->SetInputImage(movingExtractor->GetOutput());
218 
219  // Perform the registration
220  registration->Update();
221 
222  // Retrieve the final parameters
223  ParametersType finalParameters = registration->GetLastTransformParameters();
224  double value = m_Optimizer->GetValue(registration->GetLastTransformParameters());
225 
226  // Computing moving image point
227  typename FixedImageType::PointType inputPoint, outputPoint;
228  typename FixedImageType::IndexType inputIndex;
229 
230  // ensure that we have the right coord rep type
231  inputIndex[0] = static_cast<unsigned int>(p[0]);
232  inputIndex[1] = static_cast<unsigned int>(p[1]);
233 
234  fixed->TransformIndexToPhysicalPoint(inputIndex, inputPoint);
235 
236  m_Transform->SetParameters(finalParameters);
237  outputPoint = m_Transform->TransformPoint(inputPoint);
238 
239  otbMsgDevMacro(<< "Metric value: " << value);
240  otbMsgDevMacro(<< "Transform parameters: " << finalParameters);
241  otbMsgDevMacro(<< "Displacement: (" << outputPoint[0] - inputPoint[0] << ", " << outputPoint[1] - inputPoint[1] << ")");
242  otbMsgDevMacro(<< "Final parameters: " << finalParameters);
243 
244  ParametersType data(finalParameters.GetSize() + 3);
245 
246  data[0] = value;
247  data[1] = outputPoint[0] - inputPoint[0];
248  data[2] = outputPoint[1] - inputPoint[1];
249 
250  for (unsigned int i = 0; i < finalParameters.GetSize(); ++i)
251  {
252  data[i + 3] = finalParameters[i];
253  }
254 
255  // Set the parameters value in the point set data container.
256  output->SetPoint(dataId, p);
257  output->SetPointData(dataId, data);
258  // otbMsgDevMacro(<<"Point "<<dataId<<": "<<finalParameters);
259  ++pointIterator; // advance to next point
260  ++dataId;
261  }
262 }
263 template <class TFixedImage, class TMovingImage, class TPointSet>
265 {
266  Superclass::PrintSelf(os, indent);
267  os << indent << "Window size: " << m_WinSize << std::endl;
268  os << indent << "Exploration size: " << m_ExploSize << std::endl;
269 }
270 }
271 #endif
otb::DisparityMapEstimationMethod::PointSetType
TPointSet PointSetType
Definition: otbDisparityMapEstimationMethod.h:83
otb::DisparityMapEstimationMethod::ParametersType
MetricType::TransformParametersType ParametersType
Definition: otbDisparityMapEstimationMethod.h:110
otb::DisparityMapEstimationMethod::~DisparityMapEstimationMethod
~DisparityMapEstimationMethod() override
Definition: otbDisparityMapEstimationMethod.hxx:57
otb::DisparityMapEstimationMethod::SetPointSet
void SetPointSet(const TPointSet *pointset)
Definition: otbDisparityMapEstimationMethod.hxx:65
otbDisparityMapEstimationMethod.h
otb
The "otb" namespace contains all Orfeo Toolbox (OTB) classes.
Definition: otbJoinContainer.h:32
otb::DisparityMapEstimationMethod::SetFixedImage
void SetFixedImage(const TFixedImage *image)
Definition: otbDisparityMapEstimationMethod.hxx:85
otbMacro.h
otb::DisparityMapEstimationMethod::GetFixedImage
const TFixedImage * GetFixedImage(void)
Definition: otbDisparityMapEstimationMethod.hxx:95
otb::DisparityMapEstimationMethod::MovingImageType
TMovingImage MovingImageType
Definition: otbDisparityMapEstimationMethod.h:78
otb::DisparityMapEstimationMethod::DisparityMapEstimationMethod
DisparityMapEstimationMethod()
Definition: otbDisparityMapEstimationMethod.hxx:38
otbExtractROI.h
otb::DisparityMapEstimationMethod::FixedImageType
TFixedImage FixedImageType
Definition: otbDisparityMapEstimationMethod.h:69
otb::DisparityMapEstimationMethod::GenerateData
void GenerateData() override
Definition: otbDisparityMapEstimationMethod.hxx:124
otbMsgDevMacro
#define otbMsgDevMacro(x)
Definition: otbMacro.h:64
otb::DisparityMapEstimationMethod::GetPointSet
const TPointSet * GetPointSet(void)
Definition: otbDisparityMapEstimationMethod.hxx:75
otb::DisparityMapEstimationMethod::PrintSelf
void PrintSelf(std::ostream &os, itk::Indent indent) const override
Definition: otbDisparityMapEstimationMethod.hxx:264
otb::DisparityMapEstimationMethod::SetMovingImage
void SetMovingImage(const TMovingImage *image)
Definition: otbDisparityMapEstimationMethod.hxx:105
otb::ExtractROI
Extract a subset of a mono-channel image.
Definition: otbExtractROI.h:43
otb::DisparityMapEstimationMethod::GetMovingImage
const TMovingImage * GetMovingImage(void)
Definition: otbDisparityMapEstimationMethod.hxx:115