From d882042f13dd241ca396d91552931640eee4c7cf Mon Sep 17 00:00:00 2001
From: Julien Michel <julien.michel@c-s.fr>
Date: Wed, 4 Apr 2007 06:54:17 +0000
Subject: [PATCH] =?UTF-8?q?suite=20cartes=20de=20disparit=C3=A9.?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit
---
.../otbPointSetBasedResamplingFilter.txx | 182 ++++++++++++++++++
1 file changed, 182 insertions(+)
create mode 100644 Code/DisparityMap/otbPointSetBasedResamplingFilter.txx
diff --git a/Code/DisparityMap/otbPointSetBasedResamplingFilter.txx b/Code/DisparityMap/otbPointSetBasedResamplingFilter.txx
new file mode 100644
index 0000000000..305007336f
--- /dev/null
+++ b/Code/DisparityMap/otbPointSetBasedResamplingFilter.txx
@@ -0,0 +1,182 @@
+/*=========================================================================
+
+Program: ORFEO Toolbox
+Language: C++
+Date: $Date$
+Version: $Revision$
+
+
+Copyright (c) Centre National d'Etudes Spatiales. All rights reserved.
+See OTBCopyright.txt for details.
+
+
+This software is distributed WITHOUT ANY WARRANTY; without even
+the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+PURPOSE. See the above copyright notices for more information.
+
+=========================================================================*/
+#ifndef _otbPointSetBasedResamplingFilter_txx
+#define _otbPointSetBasedResamplingFilter_txx
+
+#include "otbPointSetBasedResamplingFilter.h"
+#include "itkImageRegionIteratorWithIndex.h"
+
+namespace otb
+{
+/**
+ * Constructor
+ */
+template <class TInputImage, class TPointSet, class TOutputImage>
+PointSetBasedResamplingFilter<TInputImage, TPointSet, TOutputImage>
+::PointSetBasedResamplingFilter()
+{
+ this->SetNumberOfRequiredInputs(2);
+ m_MetricThreshold = 0.0;
+ m_OutputSize.Fill(100);
+ m_OutputSpacing.Fill(1.0);
+ m_OutputOrigin.Fill(0.0);
+ m_DefaultValue = 0;
+ m_Interpolator = 0; // has to be provided by the user
+ m_Transform = 0; // has to be provided by the user
+}
+/**
+ * Set the pointset containing the disparity.
+ * \param pointset The pointset containing the disparity.
+ */
+template <class TInputImage, class TPointSet, class TOutputImage>
+void
+PointSetBasedResamplingFilter<TInputImage, TPointSet, TOutputImage>
+::SetPointSet(const TPointSet * pointset)
+{
+ this->itk::ProcessObject::SetNthInput(1,const_cast<PointSetType *>(pointset));
+}
+/**
+ * Get the pointset containing the disparity.
+ * \return The pointset containing the disparity.
+ */
+template <class TInputImage, class TPointSet, class TOutputImage>
+const TPointSet *
+PointSetBasedResamplingFilter<TInputImage, TPointSet, TOutputImage>
+::GetPointSet(void)
+{
+ return static_cast<const PointSetType *>(this->itk::ProcessObject::GetInput(1));
+}
+/** Generate output information */
+template <class TInputImage, class TPointSet, class TOutputImage>
+void
+PointSetBasedResamplingFilter<TInputImage, TPointSet, TOutputImage>
+::GenerateOutputInformation(void)
+{
+ OutputImagePointerType outputPtr = this->GetOutput();
+ typename OutputImageType::RegionType largest;
+ largest.SetSize(m_OutputSize);
+ IndexType index;
+ index.Fill(0);
+ largest.SetIndex(index);
+ outputPtr->SetLargestPossibleRegion(largest);
+ outputPtr->SetSpacing(m_OutputSpacing);
+ outputPtr->SetOrigin(m_OutputOrigin);
+}
+/** Main computation method */
+template <class TInputImage, class TPointSet, class TOutputImage>
+void
+PointSetBasedResamplingFilter<TInputImage, TPointSet, TOutputImage>
+::GenerateData(void)
+{
+ this->AllocateOutputs();
+ OutputImagePointerType outputPtr = this->GetOutput();
+ InputImageConstPointerType inputPtr = this->GetInput();
+ outputPtr->FillBuffer(0);
+ m_Interpolator->SetInputImage(inputPtr);
+
+ typedef itk::ImageRegionIteratorWithIndex<OutputImageType> IteratorType;
+ IteratorType outputIt(outputPtr,outputPtr->GetRequestedRegion());
+
+ for(outputIt.GoToBegin();!outputIt.IsAtEnd();++outputIt)
+ {
+ PointType outputPoint,inputPoint;
+ IndexType index;
+ typename TransformType::Pointer inverseTransform;
+
+ index = outputIt.GetIndex();
+ m_Transform->SetParameters(this->GetNearestPointTransformParameters(index));
+ m_Transform->GetInverse(inverseTransform);
+ outputPtr->TransformIndexToPhysicalPoint(index,outputPoint);
+ inputPoint = inverseTransform->TransformPoint(outputPoint);
+ otbMsgDevMacro(<<"back point: "<<inputPoint<<"(inverse trasform param: "<<inverseTransform->GetParameters()<<")");
+
+ if(m_Interpolator->IsInsideBuffer(inputPoint))
+ {
+ outputIt.Set(m_Interpolator->Evaluate(inputPoint));
+ }
+ else
+ {
+ outputIt.Set(m_DefaultValue);
+ }
+ }
+}
+/**
+ * \return The parameters of the transform associated with the nearest suitable point in pointset.
+ */
+template <class TInputImage, class TPointSet, class TOutputImage>
+typename PointSetBasedResamplingFilter<TInputImage, TPointSet, TOutputImage>
+::ParametersType
+PointSetBasedResamplingFilter<TInputImage, TPointSet, TOutputImage>
+::GetNearestPointTransformParameters(IndexType &index)
+{
+ typedef typename PointSetType::PointsContainer::ConstIterator PointSetIteratorType;
+ typedef typename PointSetType::PointsContainer PointsContainerType;
+ typedef typename PointSetType::PointDataContainer::ConstIterator PointDataIteratorType;
+ typedef typename PointSetType::PointDataContainer PointDataContainerType;
+
+ PointSetIteratorType it = this->GetPointSet()->GetPoints()->Begin();
+ PointDataIteratorType itData = this->GetPointSet()->GetPointData()->Begin();
+ OutputImagePointerType outputPtr = this->GetOutput();
+
+ double distance = -1.;
+ ParametersType parameters(m_Transform->GetNumberOfParameters());
+ PointType indexPoint;
+ outputPtr->TransformIndexToPhysicalPoint(index,indexPoint);
+
+ for(;it!=this->GetPointSet()->GetPoints()->End()&&itData!=this->GetPointSet()->GetPointData()->End();++it,++itData)
+ {
+ // If the point has a sufficient score
+ if(vcl_abs(itData.Value()[0])>m_MetricThreshold)
+ {
+ ParametersType tmpParameters(m_Transform->GetNumberOfParameters());
+ PointType inputPoint, outputPoint;
+
+ for(unsigned int i = 0; i<m_Transform->GetNumberOfParameters();++i)
+ {
+ tmpParameters[i]=itData.Value()[i+1];
+ }
+ inputPoint[0] = it.Value()[0];
+ inputPoint[1] = it.Value()[1];
+
+ m_Transform->SetParameters(tmpParameters);
+ outputPoint = m_Transform->TransformPoint(inputPoint);
+ // compute the distance to current point
+ double d = vcl_pow(outputPoint[0]-indexPoint[0],2) + vcl_pow(outputPoint[1]-indexPoint[1],2);
+
+ if(distance<0||distance>d)
+ {
+ distance = d;
+ parameters = tmpParameters;
+ }
+ }
+ }
+ return parameters;
+}
+
+/**
+ * PrintSelf Method
+ */
+template <class TInputImage, class TPointSet, class TOutputImage>
+void
+PointSetBasedResamplingFilter<TInputImage, TPointSet, TOutputImage>
+::PrintSelf(std::ostream& os, itk::Indent indent) const
+{
+ Superclass::PrintSelf(os, indent);
+}
+} // End namespace otb
+#endif
--
GitLab