[PROJ] Geodetic to Authalic latitude conversions

[Even Rouault]

> If we wanted to be fully precised, for values close to 1, we could 
> actually evaluate sqrt(delta * (2 - delta)) with delta = 1 - something
>
> then
>
> >>> '%.17g' % math.sqrt(1e-7 * (2 - 1e-7))
> '0.00044721358431961788'
>
> which matches quite perfectly the full precision result.
>
Well, revisiting that, the above would be correct if "delta = 1 - 
something" would be evaluated with full precision... But in reality, it 
doesn't, so:

if something = 1 - 1e-7,

delta = 1 - something = 9.9999999947364415e-08

and

'%.17g' % math.sqrt(9.9999999947364415e-08 * (2 - 9.9999999947364415e-08))

= 0.00044721358420192118

which is not better than using the naive sqrt(1 - x * x)

So my pj_cos_of_asin(x) could simplify be simplified as sqrt(1 - x * x), 
possibly with being tolerant for input values in [1, 1 + 
EPSILON_FOR_SLIGHTLY_OUT_OF_DOMAIN] or [-1 - 
EPSILON_FOR_SLIGHTLY_OUT_OF_DOMAIN, -1] as we do in aasin() (in 
src/aasincos.cpp), although I don't know if the use cases where authalic 
lat --> geographic lat is involved we can get those slightly out of 
range values, due to rounding errors in previous computations.

-- 
http://www.spatialys.com
My software is free, but my time generally not.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/proj/attachments/20240913/cc2efd2a/attachment.htm>