30 #define GEOS_USE_ONLY_R_API 35 using WKB = std::vector<uint8_t>;
37 #define MAX_GEOS_MESSAGE_LEN 200 39 static std::mutex geos_log_info_mutex;
40 static std::mutex geos_log_error_mutex;
43 static void geos_notice_handler(
const char* fmt, ...) {
44 char buffer[MAX_GEOS_MESSAGE_LEN];
47 vsnprintf(buffer, MAX_GEOS_MESSAGE_LEN, fmt, args);
50 std::lock_guard<std::mutex> guard(geos_log_info_mutex);
51 LOG(
INFO) <<
"GEOS Notice: " << std::string(buffer);
56 static void geos_error_handler(
const char* fmt, ...) {
58 char buffer[MAX_GEOS_MESSAGE_LEN];
60 vsnprintf(buffer, MAX_GEOS_MESSAGE_LEN, fmt, args);
63 std::lock_guard<std::mutex> guard(geos_log_error_mutex);
64 LOG(
ERROR) <<
"GEOS Error: " << std::string(buffer);
68 GEOSContextHandle_t create_context() {
69 #if GEOS_VERSION_MAJOR == 3 && GEOS_VERSION_MINOR < 5 70 GEOSContextHandle_t context = initGEOS_r(geos_notice_handler, geos_error_handler);
74 GEOSContextHandle_t context = GEOS_init_r();
76 GEOSContext_setNoticeHandler_r(context, geos_notice_handler);
77 GEOSContext_setErrorHandler_r(context, geos_error_handler);
82 void destroy_context(GEOSContextHandle_t context) {
84 #if GEOS_VERSION_MAJOR == 3 && GEOS_VERSION_MINOR < 5 85 finishGEOS_r(context);
87 GEOS_finish_r(context);
100 int32_t* best_planar_srid_ptr) {
103 if (static_cast<SQLTypes>(type) ==
kPOINT) {
105 if (best_planar_srid_ptr) {
108 *best_planar_srid_ptr = point.getBestPlanarSRID();
109 if (!point.transform(4326, *best_planar_srid_ptr)) {
113 return point.getWkb(wkb);
117 if (best_planar_srid_ptr) {
120 *best_planar_srid_ptr = linestring.getBestPlanarSRID();
121 if (!linestring.transform(4326, *best_planar_srid_ptr)) {
125 return linestring.getWkb(wkb);
127 std::vector<int32_t> meta1v(meta1, meta1 + meta1_size);
128 if (static_cast<SQLTypes>(type) ==
kPOLYGON) {
130 if (best_planar_srid_ptr) {
133 *best_planar_srid_ptr = poly.getBestPlanarSRID();
134 if (!poly.transform(4326, *best_planar_srid_ptr)) {
138 return poly.getWkb(wkb);
140 std::vector<int32_t> meta2v(meta2, meta2 + meta2_size);
145 if (meta1_size == 1 && meta2_size == 1) {
146 const std::vector<double> ecv = {0.0, 0.0, 0.00000012345, 0.0, 0.0, 0.00000012345};
149 return empty.getWkb(wkb);
153 if (best_planar_srid_ptr) {
156 *best_planar_srid_ptr = mpoly.getBestPlanarSRID();
157 if (!mpoly.transform(4326, *best_planar_srid_ptr)) {
161 return mpoly.getWkb(wkb);
168 bool fromWkb(WKB& wkb,
170 int8_t** result_coords,
171 int64_t* result_coords_size,
172 int32_t** result_meta1,
173 int64_t* result_meta1_size,
174 int32_t** result_meta2,
175 int64_t* result_meta2_size,
176 int32_t* best_planar_srid_ptr) {
178 if (best_planar_srid_ptr && !
result->isEmpty()) {
181 if (!
result->transform(*best_planar_srid_ptr, 4326)) {
187 std::vector<double> coords{};
188 std::vector<int32_t> ring_sizes{};
189 std::vector<int32_t> poly_rings{};
190 std::vector<double> bounds{};
197 coords = {0.0, 0.0, 0.00000012345, 0.0, 0.0, 0.00000012345};
198 ring_sizes.push_back(3);
199 poly_rings.push_back(1);
200 }
else if (
auto result_point = dynamic_cast<GeoPoint*>(
result.get())) {
201 result_point->getColumns(coords);
203 coords.push_back(coords[0] + 0.0000001);
204 coords.push_back(coords[1]);
205 coords.push_back(coords[0]);
206 coords.push_back(coords[1] + 0.0000001);
207 ring_sizes.push_back(3);
208 poly_rings.push_back(ring_sizes.size());
209 }
else if (
auto result_poly = dynamic_cast<GeoPolygon*>(
result.get())) {
210 result_poly->getColumns(coords, ring_sizes, bounds);
212 poly_rings.push_back(ring_sizes.size());
213 }
else if (
auto result_mpoly = dynamic_cast<GeoMultiPolygon*>(
result.get())) {
214 result_mpoly->getColumns(coords, ring_sizes, poly_rings, bounds);
224 *result_coords =
nullptr;
225 int64_t size = coords.size() *
sizeof(double);
228 std::memcpy(buf, coords.data(), size);
229 *result_coords =
reinterpret_cast<int8_t*
>(buf);
231 *result_coords_size = size;
233 *result_meta1 =
nullptr;
234 size = ring_sizes.size() *
sizeof(int32_t);
237 std::memcpy(buf, ring_sizes.data(), size);
238 *result_meta1 =
reinterpret_cast<int32_t*
>(buf);
240 *result_meta1_size = ring_sizes.size();
242 *result_meta2 =
nullptr;
243 size = poly_rings.size() *
sizeof(int32_t);
246 std::memcpy(buf, poly_rings.data(), size);
247 *result_meta2 =
reinterpret_cast<int32_t*
>(buf);
249 *result_meta2_size = poly_rings.size();
254 GEOSGeometry* postprocess(GEOSContextHandle_t context, GEOSGeometry* g) {
255 if (g && GEOSisEmpty_r(context, g) == 0) {
256 auto type = GEOSGeomTypeId_r(context, g);
258 if (type != GEOS_POINT && type != GEOS_POLYGON && type != GEOS_MULTIPOLYGON) {
260 double tiny_distance = 0.000000001;
261 auto ng = GEOSBuffer_r(context, g, tiny_distance, quadsegs);
262 GEOSGeom_destroy_r(context, g);
274 int64_t arg1_coords_size,
276 int64_t arg1_meta1_size,
278 int64_t arg1_meta2_size,
284 int64_t arg2_coords_size,
286 int64_t arg2_meta1_size,
288 int64_t arg2_meta2_size,
294 int8_t** result_coords,
295 int64_t* result_coords_size,
296 int32_t** result_meta1,
297 int64_t* result_meta1_size,
298 int32_t** result_meta2,
299 int64_t* result_meta2_size) {
333 auto context = create_context();
337 auto* g1 = GEOSGeomFromWKB_buf_r(context, wkb1.data(), wkb1.size());
339 auto* g2 = GEOSGeomFromWKB_buf_r(context, wkb2.data(), wkb2.size());
341 GEOSGeometry* g =
nullptr;
343 g = GEOSIntersection_r(context, g1, g2);
345 g = GEOSDifference_r(context, g1, g2);
347 g = GEOSUnion_r(context, g1, g2);
349 g = postprocess(context, g);
351 size_t wkb_size = 0ULL;
352 auto wkb_buf = GEOSGeomToWKB_buf_r(context, g, &wkb_size);
353 if (wkb_buf && wkb_size > 0ULL) {
354 WKB wkb(wkb_buf, wkb_buf + wkb_size);
356 status = fromWkb(wkb,
366 GEOSGeom_destroy_r(context, g);
368 GEOSGeom_destroy_r(context, g2);
370 GEOSGeom_destroy_r(context, g1);
372 destroy_context(context);
382 int64_t arg1_coords_size,
384 int64_t arg1_meta1_size,
386 int64_t arg1_meta2_size,
393 int8_t** result_coords,
394 int64_t* result_coords_size,
395 int32_t** result_meta1,
396 int64_t* result_meta1_size,
397 int32_t** result_meta2,
398 int64_t* result_meta2_size) {
400 int32_t best_planar_srid;
401 int32_t* best_planar_srid_ptr =
nullptr;
405 best_planar_srid_ptr = &best_planar_srid;
418 best_planar_srid_ptr)) {
423 auto context = create_context();
427 auto* g1 = GEOSGeomFromWKB_buf_r(context, wkb1.data(), wkb1.size());
429 GEOSGeometry* g =
nullptr;
432 g = GEOSBuffer_r(context, g1, arg2, quadsegs);
434 g = postprocess(context, g);
436 size_t wkb_size = 0ULL;
437 auto wkb_buf = GEOSGeomToWKB_buf_r(context, g, &wkb_size);
438 if (wkb_buf && wkb_size > 0ULL) {
439 WKB wkb(wkb_buf, wkb_buf + wkb_size);
442 status = fromWkb(wkb,
450 best_planar_srid_ptr);
452 GEOSGeom_destroy_r(context, g);
454 GEOSGeom_destroy_r(context, g1);
456 destroy_context(context);
466 int64_t arg_coords_size,
468 int64_t arg_meta1_size,
470 int64_t arg_meta2_size,
477 if (!result || !toWkb(wkb1,
491 auto context = create_context();
495 auto* g1 = GEOSGeomFromWKB_buf_r(context, wkb1.data(), wkb1.size());
498 *result = GEOSisEmpty_r(context, g1);
501 *result = GEOSisValid_r(context, g1);
504 GEOSGeom_destroy_r(context, g1);
506 destroy_context(context);
static std::unique_ptr< GeoBase > createGeoType(const std::string &wkt_or_wkb_hex)
bool Geos_Wkb(int op, int arg_type, int8_t *arg_coords, int64_t arg_coords_size, int32_t *arg_meta1, int64_t arg_meta1_size, int32_t *arg_meta2, int64_t arg_meta2_size, int32_t arg_ic, int32_t arg_srid, bool *result)
void * checked_malloc(const size_t size)
bool Geos_Wkb_Wkb(int op, int arg1_type, int8_t *arg1_coords, int64_t arg1_coords_size, int32_t *arg1_meta1, int64_t arg1_meta1_size, int32_t *arg1_meta2, int64_t arg1_meta2_size, int32_t arg1_ic, int32_t arg1_srid, int arg2_type, int8_t *arg2_coords, int64_t arg2_coords_size, int32_t *arg2_meta1, int64_t arg2_meta1_size, int32_t *arg2_meta2, int64_t arg2_meta2_size, int32_t arg2_ic, int32_t arg2_srid, int *result_type, int8_t **result_coords, int64_t *result_coords_size, int32_t **result_meta1, int64_t *result_meta1_size, int32_t **result_meta2, int64_t *result_meta2_size)
bool Geos_Wkb_double(int op, int arg1_type, int8_t *arg1_coords, int64_t arg1_coords_size, int32_t *arg1_meta1, int64_t arg1_meta1_size, int32_t *arg1_meta2, int64_t arg1_meta2_size, int32_t arg1_ic, int32_t arg1_srid, double arg2, int *result_type, int8_t **result_coords, int64_t *result_coords_size, int32_t **result_meta1, int64_t *result_meta1_size, int32_t **result_meta2, int64_t *result_meta2_size)
std::shared_ptr< std::vector< double > > decompress_coords< double, int32_t >(const int32_t &ic, const int8_t *coords, const size_t coords_sz)