[PD] haversine formula in Pd
Claude Heiland-Allen
claude at mathr.co.uk
Sun Jun 7 05:27:24 CEST 2015
On 07/06/15 03:48, Max wrote:
> In the test case the Pd implementation is 8917.74 km off the proper
> result (2887.26). However I need a precision of about 1m.
Circumference of the Earth in meters: 40,075,000
Accuracy of single precision (24bit): 16,777,216
So your input values will already be inaccurate before you even start
messing with rounding errors etc.
> So I assume the haversine formula is not implementable in Pd at all?
> (unless double precision will be there that is)
I guess so.
> Or is there a workaround?
Not easy at all, but:
Maybe use two floats at different scales 'a' (close to the true value),
'b' (the residual difference from the true value) to express each
coordinate 'c', where:
c = a + b
This would give around 48bits, which should be enough. Actually
performing the addition would give just 'a', due to limited precision.
You have to work with the unpleasant properties of floating point
numbers, like "(a + b) - b" not always being equal to "a".
There are a few blog posts out there about using it on GPU with GLSL
("emulated double", "double-single"), and there is the QD package for
"double-double" and "quad-double" in C++ and FORTRAN (maybe with C
wrappers), there's a Haskell library called "compensated" which might
give some ideas also.
You'll probably also need to do some algebraic manipulations with
trigonometric identities etc to get accurate results.
See:
<https://www.thasler.com/blog/blog/glsl-part2-emu>
<http://crd-legacy.lbl.gov/~dhbailey/mpdist/>
<http://hackage.haskell.org/package/compensated>
<https://en.wikipedia.org/wiki/List_of_trigonometric_identities>
Claude
--
http://mathr.co.uk
More information about the Pd-list
mailing list