Commit 613bfa77 authored by Cédric Traizet's avatar Cédric Traizet
Browse files

Merge branch 'refactor_geocentric_transform' into 'develop'

Refactor GeocentricTransform and remove EllipsoidAdapter

See merge request !845
parents c783367e 795199a5
Pipeline #8443 passed with stages
in 6 minutes and 56 seconds
/*
* Copyright (C) 2005-2020 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 otbEllipsoidAdapter_h
#define otbEllipsoidAdapter_h
#include "itkObject.h"
#include "itkObjectFactory.h"
#include "OTBOSSIMAdaptersExport.h"
class ossimEllipsoid;
namespace otb
{
/** \class EllipsoidAdapter
* \brief class providing the interface to the ossimEllipsoid
*
* This is a class to be used internally instead of introducing a dependence on
* ossimEllipsoid
*
* \todo{Add the support for different ellipsoid models}
*
*
* \ingroup OTBOSSIMAdapters
**/
class OTBOSSIMAdapters_EXPORT EllipsoidAdapter : public itk::Object
{
public:
/** Standard class typedefs. */
typedef EllipsoidAdapter Self;
typedef itk::Object Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
/** Method for creation through the object factory. */
itkNewMacro(Self);
/** Run-time type information (and related methods). */
itkTypeMacro(EllipsoidAdapter, itk::Object);
/** Convert a XYZ coordinate into a lon, lat, height on the ellipsoid */
void XYZToLonLatHeight(double x, double y, double z, double& lon, double& lat, double& h) const;
/** Convert a lon, lat, height on the ellipsoid into a XYZ geocentric system*/
void LonLatHeightToXYZ(double lon, double lat, double h, double& x, double& y, double& z) const;
protected:
EllipsoidAdapter();
~EllipsoidAdapter() override;
private:
EllipsoidAdapter(const Self&) = delete;
void operator=(const Self&) = delete;
ossimEllipsoid* m_Ellipsoid;
};
} // namespace otb
#endif
......@@ -23,7 +23,6 @@ set(OTBOSSIMAdapters_SRC
otbImageKeywordlist.cxx
otbRPCSolverAdapter.cxx
otbDateTimeAdapter.cxx
otbEllipsoidAdapter.cxx
)
add_library(OTBOSSIMAdapters ${OTBOSSIMAdapters_SRC})
......
/*
* Copyright (C) 2005-2020 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.
*/
#include "otbEllipsoidAdapter.h"
#if defined(__GNUC__) || defined(__clang__)
#pragma GCC diagnostic push
#pragma GCC diagnostic ignored "-Wunused-parameter"
#pragma GCC diagnostic ignored "-Woverloaded-virtual"
#pragma GCC diagnostic ignored "-Wshadow"
#include "ossim/base/ossimEllipsoid.h"
#pragma GCC diagnostic pop
#else
#include "ossim/base/ossimEllipsoid.h"
#endif
namespace otb
{
EllipsoidAdapter::EllipsoidAdapter()
{
m_Ellipsoid = new ossimEllipsoid();
}
EllipsoidAdapter::~EllipsoidAdapter()
{
if (m_Ellipsoid != nullptr)
{
delete m_Ellipsoid;
}
}
void EllipsoidAdapter::XYZToLonLatHeight(double x, double y, double z, double& lon, double& lat, double& h) const
{
// Note the lat/lon convension for ossim vs lon/lat for OTB
m_Ellipsoid->XYZToLatLonHeight(x, y, z, lat, lon, h);
}
void EllipsoidAdapter::LonLatHeightToXYZ(double lon, double lat, double h, double& x, double& y, double& z) const
{
// Note the lat/lon convension for ossim vs lon/lat for OTB
m_Ellipsoid->latLonHeightToXYZ(lat, lon, h, x, y, z);
}
} // namespace otb
......@@ -21,8 +21,8 @@
#ifndef otbGeocentricTransform_h
#define otbGeocentricTransform_h
#include "otbGenericMapProjection.h"
#include "otbEllipsoidAdapter.h"
#include "otbCoordinateTransformation.h"
#include "otbTransform.h"
namespace otb
{
......@@ -33,22 +33,22 @@ namespace otb
*
* \ingroup OTBTransform
*/
template <TransformDirection TDirectionOfMapping, class TScalarType = double, unsigned int NInputDimensions = 3,
unsigned int NOutputDimensions = 3>
template <TransformDirection TDirectionOfMapping, class TScalarType = double>
class ITK_EXPORT GeocentricTransform : public Transform<TScalarType, // Data type for scalars
NInputDimensions, // Number of dimensions in the input space
NOutputDimensions> // Number of dimensions in the output space
3, // Number of dimensions in the input space
3> // Number of dimensions in the output space
{
public:
/** Standard class typedefs. */
typedef GeocentricTransform Self;
typedef Transform<TScalarType, NInputDimensions, NOutputDimensions> Superclass;
typedef Transform<TScalarType, 3, 3> Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
typedef typename Superclass::ScalarType ScalarType;
typedef itk::Point<ScalarType, NInputDimensions> InputPointType;
typedef itk::Point<ScalarType, NOutputDimensions> OutputPointType;
typedef itk::Point<ScalarType, 3> InputPointType;
typedef itk::Point<ScalarType, 3> OutputPointType;
/** Method for creation through the object factory. */
itkNewMacro(Self);
......@@ -58,21 +58,22 @@ public:
static const TransformDirection DirectionOfMapping = TDirectionOfMapping;
itkStaticConstMacro(InputSpaceDimension, unsigned int, NInputDimensions);
itkStaticConstMacro(OutputSpaceDimension, unsigned int, NOutputDimensions);
itkStaticConstMacro(SpaceDimension, unsigned int, NInputDimensions);
itkStaticConstMacro(ParametersDimension, unsigned int, NInputDimensions*(NInputDimensions + 1));
itkStaticConstMacro(InputSpaceDimension, unsigned int, 3);
itkStaticConstMacro(OutputSpaceDimension, unsigned int, 3);
itkStaticConstMacro(SpaceDimension, unsigned int, 3);
itkStaticConstMacro(ParametersDimension, unsigned int, 3*(3 + 1));
OutputPointType TransformPoint(const InputPointType& point) const override;
protected:
GeocentricTransform();
~GeocentricTransform() override;
EllipsoidAdapter::Pointer m_Ellipsoid;
~GeocentricTransform() override = default;
private:
GeocentricTransform(const Self&) = delete;
void operator=(const Self&) = delete;
std::unique_ptr<CoordinateTransformation> m_MapProjection;
};
} // namespace otb
......
......@@ -26,32 +26,42 @@
namespace otb
{
template <TransformDirection TDirectionOfMapping, class TScalarType, unsigned int NInputDimensions, unsigned int NOutputDimensions>
GeocentricTransform<TDirectionOfMapping, TScalarType, NInputDimensions, NOutputDimensions>::GeocentricTransform() : Superclass(ParametersDimension)
template <TransformDirection TDirectionOfMapping, class TScalarType>
GeocentricTransform<TDirectionOfMapping, TScalarType>::GeocentricTransform() : Superclass(ParametersDimension)
{
m_Ellipsoid = EllipsoidAdapter::New();
}
template <TransformDirection TDirectionOfMapping, class TScalarType, unsigned int NInputDimensions, unsigned int NOutputDimensions>
GeocentricTransform<TDirectionOfMapping, TScalarType, NInputDimensions, NOutputDimensions>::~GeocentricTransform()
{
}
SpatialReference wgs84 = SpatialReference::FromWGS84();
SpatialReference ecefSpatialReference = SpatialReference::FromDescription("EPSG:4978");
template <TransformDirection TDirectionOfMapping, class TScalarType, unsigned int NInputDimensions, unsigned int NOutputDimensions>
typename GeocentricTransform<TDirectionOfMapping, TScalarType, NInputDimensions, NOutputDimensions>::OutputPointType
GeocentricTransform<TDirectionOfMapping, TScalarType, NInputDimensions, NOutputDimensions>::TransformPoint(const InputPointType& point) const
{
OutputPointType outputPoint;
#if GDAL_VERSION_NUM >= 3000000
wgs84.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
ecefSpatialReference.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
#endif
if (DirectionOfMapping == TransformDirection::INVERSE)
{
m_Ellipsoid->XYZToLonLatHeight(point[0], point[1], point[2], outputPoint[0], outputPoint[1], outputPoint[2]);
std::unique_ptr<CoordinateTransformation> newMapProjection(new CoordinateTransformation(ecefSpatialReference, wgs84));
if (newMapProjection)
m_MapProjection.swap(newMapProjection);
}
if (DirectionOfMapping == TransformDirection::FORWARD)
else
{
m_Ellipsoid->LonLatHeightToXYZ(point[0], point[1], point[2], outputPoint[0], outputPoint[1], outputPoint[2]);
std::unique_ptr<CoordinateTransformation> newMapProjection(new CoordinateTransformation(wgs84, ecefSpatialReference));
if (newMapProjection)
m_MapProjection.swap(newMapProjection);
}
// To be completed
}
template <TransformDirection TDirectionOfMapping, class TScalarType>
typename GeocentricTransform<TDirectionOfMapping, TScalarType>::OutputPointType
GeocentricTransform<TDirectionOfMapping, TScalarType>::TransformPoint(const InputPointType& point) const
{
OutputPointType outputPoint;
std::tie(outputPoint[0], outputPoint[1], outputPoint[2]) = m_MapProjection->Transform(std::make_tuple(point[0], point[1], point[2]));
return outputPoint;
}
......
......@@ -20,8 +20,7 @@
#include "otbSarSensorModel.h"
#include "otbGenericMapProjection.h"
#include "otbGeocentricTransform.h"
#include "otbDEMHandler.h"
#include <numeric>
......@@ -30,18 +29,14 @@ namespace otb
{
itk::Point<double, 3> EcefToWorld(const itk::Point<double, 3> & ecefPoint)
{
auto transform = otb::GenericMapProjection<otb::TransformDirection::INVERSE, double,3, 3>::New();
transform->SetWkt("EPSG:4978");
auto transform = otb::GeocentricTransform<otb::TransformDirection::INVERSE, double>::New();
return transform->TransformPoint(ecefPoint);
}
itk::Point<double, 3> WorldToEcef(const itk::Point<double, 3> & worldPoint)
{
auto transform = otb::GenericMapProjection<otb::TransformDirection::FORWARD, double,3, 3>::New();
transform->SetWkt("EPSG:4978");
auto transform = otb::GeocentricTransform<otb::TransformDirection::FORWARD, double>::New();
return transform->TransformPoint(worldPoint);
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment