Writing the Inverse Projection
Before writing the inverse function for the Mercator projection, you should be familiar with the IlvProjection::inverse function.
The IlvProjection::inverse Function
The inverse function prepares the data for inversion and processes the data for appropriate offset.
This function:
suppresses the offset produced by the Cartesian coordinates and converts these coordinates to meters.
reverts the coordinates to their geographic values and applies them to a standard ellipsoid with a semi-major axis of value 1.
calls the function
sInverse or
eInverse depending on whether the ellipsoid is a sphere or not.
adds the value of the central meridian to the longitude and adjusts the longitude to the range
[-PI;PI] if longitude reduction is used (the default value).
converts the latitude if the coordinates are geocentric.
Inverse Projection onto a Sphere
The inverse projection onto a sphere is performed via the sInverse function.
It is not necessary for you to implement the sInverse function. If you call the IlvProjection::inverse function for a projection that does not support the inverse function, the error code IlvProjection::UnsupportedFeatureError() is returned.
As we saw earlier with the sForward function, the projection does not modify the x value. Therefore, the inverse equation is applied only to the y value.
IlvMapsError
MercatorProjection::sInverse(IlvCoordinate& xy) const
{
xy.setY(IlvMaps::Pi()/2. - 2.*atan(exp(-xy.y())));
return IlvMaps::NoError();
}
Inverse Projection onto an Ellipsoid
The inverse projection onto an ellipsoid is performed via the eInverse function. This function assumes that the value of the semi-major axis of the ellipsoid is 1.
In the particular case of the Mercator projection, the implementation of this function is more complex for an ellipsoid than for a sphere. It requires iterations and might fail since there is no simple analytical inverse equation for the Mercator projection from a nonspherical ellipsoid.
IlvMapsError
Mercator::eInverse(IlvCoordinate& xy) const
{
IlvDouble ts = exp(-xy.y());
IlvDouble e = sqrt(getEllipsoid()->getES());
IlvDouble eccnth = 0.5 * e;
IlvDouble Phi = IlvMaps::Pi() / 2. -2. * atan(ts);
int i = 15;
IlvDoubledphi;
do {
IlvDouble con = e * sin(Phi);
dphi = IlvMaps::Pi()/2. -
2. * atan(ts * pow((1 - con)/(1 + con), eccnth)) - Phi;
Phi += dphi;
} while(fabs(dphi) > 1.e-10 && --i != 0);
if(i <= 0)
return ToleranceConditionError();
xy.setY(Phi);
return IlvMaps::NoError();
}
Published date: 05/24/2022
Last modified date: 02/24/2022