OTB  6.7.0
Orfeo Toolbox
otbDisparityTranslateFilter.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 otbDisparityTranslateFilter_hxx
22 #define otbDisparityTranslateFilter_hxx
23 
25 
27 #include "itkImageRegionIterator.h"
28 
29 namespace otb
30 {
31 
32 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
35 {
36  m_NoDataValue = -32768;
37  // Set the number of inputs (1 moving image by default -> 3 inputs)
38  this->SetNumberOfRequiredInputs(6);
39  this->SetNumberOfRequiredInputs(1);
40 
41  // Set the outputs
42  this->SetNumberOfRequiredOutputs(2);
43  this->SetNthOutput(0,TDisparityImage::New());
44  this->SetNthOutput(1,TDisparityImage::New());
45 }
46 
47 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
48 void
50 ::SetHorizontalDisparityMapInput( const TDisparityImage * hmap )
51 {
52  // Process object is not const-correct so the const casting is required.
53  this->SetNthInput(0, const_cast<TDisparityImage *>( hmap ));
54 }
55 
56 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
57 void
59 ::SetVerticalDisparityMapInput( const TDisparityImage * vmap )
60 {
61  // Process object is not const-correct so the const casting is required.
62  this->SetNthInput(1, const_cast<TDisparityImage *>( vmap ));
63 }
64 
65 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
66 void
68 ::SetInverseEpipolarLeftGrid( const TGridImage * lgrid )
69 {
70  // Process object is not const-correct so the const casting is required.
71  this->SetNthInput(2, const_cast<TGridImage *>( lgrid ));
72 }
73 
74 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
75 void
77 ::SetDirectEpipolarRightGrid( const TGridImage * rgrid )
78 {
79  // Process object is not const-correct so the const casting is required.
80  this->SetNthInput(3, const_cast<TGridImage *>( rgrid ));
81 }
82 
83 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
84 void
86 ::SetDisparityMaskInput( const TMaskImage * mask )
87 {
88  // Process object is not const-correct so the const casting is required.
89  this->SetNthInput(4, const_cast<TMaskImage *>( mask ));
90 }
91 
92 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
93 void
95 ::SetLeftSensorImageInput(const TSensorImage * left)
96 {
97  // Process object is not const-correct so the const casting is required.
98  this->SetNthInput(5, const_cast<TSensorImage *>( left ));
99 }
100 
101 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
102 const TDisparityImage *
105 {
106  if (this->GetNumberOfInputs()<1)
107  {
108  return nullptr;
109  }
110  return static_cast<const TDisparityImage *>(this->itk::ProcessObject::GetInput(0));
111 }
112 
113 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
114 const TDisparityImage *
117 {
118  if (this->GetNumberOfInputs()<2)
119  {
120  return nullptr;
121  }
122  return static_cast<const TDisparityImage *>(this->itk::ProcessObject::GetInput(1));
123 }
124 
125 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
126 const TGridImage *
129 {
130  if (this->GetNumberOfInputs()<3)
131  {
132  return nullptr;
133  }
134  return static_cast<const TGridImage *>(this->itk::ProcessObject::GetInput(2));
135 }
136 
137 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
138 const TGridImage *
141 {
142  if (this->GetNumberOfInputs()<4)
143  {
144  return nullptr;
145  }
146  return static_cast<const TGridImage *>(this->itk::ProcessObject::GetInput(3));
147 }
148 
149 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
150 const TMaskImage *
153 {
154  if (this->GetNumberOfInputs()<5)
155  {
156  return nullptr;
157  }
158  return static_cast<const TMaskImage *>(this->itk::ProcessObject::GetInput(4));
159 }
160 
161 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
162 const TSensorImage *
165 {
166  if (this->GetNumberOfInputs()<6)
167  {
168  return nullptr;
169  }
170  return static_cast<const TSensorImage *>(this->itk::ProcessObject::GetInput(5));
171 }
172 
173 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
174 TDisparityImage *
177 {
178  if (this->GetNumberOfOutputs()<1)
179  {
180  return nullptr;
181  }
182  return static_cast<TDisparityImage *>(this->itk::ProcessObject::GetOutput(0));
183 }
184 
185 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
186 TDisparityImage *
189 {
190  if (this->GetNumberOfOutputs()<2)
191  {
192  return nullptr;
193  }
194  return static_cast<TDisparityImage *>(this->itk::ProcessObject::GetOutput(1));
195 }
196 
197 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
198 void
201 {
202 // this->Superclass::GenerateOutputInformation();
203 
204  TDisparityImage * horizOut = this->GetHorizontalDisparityMapOutput();
205  TDisparityImage * vertiOut = this->GetVerticalDisparityMapOutput();
206 
207  const TSensorImage * leftIn = this->GetLeftSensorImageInput();
208 
209  horizOut->CopyInformation(leftIn);
210  vertiOut->CopyInformation(leftIn);
211 
212  // Set the NoData value
213  std::vector<bool> noDataValueAvailable;
214  noDataValueAvailable.push_back(true);
215  std::vector<double> noDataValue;
216  noDataValue.push_back(m_NoDataValue);
217  itk::MetaDataDictionary& dict = horizOut->GetMetaDataDictionary();
218  itk::EncapsulateMetaData<std::vector<bool> >(dict,MetaDataKey::NoDataValueAvailable,noDataValueAvailable);
219  itk::EncapsulateMetaData<std::vector<double> >(dict,MetaDataKey::NoDataValue,noDataValue);
220  dict = vertiOut->GetMetaDataDictionary();
221  itk::EncapsulateMetaData<std::vector<bool> >(dict,MetaDataKey::NoDataValueAvailable,noDataValueAvailable);
222  itk::EncapsulateMetaData<std::vector<double> >(dict,MetaDataKey::NoDataValue,noDataValue);
223 }
224 
225 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
226 void
229 {
230  this->Superclass::GenerateInputRequestedRegion();
231 
232  TGridImage * leftGrid = const_cast<TGridImage*>(this->GetInverseEpipolarLeftGrid());
233  TGridImage * rightGrid = const_cast<TGridImage*>(this->GetDirectEpipolarRightGrid());
234 
235  leftGrid->SetRequestedRegionToLargestPossibleRegion();
236  rightGrid->SetRequestedRegionToLargestPossibleRegion();
237 
238  leftGrid->UpdateOutputData();
239  rightGrid->UpdateOutputData();
240 
241  TSensorImage * leftIn = const_cast<TSensorImage*>(this->GetLeftSensorImageInput());
242  RegionType emptyRegion = leftIn->GetLargestPossibleRegion();
243  emptyRegion.SetSize(0,0);
244  emptyRegion.SetSize(1,0);
245  leftIn->SetRequestedRegion(emptyRegion);
246 
247  TDisparityImage * horizOut = this->GetHorizontalDisparityMapOutput();
248 
249  TDisparityImage * horizIn = const_cast<TDisparityImage*>(this->GetHorizontalDisparityMapInput());
250  TDisparityImage * vertiIn = const_cast<TDisparityImage*>(this->GetVerticalDisparityMapInput());
251 
252  TMaskImage * maskIn = const_cast<TMaskImage*>(this->GetDisparityMaskInput());
253 
254  RegionType requested = this->GetHorizontalDisparityMapOutput()->GetRequestedRegion();
255  RegionType inputlargest = this->GetHorizontalDisparityMapInput()->GetLargestPossibleRegion();
256  RegionType inputRequested;
257  GridRegionType gridLargest = leftGrid->GetLargestPossibleRegion();
258 
259  IndexType minIndex,maxIndex;
260  // -Wmaybe-uninitialized
263 
264  IndexType corners[4];
265  for(int i = 0; i < 4; i++)
267 
268  corners[0] = requested.GetIndex();
269  corners[1] = requested.GetIndex();
270  corners[1][0] += static_cast<IndexValueType>(requested.GetSize()[0]) - 1;
271  corners[2] = requested.GetIndex();
272  corners[2][1] += static_cast<IndexValueType>(requested.GetSize()[1]) - 1;
273  corners[3] = requested.GetIndex();
274  corners[3][0] += static_cast<IndexValueType>(requested.GetSize()[0]) - 1;
275  corners[3][1] += static_cast<IndexValueType>(requested.GetSize()[1]) - 1;
276  for (unsigned int k=0; k<4; ++k)
277  {
278  PointType pointSensor;
279  horizOut->TransformIndexToPhysicalPoint(corners[k],pointSensor);
281  leftGrid->TransformPhysicalPointToContinuousIndex(pointSensor,indexGrid);
282  IndexType ul;
283  ul[0] = static_cast<long>(std::floor(indexGrid[0]));
284  ul[1] = static_cast<long>(std::floor(indexGrid[1]));
285  if (ul[0] < gridLargest.GetIndex()[0]) ul[0] = gridLargest.GetIndex()[0];
286  if (ul[1] < gridLargest.GetIndex()[1]) ul[1] = gridLargest.GetIndex()[1];
287  if (ul[0] > static_cast<IndexValueType>(gridLargest.GetIndex()[0] + gridLargest.GetSize()[0]-2))
288  ul[0] = (gridLargest.GetIndex()[0] + gridLargest.GetSize()[0]-2);
289  if (ul[1] > static_cast<IndexValueType>(gridLargest.GetIndex()[1] + gridLargest.GetSize()[1]-2))
290  ul[1] = (gridLargest.GetIndex()[1] + gridLargest.GetSize()[1]-2);
291 
292  IndexType ur = ul;
293  ur[0] += 1;
294  IndexType ll = ul;
295  ll[1] += 1;
296  IndexType lr = ul;
297  lr[0] += 1;
298  lr[1] += 1;
299 
300  double rx = indexGrid[0] - static_cast<double>(ul[0]);
301  double ry = indexGrid[1] - static_cast<double>(ul[1]);
302  PointType pointEpi = pointSensor;
303 
304  pointEpi[0] += (1. - ry) * ((1. - rx) * leftGrid->GetPixel(ul)[0] + rx * leftGrid->GetPixel(ur)[0]) +
305  ry * ((1. - rx) * leftGrid->GetPixel(ll)[0] + rx * leftGrid->GetPixel(lr)[0]);
306  pointEpi[1] += (1. - ry) * ((1. - rx) * leftGrid->GetPixel(ul)[1] + rx * leftGrid->GetPixel(ur)[1]) +
307  ry * ((1. - rx) * leftGrid->GetPixel(ll)[1] + rx * leftGrid->GetPixel(lr)[1]);
309 
310  horizIn->TransformPhysicalPointToContinuousIndex(pointEpi,indexEpi);
311  if (k == 0)
312  {
313  minIndex[0] = static_cast<long>(std::floor(indexEpi[0]));
314  minIndex[1] = static_cast<long>(std::floor(indexEpi[1]));
315  maxIndex[0] = static_cast<long>(std::ceil(indexEpi[0]));
316  maxIndex[1] = static_cast<long>(std::ceil(indexEpi[1]));
317  }
318  else
319  {
320  if (minIndex[0]>static_cast<long>(std::floor(indexEpi[0]))) minIndex[0]=static_cast<long>(std::floor(indexEpi[0]));
321  if (minIndex[1]>static_cast<long>(std::floor(indexEpi[1]))) minIndex[1]=static_cast<long>(std::floor(indexEpi[1]));
322  if (maxIndex[0]<static_cast<long>(std::ceil(indexEpi[0]))) maxIndex[0]=static_cast<long>(std::ceil(indexEpi[0]));
323  if (maxIndex[1]<static_cast<long>(std::ceil(indexEpi[1]))) maxIndex[1]=static_cast<long>(std::ceil(indexEpi[1]));
324  }
325  }
326 
327  inputRequested.SetIndex(minIndex);
328  inputRequested.SetSize(0,static_cast<unsigned long>(maxIndex[0]-minIndex[0]));
329  inputRequested.SetSize(1,static_cast<unsigned long>(maxIndex[1]-minIndex[1]));
330 
331  if (!inputRequested.Crop(inputlargest))
332  {
333  inputRequested.SetSize(0,0);
334  inputRequested.SetSize(1,0);
335  inputRequested.SetIndex(inputlargest.GetIndex());
336  }
337 
338  horizIn->SetRequestedRegion(inputRequested);
339  if (vertiIn) vertiIn->SetRequestedRegion(inputRequested);
340  if (maskIn) maskIn->SetRequestedRegion(inputRequested);
341 }
342 
343 template <class TDisparityImage, class TGridImage, class TSensorImage, class TMaskImage>
344 void
346 ::ThreadedGenerateData(const RegionType & outputRegionForThread, itk::ThreadIdType itkNotUsed(threadId))
347 {
348  const TGridImage * leftGrid = this->GetInverseEpipolarLeftGrid();
349  const TGridImage * rightGrid = this->GetDirectEpipolarRightGrid();
350 
351  TDisparityImage * horizOut = this->GetHorizontalDisparityMapOutput();
352  TDisparityImage * vertiOut = this->GetVerticalDisparityMapOutput();
353 
354  const TDisparityImage * horizIn = this->GetHorizontalDisparityMapInput();
355  const TDisparityImage * vertiIn = this->GetVerticalDisparityMapInput();
356 
357  const TMaskImage * maskIn = this->GetDisparityMaskInput();
358 
359  GridRegionType leftLargest = leftGrid->GetLargestPossibleRegion();
360  GridRegionType rightLargest = rightGrid->GetLargestPossibleRegion();
361  RegionType buffered = horizIn->GetBufferedRegion();
362 
363  const bool emptyInputRegion=buffered.GetNumberOfPixels() == 0;
364 
366  DispIterator horizIter(horizOut, outputRegionForThread);
367  DispIterator vertiIter(vertiOut, outputRegionForThread);
368 
369  horizIter.GoToBegin();
370  vertiIter.GoToBegin();
371 
372  while (!horizIter.IsAtEnd() && !vertiIter.IsAtEnd())
373  {
374 
375  if(!emptyInputRegion)
376  {
377  PointType pointSensor;
378  horizOut->TransformIndexToPhysicalPoint(horizIter.GetIndex(), pointSensor);
379 
381  leftGrid->TransformPhysicalPointToContinuousIndex(pointSensor, indexGrid);
382 
383  // Interpolate in left grid
384  IndexType ul;
385  ul[0] = static_cast<long> (std::floor(indexGrid[0]));
386  ul[1] = static_cast<long> (std::floor(indexGrid[1]));
387  if (ul[0] < leftLargest.GetIndex()[0]) ul[0] = leftLargest.GetIndex()[0];
388  if (ul[1] < leftLargest.GetIndex()[1]) ul[1] = leftLargest.GetIndex()[1];
389  if (ul[0] > static_cast<IndexValueType>(leftLargest.GetIndex()[0] + leftLargest.GetSize()[0] - 2))
390  ul[0] = (leftLargest.GetIndex()[0] + leftLargest.GetSize()[0] - 2);
391  if (ul[1] > static_cast<IndexValueType>(leftLargest.GetIndex()[1] + leftLargest.GetSize()[1] - 2))
392  ul[1] = (leftLargest.GetIndex()[1] + leftLargest.GetSize()[1] - 2);
393 
394  IndexType ur = ul;
395  ur[0] += 1;
396  IndexType ll = ul;
397  ll[1] += 1;
398  IndexType lr = ul;
399  lr[0] += 1;
400  lr[1] += 1;
401 
402  double rx = indexGrid[0] - static_cast<double> (ul[0]);
403  double ry = indexGrid[1] - static_cast<double> (ul[1]);
404  PointType pointEpi = pointSensor;
405 
406  pointEpi[0] += (1. - ry) * ((1. - rx) * leftGrid->GetPixel(ul)[0] + rx * leftGrid->GetPixel(ur)[0]) + ry * ((1.
407  - rx) * leftGrid->GetPixel(ll)[0] + rx * leftGrid->GetPixel(lr)[0]);
408  pointEpi[1] += (1. - ry) * ((1. - rx) * leftGrid->GetPixel(ul)[1] + rx * leftGrid->GetPixel(ur)[1]) + ry * ((1.
409  - rx) * leftGrid->GetPixel(ll)[1] + rx * leftGrid->GetPixel(lr)[1]);
410 
412  horizIn->TransformPhysicalPointToContinuousIndex(pointEpi, indexEpi);
413 
414  // Interpolate in disparity map
415  ul[0] = static_cast<long> (std::floor(indexEpi[0]));
416  ul[1] = static_cast<long> (std::floor(indexEpi[1]));
417  if (ul[0] < buffered.GetIndex()[0]) ul[0] = buffered.GetIndex()[0];
418  if (ul[1] < buffered.GetIndex()[1]) ul[1] = buffered.GetIndex()[1];
419  if (ul[0] > static_cast<IndexValueType>(buffered.GetIndex()[0] + buffered.GetSize()[0] - 2))
420  ul[0] = (buffered.GetIndex()[0] + buffered.GetSize()[0] - 2);
421  if (ul[1] > static_cast<IndexValueType>(buffered.GetIndex()[1] + buffered.GetSize()[1] - 2))
422  ul[1] = (buffered.GetIndex()[1] + buffered.GetSize()[1] - 2);
423 
424  ur = ul;
425  ur[0] += 1;
426  ll = ul;
427  ll[1] += 1;
428  lr = ul;
429  lr[0] += 1;
430  lr[1] += 1;
431 
432  // check if all corners are valid
433  if (!maskIn || (maskIn && maskIn->GetPixel(ul) > 0 && maskIn->GetPixel(ur) > 0 && maskIn->GetPixel(ll) > 0 &&
434  maskIn->GetPixel(lr) > 0))
435  {
436  rx = indexEpi[0] - static_cast<double> (ul[0]);
437  ry = indexEpi[1] - static_cast<double> (ul[1]);
438 
439  itk::ContinuousIndex<double, 2> indexRight(indexEpi);
440 
441  indexRight[0] += (1. - ry) * ((1. - rx) * horizIn->GetPixel(ul) + rx * horizIn->GetPixel(ur)) + ry * ((1. - rx)
442  * horizIn->GetPixel(ll) + rx * horizIn->GetPixel(lr));
443  if (vertiIn)
444  {
445  indexRight[1] += (1. - ry) * ((1. - rx) * vertiIn->GetPixel(ul) + rx * vertiIn->GetPixel(ur)) + ry * ((1.
446  - rx) * vertiIn->GetPixel(ll) + rx * vertiIn->GetPixel(lr));
447  }
448 
449  PointType pointRight;
450  horizIn->TransformContinuousIndexToPhysicalPoint(indexRight, pointRight);
451 
452  itk::ContinuousIndex<double, 2> indexGridRight;
453  rightGrid->TransformPhysicalPointToContinuousIndex(pointRight, indexGridRight);
454 
455  // Interpolate in right grid
456  ul[0] = static_cast<long> (std::floor(indexGridRight[0]));
457  ul[1] = static_cast<long> (std::floor(indexGridRight[1]));
458  if (ul[0] < rightLargest.GetIndex()[0]) ul[0] = rightLargest.GetIndex()[0];
459  if (ul[1] < rightLargest.GetIndex()[1]) ul[1] = rightLargest.GetIndex()[1];
460  if (ul[0] > static_cast<IndexValueType>(rightLargest.GetIndex()[0] + rightLargest.GetSize()[0] - 2))
461  ul[0] = (rightLargest.GetIndex()[0] + rightLargest.GetSize()[0] - 2);
462  if (ul[1] > static_cast<IndexValueType>(rightLargest.GetIndex()[1] + rightLargest.GetSize()[1] - 2))
463  ul[1] = (rightLargest.GetIndex()[1] + rightLargest.GetSize()[1] - 2);
464 
465  ur = ul;
466  ur[0] += 1;
467  ll = ul;
468  ll[1] += 1;
469  lr = ul;
470  lr[0] += 1;
471  lr[1] += 1;
472 
473  rx = indexGridRight[0] - static_cast<double> (ul[0]);
474  ry = indexGridRight[1] - static_cast<double> (ul[1]);
475 
476  PointType pointSensorRight = pointRight;
477 
478  pointSensorRight[0] += (1. - ry) * ((1. - rx) * rightGrid->GetPixel(ul)[0] + rx * rightGrid->GetPixel(ur)[0]) +
479  ry * ((1. - rx) * rightGrid->GetPixel(ll)[0] + rx * rightGrid->GetPixel(lr)[0]);
480  pointSensorRight[1] += (1. - ry) * ((1. - rx) * rightGrid->GetPixel(ul)[1] + rx * rightGrid->GetPixel(ur)[1]) +
481  ry * ((1. - rx) * rightGrid->GetPixel(ll)[1] + rx * rightGrid->GetPixel(lr)[1]);
482 
483  horizIter.Set(pointSensorRight[0] - pointSensor[0]);
484  vertiIter.Set(pointSensorRight[1] - pointSensor[1]);
485  }
486  else
487  {
488  horizIter.Set(m_NoDataValue);
489  vertiIter.Set(m_NoDataValue);
490  }
491  }
492  else
493  {
494  horizIter.Set(m_NoDataValue);
495  vertiIter.Set(m_NoDataValue);
496  }
497  ++horizIter;
498  ++vertiIter;
499  }
500 
501 
502 }
503 
504 
505 }
506 
507 #endif
const TMaskImage * GetDisparityMaskInput() const
void SetLeftSensorImageInput(const TSensorImage *left)
void SetInverseEpipolarLeftGrid(const TGridImage *lgrid)
const TGridImage * GetDirectEpipolarRightGrid() const
void ThreadedGenerateData(const RegionType &outputRegionForThread, itk::ThreadIdType threadId) override
const TSensorImage * GetLeftSensorImageInput() const
OTBOSSIMAdapters_EXPORT char const * NoDataValue
const TDisparityImage * GetHorizontalDisparityMapInput() const
void SetHorizontalDisparityMapInput(const TDisparityImage *hmap)
void SetDirectEpipolarRightGrid(const TGridImage *rgrid)
OTBOSSIMAdapters_EXPORT char const * NoDataValueAvailable
void SetDisparityMaskInput(const TMaskImage *mask)
unsigned int ThreadIdType
const TDisparityImage * GetVerticalDisparityMapInput() const
DataObject * GetInput(const DataObjectIdentifierType &key)
DispMapType::IndexValueType IndexValueType
void SetVerticalDisparityMapInput(const TDisparityImage *vmap)
const TGridImage * GetInverseEpipolarLeftGrid() const
DataObject * GetOutput(const DataObjectIdentifierType &key)