/* * Copyright (C) 2005-2019 Centre National d'Etudes Spatiales (CNES) * * This file is part of Orfeo Toolbox * * https://www.orfeo-toolbox.org/ * * Licensed under the Apache License, Version 2.0 (the "License"); * you may not use this file except in compliance with the License. * You may obtain a copy of the License at * * http://www.apache.org/licenses/LICENSE-2.0 * * Unless required by applicable law or agreed to in writing, software * distributed under the License is distributed on an "AS IS" BASIS, * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. * See the License for the specific language governing permissions and * limitations under the License. */ #ifndef otbStereorectificationDisplacementFieldSource_h #define otbStereorectificationDisplacementFieldSource_h #include "itkImageSource.h" #include "otbGenericRSTransform.h" #include "itkImageRegionIteratorWithIndex.h" #include "itkVector.h" #include "otbMacro.h" namespace otb { /** \class StereorectificationDisplacementFieldSource * \brief Compute the deformation fields for stereo-rectification * * The geometry of acqusition related to stereo pairs is such that * the displacement related to the elevation between the two images * of the pair always occurs along lines that are called epipolar * lines. * * In the case of conic acqusition, these lines are perfectly * parallel by definition, while in push-broom geometry, these lines * might not be exactly parallel due to perturbations during * acquisition, especially when considering a large field of view. * * The purpose of stereo-rectification is to warp both images of the * pairs so that the displacement related to the elevation only * occurs in the horizontal direction (i.e. epipolar lines are * horizontal). This operation is useful for mainly two reasons: it * allows searching for disparities in one direction only, and it * allows deriving anaglyph for 3D viewing with 3D glasses. * * This filter allows you to compute the deformation fields up to the * sensor model precision needed to warp a pair of stereo images into * epipolar geometry. Warping can be done using the * otb::StreamingWarpImageFilter. * * Since lines might not be perfectly regular, the algorithm * performed by this filter uses the otb::GenericRSTransform * capabilities to compute the local epipolar lines, and iteratively * build a resampling grid by propagating along these locally * estimated lines. * * Epipolar images will have a null origin and a size as given by the * GetRectifiedImageSize() method. The deformation fields and size * are derived to produce epipolar images covering the whole extent * of the left image. * * The SetAverageElevation() method allows you to set the elevation * hypothesis on which the epipolar geometry is built. It means that * any pair of pixels in the stereo pair whose elevation is exactly * equal to the average elevation will have a null disparity (no * displacement). The SetElevationOffset() method allows tuning the * elevation offset which is only used for local epipolar lines * estimation. The default value of 50 meters should do. * * Additionnaly, the SetScale() method allows deriving deformation * fields and images size at a coarser (scale > 1) or finer (scale < * 1) resolution. The SetGridStep() allows tuning the step of the * resampling grid. Please keep in mind that the whole grid is loaded * into memory, and that the epipolar lines direction may only vary * smoothly. When working with large images, a coarse grid-step will * generally be accurate enough and will preserve the memory resources. * * \sa StreamingWarpImageFilter * \sa StereoSensorModelToElevationMapFilter * * \ingroup OTBStereo */ template class ITK_EXPORT StereorectificationDisplacementFieldSource : public itk::ImageSource { public: /** Standard class typedefs */ typedef StereorectificationDisplacementFieldSource Self; typedef itk::ImageSource Superclass; typedef itk::SmartPointer Pointer; typedef itk::SmartPointer ConstPointer; // Inputs and outputs relative typedefs typedef TInputImage InputImageType; typedef typename InputImageType::Pointer InputImagePointerType; typedef TOutputImage OutputImageType; typedef typename OutputImageType::Pointer OutputImagePointerType; // Image related typedefs typedef typename OutputImageType::SizeType SizeType; typedef typename OutputImageType::PointType PointType; typedef typename OutputImageType::SpacingType SpacingType; typedef typename OutputImageType::RegionType RegionType; // 3D RS transform // TODO: Allow tuning precision (i.e. double or float) typedef otb::GenericRSTransform RSTransformType; typedef typename RSTransformType::Pointer RSTransformPointerType; // 3D points typedef typename RSTransformType::InputPointType TDPointType; /** Method for creation through the object factory. */ itkNewMacro(Self); /** Run-time type information (and related methods). */ itkTypeMacro(StereorectificationDisplacementFieldSource, ImageSource); /** Get/Set the scale */ itkSetMacro(Scale, double); itkGetConstReferenceMacro(Scale, double); /** Get/Set the grid scale */ itkSetMacro(GridStep, double); itkGetMacro(GridStep, double); /** Get/Set left image for stereo-rectification */ itkSetObjectMacro(LeftImage, InputImageType); itkGetObjectMacro(LeftImage, InputImageType); /** Get/Set right image for stereo-rectification */ itkSetObjectMacro(RightImage, InputImageType); itkGetObjectMacro(RightImage, InputImageType); /** Get the size of the rectified image */ itkGetConstReferenceMacro(RectifiedImageSize, SizeType); /** Get the estimated mean baseline ratio */ itkGetConstReferenceMacro(MeanBaselineRatio, double); /** Return the left deformation field (const version) */ const OutputImageType* GetLeftDisplacementFieldOutput() const; /** Return the left deformation field */ OutputImageType* GetLeftDisplacementFieldOutput(); /** Return the left deformation field (const version) */ const OutputImageType* GetRightDisplacementFieldOutput() const; /** Return the left deformation field */ OutputImageType* GetRightDisplacementFieldOutput(); itkSetMacro(UseDEM, bool); itkGetMacro(UseDEM, bool); itkBooleanMacro(UseDEM); protected: /** Constructor */ StereorectificationDisplacementFieldSource(void); /** Destructor */ ~StereorectificationDisplacementFieldSource(void) override{}; /** Generate output images information */ void GenerateOutputInformation() override; /** Enlarge output requested region (no streaming) */ void EnlargeOutputRequestedRegion(itk::DataObject* itkNotUsed(output)) override; /** Compute the deformation field */ void GenerateData() override; /** PrintSelf method */ void PrintSelf(std::ostream& os, itk::Indent indent) const override; private: StereorectificationDisplacementFieldSource(const Self&); // purposely // not // implemented void operator=(const Self&) = delete; /** This elevation offset is used to compute the epipolar direction */ double m_ElevationOffset; /** A scale greater than 1 will lead to zoomed stereo-rectified * pairs */ double m_Scale; /** Controls the step of the resampling grid (in pixels). A finer * step will lead to more memory consumption. */ double m_GridStep; /** Left image */ InputImagePointerType m_LeftImage; /** Right image */ InputImagePointerType m_RightImage; /** Left to right transform */ RSTransformPointerType m_LeftToRightTransform; /** Right to left transform */ RSTransformPointerType m_RightToLeftTransform; /** Size of the rectified images */ SizeType m_RectifiedImageSize; /** Output origin in left image (internal use) */ TDPointType m_OutputOriginInLeftImage; /** This variable contains the estimated mean baseline ratio over * the image */ double m_MeanBaselineRatio; /** If set to true, elevation is retrieved through * DEMHandler::GetHeightAboveEllipsoid(). If false, elevation is * retrieved from DEMHandler::GetDefaultHeightAboveEllipsoid() */ bool m_UseDEM; }; } // End namespace otb #ifndef OTB_MANUAL_INSTANTIATION #include "otbStereorectificationDisplacementFieldSource.hxx" #endif #endif