linrrzqqq commented on code in PR #60170: URL: https://github.com/apache/doris/pull/60170#discussion_r2726287939
########## be/src/geo/st_distance.cpp: ########## @@ -0,0 +1,320 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +#include <s2/s2cap.h> +#include <s2/s2earth.h> +#include <s2/s2loop.h> +#include <s2/s2polygon.h> +#include <s2/s2polyline.h> + +#include <cmath> +#include <limits> + +#include "geo/geo_types.h" + +namespace doris { + +// Helper function to compute distance from a point to a line segment +double distance_point_to_segment(const S2Point& point, const S2Point& line_start, + const S2Point& line_end) { + S2LatLng point_ll = S2LatLng(point); + S2LatLng start_ll = S2LatLng(line_start); + S2LatLng end_ll = S2LatLng(line_end); + + double px = point_ll.lng().degrees(); + double py = point_ll.lat().degrees(); + double x1 = start_ll.lng().degrees(); + double y1 = start_ll.lat().degrees(); + double x2 = end_ll.lng().degrees(); + double y2 = end_ll.lat().degrees(); + + double dx = x2 - x1; + double dy = y2 - y1; + + if (dx == 0 && dy == 0) { + return S2Earth::GetDistanceMeters(point_ll, start_ll); + } + + double t = ((px - x1) * dx + (py - y1) * dy) / (dx * dx + dy * dy); + t = std::max(0.0, std::min(1.0, t)); + + double closest_x = x1 + t * dx; + double closest_y = y1 + t * dy; + + S2LatLng closest_point = S2LatLng::FromDegrees(closest_y, closest_x); + return S2Earth::GetDistanceMeters(point_ll, closest_point); +} + +// Helper function to compute distance from a point to a polyline +double distance_point_to_polyline(const S2Point& point, const S2Polyline* polyline) { + double min_distance = std::numeric_limits<double>::max(); + + for (int i = 0; i < polyline->num_vertices() - 1; ++i) { + const S2Point& p1 = polyline->vertex(i); + const S2Point& p2 = polyline->vertex(i + 1); + + double dist = distance_point_to_segment(point, p1, p2); + min_distance = std::min(min_distance, dist); + } + + return min_distance; +} + +// Helper function to compute distance from a point to a polygon +double distance_point_to_polygon(const S2Point& point, const S2Polygon* polygon) { + // Check if point is inside polygon + if (polygon->Contains(point)) { + return 0.0; + } + + // Find minimum distance to polygon boundary + double min_distance = std::numeric_limits<double>::max(); + + for (int i = 0; i < polygon->num_loops(); ++i) { + const S2Loop* loop = polygon->loop(i); + + for (int j = 0; j < loop->num_vertices(); ++j) { + const S2Point& p1 = loop->vertex(j); + const S2Point& p2 = loop->vertex((j + 1) % loop->num_vertices()); + + double dist = distance_point_to_segment(point, p1, p2); + min_distance = std::min(min_distance, dist); + } + } + + return min_distance; +} + +double GeoPoint::Distance(const GeoShape* rhs) const { + // rhs is guaranteed to be valid by StDualGeoDoubleFunction (functions_geo.cpp) which checks shapes[0] && shapes[1] + switch (rhs->type()) { + case GEO_SHAPE_POINT: { + const GeoPoint* point = static_cast<const GeoPoint*>(rhs); + S2LatLng this_ll = S2LatLng(*_point); + S2LatLng other_ll = S2LatLng(*point->point()); + return S2Earth::GetDistanceMeters(this_ll, other_ll); + } + case GEO_SHAPE_LINE_STRING: { + const GeoLine* line = static_cast<const GeoLine*>(rhs); + return distance_point_to_polyline(*_point, line->polyline()); + } + case GEO_SHAPE_POLYGON: { + const GeoPolygon* polygon = static_cast<const GeoPolygon*>(rhs); + return distance_point_to_polygon(*_point, polygon->polygon()); + } + case GEO_SHAPE_CIRCLE: { + const GeoCircle* circle = static_cast<const GeoCircle*>(rhs); + S2LatLng this_ll = S2LatLng(*_point); + S2LatLng center_ll = S2LatLng(circle->circle()->center()); + double dist_to_center = S2Earth::GetDistanceMeters(this_ll, center_ll); + double circle_radius = S2Earth::ToMeters(circle->circle()->radius()); + + // Distance from point to circle is distance to center minus radius + return std::max(0.0, dist_to_center - circle_radius); + } + case GEO_SHAPE_MULTI_POLYGON: { + return rhs->Distance(this); // Delegate to MultiPolygon's implementation + } + default: + return -1.0; + } +} + +double GeoLine::Distance(const GeoShape* rhs) const { + // rhs is guaranteed to be valid by StDualGeoDoubleFunction (functions_geo.cpp) which checks shapes[0] && shapes[1] + switch (rhs->type()) { + case GEO_SHAPE_POINT: { + const GeoPoint* point = static_cast<const GeoPoint*>(rhs); + return distance_point_to_polyline(*point->point(), _polyline.get()); + } + case GEO_SHAPE_LINE_STRING: { Review Comment: plygon-polygon and circle-line seem to have this issue as well -- This is an automated message from the Apache Git Service. To respond to the message, please log on to GitHub and use the URL above to go to the specific comment. To unsubscribe, e-mail: [email protected] For queries about this service, please contact Infrastructure at: [email protected] --------------------------------------------------------------------- To unsubscribe, e-mail: [email protected] For additional commands, e-mail: [email protected]
