diff --git a/src/Functions/greatCircleDistance.cpp b/src/Functions/greatCircleDistance.cpp index 94f44ba37ba..8b9df91e117 100644 --- a/src/Functions/greatCircleDistance.cpp +++ b/src/Functions/greatCircleDistance.cpp @@ -99,6 +99,12 @@ void geodistInit() } } +inline NO_SANITIZE_UNDEFINED size_t floatToIndex(float x) +{ + /// Implementation specific behaviour on overflow or infinite value. + return static_cast(x); +} + inline float geodistDegDiff(float f) { f = fabsf(f); @@ -110,7 +116,7 @@ inline float geodistDegDiff(float f) inline float geodistFastCos(float x) { float y = fabsf(x) * (COS_LUT_SIZE / PI / 2); - size_t i = static_cast(y); + size_t i = floatToIndex(y); y -= i; i &= (COS_LUT_SIZE - 1); return cos_lut[i] + (cos_lut[i + 1] - cos_lut[i]) * y; @@ -119,7 +125,7 @@ inline float geodistFastCos(float x) inline float geodistFastSin(float x) { float y = fabsf(x) * (COS_LUT_SIZE / PI / 2); - size_t i = static_cast(y); + size_t i = floatToIndex(y); y -= i; i = (i - COS_LUT_SIZE / 4) & (COS_LUT_SIZE - 1); // cos(x - pi / 2) = sin(x), costable / 4 = pi / 2 return cos_lut[i] + (cos_lut[i + 1] - cos_lut[i]) * y; @@ -139,7 +145,7 @@ inline float geodistFastAsinSqrt(float x) { // distance under 17083 km, 512-entry LUT error under 0.00072% x *= ASIN_SQRT_LUT_SIZE; - size_t i = static_cast(x); + size_t i = floatToIndex(x); return asin_sqrt_lut[i] + (asin_sqrt_lut[i + 1] - asin_sqrt_lut[i]) * (x - i); } return asinf(sqrtf(x)); // distance over 17083 km, just compute exact @@ -160,12 +166,6 @@ DECLARE_MULTITARGET_CODE( namespace { -inline NO_SANITIZE_UNDEFINED size_t metricLUTIndex(float latitude_midpoint) -{ - /// Implementation specific behaviour on overflow or infinite value. - return static_cast(latitude_midpoint) & (METRIC_LUT_SIZE - 1); -} - template float distance(float lon1deg, float lat1deg, float lon2deg, float lat2deg) { @@ -183,7 +183,7 @@ float distance(float lon1deg, float lat1deg, float lon2deg, float lat2deg) /// But if longitude is close but latitude is different enough, there is no difference between meridian and great circle line. float latitude_midpoint = (lat1deg + lat2deg + 180) * METRIC_LUT_SIZE / 360; // [-90, 90] degrees -> [0, KTABLE] indexes - size_t latitude_midpoint_index = metricLUTIndex(latitude_midpoint); + size_t latitude_midpoint_index = floatToIndex(latitude_midpoint) & (METRIC_LUT_SIZE - 1); /// This is linear interpolation between two table items at index "latitude_midpoint_index" and "latitude_midpoint_index + 1".