gvsafronov/futriix
1
0
Fork 0
Code Issues Pull Requests Actions 4 Packages Projects Releases Wiki Activity

GEOSEARCH BYBOX: Simplified haversine distance formula when longitude diff is 0 (#11579)

Browse Source 
This is take 2 of `GEOSEARCH BYBOX` optimizations based on haversine
distance formula when longitude diff is 0.
The first one was in #11535 . 

- Given longitude diff is 0 the asin(sqrt(a)) on the haversine is asin(sin(abs(u))).
- arcsin(sin(x)) equal to x when x ∈[−𝜋/2,𝜋/2]. 
- Given latitude is between [−𝜋/2,𝜋/2] we can simplifiy arcsin(sin(x)) to x.

On the sample dataset with 60M datapoints, we've measured 55% increase
in the achievable ops/sec.
This commit is contained in:
filipe oliveira 2022-12-05 13:45:04 +00:00 committed by GitHub
parent 2d80cd7840
commit e48ac075c0
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
1 changed files with 18 additions and 10 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

28
src/geohash_helper.c
View File

@ -216,20 +216,28 @@ GeoHashFix52Bits geohashAlign52Bits(const GeoHashBits hash) {
return bits;
}
/* Calculate distance using simplified haversine great circle distance formula.
* Given longitude diff is 0 the asin(sqrt(a)) on the haversine is asin(sin(abs(u))).
* arcsin(sin(x)) equal to x when x ∈[−𝜋/2,𝜋/2]. Given latitude is between [−𝜋/2,𝜋/2]
* we can simplify arcsin(sin(x)) to x.
*/
double geohashGetLatDistance(double lat1d, double lat2d) {
return EARTH_RADIUS_IN_METERS * fabs(deg_rad(lat2d) - deg_rad(lat1d));
}
/* Calculate distance using haversine great circle distance formula. */
double geohashGetDistance(double lon1d, double lat1d, double lon2d, double lat2d) {
double lat1r, lon1r, lat2r, lon2r, u, v;
lat1r = deg_rad(lat1d);
double lat1r, lon1r, lat2r, lon2r, u, v, a;
lon1r = deg_rad(lon1d);
lat2r = deg_rad(lat2d);
lon2r = deg_rad(lon2d);
u = sin((lat2r - lat1r) / 2);
v = sin((lon2r - lon1r) / 2);
double a = u * u;
/* if v == 0 we can avoid doing expensive math */
if (v != 0.0){
a += cos(lat1r) * cos(lat2r) * v * v;
}
/* if v == 0 we can avoid doing expensive math when lons are practically the same */
if (v == 0.0)
return geohashGetLatDistance(lat1d, lat2d);
lat1r = deg_rad(lat1d);
lat2r = deg_rad(lat2d);
u = sin((lat2r - lat1r) / 2);
a = u * u + cos(lat1r) * cos(lat2r) * v * v;
return 2.0 * EARTH_RADIUS_IN_METERS * asin(sqrt(a));
}
@ -259,7 +267,7 @@ int geohashGetDistanceIfInRectangle(double width_m, double height_m, double x1,
double x2, double y2, double *distance) {
/* latitude distance is less expensive to compute than longitude distance
* so we check first for the latitude condition */
double lat_distance = geohashGetDistance(x2, y2, x2, y1);
double lat_distance = geohashGetLatDistance(y2, y1);
if (lat_distance > height_m/2) {
return 0;
}