1 // Copyright (c) 2017, 2020, Oracle and/or its affiliates. All rights reserved.
2 //
3 // This program is free software; you can redistribute it and/or modify
4 // it under the terms of the GNU General Public License, version 2.0,
5 // as published by the Free Software Foundation.
6 //
7 // This program is also distributed with certain software (including
8 // but not limited to OpenSSL) that is licensed under separate terms,
9 // as designated in a particular file or component or in included license
10 // documentation.  The authors of MySQL hereby grant you an additional
11 // permission to link the program and your derivative works with the
12 // separately licensed software that they have included with MySQL.
13 //
14 // This program is distributed in the hope that it will be useful,
15 // but WITHOUT ANY WARRANTY; without even the implied warranty of
16 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17 // GNU General Public License, version 2.0, for more details.
18 //
19 // You should have received a copy of the GNU General Public License
20 // along with this program; if not, write to the Free Software
21 // Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301  USA
22 
23 /// @file
24 ///
25 /// This file implements the set of functions that storage engines can call to
26 /// do geometrical operations.
27 
28 #include "sql/gis/rtree_support.h"
29 
30 #include <algorithm>  // std::min, std::max
31 #include <cmath>      // std::isfinite, std::isinf, std::isnan
32 #include <limits>
33 
34 #include <boost/geometry.hpp>
35 
36 #include "my_byteorder.h"  // doubleget, float8get
37 #include "my_inttypes.h"   // uchar
38 #include "sql/current_thd.h"
39 #include "sql/dd/cache/dictionary_client.h"
40 #include "sql/dd/types/spatial_reference_system.h"
41 #include "sql/gis/box.h"
42 #include "sql/gis/box_traits.h"
43 #include "sql/gis/covered_by_functor.h"
44 #include "sql/gis/disjoint_functor.h"
45 #include "sql/gis/equals_functor.h"
46 #include "sql/gis/geometries.h"
47 #include "sql/gis/geometries_cs.h"
48 #include "sql/gis/intersects_functor.h"
49 #include "sql/gis/mbr_utils.h"
50 #include "sql/gis/srid.h"
51 #include "sql/gis/wkb.h"
52 #include "sql/spatial.h"    // SRID_SIZE
53 #include "sql/sql_class.h"  // THD
54 #include "sql/srs_fetcher.h"
55 #include "template_utils.h"  // pointer_cast
56 
57 namespace bg = boost::geometry;
58 
fetch_srs(gis::srid_t srid)59 dd::Spatial_reference_system *fetch_srs(gis::srid_t srid) {
60   const dd::Spatial_reference_system *srs = nullptr;
61   dd::cache::Dictionary_client::Auto_releaser m_releaser(
62       current_thd->dd_client());
63   Srs_fetcher fetcher(current_thd);
64   if (srid != 0 && fetcher.acquire(srid, &srs)) return nullptr;
65 
66   if (srs)
67     return srs->clone();
68   else
69     return nullptr;
70 }
71 
mbr_contain_cmp(const dd::Spatial_reference_system * srs,rtr_mbr_t * a,rtr_mbr_t * b)72 bool mbr_contain_cmp(const dd::Spatial_reference_system *srs, rtr_mbr_t *a,
73                      rtr_mbr_t *b) {
74   DBUG_ASSERT(a->xmin <= a->xmax && a->ymin <= a->ymax);
75   DBUG_ASSERT(b->xmin <= b->xmax && b->ymin <= b->ymax);
76 
77   bool result = false;
78   try {
79     gis::Covered_by covered_by(srs ? srs->semi_major_axis() : 0.0,
80                                srs ? srs->semi_minor_axis() : 0.0);
81     if (srs == nullptr || srs->is_cartesian()) {
82       gis::Cartesian_box a_box(gis::Cartesian_point(a->xmin, a->ymin),
83                                gis::Cartesian_point(a->xmax, a->ymax));
84       gis::Cartesian_box b_box(gis::Cartesian_point(b->xmin, b->ymin),
85                                gis::Cartesian_point(b->xmax, b->ymax));
86       result = covered_by(&b_box, &a_box);
87     } else {
88       DBUG_ASSERT(srs->is_geographic());
89       gis::Geographic_box a_box(
90           gis::Geographic_point(srs->to_radians(a->xmin),
91                                 srs->to_radians(a->ymin)),
92           gis::Geographic_point(srs->to_radians(a->xmax),
93                                 srs->to_radians(a->ymax)));
94       gis::Geographic_box b_box(
95           gis::Geographic_point(srs->to_radians(b->xmin),
96                                 srs->to_radians(b->ymin)),
97           gis::Geographic_point(srs->to_radians(b->xmax),
98                                 srs->to_radians(b->ymax)));
99       result = covered_by(&b_box, &a_box);
100     }
101   } catch (...) {
102     DBUG_ASSERT(false); /* purecov: inspected */
103   }
104 
105   return result;
106 }
107 
mbr_equal_cmp(const dd::Spatial_reference_system * srs,rtr_mbr_t * a,rtr_mbr_t * b)108 bool mbr_equal_cmp(const dd::Spatial_reference_system *srs, rtr_mbr_t *a,
109                    rtr_mbr_t *b) {
110   // These points should not have initialized values at this point,
111   // which are min == DBL_MAX and max == -DBL_MAX.
112   DBUG_ASSERT(a->xmin <= a->xmax && a->ymin <= a->ymax);
113   DBUG_ASSERT(b->xmin <= b->xmax && b->ymin <= b->ymax);
114 
115   bool result = false;
116   try {
117     gis::Equals equals(srs ? srs->semi_major_axis() : 0.0,
118                        srs ? srs->semi_minor_axis() : 0.0);
119     if (srs == nullptr || srs->is_cartesian()) {
120       gis::Cartesian_box a_box(gis::Cartesian_point(a->xmin, a->ymin),
121                                gis::Cartesian_point(a->xmax, a->ymax));
122       gis::Cartesian_box b_box(gis::Cartesian_point(b->xmin, b->ymin),
123                                gis::Cartesian_point(b->xmax, b->ymax));
124       result = equals(&a_box, &b_box);
125     } else {
126       DBUG_ASSERT(srs->is_geographic());
127       gis::Geographic_box a_box(
128           gis::Geographic_point(srs->to_radians(a->xmin),
129                                 srs->to_radians(a->ymin)),
130           gis::Geographic_point(srs->to_radians(a->xmax),
131                                 srs->to_radians(a->ymax)));
132       gis::Geographic_box b_box(
133           gis::Geographic_point(srs->to_radians(b->xmin),
134                                 srs->to_radians(b->ymin)),
135           gis::Geographic_point(srs->to_radians(b->xmax),
136                                 srs->to_radians(b->ymax)));
137       result = equals(&a_box, &b_box);
138     }
139   } catch (...) {
140     DBUG_ASSERT(false); /* purecov: inspected */
141   }
142 
143   return result;
144 }
145 
mbr_intersect_cmp(const dd::Spatial_reference_system * srs,rtr_mbr_t * a,rtr_mbr_t * b)146 bool mbr_intersect_cmp(const dd::Spatial_reference_system *srs, rtr_mbr_t *a,
147                        rtr_mbr_t *b) {
148   try {
149     gis::Intersects intersects(srs ? srs->semi_major_axis() : 0.0,
150                                srs ? srs->semi_minor_axis() : 0.0);
151     if (srs == nullptr || srs->is_cartesian()) {
152       gis::Cartesian_box a_box(gis::Cartesian_point(a->xmin, a->ymin),
153                                gis::Cartesian_point(a->xmax, a->ymax));
154       gis::Cartesian_box b_box(gis::Cartesian_point(b->xmin, b->ymin),
155                                gis::Cartesian_point(b->xmax, b->ymax));
156       return intersects(&a_box, &b_box);
157     } else {
158       DBUG_ASSERT(srs->is_geographic());
159       gis::Geographic_box a_box(
160           gis::Geographic_point(srs->to_radians(a->xmin),
161                                 srs->to_radians(a->ymin)),
162           gis::Geographic_point(srs->to_radians(a->xmax),
163                                 srs->to_radians(a->ymax)));
164       gis::Geographic_box b_box(
165           gis::Geographic_point(srs->to_radians(b->xmin),
166                                 srs->to_radians(b->ymin)),
167           gis::Geographic_point(srs->to_radians(b->xmax),
168                                 srs->to_radians(b->ymax)));
169       return intersects(&a_box, &b_box);
170     }
171   } catch (...) {
172     assert(false); /* purecov: inspected */
173   }
174   return false; /* purecov: dead code */
175 }
176 
mbr_disjoint_cmp(const dd::Spatial_reference_system * srs,rtr_mbr_t * a,rtr_mbr_t * b)177 bool mbr_disjoint_cmp(const dd::Spatial_reference_system *srs, rtr_mbr_t *a,
178                       rtr_mbr_t *b) {
179   try {
180     gis::Disjoint disjoint(srs ? srs->semi_major_axis() : 0.0,
181                            srs ? srs->semi_minor_axis() : 0.0);
182     if (srs == nullptr || srs->is_cartesian()) {
183       gis::Cartesian_box a_box(gis::Cartesian_point(a->xmin, a->ymin),
184                                gis::Cartesian_point(a->xmax, a->ymax));
185       gis::Cartesian_box b_box(gis::Cartesian_point(b->xmin, b->ymin),
186                                gis::Cartesian_point(b->xmax, b->ymax));
187       return disjoint(&a_box, &b_box);
188     } else {
189       DBUG_ASSERT(srs->is_geographic());
190       gis::Geographic_box a_box(
191           gis::Geographic_point(srs->to_radians(a->xmin),
192                                 srs->to_radians(a->ymin)),
193           gis::Geographic_point(srs->to_radians(a->xmax),
194                                 srs->to_radians(a->ymax)));
195       gis::Geographic_box b_box(
196           gis::Geographic_point(srs->to_radians(b->xmin),
197                                 srs->to_radians(b->ymin)),
198           gis::Geographic_point(srs->to_radians(b->xmax),
199                                 srs->to_radians(b->ymax)));
200       return disjoint(&a_box, &b_box);
201     }
202   } catch (...) {
203     assert(false); /* purecov: inspected */
204   }
205   return false; /* purecov: dead code */
206 }
207 
mbr_within_cmp(const dd::Spatial_reference_system * srs,rtr_mbr_t * a,rtr_mbr_t * b)208 bool mbr_within_cmp(const dd::Spatial_reference_system *srs, rtr_mbr_t *a,
209                     rtr_mbr_t *b) {
210   bool result = false;
211   try {
212     // If min and max coordinates have been reversed, InnoDB expects the result
213     // to be inverse too. But not if a and b have the exact same coordinates.
214     bool invert = false;
215     if (a->xmin > a->xmax && a->ymin > a->ymax &&
216         !(a->xmin == b->xmin && a->ymin == b->ymin && a->xmax == b->xmax &&
217           a->ymax == b->ymax)) {
218       invert = true;
219     }
220 
221     // Correct the min and max corners to generate proper boxes.
222     double a_xmin = std::min(a->xmin, a->xmax);
223     double a_ymin = std::min(a->ymin, a->ymax);
224     double a_xmax = std::max(a->xmin, a->xmax);
225     double a_ymax = std::max(a->ymin, a->ymax);
226     double b_xmin = std::min(b->xmin, b->xmax);
227     double b_ymin = std::min(b->ymin, b->ymax);
228     double b_xmax = std::max(b->xmin, b->xmax);
229     double b_ymax = std::max(b->ymin, b->ymax);
230 
231     gis::Covered_by covered_by(srs ? srs->semi_major_axis() : 0.0,
232                                srs ? srs->semi_minor_axis() : 0.0);
233     if (srs == nullptr || srs->is_cartesian()) {
234       gis::Cartesian_box a_box(gis::Cartesian_point(a_xmin, a_ymin),
235                                gis::Cartesian_point(a_xmax, a_ymax));
236       gis::Cartesian_box b_box(gis::Cartesian_point(b_xmin, b_ymin),
237                                gis::Cartesian_point(b_xmax, b_ymax));
238       result = covered_by(&a_box, &b_box);
239     } else {
240       DBUG_ASSERT(srs->is_geographic());
241       gis::Geographic_box a_box(gis::Geographic_point(srs->to_radians(a_xmin),
242                                                       srs->to_radians(a_ymin)),
243                                 gis::Geographic_point(srs->to_radians(a_xmax),
244                                                       srs->to_radians(a_ymax)));
245       gis::Geographic_box b_box(gis::Geographic_point(srs->to_radians(b_xmin),
246                                                       srs->to_radians(b_ymin)),
247                                 gis::Geographic_point(srs->to_radians(b_xmax),
248                                                       srs->to_radians(b_ymax)));
249       result = covered_by(&a_box, &b_box);
250     }
251     if (invert) result = !result;
252   } catch (...) {
253     DBUG_ASSERT(false); /* purecov: inspected */
254   }
255 
256   return result;
257 }
258 
mbr_join(const dd::Spatial_reference_system * srs,double * a,const double * b,int n_dim MY_ATTRIBUTE ((unused)))259 void mbr_join(const dd::Spatial_reference_system *srs, double *a,
260               const double *b, int n_dim MY_ATTRIBUTE((unused))) {
261   DBUG_ASSERT(n_dim == 2);
262 
263   try {
264     if (srs == nullptr || srs->is_cartesian()) {
265       gis::Cartesian_box a_box(gis::Cartesian_point(a[0], a[2]),
266                                gis::Cartesian_point(a[1], a[3]));
267       gis::Cartesian_box b_box(gis::Cartesian_point(b[0], b[2]),
268                                gis::Cartesian_point(b[1], b[3]));
269       bg::expand(a_box, b_box);
270       a[0] = a_box.min_corner().x();
271       a[1] = a_box.max_corner().x();
272       a[2] = a_box.min_corner().y();
273       a[3] = a_box.max_corner().y();
274     } else {
275       DBUG_ASSERT(srs->is_geographic());
276       gis::Geographic_box a_box(
277           gis::Geographic_point(srs->to_radians(a[0]), srs->to_radians(a[2])),
278           gis::Geographic_point(srs->to_radians(a[1]), srs->to_radians(a[3])));
279       gis::Geographic_box b_box(
280           gis::Geographic_point(srs->to_radians(b[0]), srs->to_radians(b[2])),
281           gis::Geographic_point(srs->to_radians(b[1]), srs->to_radians(b[3])));
282       bg::expand(a_box, b_box);
283       a[0] = srs->from_radians(a_box.min_corner().x());
284       a[1] = srs->from_radians(a_box.max_corner().x());
285       a[2] = srs->from_radians(a_box.min_corner().y());
286       a[3] = srs->from_radians(a_box.max_corner().y());
287     }
288   } catch (...) {
289     DBUG_ASSERT(false); /* purecov: inspected */
290   }
291 }
292 
mbr_join_area(const dd::Spatial_reference_system * srs,const double * a,const double * b,int n_dim MY_ATTRIBUTE ((unused)))293 double mbr_join_area(const dd::Spatial_reference_system *srs, const double *a,
294                      const double *b, int n_dim MY_ATTRIBUTE((unused))) {
295   DBUG_ASSERT(n_dim == 2);
296 
297   double area = 0.0;
298   try {
299     if (srs == nullptr || srs->is_cartesian()) {
300       gis::Cartesian_box a_box(gis::Cartesian_point(a[0], a[2]),
301                                gis::Cartesian_point(a[1], a[3]));
302       gis::Cartesian_box b_box(gis::Cartesian_point(b[0], b[2]),
303                                gis::Cartesian_point(b[1], b[3]));
304       bg::expand(a_box, b_box);
305       area = bg::area(a_box);
306     } else {
307       DBUG_ASSERT(srs->is_geographic());
308       gis::Geographic_box a_box(
309           gis::Geographic_point(srs->to_radians(a[0]), srs->to_radians(a[2])),
310           gis::Geographic_point(srs->to_radians(a[1]), srs->to_radians(a[3])));
311       gis::Geographic_box b_box(
312           gis::Geographic_point(srs->to_radians(b[0]), srs->to_radians(b[2])),
313           gis::Geographic_point(srs->to_radians(b[1]), srs->to_radians(b[3])));
314       bg::expand(a_box, b_box);
315       area = bg::area(
316           a_box, bg::strategy::area::geographic<>(bg::srs::spheroid<double>(
317                      srs->semi_major_axis(), srs->semi_minor_axis())));
318     }
319   } catch (...) {
320     DBUG_ASSERT(false); /* purecov: inspected */
321   }
322 
323   if (!std::isfinite(area)) area = std::numeric_limits<double>::max();
324   return area;
325 }
326 
compute_area(const dd::Spatial_reference_system * srs,const double * a,int n_dim MY_ATTRIBUTE ((unused)))327 double compute_area(const dd::Spatial_reference_system *srs, const double *a,
328                     int n_dim MY_ATTRIBUTE((unused))) {
329   DBUG_ASSERT(n_dim == 2);
330 
331   double area = 0.0;
332   try {
333     if (srs == nullptr || srs->is_cartesian()) {
334       gis::Cartesian_box a_box(gis::Cartesian_point(a[0], a[2]),
335                                gis::Cartesian_point(a[1], a[3]));
336       area = bg::area(a_box);
337     } else {
338       DBUG_ASSERT(srs->is_geographic());
339       gis::Geographic_box a_box(
340           gis::Geographic_point(srs->to_radians(a[0]), srs->to_radians(a[2])),
341           gis::Geographic_point(srs->to_radians(a[1]), srs->to_radians(a[3])));
342       area = bg::area(
343           a_box, bg::strategy::area::geographic<>(bg::srs::spheroid<double>(
344                      srs->semi_major_axis(), srs->semi_minor_axis())));
345     }
346   } catch (...) {
347     DBUG_ASSERT(false); /* purecov: inspected */
348   }
349 
350   return area;
351 }
352 
get_mbr_from_store(const dd::Spatial_reference_system * srs,const uchar * store,uint size,uint n_dims MY_ATTRIBUTE ((unused)),double * mbr,gis::srid_t * srid)353 int get_mbr_from_store(const dd::Spatial_reference_system *srs,
354                        const uchar *store, uint size,
355                        uint n_dims MY_ATTRIBUTE((unused)), double *mbr,
356                        gis::srid_t *srid) {
357   DBUG_ASSERT(n_dims == 2);
358   // The SRS should match the SRID of the geometry, with one exception: For
359   // backwards compatibility it is allowed to create indexes with mixed
360   // SRIDs. Although these indexes can never be used to optimize queries, the
361   // user is allowed to create them. These indexes will call get_mbr_from_store
362   // with srs == nullptr. There is, unfortunately, no way to differentiate mixed
363   // SRID indexes from SRID 0 indexes here, so the assertion is not perfect.
364   DBUG_ASSERT(srs == nullptr || (srs->id() == uint4korr(store)));
365 
366   if (srid != nullptr) *srid = uint4korr(store);
367 
368   try {
369     // Note: current_thd may be nullptr here if this function was called from an
370     // internal InnoDB thread. In that case, we won't get any stack size check
371     // in gis::parse_wkb, but the geometry has been parsed before with the stack
372     // size check enabled. We assume we have at least the same amount of stack
373     // when called from an internal thread as when called from a MySQL thread.
374     std::unique_ptr<gis::Geometry> g =
375         gis::parse_wkb(current_thd, srs,
376                        pointer_cast<const char *>(store) + sizeof(gis::srid_t),
377                        size - sizeof(gis::srid_t), true);
378     if (g.get() == nullptr) {
379       return -1; /* purecov: inspected */
380     }
381     if (srs == nullptr || srs->is_cartesian()) {
382       gis::Cartesian_box box;
383       gis::box_envelope(g.get(), srs, &box);
384       mbr[0] = box.min_corner().x();
385       mbr[1] = box.max_corner().x();
386       mbr[2] = box.min_corner().y();
387       mbr[3] = box.max_corner().y();
388     } else {
389       DBUG_ASSERT(srs->is_geographic());
390       gis::Geographic_box box;
391       gis::box_envelope(g.get(), srs, &box);
392       mbr[0] = srs->from_radians(box.min_corner().x());
393       mbr[1] = srs->from_radians(box.max_corner().x());
394       mbr[2] = srs->from_radians(box.min_corner().y());
395       mbr[3] = srs->from_radians(box.max_corner().y());
396     }
397   } catch (...) {
398     DBUG_ASSERT(false); /* purecov: inspected */
399     return -1;
400   }
401 
402   if (std::isnan(mbr[0])) {
403     /* purecov: begin inspected */
404     DBUG_ASSERT(std::isnan(mbr[1]) && std::isnan(mbr[2]) && std::isnan(mbr[3]));
405     // The geometry is empty, so there is no bounding box. Return a box that
406     // covers the entire domain.
407     mbr[0] = std::numeric_limits<double>::lowest();
408     mbr[1] = std::numeric_limits<double>::max();
409     mbr[2] = std::numeric_limits<double>::lowest();
410     mbr[3] = std::numeric_limits<double>::max();
411     /* purecov: end inspected */
412   }
413 
414   // xmin <= xmax && ymin <= ymax
415   DBUG_ASSERT(mbr[0] <= mbr[1] && mbr[2] <= mbr[3]);
416 
417   return 0;
418 }
419 
rtree_area_increase(const dd::Spatial_reference_system * srs,const uchar * mbr_a,const uchar * mbr_b,int mbr_len MY_ATTRIBUTE ((unused)),double * ab_area)420 double rtree_area_increase(const dd::Spatial_reference_system *srs,
421                            const uchar *mbr_a, const uchar *mbr_b,
422                            int mbr_len MY_ATTRIBUTE((unused)),
423                            double *ab_area) {
424   DBUG_ASSERT(mbr_len == sizeof(double) * 4);
425 
426   double a_xmin = float8get(mbr_a);
427   double a_xmax = float8get(mbr_a + sizeof(double));
428   double a_ymin = float8get(mbr_a + sizeof(double) * 2);
429   double a_ymax = float8get(mbr_a + sizeof(double) * 3);
430   double b_xmin = float8get(mbr_b);
431   double b_xmax = float8get(mbr_b + sizeof(double));
432   double b_ymin = float8get(mbr_b + sizeof(double) * 2);
433   double b_ymax = float8get(mbr_b + sizeof(double) * 3);
434 
435   DBUG_ASSERT(a_xmin <= a_xmax && a_ymin <= a_ymax);
436   DBUG_ASSERT(b_xmin <= b_xmax && b_ymin <= b_ymax);
437 
438   double a_area = 0.0;
439   try {
440     if (srs == nullptr || srs->is_cartesian()) {
441       gis::Cartesian_box a_box(gis::Cartesian_point(a_xmin, a_ymin),
442                                gis::Cartesian_point(a_xmax, a_ymax));
443       gis::Cartesian_box b_box(gis::Cartesian_point(b_xmin, b_ymin),
444                                gis::Cartesian_point(b_xmax, b_ymax));
445       a_area = bg::area(a_box);
446       if (a_area == 0.0) a_area = 0.001 * 0.001;
447       bg::expand(a_box, b_box);
448       *ab_area = bg::area(a_box);
449     } else {
450       DBUG_ASSERT(srs->is_geographic());
451       gis::Geographic_box a_box(gis::Geographic_point(srs->to_radians(a_xmin),
452                                                       srs->to_radians(a_ymin)),
453                                 gis::Geographic_point(srs->to_radians(a_xmax),
454                                                       srs->to_radians(a_ymax)));
455       gis::Geographic_box b_box(gis::Geographic_point(srs->to_radians(b_xmin),
456                                                       srs->to_radians(b_ymin)),
457                                 gis::Geographic_point(srs->to_radians(b_xmax),
458                                                       srs->to_radians(b_ymax)));
459       a_area = bg::area(
460           a_box, bg::strategy::area::geographic<>(bg::srs::spheroid<double>(
461                      srs->semi_major_axis(), srs->semi_minor_axis())));
462       bg::expand(a_box, b_box);
463       *ab_area = bg::area(
464           a_box, bg::strategy::area::geographic<>(bg::srs::spheroid<double>(
465                      srs->semi_major_axis(), srs->semi_minor_axis())));
466     }
467     if (std::isinf(a_area)) a_area = std::numeric_limits<double>::max();
468     if (std::isinf(*ab_area)) *ab_area = std::numeric_limits<double>::max();
469   } catch (...) {
470     DBUG_ASSERT(false); /* purecov: inspected */
471   }
472 
473   DBUG_ASSERT(std::isfinite(*ab_area - a_area));
474   return *ab_area - a_area;
475 }
476 
rtree_area_overlapping(const dd::Spatial_reference_system * srs,const uchar * mbr_a,const uchar * mbr_b,int mbr_len MY_ATTRIBUTE ((unused)))477 double rtree_area_overlapping(const dd::Spatial_reference_system *srs,
478                               const uchar *mbr_a, const uchar *mbr_b,
479                               int mbr_len MY_ATTRIBUTE((unused))) {
480   DBUG_ASSERT(mbr_len == sizeof(double) * 4);
481 
482   double a_xmin = float8get(mbr_a);
483   double a_xmax = float8get(mbr_a + sizeof(double));
484   double a_ymin = float8get(mbr_a + sizeof(double) * 2);
485   double a_ymax = float8get(mbr_a + sizeof(double) * 3);
486   double b_xmin = float8get(mbr_b);
487   double b_xmax = float8get(mbr_b + sizeof(double));
488   double b_ymin = float8get(mbr_b + sizeof(double) * 2);
489   double b_ymax = float8get(mbr_b + sizeof(double) * 3);
490 
491   DBUG_ASSERT(a_xmin <= a_xmax && a_ymin <= a_ymax);
492   DBUG_ASSERT(b_xmin <= b_xmax && b_ymin <= b_ymax);
493 
494   double area = 0.0;
495   try {
496     if (srs == nullptr || srs->is_cartesian()) {
497       gis::Cartesian_box a_box(gis::Cartesian_point(a_xmin, a_ymin),
498                                gis::Cartesian_point(a_xmax, a_ymax));
499       gis::Cartesian_box b_box(gis::Cartesian_point(b_xmin, b_ymin),
500                                gis::Cartesian_point(b_xmax, b_ymax));
501       gis::Cartesian_box overlapping_box;
502       bg::intersection(a_box, b_box, overlapping_box);
503       area = bg::area(overlapping_box);
504     } else {
505       DBUG_ASSERT(srs->is_geographic());
506       gis::Geographic_box a_box(gis::Geographic_point(srs->to_radians(a_xmin),
507                                                       srs->to_radians(a_ymin)),
508                                 gis::Geographic_point(srs->to_radians(a_xmax),
509                                                       srs->to_radians(a_ymax)));
510       gis::Geographic_box b_box(gis::Geographic_point(srs->to_radians(b_xmin),
511                                                       srs->to_radians(b_ymin)),
512                                 gis::Geographic_point(srs->to_radians(b_xmax),
513                                                       srs->to_radians(b_ymax)));
514       gis::Geographic_box overlapping_box;
515       bg::intersection(a_box, b_box, overlapping_box,
516                        bg::strategy::intersection::geographic_segments<>(
517                            bg::srs::spheroid<double>(srs->semi_major_axis(),
518                                                      srs->semi_minor_axis())));
519       area =
520           bg::area(overlapping_box,
521                    bg::strategy::area::geographic<>(bg::srs::spheroid<double>(
522                        srs->semi_major_axis(), srs->semi_minor_axis())));
523     }
524   } catch (...) {
525     DBUG_ASSERT(false); /* purecov: inspected */
526   }
527 
528   if (std::isnan(area)) area = 0.0;
529   return area;
530 }
531