feat(query): implement HAVERSINE

src/query/functions/src/scalars/geometry.rs

@@ -23,6 +23,7 @@ use databend_common_expression::types::VariantType;
 use databend_common_expression::types::F64;
 use databend_common_expression::vectorize_with_builder_1_arg;
 use databend_common_expression::vectorize_with_builder_2_arg;
+use databend_common_expression::vectorize_with_builder_4_arg;
 use databend_common_expression::FunctionDomain;
 use databend_common_expression::FunctionRegistry;
 use databend_common_io::geometry_format;
@@ -34,6 +35,7 @@ use databend_common_io::GeometryDataType;
 use geo::dimensions::Dimensions;
 use geo::BoundingRect;
 use geo::HasDimensions;
+use geo::HaversineDistance;
 use geo::Point;
 use geo_types::Polygon;
 use geohash::decode_bbox;
@@ -55,6 +57,7 @@ use geozero::ToWkt;
 use jsonb::parse_value;
 use jsonb::to_string;
 use num_traits::AsPrimitive;
+use ordered_float::OrderedFloat;
 
 pub fn register(registry: &mut FunctionRegistry) {
     // aliases
@@ -77,6 +80,16 @@ pub fn register(registry: &mut FunctionRegistry) {
     ]);
 
     // functions
+    registry.register_passthrough_nullable_4_arg::<NumberType<F64>, NumberType<F64>, NumberType<F64>, NumberType<F64>, NumberType<F64>, _, _>(
+        "haversine",
+        |_, _, _, _, _| FunctionDomain::Full,
+        vectorize_with_builder_4_arg::<NumberType<F64>, NumberType<F64>, NumberType<F64>, NumberType<F64>, NumberType<F64>,>(|lat1, lon1, lat2, lon2, builder, _| {
+            let p1 = Point::new(lon1, lat1);
+            let p2 = Point::new(lon2, lat2);
+            builder.push(p1.haversine_distance(&p2) * 0.001);
+        }),
+    );
+
     registry.register_passthrough_nullable_1_arg::<GeometryType, VariantType, _, _>(
         "st_asgeojson",
         |_, _| FunctionDomain::MayThrow,

src/query/functions/tests/it/scalars/geometry.rs

@@ -26,6 +26,7 @@ use crate::scalars::run_ast;
 fn test_geometry() {
     let mut mint = Mint::new("tests/it/scalars/testdata");
     let file = &mut mint.new_goldenfile("geometry.txt").unwrap();
+    test_haversine(file);
     test_st_asewkb(file);
     test_st_aswkb(file);
     test_st_asewkt(file);
@@ -56,6 +57,11 @@ fn test_geometry() {
     // test_st_transform(file);
 }
 
+fn test_haversine(file: &mut impl Write) {
+    run_ast(file, "haversine(40.7127, -74.0059, 34.0500, -118.2500)", &[
+    ]);
+}
+
 fn test_st_asewkb(file: &mut impl Write) {
     run_ast(
         file,

src/query/functions/tests/it/scalars/testdata/function_list.txt

@@ -2187,6 +2187,8 @@ Functions overloads:
 1 h3_to_string(UInt64 NULL) :: String NULL
 0 h3_unidirectional_edge_is_valid(UInt64) :: Boolean
 1 h3_unidirectional_edge_is_valid(UInt64 NULL) :: Boolean NULL
+0 haversine(Float64, Float64, Float64, Float64) :: Float64
+1 haversine(Float64 NULL, Float64 NULL, Float64 NULL, Float64 NULL) :: Float64 NULL
 0 humanize_number(Float64) :: String
 1 humanize_number(Float64 NULL) :: String NULL
 0 humanize_size(Float64) :: String

src/query/functions/tests/it/scalars/testdata/geometry.txt

@@ -1,3 +1,12 @@
+ast            : haversine(40.7127, -74.0059, 34.0500, -118.2500)
+raw expr       : haversine(40.7127, minus(74.0059), 34.0500, minus(118.2500))
+checked expr   : haversine<Float64, Float64, Float64, Float64>(to_float64<Decimal(6, 4)>(40.7127_d128(6,4)), to_float64<Decimal(6, 4)>(minus<Decimal(6, 4)>(74.0059_d128(6,4))), to_float64<Decimal(6, 4)>(34.0500_d128(6,4)), to_float64<Decimal(7, 4)>(minus<Decimal(7, 4)>(118.2500_d128(7,4))))
+optimized expr : 4911678.33444826_f64
+output type    : Float64
+output domain  : {4911678.33444826..=4911678.33444826}
+output         : 4911678.33444826
+
+
 ast            : st_asewkb(to_geometry('SRID=4326;POINT(-122.35 37.55)'))
 raw expr       : st_asewkb(to_geometry('SRID=4326;POINT(-122.35 37.55)'))
 checked expr   : st_asewkb<Geometry>(to_geometry<String>("SRID=4326;POINT(-122.35 37.55)"))

tests/sqllogictests/suites/query/functions/02_0060_function_geometry.test

@@ -1,6 +1,11 @@
 statement ok
 DROP TABLE IF EXISTS t1
 
+query T
+SELECT haversine(40.7127, -74.0059, 34.0500, -118.2500);
+----
+3936.3905335559616
+
 statement ok
 SET geometry_output_format='EWKT'