diff options
author | Doug Junkins <junkins@foghead.com> | 2019-04-30 06:53:31 -0700 |
---|---|---|
committer | bstoeger <32835590+bstoeger@users.noreply.github.com> | 2019-04-30 23:32:50 +0200 |
commit | 14a763a6a0148a42fad06717bd2532181a54c626 (patch) | |
tree | d8b71f0c9a2d33d40d60aeccb44e12111abcb2e0 | |
parent | 154792ffd1c217e5c1c9774b5593a753ae553601 (diff) | |
download | subsurface-14a763a6a0148a42fad06717bd2532181a54c626.tar.gz |
Bugfix for algorithm in get_distance()
Fixed bug in the Haversine function in get_distance() based on algorithm
at https://www.movable-type.co.uk/scripts/latlong.html and added bounds
to the 'a' term to avoid floating point errors for antipodal points.
Signed-off-by: Doug Junkins <junkins@foghead.com>
-rw-r--r-- | CHANGELOG.md | 1 | ||||
-rw-r--r-- | core/divesite.c | 7 |
2 files changed, 6 insertions, 2 deletions
diff --git a/CHANGELOG.md b/CHANGELOG.md index cce954eb1..6b3930ac7 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,3 +1,4 @@ +- Core: fix bug in get_distance() to correctly compute spherical distance - Desktop: For videos, add save data export as subtitle file - Desktop: make dive sites 1st class citizens with their own dive site table - Desktop: only show dives at the dive sites selected in dive site table diff --git a/core/divesite.c b/core/divesite.c index b5ce99a8d..45601537c 100644 --- a/core/divesite.c +++ b/core/divesite.c @@ -77,15 +77,18 @@ struct dive_site *get_dive_site_by_gps_and_name(char *name, const location_t *lo // Calculate the distance in meters between two coordinates. unsigned int get_distance(const location_t *loc1, const location_t *loc2) { + double lat1_r = udeg_to_radians(loc1->lat.udeg); double lat2_r = udeg_to_radians(loc2->lat.udeg); double lat_d_r = udeg_to_radians(loc2->lat.udeg - loc1->lat.udeg); double lon_d_r = udeg_to_radians(loc2->lon.udeg - loc1->lon.udeg); double a = sin(lat_d_r/2) * sin(lat_d_r/2) + - cos(lat2_r) * cos(lat2_r) * sin(lon_d_r/2) * sin(lon_d_r/2); + cos(lat1_r) * cos(lat2_r) * sin(lon_d_r/2) * sin(lon_d_r/2); + if (a < 0.0) a = 0.0; + if (a > 1.0) a = 1.0; double c = 2 * atan2(sqrt(a), sqrt(1.0 - a)); - // Earth radious in metres + // Earth radius in metres return lrint(6371000 * c); } |