pgLatLon
changeset 15:95f185a648a4
Added files for next version 0.4 (latlon-v0004.c)
author | jbe |
---|---|
date | Sat Sep 03 16:00:22 2016 +0200 (2016-09-03) |
parents | 7df8eed7539b |
children | e319679cefbd |
files | GNUmakefile latlon--0.4.sql latlon-v0004.c |
line diff
1.1 --- a/GNUmakefile Fri Sep 02 14:32:38 2016 +0200 1.2 +++ b/GNUmakefile Sat Sep 03 16:00:22 2016 +0200 1.3 @@ -1,6 +1,6 @@ 1.4 EXTENSION = latlon 1.5 -DATA = latlon--0.1--0.2.sql latlon--0.2--0.3.sql latlon--0.3.sql 1.6 -MODULES = latlon-v0003 1.7 +DATA = latlon--0.1--0.2.sql latlon--0.2--0.3.sql latlon--0.3.sql latlon--0.4.sql 1.8 +MODULES = latlon-v0003 latlon-v0004 1.9 1.10 PG_CONFIG = pg_config 1.11 PGXS := $(shell $(PG_CONFIG) --pgxs)
2.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000 2.2 +++ b/latlon--0.4.sql Sat Sep 03 16:00:22 2016 +0200 2.3 @@ -0,0 +1,1335 @@ 2.4 + 2.5 +---------------------------------------- 2.6 +-- forward declarations (shell types) -- 2.7 +---------------------------------------- 2.8 + 2.9 +CREATE TYPE epoint; 2.10 +CREATE TYPE ebox; 2.11 +CREATE TYPE ecircle; 2.12 +CREATE TYPE ecluster; 2.13 + 2.14 + 2.15 +------------------------------------------------------------ 2.16 +-- dummy input/output functions for dummy index key types -- 2.17 +------------------------------------------------------------ 2.18 + 2.19 +CREATE FUNCTION ekey_point_in_dummy(cstring) 2.20 + RETURNS ekey_point 2.21 + LANGUAGE C IMMUTABLE STRICT 2.22 + AS '$libdir/latlon-v0004', 'pgl_notimpl'; 2.23 + 2.24 +CREATE FUNCTION ekey_point_out_dummy(ekey_point) 2.25 + RETURNS cstring 2.26 + LANGUAGE C IMMUTABLE STRICT 2.27 + AS '$libdir/latlon-v0004', 'pgl_notimpl'; 2.28 + 2.29 +CREATE FUNCTION ekey_area_in_dummy(cstring) 2.30 + RETURNS ekey_area 2.31 + LANGUAGE C IMMUTABLE STRICT 2.32 + AS '$libdir/latlon-v0004', 'pgl_notimpl'; 2.33 + 2.34 +CREATE FUNCTION ekey_area_out_dummy(ekey_area) 2.35 + RETURNS cstring 2.36 + LANGUAGE C IMMUTABLE STRICT 2.37 + AS '$libdir/latlon-v0004', 'pgl_notimpl'; 2.38 + 2.39 + 2.40 +-------------------------- 2.41 +-- text input functions -- 2.42 +-------------------------- 2.43 + 2.44 +CREATE FUNCTION epoint_in(cstring) 2.45 + RETURNS epoint 2.46 + LANGUAGE C IMMUTABLE STRICT 2.47 + AS '$libdir/latlon-v0004', 'pgl_epoint_in'; 2.48 + 2.49 +CREATE FUNCTION ebox_in(cstring) 2.50 + RETURNS ebox 2.51 + LANGUAGE C IMMUTABLE STRICT 2.52 + AS '$libdir/latlon-v0004', 'pgl_ebox_in'; 2.53 + 2.54 +CREATE FUNCTION ecircle_in(cstring) 2.55 + RETURNS ecircle 2.56 + LANGUAGE C IMMUTABLE STRICT 2.57 + AS '$libdir/latlon-v0004', 'pgl_ecircle_in'; 2.58 + 2.59 +CREATE FUNCTION ecluster_in(cstring) 2.60 + RETURNS ecluster 2.61 + LANGUAGE C IMMUTABLE STRICT 2.62 + AS '$libdir/latlon-v0004', 'pgl_ecluster_in'; 2.63 + 2.64 + 2.65 +--------------------------- 2.66 +-- text output functions -- 2.67 +--------------------------- 2.68 + 2.69 +CREATE FUNCTION epoint_out(epoint) 2.70 + RETURNS cstring 2.71 + LANGUAGE C IMMUTABLE STRICT 2.72 + AS '$libdir/latlon-v0004', 'pgl_epoint_out'; 2.73 + 2.74 +CREATE FUNCTION ebox_out(ebox) 2.75 + RETURNS cstring 2.76 + LANGUAGE C IMMUTABLE STRICT 2.77 + AS '$libdir/latlon-v0004', 'pgl_ebox_out'; 2.78 + 2.79 +CREATE FUNCTION ecircle_out(ecircle) 2.80 + RETURNS cstring 2.81 + LANGUAGE C IMMUTABLE STRICT 2.82 + AS '$libdir/latlon-v0004', 'pgl_ecircle_out'; 2.83 + 2.84 +CREATE FUNCTION ecluster_out(ecluster) 2.85 + RETURNS cstring 2.86 + LANGUAGE C IMMUTABLE STRICT 2.87 + AS '$libdir/latlon-v0004', 'pgl_ecluster_out'; 2.88 + 2.89 + 2.90 +-------------------------- 2.91 +-- binary I/O functions -- 2.92 +-------------------------- 2.93 + 2.94 +CREATE FUNCTION epoint_recv(internal) 2.95 + RETURNS epoint 2.96 + LANGUAGE C IMMUTABLE STRICT 2.97 + AS '$libdir/latlon-v0004', 'pgl_epoint_recv'; 2.98 + 2.99 +CREATE FUNCTION ebox_recv(internal) 2.100 + RETURNS ebox 2.101 + LANGUAGE C IMMUTABLE STRICT 2.102 + AS '$libdir/latlon-v0004', 'pgl_ebox_recv'; 2.103 + 2.104 +CREATE FUNCTION ecircle_recv(internal) 2.105 + RETURNS ecircle 2.106 + LANGUAGE C IMMUTABLE STRICT 2.107 + AS '$libdir/latlon-v0004', 'pgl_ecircle_recv'; 2.108 + 2.109 +CREATE FUNCTION epoint_send(epoint) 2.110 + RETURNS bytea 2.111 + LANGUAGE C IMMUTABLE STRICT 2.112 + AS '$libdir/latlon-v0004', 'pgl_epoint_send'; 2.113 + 2.114 +CREATE FUNCTION ebox_send(ebox) 2.115 + RETURNS bytea 2.116 + LANGUAGE C IMMUTABLE STRICT 2.117 + AS '$libdir/latlon-v0004', 'pgl_ebox_send'; 2.118 + 2.119 +CREATE FUNCTION ecircle_send(ecircle) 2.120 + RETURNS bytea 2.121 + LANGUAGE C IMMUTABLE STRICT 2.122 + AS '$libdir/latlon-v0004', 'pgl_ecircle_send'; 2.123 + 2.124 + 2.125 +----------------------------------------------- 2.126 +-- type definitions of dummy index key types -- 2.127 +----------------------------------------------- 2.128 + 2.129 +CREATE TYPE ekey_point ( 2.130 + internallength = 8, 2.131 + input = ekey_point_in_dummy, 2.132 + output = ekey_point_out_dummy, 2.133 + alignment = char ); 2.134 + 2.135 +CREATE TYPE ekey_area ( 2.136 + internallength = 9, 2.137 + input = ekey_area_in_dummy, 2.138 + output = ekey_area_out_dummy, 2.139 + alignment = char ); 2.140 + 2.141 + 2.142 +------------------------------------------ 2.143 +-- definitions of geographic data types -- 2.144 +------------------------------------------ 2.145 + 2.146 +CREATE TYPE epoint ( 2.147 + internallength = 16, 2.148 + input = epoint_in, 2.149 + output = epoint_out, 2.150 + receive = epoint_recv, 2.151 + send = epoint_send, 2.152 + alignment = double ); 2.153 + 2.154 +CREATE TYPE ebox ( 2.155 + internallength = 32, 2.156 + input = ebox_in, 2.157 + output = ebox_out, 2.158 + receive = ebox_recv, 2.159 + send = ebox_send, 2.160 + alignment = double ); 2.161 + 2.162 +CREATE TYPE ecircle ( 2.163 + internallength = 24, 2.164 + input = ecircle_in, 2.165 + output = ecircle_out, 2.166 + receive = ecircle_recv, 2.167 + send = ecircle_send, 2.168 + alignment = double ); 2.169 + 2.170 +CREATE TYPE ecluster ( 2.171 + internallength = VARIABLE, 2.172 + input = ecluster_in, 2.173 + output = ecluster_out, 2.174 + alignment = double, 2.175 + storage = external ); 2.176 + 2.177 + 2.178 +-------------------- 2.179 +-- B-tree support -- 2.180 +-------------------- 2.181 + 2.182 +-- begin of B-tree support for epoint 2.183 + 2.184 +CREATE FUNCTION epoint_btree_lt(epoint, epoint) 2.185 + RETURNS boolean 2.186 + LANGUAGE C IMMUTABLE STRICT 2.187 + AS '$libdir/latlon-v0004', 'pgl_btree_epoint_lt'; 2.188 + 2.189 +CREATE FUNCTION epoint_btree_le(epoint, epoint) 2.190 + RETURNS boolean 2.191 + LANGUAGE C IMMUTABLE STRICT 2.192 + AS '$libdir/latlon-v0004', 'pgl_btree_epoint_le'; 2.193 + 2.194 +CREATE FUNCTION epoint_btree_eq(epoint, epoint) 2.195 + RETURNS boolean 2.196 + LANGUAGE C IMMUTABLE STRICT 2.197 + AS '$libdir/latlon-v0004', 'pgl_btree_epoint_eq'; 2.198 + 2.199 +CREATE FUNCTION epoint_btree_ne(epoint, epoint) 2.200 + RETURNS boolean 2.201 + LANGUAGE C IMMUTABLE STRICT 2.202 + AS '$libdir/latlon-v0004', 'pgl_btree_epoint_ne'; 2.203 + 2.204 +CREATE FUNCTION epoint_btree_ge(epoint, epoint) 2.205 + RETURNS boolean 2.206 + LANGUAGE C IMMUTABLE STRICT 2.207 + AS '$libdir/latlon-v0004', 'pgl_btree_epoint_ge'; 2.208 + 2.209 +CREATE FUNCTION epoint_btree_gt(epoint, epoint) 2.210 + RETURNS boolean 2.211 + LANGUAGE C IMMUTABLE STRICT 2.212 + AS '$libdir/latlon-v0004', 'pgl_btree_epoint_gt'; 2.213 + 2.214 +CREATE OPERATOR <<< ( 2.215 + leftarg = epoint, 2.216 + rightarg = epoint, 2.217 + procedure = epoint_btree_lt, 2.218 + commutator = >>>, 2.219 + negator = >>>=, 2.220 + restrict = scalarltsel, 2.221 + join = scalarltjoinsel 2.222 +); 2.223 + 2.224 +CREATE OPERATOR <<<= ( 2.225 + leftarg = epoint, 2.226 + rightarg = epoint, 2.227 + procedure = epoint_btree_le, 2.228 + commutator = >>>=, 2.229 + negator = >>>, 2.230 + restrict = scalarltsel, 2.231 + join = scalarltjoinsel 2.232 +); 2.233 + 2.234 +CREATE OPERATOR = ( 2.235 + leftarg = epoint, 2.236 + rightarg = epoint, 2.237 + procedure = epoint_btree_eq, 2.238 + commutator = =, 2.239 + negator = <>, 2.240 + restrict = eqsel, 2.241 + join = eqjoinsel, 2.242 + merges 2.243 +); 2.244 + 2.245 +CREATE OPERATOR <> ( 2.246 + leftarg = epoint, 2.247 + rightarg = epoint, 2.248 + procedure = epoint_btree_eq, 2.249 + commutator = <>, 2.250 + negator = =, 2.251 + restrict = neqsel, 2.252 + join = neqjoinsel 2.253 +); 2.254 + 2.255 +CREATE OPERATOR >>>= ( 2.256 + leftarg = epoint, 2.257 + rightarg = epoint, 2.258 + procedure = epoint_btree_ge, 2.259 + commutator = <<<=, 2.260 + negator = <<<, 2.261 + restrict = scalargtsel, 2.262 + join = scalargtjoinsel 2.263 +); 2.264 + 2.265 +CREATE OPERATOR >>> ( 2.266 + leftarg = epoint, 2.267 + rightarg = epoint, 2.268 + procedure = epoint_btree_gt, 2.269 + commutator = <<<, 2.270 + negator = <<<=, 2.271 + restrict = scalargtsel, 2.272 + join = scalargtjoinsel 2.273 +); 2.274 + 2.275 +CREATE FUNCTION epoint_btree_cmp(epoint, epoint) 2.276 + RETURNS int4 2.277 + LANGUAGE C IMMUTABLE STRICT 2.278 + AS '$libdir/latlon-v0004', 'pgl_btree_epoint_cmp'; 2.279 + 2.280 +CREATE OPERATOR CLASS epoint_btree_ops 2.281 + DEFAULT FOR TYPE epoint USING btree AS 2.282 + OPERATOR 1 <<< , 2.283 + OPERATOR 2 <<<= , 2.284 + OPERATOR 3 = , 2.285 + OPERATOR 4 >>>= , 2.286 + OPERATOR 5 >>> , 2.287 + FUNCTION 1 epoint_btree_cmp(epoint, epoint); 2.288 + 2.289 +-- end of B-tree support for epoint 2.290 + 2.291 +-- begin of B-tree support for ebox 2.292 + 2.293 +CREATE FUNCTION ebox_btree_lt(ebox, ebox) 2.294 + RETURNS boolean 2.295 + LANGUAGE C IMMUTABLE STRICT 2.296 + AS '$libdir/latlon-v0004', 'pgl_btree_ebox_lt'; 2.297 + 2.298 +CREATE FUNCTION ebox_btree_le(ebox, ebox) 2.299 + RETURNS boolean 2.300 + LANGUAGE C IMMUTABLE STRICT 2.301 + AS '$libdir/latlon-v0004', 'pgl_btree_ebox_le'; 2.302 + 2.303 +CREATE FUNCTION ebox_btree_eq(ebox, ebox) 2.304 + RETURNS boolean 2.305 + LANGUAGE C IMMUTABLE STRICT 2.306 + AS '$libdir/latlon-v0004', 'pgl_btree_ebox_eq'; 2.307 + 2.308 +CREATE FUNCTION ebox_btree_ne(ebox, ebox) 2.309 + RETURNS boolean 2.310 + LANGUAGE C IMMUTABLE STRICT 2.311 + AS '$libdir/latlon-v0004', 'pgl_btree_ebox_ne'; 2.312 + 2.313 +CREATE FUNCTION ebox_btree_ge(ebox, ebox) 2.314 + RETURNS boolean 2.315 + LANGUAGE C IMMUTABLE STRICT 2.316 + AS '$libdir/latlon-v0004', 'pgl_btree_ebox_ge'; 2.317 + 2.318 +CREATE FUNCTION ebox_btree_gt(ebox, ebox) 2.319 + RETURNS boolean 2.320 + LANGUAGE C IMMUTABLE STRICT 2.321 + AS '$libdir/latlon-v0004', 'pgl_btree_ebox_gt'; 2.322 + 2.323 +CREATE OPERATOR <<< ( 2.324 + leftarg = ebox, 2.325 + rightarg = ebox, 2.326 + procedure = ebox_btree_lt, 2.327 + commutator = >>>, 2.328 + negator = >>>=, 2.329 + restrict = scalarltsel, 2.330 + join = scalarltjoinsel 2.331 +); 2.332 + 2.333 +CREATE OPERATOR <<<= ( 2.334 + leftarg = ebox, 2.335 + rightarg = ebox, 2.336 + procedure = ebox_btree_le, 2.337 + commutator = >>>=, 2.338 + negator = >>>, 2.339 + restrict = scalarltsel, 2.340 + join = scalarltjoinsel 2.341 +); 2.342 + 2.343 +CREATE OPERATOR = ( 2.344 + leftarg = ebox, 2.345 + rightarg = ebox, 2.346 + procedure = ebox_btree_eq, 2.347 + commutator = =, 2.348 + negator = <>, 2.349 + restrict = eqsel, 2.350 + join = eqjoinsel, 2.351 + merges 2.352 +); 2.353 + 2.354 +CREATE OPERATOR <> ( 2.355 + leftarg = ebox, 2.356 + rightarg = ebox, 2.357 + procedure = ebox_btree_eq, 2.358 + commutator = <>, 2.359 + negator = =, 2.360 + restrict = neqsel, 2.361 + join = neqjoinsel 2.362 +); 2.363 + 2.364 +CREATE OPERATOR >>>= ( 2.365 + leftarg = ebox, 2.366 + rightarg = ebox, 2.367 + procedure = ebox_btree_ge, 2.368 + commutator = <<<=, 2.369 + negator = <<<, 2.370 + restrict = scalargtsel, 2.371 + join = scalargtjoinsel 2.372 +); 2.373 + 2.374 +CREATE OPERATOR >>> ( 2.375 + leftarg = ebox, 2.376 + rightarg = ebox, 2.377 + procedure = ebox_btree_gt, 2.378 + commutator = <<<, 2.379 + negator = <<<=, 2.380 + restrict = scalargtsel, 2.381 + join = scalargtjoinsel 2.382 +); 2.383 + 2.384 +CREATE FUNCTION ebox_btree_cmp(ebox, ebox) 2.385 + RETURNS int4 2.386 + LANGUAGE C IMMUTABLE STRICT 2.387 + AS '$libdir/latlon-v0004', 'pgl_btree_ebox_cmp'; 2.388 + 2.389 +CREATE OPERATOR CLASS ebox_btree_ops 2.390 + DEFAULT FOR TYPE ebox USING btree AS 2.391 + OPERATOR 1 <<< , 2.392 + OPERATOR 2 <<<= , 2.393 + OPERATOR 3 = , 2.394 + OPERATOR 4 >>>= , 2.395 + OPERATOR 5 >>> , 2.396 + FUNCTION 1 ebox_btree_cmp(ebox, ebox); 2.397 + 2.398 +-- end of B-tree support for ebox 2.399 + 2.400 +-- begin of B-tree support for ecircle 2.401 + 2.402 +CREATE FUNCTION ecircle_btree_lt(ecircle, ecircle) 2.403 + RETURNS boolean 2.404 + LANGUAGE C IMMUTABLE STRICT 2.405 + AS '$libdir/latlon-v0004', 'pgl_btree_ecircle_lt'; 2.406 + 2.407 +CREATE FUNCTION ecircle_btree_le(ecircle, ecircle) 2.408 + RETURNS boolean 2.409 + LANGUAGE C IMMUTABLE STRICT 2.410 + AS '$libdir/latlon-v0004', 'pgl_btree_ecircle_le'; 2.411 + 2.412 +CREATE FUNCTION ecircle_btree_eq(ecircle, ecircle) 2.413 + RETURNS boolean 2.414 + LANGUAGE C IMMUTABLE STRICT 2.415 + AS '$libdir/latlon-v0004', 'pgl_btree_ecircle_eq'; 2.416 + 2.417 +CREATE FUNCTION ecircle_btree_ne(ecircle, ecircle) 2.418 + RETURNS boolean 2.419 + LANGUAGE C IMMUTABLE STRICT 2.420 + AS '$libdir/latlon-v0004', 'pgl_btree_ecircle_ne'; 2.421 + 2.422 +CREATE FUNCTION ecircle_btree_ge(ecircle, ecircle) 2.423 + RETURNS boolean 2.424 + LANGUAGE C IMMUTABLE STRICT 2.425 + AS '$libdir/latlon-v0004', 'pgl_btree_ecircle_ge'; 2.426 + 2.427 +CREATE FUNCTION ecircle_btree_gt(ecircle, ecircle) 2.428 + RETURNS boolean 2.429 + LANGUAGE C IMMUTABLE STRICT 2.430 + AS '$libdir/latlon-v0004', 'pgl_btree_ecircle_gt'; 2.431 + 2.432 +CREATE OPERATOR <<< ( 2.433 + leftarg = ecircle, 2.434 + rightarg = ecircle, 2.435 + procedure = ecircle_btree_lt, 2.436 + commutator = >>>, 2.437 + negator = >>>=, 2.438 + restrict = scalarltsel, 2.439 + join = scalarltjoinsel 2.440 +); 2.441 + 2.442 +CREATE OPERATOR <<<= ( 2.443 + leftarg = ecircle, 2.444 + rightarg = ecircle, 2.445 + procedure = ecircle_btree_le, 2.446 + commutator = >>>=, 2.447 + negator = >>>, 2.448 + restrict = scalarltsel, 2.449 + join = scalarltjoinsel 2.450 +); 2.451 + 2.452 +CREATE OPERATOR = ( 2.453 + leftarg = ecircle, 2.454 + rightarg = ecircle, 2.455 + procedure = ecircle_btree_eq, 2.456 + commutator = =, 2.457 + negator = <>, 2.458 + restrict = eqsel, 2.459 + join = eqjoinsel, 2.460 + merges 2.461 +); 2.462 + 2.463 +CREATE OPERATOR <> ( 2.464 + leftarg = ecircle, 2.465 + rightarg = ecircle, 2.466 + procedure = ecircle_btree_eq, 2.467 + commutator = <>, 2.468 + negator = =, 2.469 + restrict = neqsel, 2.470 + join = neqjoinsel 2.471 +); 2.472 + 2.473 +CREATE OPERATOR >>>= ( 2.474 + leftarg = ecircle, 2.475 + rightarg = ecircle, 2.476 + procedure = ecircle_btree_ge, 2.477 + commutator = <<<=, 2.478 + negator = <<<, 2.479 + restrict = scalargtsel, 2.480 + join = scalargtjoinsel 2.481 +); 2.482 + 2.483 +CREATE OPERATOR >>> ( 2.484 + leftarg = ecircle, 2.485 + rightarg = ecircle, 2.486 + procedure = ecircle_btree_gt, 2.487 + commutator = <<<, 2.488 + negator = <<<=, 2.489 + restrict = scalargtsel, 2.490 + join = scalargtjoinsel 2.491 +); 2.492 + 2.493 +CREATE FUNCTION ecircle_btree_cmp(ecircle, ecircle) 2.494 + RETURNS int4 2.495 + LANGUAGE C IMMUTABLE STRICT 2.496 + AS '$libdir/latlon-v0004', 'pgl_btree_ecircle_cmp'; 2.497 + 2.498 +CREATE OPERATOR CLASS ecircle_btree_ops 2.499 + DEFAULT FOR TYPE ecircle USING btree AS 2.500 + OPERATOR 1 <<< , 2.501 + OPERATOR 2 <<<= , 2.502 + OPERATOR 3 = , 2.503 + OPERATOR 4 >>>= , 2.504 + OPERATOR 5 >>> , 2.505 + FUNCTION 1 ecircle_btree_cmp(ecircle, ecircle); 2.506 + 2.507 +-- end of B-tree support for ecircle 2.508 + 2.509 + 2.510 +---------------- 2.511 +-- type casts -- 2.512 +---------------- 2.513 + 2.514 +CREATE FUNCTION cast_epoint_to_ebox(epoint) 2.515 + RETURNS ebox 2.516 + LANGUAGE C IMMUTABLE STRICT 2.517 + AS '$libdir/latlon-v0004', 'pgl_epoint_to_ebox'; 2.518 + 2.519 +CREATE CAST (epoint AS ebox) WITH FUNCTION cast_epoint_to_ebox(epoint); 2.520 + 2.521 +CREATE FUNCTION cast_epoint_to_ecircle(epoint) 2.522 + RETURNS ecircle 2.523 + LANGUAGE C IMMUTABLE STRICT 2.524 + AS '$libdir/latlon-v0004', 'pgl_epoint_to_ecircle'; 2.525 + 2.526 +CREATE CAST (epoint AS ecircle) WITH FUNCTION cast_epoint_to_ecircle(epoint); 2.527 + 2.528 +CREATE FUNCTION cast_epoint_to_ecluster(epoint) 2.529 + RETURNS ecluster 2.530 + LANGUAGE C IMMUTABLE STRICT 2.531 + AS '$libdir/latlon-v0004', 'pgl_epoint_to_ecluster'; 2.532 + 2.533 +CREATE CAST (epoint AS ecluster) WITH FUNCTION cast_epoint_to_ecluster(epoint); 2.534 + 2.535 +CREATE FUNCTION cast_ebox_to_ecluster(ebox) 2.536 + RETURNS ecluster 2.537 + LANGUAGE C IMMUTABLE STRICT 2.538 + AS '$libdir/latlon-v0004', 'pgl_ebox_to_ecluster'; 2.539 + 2.540 +CREATE CAST (ebox AS ecluster) WITH FUNCTION cast_ebox_to_ecluster(ebox); 2.541 + 2.542 + 2.543 +--------------------------- 2.544 +-- constructor functions -- 2.545 +--------------------------- 2.546 + 2.547 +CREATE FUNCTION epoint(float8, float8) 2.548 + RETURNS epoint 2.549 + LANGUAGE C IMMUTABLE STRICT 2.550 + AS '$libdir/latlon-v0004', 'pgl_create_epoint'; 2.551 + 2.552 +CREATE FUNCTION epoint_latlon(float8, float8) 2.553 + RETURNS epoint 2.554 + LANGUAGE SQL IMMUTABLE STRICT AS $$ 2.555 + SELECT epoint($1, $2) 2.556 + $$; 2.557 + 2.558 +CREATE FUNCTION epoint_lonlat(float8, float8) 2.559 + RETURNS epoint 2.560 + LANGUAGE SQL IMMUTABLE STRICT AS $$ 2.561 + SELECT epoint($2, $1) 2.562 + $$; 2.563 + 2.564 +CREATE FUNCTION empty_ebox() 2.565 + RETURNS ebox 2.566 + LANGUAGE C IMMUTABLE STRICT 2.567 + AS '$libdir/latlon-v0004', 'pgl_create_empty_ebox'; 2.568 + 2.569 +CREATE FUNCTION ebox(float8, float8, float8, float8) 2.570 + RETURNS ebox 2.571 + LANGUAGE C IMMUTABLE STRICT 2.572 + AS '$libdir/latlon-v0004', 'pgl_create_ebox'; 2.573 + 2.574 +CREATE FUNCTION ebox(epoint, epoint) 2.575 + RETURNS ebox 2.576 + LANGUAGE C IMMUTABLE STRICT 2.577 + AS '$libdir/latlon-v0004', 'pgl_create_ebox_from_epoints'; 2.578 + 2.579 +CREATE FUNCTION ecircle(float8, float8, float8) 2.580 + RETURNS ecircle 2.581 + LANGUAGE C IMMUTABLE STRICT 2.582 + AS '$libdir/latlon-v0004', 'pgl_create_ecircle'; 2.583 + 2.584 +CREATE FUNCTION ecircle(epoint, float8) 2.585 + RETURNS ecircle 2.586 + LANGUAGE C IMMUTABLE STRICT 2.587 + AS '$libdir/latlon-v0004', 'pgl_create_ecircle_from_epoint'; 2.588 + 2.589 +CREATE FUNCTION ecluster_concat(ecluster[]) 2.590 + RETURNS ecluster 2.591 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.592 + SELECT array_to_string($1, ' ')::ecluster 2.593 + $$; 2.594 + 2.595 +CREATE FUNCTION ecluster_concat(ecluster, ecluster) 2.596 + RETURNS ecluster 2.597 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.598 + SELECT ($1::text || ' ' || $2::text)::ecluster 2.599 + $$; 2.600 + 2.601 +CREATE FUNCTION ecluster_create_multipoint(epoint[]) 2.602 + RETURNS ecluster 2.603 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.604 + SELECT 2.605 + array_to_string(array_agg('point (' || unnest || ')'), ' ')::ecluster 2.606 + FROM unnest($1) 2.607 + $$; 2.608 + 2.609 +CREATE FUNCTION ecluster_create_path(epoint[]) 2.610 + RETURNS ecluster 2.611 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.612 + SELECT CASE WHEN "str" = '' THEN 'empty'::ecluster ELSE 2.613 + ('path (' || array_to_string($1, ' ') || ')')::ecluster 2.614 + END 2.615 + FROM array_to_string($1, ' ') AS "str" 2.616 + $$; 2.617 + 2.618 +CREATE FUNCTION ecluster_create_outline(epoint[]) 2.619 + RETURNS ecluster 2.620 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.621 + SELECT CASE WHEN "str" = '' THEN 'empty'::ecluster ELSE 2.622 + ('outline (' || array_to_string($1, ' ') || ')')::ecluster 2.623 + END 2.624 + FROM array_to_string($1, ' ') AS "str" 2.625 + $$; 2.626 + 2.627 +CREATE FUNCTION ecluster_create_polygon(epoint[]) 2.628 + RETURNS ecluster 2.629 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.630 + SELECT CASE WHEN "str" = '' THEN 'empty'::ecluster ELSE 2.631 + ('polygon (' || array_to_string($1, ' ') || ')')::ecluster 2.632 + END 2.633 + FROM array_to_string($1, ' ') AS "str" 2.634 + $$; 2.635 + 2.636 + 2.637 +---------------------- 2.638 +-- getter functions -- 2.639 +---------------------- 2.640 + 2.641 +CREATE FUNCTION latitude(epoint) 2.642 + RETURNS float8 2.643 + LANGUAGE C IMMUTABLE STRICT 2.644 + AS '$libdir/latlon-v0004', 'pgl_epoint_lat'; 2.645 + 2.646 +CREATE FUNCTION longitude(epoint) 2.647 + RETURNS float8 2.648 + LANGUAGE C IMMUTABLE STRICT 2.649 + AS '$libdir/latlon-v0004', 'pgl_epoint_lon'; 2.650 + 2.651 +CREATE FUNCTION min_latitude(ebox) 2.652 + RETURNS float8 2.653 + LANGUAGE C IMMUTABLE STRICT 2.654 + AS '$libdir/latlon-v0004', 'pgl_ebox_lat_min'; 2.655 + 2.656 +CREATE FUNCTION max_latitude(ebox) 2.657 + RETURNS float8 2.658 + LANGUAGE C IMMUTABLE STRICT 2.659 + AS '$libdir/latlon-v0004', 'pgl_ebox_lat_max'; 2.660 + 2.661 +CREATE FUNCTION min_longitude(ebox) 2.662 + RETURNS float8 2.663 + LANGUAGE C IMMUTABLE STRICT 2.664 + AS '$libdir/latlon-v0004', 'pgl_ebox_lon_min'; 2.665 + 2.666 +CREATE FUNCTION max_longitude(ebox) 2.667 + RETURNS float8 2.668 + LANGUAGE C IMMUTABLE STRICT 2.669 + AS '$libdir/latlon-v0004', 'pgl_ebox_lon_max'; 2.670 + 2.671 +CREATE FUNCTION center(ecircle) 2.672 + RETURNS epoint 2.673 + LANGUAGE C IMMUTABLE STRICT 2.674 + AS '$libdir/latlon-v0004', 'pgl_ecircle_center'; 2.675 + 2.676 +CREATE FUNCTION radius(ecircle) 2.677 + RETURNS float8 2.678 + LANGUAGE C IMMUTABLE STRICT 2.679 + AS '$libdir/latlon-v0004', 'pgl_ecircle_radius'; 2.680 + 2.681 +CREATE FUNCTION ecluster_extract_points(ecluster) 2.682 + RETURNS SETOF epoint 2.683 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.684 + SELECT "match"[2]::epoint 2.685 + FROM regexp_matches($1::text, e'(^| )point \\(([^)]+)\\)', 'g') AS "match" 2.686 + $$; 2.687 + 2.688 +CREATE FUNCTION ecluster_extract_paths(ecluster) 2.689 + RETURNS SETOF epoint[] 2.690 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.691 + SELECT ( 2.692 + SELECT array_agg("m2"[1]::epoint) 2.693 + FROM regexp_matches("m1"[2], e'[^ ]+ [^ ]+', 'g') AS "m2" 2.694 + ) 2.695 + FROM regexp_matches($1::text, e'(^| )path \\(([^)]+)\\)', 'g') AS "m1" 2.696 + $$; 2.697 + 2.698 +CREATE FUNCTION ecluster_extract_outlines(ecluster) 2.699 + RETURNS SETOF epoint[] 2.700 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.701 + SELECT ( 2.702 + SELECT array_agg("m2"[1]::epoint) 2.703 + FROM regexp_matches("m1"[2], e'[^ ]+ [^ ]+', 'g') AS "m2" 2.704 + ) 2.705 + FROM regexp_matches($1::text, e'(^| )outline \\(([^)]+)\\)', 'g') AS "m1" 2.706 + $$; 2.707 + 2.708 +CREATE FUNCTION ecluster_extract_polygons(ecluster) 2.709 + RETURNS SETOF epoint[] 2.710 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.711 + SELECT ( 2.712 + SELECT array_agg("m2"[1]::epoint) 2.713 + FROM regexp_matches("m1"[2], e'[^ ]+ [^ ]+', 'g') AS "m2" 2.714 + ) 2.715 + FROM regexp_matches($1::text, e'(^| )polygon \\(([^)]+)\\)', 'g') AS "m1" 2.716 + $$; 2.717 + 2.718 + 2.719 +--------------- 2.720 +-- operators -- 2.721 +--------------- 2.722 + 2.723 +CREATE FUNCTION epoint_ebox_overlap_proc(epoint, ebox) 2.724 + RETURNS boolean 2.725 + LANGUAGE C IMMUTABLE STRICT 2.726 + AS '$libdir/latlon-v0004', 'pgl_epoint_ebox_overlap'; 2.727 + 2.728 +CREATE FUNCTION epoint_ecircle_overlap_proc(epoint, ecircle) 2.729 + RETURNS boolean 2.730 + LANGUAGE C IMMUTABLE STRICT 2.731 + AS '$libdir/latlon-v0004', 'pgl_epoint_ecircle_overlap'; 2.732 + 2.733 +CREATE FUNCTION epoint_ecluster_overlap_proc(epoint, ecluster) 2.734 + RETURNS boolean 2.735 + LANGUAGE C IMMUTABLE STRICT 2.736 + AS '$libdir/latlon-v0004', 'pgl_epoint_ecluster_overlap'; 2.737 + 2.738 +CREATE FUNCTION epoint_ecluster_may_overlap_proc(epoint, ecluster) 2.739 + RETURNS boolean 2.740 + LANGUAGE C IMMUTABLE STRICT 2.741 + AS '$libdir/latlon-v0004', 'pgl_epoint_ecluster_may_overlap'; 2.742 + 2.743 +CREATE FUNCTION ebox_overlap_proc(ebox, ebox) 2.744 + RETURNS boolean 2.745 + LANGUAGE C IMMUTABLE STRICT 2.746 + AS '$libdir/latlon-v0004', 'pgl_ebox_overlap'; 2.747 + 2.748 +CREATE FUNCTION ebox_ecircle_may_overlap_proc(ebox, ecircle) 2.749 + RETURNS boolean 2.750 + LANGUAGE C IMMUTABLE STRICT 2.751 + AS '$libdir/latlon-v0004', 'pgl_ebox_ecircle_may_overlap'; 2.752 + 2.753 +CREATE FUNCTION ebox_ecluster_may_overlap_proc(ebox, ecluster) 2.754 + RETURNS boolean 2.755 + LANGUAGE C IMMUTABLE STRICT 2.756 + AS '$libdir/latlon-v0004', 'pgl_ebox_ecluster_may_overlap'; 2.757 + 2.758 +CREATE FUNCTION ecircle_overlap_proc(ecircle, ecircle) 2.759 + RETURNS boolean 2.760 + LANGUAGE C IMMUTABLE STRICT 2.761 + AS '$libdir/latlon-v0004', 'pgl_ecircle_overlap'; 2.762 + 2.763 +CREATE FUNCTION ecircle_ecluster_overlap_proc(ecircle, ecluster) 2.764 + RETURNS boolean 2.765 + LANGUAGE C IMMUTABLE STRICT 2.766 + AS '$libdir/latlon-v0004', 'pgl_ecircle_ecluster_overlap'; 2.767 + 2.768 +CREATE FUNCTION ecircle_ecluster_may_overlap_proc(ecircle, ecluster) 2.769 + RETURNS boolean 2.770 + LANGUAGE C IMMUTABLE STRICT 2.771 + AS '$libdir/latlon-v0004', 'pgl_ecircle_ecluster_may_overlap'; 2.772 + 2.773 +CREATE FUNCTION ecluster_may_overlap_proc(ecluster, ecluster) 2.774 + RETURNS boolean 2.775 + LANGUAGE C IMMUTABLE STRICT 2.776 + AS '$libdir/latlon-v0004', 'pgl_ecluster_may_overlap'; 2.777 + 2.778 +CREATE FUNCTION epoint_distance_proc(epoint, epoint) 2.779 + RETURNS float8 2.780 + LANGUAGE C IMMUTABLE STRICT 2.781 + AS '$libdir/latlon-v0004', 'pgl_epoint_distance'; 2.782 + 2.783 +CREATE FUNCTION epoint_ecircle_distance_proc(epoint, ecircle) 2.784 + RETURNS float8 2.785 + LANGUAGE C IMMUTABLE STRICT 2.786 + AS '$libdir/latlon-v0004', 'pgl_epoint_ecircle_distance'; 2.787 + 2.788 +CREATE FUNCTION epoint_ecluster_distance_proc(epoint, ecluster) 2.789 + RETURNS float8 2.790 + LANGUAGE C IMMUTABLE STRICT 2.791 + AS '$libdir/latlon-v0004', 'pgl_epoint_ecluster_distance'; 2.792 + 2.793 +CREATE FUNCTION ecircle_distance_proc(ecircle, ecircle) 2.794 + RETURNS float8 2.795 + LANGUAGE C IMMUTABLE STRICT 2.796 + AS '$libdir/latlon-v0004', 'pgl_ecircle_distance'; 2.797 + 2.798 +CREATE FUNCTION ecircle_ecluster_distance_proc(ecircle, ecluster) 2.799 + RETURNS float8 2.800 + LANGUAGE C IMMUTABLE STRICT 2.801 + AS '$libdir/latlon-v0004', 'pgl_ecircle_ecluster_distance'; 2.802 + 2.803 +CREATE OPERATOR && ( 2.804 + leftarg = epoint, 2.805 + rightarg = ebox, 2.806 + procedure = epoint_ebox_overlap_proc, 2.807 + commutator = &&, 2.808 + restrict = areasel, 2.809 + join = areajoinsel 2.810 +); 2.811 + 2.812 +CREATE FUNCTION epoint_ebox_overlap_commutator(ebox, epoint) 2.813 + RETURNS boolean 2.814 + LANGUAGE sql IMMUTABLE AS 'SELECT $2 && $1'; 2.815 + 2.816 +CREATE OPERATOR && ( 2.817 + leftarg = ebox, 2.818 + rightarg = epoint, 2.819 + procedure = epoint_ebox_overlap_commutator, 2.820 + commutator = &&, 2.821 + restrict = areasel, 2.822 + join = areajoinsel 2.823 +); 2.824 + 2.825 +CREATE OPERATOR && ( 2.826 + leftarg = epoint, 2.827 + rightarg = ecircle, 2.828 + procedure = epoint_ecircle_overlap_proc, 2.829 + commutator = &&, 2.830 + restrict = areasel, 2.831 + join = areajoinsel 2.832 +); 2.833 + 2.834 +CREATE FUNCTION epoint_ecircle_overlap_commutator(ecircle, epoint) 2.835 + RETURNS boolean 2.836 + LANGUAGE sql IMMUTABLE AS 'SELECT $2 && $1'; 2.837 + 2.838 +CREATE OPERATOR && ( 2.839 + leftarg = ecircle, 2.840 + rightarg = epoint, 2.841 + procedure = epoint_ecircle_overlap_commutator, 2.842 + commutator = &&, 2.843 + restrict = areasel, 2.844 + join = areajoinsel 2.845 +); 2.846 + 2.847 +CREATE OPERATOR && ( 2.848 + leftarg = epoint, 2.849 + rightarg = ecluster, 2.850 + procedure = epoint_ecluster_overlap_proc, 2.851 + commutator = &&, 2.852 + restrict = areasel, 2.853 + join = areajoinsel 2.854 +); 2.855 + 2.856 +CREATE FUNCTION epoint_ecluster_overlap_commutator(ecluster, epoint) 2.857 + RETURNS boolean 2.858 + LANGUAGE sql IMMUTABLE AS 'SELECT $2 && $1'; 2.859 + 2.860 +CREATE OPERATOR && ( 2.861 + leftarg = ecluster, 2.862 + rightarg = epoint, 2.863 + procedure = epoint_ecluster_overlap_commutator, 2.864 + commutator = &&, 2.865 + restrict = areasel, 2.866 + join = areajoinsel 2.867 +); 2.868 + 2.869 +CREATE OPERATOR && ( 2.870 + leftarg = ebox, 2.871 + rightarg = ebox, 2.872 + procedure = ebox_overlap_proc, 2.873 + commutator = &&, 2.874 + restrict = areasel, 2.875 + join = areajoinsel 2.876 +); 2.877 + 2.878 +CREATE OPERATOR && ( 2.879 + leftarg = ecircle, 2.880 + rightarg = ecircle, 2.881 + procedure = ecircle_overlap_proc, 2.882 + commutator = &&, 2.883 + restrict = areasel, 2.884 + join = areajoinsel 2.885 +); 2.886 + 2.887 +CREATE OPERATOR && ( 2.888 + leftarg = ecircle, 2.889 + rightarg = ecluster, 2.890 + procedure = ecircle_ecluster_overlap_proc, 2.891 + commutator = &&, 2.892 + restrict = areasel, 2.893 + join = areajoinsel 2.894 +); 2.895 + 2.896 +CREATE FUNCTION ecircle_ecluster_overlap_commutator(ecluster, ecircle) 2.897 + RETURNS boolean 2.898 + LANGUAGE sql IMMUTABLE AS 'SELECT $2 && $1'; 2.899 + 2.900 +CREATE OPERATOR && ( 2.901 + leftarg = ecluster, 2.902 + rightarg = ecircle, 2.903 + procedure = ecircle_ecluster_overlap_commutator, 2.904 + commutator = &&, 2.905 + restrict = areasel, 2.906 + join = areajoinsel 2.907 +); 2.908 + 2.909 +CREATE OPERATOR &&+ ( 2.910 + leftarg = epoint, 2.911 + rightarg = ecluster, 2.912 + procedure = epoint_ecluster_may_overlap_proc, 2.913 + commutator = &&+, 2.914 + restrict = areasel, 2.915 + join = areajoinsel 2.916 +); 2.917 + 2.918 +CREATE FUNCTION epoint_ecluster_may_overlap_commutator(ecluster, epoint) 2.919 + RETURNS boolean 2.920 + LANGUAGE sql IMMUTABLE AS 'SELECT $2 &&+ $1'; 2.921 + 2.922 +CREATE OPERATOR &&+ ( 2.923 + leftarg = ecluster, 2.924 + rightarg = epoint, 2.925 + procedure = epoint_ecluster_may_overlap_commutator, 2.926 + commutator = &&+, 2.927 + restrict = areasel, 2.928 + join = areajoinsel 2.929 +); 2.930 + 2.931 +CREATE OPERATOR &&+ ( 2.932 + leftarg = ebox, 2.933 + rightarg = ecircle, 2.934 + procedure = ebox_ecircle_may_overlap_proc, 2.935 + commutator = &&+, 2.936 + restrict = areasel, 2.937 + join = areajoinsel 2.938 +); 2.939 + 2.940 +CREATE FUNCTION ebox_ecircle_may_overlap_commutator(ecircle, ebox) 2.941 + RETURNS boolean 2.942 + LANGUAGE sql IMMUTABLE AS 'SELECT $2 &&+ $1'; 2.943 + 2.944 +CREATE OPERATOR &&+ ( 2.945 + leftarg = ecircle, 2.946 + rightarg = ebox, 2.947 + procedure = ebox_ecircle_may_overlap_commutator, 2.948 + commutator = &&+, 2.949 + restrict = areasel, 2.950 + join = areajoinsel 2.951 +); 2.952 + 2.953 +CREATE OPERATOR &&+ ( 2.954 + leftarg = ebox, 2.955 + rightarg = ecluster, 2.956 + procedure = ebox_ecluster_may_overlap_proc, 2.957 + commutator = &&+, 2.958 + restrict = areasel, 2.959 + join = areajoinsel 2.960 +); 2.961 + 2.962 +CREATE FUNCTION ebox_ecluster_may_overlap_commutator(ecluster, ebox) 2.963 + RETURNS boolean 2.964 + LANGUAGE sql IMMUTABLE AS 'SELECT $2 &&+ $1'; 2.965 + 2.966 +CREATE OPERATOR &&+ ( 2.967 + leftarg = ecluster, 2.968 + rightarg = ebox, 2.969 + procedure = ebox_ecluster_may_overlap_commutator, 2.970 + commutator = &&+, 2.971 + restrict = areasel, 2.972 + join = areajoinsel 2.973 +); 2.974 + 2.975 +CREATE OPERATOR &&+ ( 2.976 + leftarg = ecircle, 2.977 + rightarg = ecluster, 2.978 + procedure = ecircle_ecluster_may_overlap_proc, 2.979 + commutator = &&+, 2.980 + restrict = areasel, 2.981 + join = areajoinsel 2.982 +); 2.983 + 2.984 +CREATE FUNCTION ecircle_ecluster_may_overlap_commutator(ecluster, ecircle) 2.985 + RETURNS boolean 2.986 + LANGUAGE sql IMMUTABLE AS 'SELECT $2 &&+ $1'; 2.987 + 2.988 +CREATE OPERATOR &&+ ( 2.989 + leftarg = ecluster, 2.990 + rightarg = ecircle, 2.991 + procedure = ecircle_ecluster_may_overlap_commutator, 2.992 + commutator = &&+, 2.993 + restrict = areasel, 2.994 + join = areajoinsel 2.995 +); 2.996 + 2.997 +CREATE OPERATOR &&+ ( 2.998 + leftarg = ecluster, 2.999 + rightarg = ecluster, 2.1000 + procedure = ecluster_may_overlap_proc, 2.1001 + commutator = &&+, 2.1002 + restrict = areasel, 2.1003 + join = areajoinsel 2.1004 +); 2.1005 + 2.1006 +CREATE OPERATOR <-> ( 2.1007 + leftarg = epoint, 2.1008 + rightarg = epoint, 2.1009 + procedure = epoint_distance_proc, 2.1010 + commutator = <-> 2.1011 +); 2.1012 + 2.1013 +CREATE OPERATOR <-> ( 2.1014 + leftarg = epoint, 2.1015 + rightarg = ecircle, 2.1016 + procedure = epoint_ecircle_distance_proc, 2.1017 + commutator = <-> 2.1018 +); 2.1019 + 2.1020 +CREATE FUNCTION epoint_ecircle_distance_commutator(ecircle, epoint) 2.1021 + RETURNS float8 2.1022 + LANGUAGE sql IMMUTABLE AS 'SELECT $2 <-> $1'; 2.1023 + 2.1024 +CREATE OPERATOR <-> ( 2.1025 + leftarg = ecircle, 2.1026 + rightarg = epoint, 2.1027 + procedure = epoint_ecircle_distance_commutator, 2.1028 + commutator = <-> 2.1029 +); 2.1030 + 2.1031 +CREATE OPERATOR <-> ( 2.1032 + leftarg = epoint, 2.1033 + rightarg = ecluster, 2.1034 + procedure = epoint_ecluster_distance_proc, 2.1035 + commutator = <-> 2.1036 +); 2.1037 + 2.1038 +CREATE FUNCTION epoint_ecluster_distance_commutator(ecluster, epoint) 2.1039 + RETURNS float8 2.1040 + LANGUAGE sql IMMUTABLE AS 'SELECT $2 <-> $1'; 2.1041 + 2.1042 +CREATE OPERATOR <-> ( 2.1043 + leftarg = ecluster, 2.1044 + rightarg = epoint, 2.1045 + procedure = epoint_ecluster_distance_commutator, 2.1046 + commutator = <-> 2.1047 +); 2.1048 + 2.1049 +CREATE OPERATOR <-> ( 2.1050 + leftarg = ecircle, 2.1051 + rightarg = ecircle, 2.1052 + procedure = ecircle_distance_proc, 2.1053 + commutator = <-> 2.1054 +); 2.1055 + 2.1056 +CREATE OPERATOR <-> ( 2.1057 + leftarg = ecircle, 2.1058 + rightarg = ecluster, 2.1059 + procedure = ecircle_ecluster_distance_proc, 2.1060 + commutator = <-> 2.1061 +); 2.1062 + 2.1063 +CREATE FUNCTION ecircle_ecluster_distance_commutator(ecluster, ecircle) 2.1064 + RETURNS float8 2.1065 + LANGUAGE sql IMMUTABLE AS 'SELECT $2 <-> $1'; 2.1066 + 2.1067 +CREATE OPERATOR <-> ( 2.1068 + leftarg = ecluster, 2.1069 + rightarg = ecircle, 2.1070 + procedure = ecircle_ecluster_distance_commutator, 2.1071 + commutator = <-> 2.1072 +); 2.1073 + 2.1074 + 2.1075 +---------------- 2.1076 +-- GiST index -- 2.1077 +---------------- 2.1078 + 2.1079 +CREATE FUNCTION pgl_gist_consistent(internal, internal, smallint, oid, internal) 2.1080 + RETURNS boolean 2.1081 + LANGUAGE C STRICT 2.1082 + AS '$libdir/latlon-v0004', 'pgl_gist_consistent'; 2.1083 + 2.1084 +CREATE FUNCTION pgl_gist_union(internal, internal) 2.1085 + RETURNS internal 2.1086 + LANGUAGE C STRICT 2.1087 + AS '$libdir/latlon-v0004', 'pgl_gist_union'; 2.1088 + 2.1089 +CREATE FUNCTION pgl_gist_compress_epoint(internal) 2.1090 + RETURNS internal 2.1091 + LANGUAGE C STRICT 2.1092 + AS '$libdir/latlon-v0004', 'pgl_gist_compress_epoint'; 2.1093 + 2.1094 +CREATE FUNCTION pgl_gist_compress_ecircle(internal) 2.1095 + RETURNS internal 2.1096 + LANGUAGE C STRICT 2.1097 + AS '$libdir/latlon-v0004', 'pgl_gist_compress_ecircle'; 2.1098 + 2.1099 +CREATE FUNCTION pgl_gist_compress_ecluster(internal) 2.1100 + RETURNS internal 2.1101 + LANGUAGE C STRICT 2.1102 + AS '$libdir/latlon-v0004', 'pgl_gist_compress_ecluster'; 2.1103 + 2.1104 +CREATE FUNCTION pgl_gist_decompress(internal) 2.1105 + RETURNS internal 2.1106 + LANGUAGE C STRICT 2.1107 + AS '$libdir/latlon-v0004', 'pgl_gist_decompress'; 2.1108 + 2.1109 +CREATE FUNCTION pgl_gist_penalty(internal, internal, internal) 2.1110 + RETURNS internal 2.1111 + LANGUAGE C STRICT 2.1112 + AS '$libdir/latlon-v0004', 'pgl_gist_penalty'; 2.1113 + 2.1114 +CREATE FUNCTION pgl_gist_picksplit(internal, internal) 2.1115 + RETURNS internal 2.1116 + LANGUAGE C STRICT 2.1117 + AS '$libdir/latlon-v0004', 'pgl_gist_picksplit'; 2.1118 + 2.1119 +CREATE FUNCTION pgl_gist_same(internal, internal, internal) 2.1120 + RETURNS internal 2.1121 + LANGUAGE C STRICT 2.1122 + AS '$libdir/latlon-v0004', 'pgl_gist_same'; 2.1123 + 2.1124 +CREATE FUNCTION pgl_gist_distance(internal, internal, smallint, oid) 2.1125 + RETURNS internal 2.1126 + LANGUAGE C STRICT 2.1127 + AS '$libdir/latlon-v0004', 'pgl_gist_distance'; 2.1128 + 2.1129 +CREATE OPERATOR CLASS epoint_ops 2.1130 + DEFAULT FOR TYPE epoint USING gist AS 2.1131 + OPERATOR 11 = , 2.1132 + OPERATOR 22 && (epoint, ebox), 2.1133 + OPERATOR 23 && (epoint, ecircle), 2.1134 + OPERATOR 24 && (epoint, ecluster), 2.1135 + OPERATOR 124 &&+ (epoint, ecluster), 2.1136 + OPERATOR 31 <-> (epoint, epoint) FOR ORDER BY float_ops, 2.1137 + OPERATOR 33 <-> (epoint, ecircle) FOR ORDER BY float_ops, 2.1138 + OPERATOR 34 <-> (epoint, ecluster) FOR ORDER BY float_ops, 2.1139 + FUNCTION 1 pgl_gist_consistent(internal, internal, smallint, oid, internal), 2.1140 + FUNCTION 2 pgl_gist_union(internal, internal), 2.1141 + FUNCTION 3 pgl_gist_compress_epoint(internal), 2.1142 + FUNCTION 4 pgl_gist_decompress(internal), 2.1143 + FUNCTION 5 pgl_gist_penalty(internal, internal, internal), 2.1144 + FUNCTION 6 pgl_gist_picksplit(internal, internal), 2.1145 + FUNCTION 7 pgl_gist_same(internal, internal, internal), 2.1146 + FUNCTION 8 pgl_gist_distance(internal, internal, smallint, oid), 2.1147 + STORAGE ekey_point; 2.1148 + 2.1149 +CREATE OPERATOR CLASS ecircle_ops 2.1150 + DEFAULT FOR TYPE ecircle USING gist AS 2.1151 + OPERATOR 13 = , 2.1152 + OPERATOR 21 && (ecircle, epoint), 2.1153 + OPERATOR 122 &&+ (ecircle, ebox), 2.1154 + OPERATOR 23 && (ecircle, ecircle), 2.1155 + OPERATOR 24 && (ecircle, ecluster), 2.1156 + OPERATOR 124 &&+ (ecircle, ecluster), 2.1157 + OPERATOR 31 <-> (ecircle, epoint) FOR ORDER BY float_ops, 2.1158 + OPERATOR 33 <-> (ecircle, ecircle) FOR ORDER BY float_ops, 2.1159 + OPERATOR 34 <-> (ecircle, ecluster) FOR ORDER BY float_ops, 2.1160 + FUNCTION 1 pgl_gist_consistent(internal, internal, smallint, oid, internal), 2.1161 + FUNCTION 2 pgl_gist_union(internal, internal), 2.1162 + FUNCTION 3 pgl_gist_compress_ecircle(internal), 2.1163 + FUNCTION 4 pgl_gist_decompress(internal), 2.1164 + FUNCTION 5 pgl_gist_penalty(internal, internal, internal), 2.1165 + FUNCTION 6 pgl_gist_picksplit(internal, internal), 2.1166 + FUNCTION 7 pgl_gist_same(internal, internal, internal), 2.1167 + FUNCTION 8 pgl_gist_distance(internal, internal, smallint, oid), 2.1168 + STORAGE ekey_area; 2.1169 + 2.1170 +CREATE OPERATOR CLASS ecluster_ops 2.1171 + DEFAULT FOR TYPE ecluster USING gist AS 2.1172 + OPERATOR 21 && (ecluster, epoint), 2.1173 + OPERATOR 121 &&+ (ecluster, epoint), 2.1174 + OPERATOR 122 &&+ (ecluster, ebox), 2.1175 + OPERATOR 23 && (ecluster, ecircle), 2.1176 + OPERATOR 123 &&+ (ecluster, ecircle), 2.1177 + OPERATOR 124 &&+ (ecluster, ecluster), 2.1178 + FUNCTION 1 pgl_gist_consistent(internal, internal, smallint, oid, internal), 2.1179 + FUNCTION 2 pgl_gist_union(internal, internal), 2.1180 + FUNCTION 3 pgl_gist_compress_ecluster(internal), 2.1181 + FUNCTION 4 pgl_gist_decompress(internal), 2.1182 + FUNCTION 5 pgl_gist_penalty(internal, internal, internal), 2.1183 + FUNCTION 6 pgl_gist_picksplit(internal, internal), 2.1184 + FUNCTION 7 pgl_gist_same(internal, internal, internal), 2.1185 + FUNCTION 8 pgl_gist_distance(internal, internal, smallint, oid), 2.1186 + STORAGE ekey_area; 2.1187 + 2.1188 + 2.1189 +--------------------- 2.1190 +-- alias functions -- 2.1191 +--------------------- 2.1192 + 2.1193 +CREATE FUNCTION distance(epoint, epoint) 2.1194 + RETURNS float8 2.1195 + LANGUAGE sql IMMUTABLE AS 'SELECT $1 <-> $2'; 2.1196 + 2.1197 +CREATE FUNCTION distance(ecluster, epoint) 2.1198 + RETURNS float8 2.1199 + LANGUAGE sql IMMUTABLE AS 'SELECT $1 <-> $2'; 2.1200 + 2.1201 +CREATE FUNCTION distance_within(epoint, epoint, float8) 2.1202 + RETURNS boolean 2.1203 + LANGUAGE sql IMMUTABLE AS 'SELECT $1 && ecircle($2, $3)'; 2.1204 + 2.1205 +CREATE FUNCTION distance_within(ecluster, epoint, float8) 2.1206 + RETURNS boolean 2.1207 + LANGUAGE sql IMMUTABLE AS 'SELECT $1 && ecircle($2, $3)'; 2.1208 + 2.1209 + 2.1210 +-------------------------------- 2.1211 +-- other data storage formats -- 2.1212 +-------------------------------- 2.1213 + 2.1214 +CREATE FUNCTION coords_to_epoint(float8, float8, text = 'epoint_lonlat') 2.1215 + RETURNS epoint 2.1216 + LANGUAGE plpgsql IMMUTABLE STRICT AS $$ 2.1217 + DECLARE 2.1218 + "result" epoint; 2.1219 + BEGIN 2.1220 + IF $3 = 'epoint_lonlat' THEN 2.1221 + -- avoid dynamic command execution for better performance 2.1222 + RETURN epoint($2, $1); 2.1223 + END IF; 2.1224 + IF $3 = 'epoint' OR $3 = 'epoint_latlon' THEN 2.1225 + -- avoid dynamic command execution for better performance 2.1226 + RETURN epoint($1, $2); 2.1227 + END IF; 2.1228 + EXECUTE 'SELECT ' || $3 || '($1, $2)' INTO STRICT "result" USING $1, $2; 2.1229 + RETURN "result"; 2.1230 + END; 2.1231 + $$; 2.1232 + 2.1233 +CREATE FUNCTION GeoJSON_to_epoint(jsonb, text = 'epoint_lonlat') 2.1234 + RETURNS epoint 2.1235 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.1236 + SELECT CASE 2.1237 + WHEN $1->>'type' = 'Point' THEN 2.1238 + coords_to_epoint( 2.1239 + ($1->'coordinates'->>1)::float8, 2.1240 + ($1->'coordinates'->>0)::float8, 2.1241 + $2 2.1242 + ) 2.1243 + WHEN $1->>'type' = 'Feature' THEN 2.1244 + GeoJSON_to_epoint($1->'geometry', $2) 2.1245 + ELSE 2.1246 + NULL 2.1247 + END 2.1248 + $$; 2.1249 + 2.1250 +CREATE FUNCTION GeoJSON_to_ecluster(jsonb, text = 'epoint_lonlat') 2.1251 + RETURNS ecluster 2.1252 + LANGUAGE sql IMMUTABLE STRICT AS $$ 2.1253 + SELECT CASE $1->>'type' 2.1254 + WHEN 'Point' THEN 2.1255 + coords_to_epoint( 2.1256 + ($1->'coordinates'->>1)::float8, 2.1257 + ($1->'coordinates'->>0)::float8, 2.1258 + $2 2.1259 + )::ecluster 2.1260 + WHEN 'MultiPoint' THEN 2.1261 + ( SELECT ecluster_create_multipoint(array_agg( 2.1262 + coords_to_epoint( 2.1263 + ("coord"->>1)::float8, 2.1264 + ("coord"->>0)::float8, 2.1265 + $2 2.1266 + ) 2.1267 + )) 2.1268 + FROM jsonb_array_elements($1->'coordinates') AS "coord" 2.1269 + ) 2.1270 + WHEN 'LineString' THEN 2.1271 + ( SELECT ecluster_create_path(array_agg( 2.1272 + coords_to_epoint( 2.1273 + ("coord"->>1)::float8, 2.1274 + ("coord"->>0)::float8, 2.1275 + $2 2.1276 + ) 2.1277 + )) 2.1278 + FROM jsonb_array_elements($1->'coordinates') AS "coord" 2.1279 + ) 2.1280 + WHEN 'MultiLineString' THEN 2.1281 + ( SELECT ecluster_concat(array_agg( 2.1282 + ( SELECT ecluster_create_path(array_agg( 2.1283 + coords_to_epoint( 2.1284 + ("coord"->>1)::float8, 2.1285 + ("coord"->>0)::float8, 2.1286 + $2 2.1287 + ) 2.1288 + )) 2.1289 + FROM jsonb_array_elements("coord_array") AS "coord" 2.1290 + ) 2.1291 + )) 2.1292 + FROM jsonb_array_elements($1->'coordinates') AS "coord_array" 2.1293 + ) 2.1294 + WHEN 'Polygon' THEN 2.1295 + ( SELECT ecluster_concat(array_agg( 2.1296 + ( SELECT ecluster_create_polygon(array_agg( 2.1297 + coords_to_epoint( 2.1298 + ("coord"->>1)::float8, 2.1299 + ("coord"->>0)::float8, 2.1300 + $2 2.1301 + ) 2.1302 + )) 2.1303 + FROM jsonb_array_elements("coord_array") AS "coord" 2.1304 + ) 2.1305 + )) 2.1306 + FROM jsonb_array_elements($1->'coordinates') AS "coord_array" 2.1307 + ) 2.1308 + WHEN 'MultiPolygon' THEN 2.1309 + ( SELECT ecluster_concat(array_agg( 2.1310 + ( SELECT ecluster_concat(array_agg( 2.1311 + ( SELECT ecluster_create_polygon(array_agg( 2.1312 + coords_to_epoint( 2.1313 + ("coord"->>1)::float8, 2.1314 + ("coord"->>0)::float8, 2.1315 + $2 2.1316 + ) 2.1317 + )) 2.1318 + FROM jsonb_array_elements("coord_array") AS "coord" 2.1319 + ) 2.1320 + )) 2.1321 + FROM jsonb_array_elements("coord_array_array") AS "coord_array" 2.1322 + ) 2.1323 + )) 2.1324 + FROM jsonb_array_elements($1->'coordinates') AS "coord_array_array" 2.1325 + ) 2.1326 + WHEN 'Feature' THEN 2.1327 + GeoJSON_to_ecluster($1->'geometry', $2) 2.1328 + WHEN 'FeatureCollection' THEN 2.1329 + ( SELECT ecluster_concat(array_agg( 2.1330 + GeoJSON_to_ecluster("feature", $2) 2.1331 + )) 2.1332 + FROM jsonb_array_elements($1->'features') AS "feature" 2.1333 + ) 2.1334 + ELSE 2.1335 + NULL 2.1336 + END 2.1337 + $$; 2.1338 +
3.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000 3.2 +++ b/latlon-v0004.c Sat Sep 03 16:00:22 2016 +0200 3.3 @@ -0,0 +1,2768 @@ 3.4 + 3.5 +/*-------------* 3.6 + * C prelude * 3.7 + *-------------*/ 3.8 + 3.9 +#include "postgres.h" 3.10 +#include "fmgr.h" 3.11 +#include "libpq/pqformat.h" 3.12 +#include "access/gist.h" 3.13 +#include "access/stratnum.h" 3.14 +#include "utils/array.h" 3.15 +#include <math.h> 3.16 + 3.17 +#ifdef PG_MODULE_MAGIC 3.18 +PG_MODULE_MAGIC; 3.19 +#endif 3.20 + 3.21 +#if INT_MAX < 2147483647 3.22 +#error Expected int type to be at least 32 bit wide 3.23 +#endif 3.24 + 3.25 + 3.26 +/*---------------------------------* 3.27 + * distance calculation on earth * 3.28 + * (using WGS-84 spheroid) * 3.29 + *---------------------------------*/ 3.30 + 3.31 +/* WGS-84 spheroid with following parameters: 3.32 + semi-major axis a = 6378137 3.33 + semi-minor axis b = a * (1 - 1/298.257223563) 3.34 + estimated diameter = 2 * (2*a+b)/3 3.35 +*/ 3.36 +#define PGL_SPHEROID_A 6378137.0 /* semi major axis */ 3.37 +#define PGL_SPHEROID_F (1.0/298.257223563) /* flattening */ 3.38 +#define PGL_SPHEROID_B (PGL_SPHEROID_A * (1.0-PGL_SPHEROID_F)) 3.39 +#define PGL_EPS2 ( ( PGL_SPHEROID_A * PGL_SPHEROID_A - \ 3.40 + PGL_SPHEROID_B * PGL_SPHEROID_B ) / \ 3.41 + ( PGL_SPHEROID_A * PGL_SPHEROID_A ) ) 3.42 +#define PGL_SUBEPS2 (1.0-PGL_EPS2) 3.43 +#define PGL_DIAMETER ((4.0*PGL_SPHEROID_A + 2.0*PGL_SPHEROID_B) / 3.0) 3.44 +#define PGL_SCALE (PGL_SPHEROID_A / PGL_DIAMETER) /* semi-major ref. */ 3.45 +#define PGL_FADELIMIT (PGL_DIAMETER * M_PI / 6.0) /* 1/6 circumference */ 3.46 +#define PGL_MAXDIST (PGL_DIAMETER * M_PI / 2.0) /* maximum distance */ 3.47 + 3.48 +/* calculate distance between two points on earth (given in degrees) */ 3.49 +static inline double pgl_distance( 3.50 + double lat1, double lon1, double lat2, double lon2 3.51 +) { 3.52 + float8 lat1cos, lat1sin, lat2cos, lat2sin, lon2cos, lon2sin; 3.53 + float8 nphi1, nphi2, x1, z1, x2, y2, z2, g, s, t; 3.54 + /* normalize delta longitude (lon2 > 0 && lon1 = 0) */ 3.55 + /* lon1 = 0 (not used anymore) */ 3.56 + lon2 = fabs(lon2-lon1); 3.57 + /* convert to radians (first divide, then multiply) */ 3.58 + lat1 = (lat1 / 180.0) * M_PI; 3.59 + lat2 = (lat2 / 180.0) * M_PI; 3.60 + lon2 = (lon2 / 180.0) * M_PI; 3.61 + /* make lat2 >= lat1 to ensure reversal-symmetry despite floating point 3.62 + operations (lon2 >= lon1 is already ensured in a previous step) */ 3.63 + if (lat2 < lat1) { float8 swap = lat1; lat1 = lat2; lat2 = swap; } 3.64 + /* calculate 3d coordinates on scaled ellipsoid which has an average diameter 3.65 + of 1.0 */ 3.66 + lat1cos = cos(lat1); lat1sin = sin(lat1); 3.67 + lat2cos = cos(lat2); lat2sin = sin(lat2); 3.68 + lon2cos = cos(lon2); lon2sin = sin(lon2); 3.69 + nphi1 = PGL_SCALE / sqrt(1 - PGL_EPS2 * lat1sin * lat1sin); 3.70 + nphi2 = PGL_SCALE / sqrt(1 - PGL_EPS2 * lat2sin * lat2sin); 3.71 + x1 = nphi1 * lat1cos; 3.72 + z1 = nphi1 * PGL_SUBEPS2 * lat1sin; 3.73 + x2 = nphi2 * lat2cos * lon2cos; 3.74 + y2 = nphi2 * lat2cos * lon2sin; 3.75 + z2 = nphi2 * PGL_SUBEPS2 * lat2sin; 3.76 + /* calculate tunnel distance through scaled (diameter 1.0) ellipsoid */ 3.77 + g = sqrt((x2-x1)*(x2-x1) + y2*y2 + (z2-z1)*(z2-z1)); 3.78 + /* convert tunnel distance through scaled ellipsoid to approximated surface 3.79 + distance on original ellipsoid */ 3.80 + if (g > 1.0) g = 1.0; 3.81 + s = PGL_DIAMETER * asin(g); 3.82 + /* return result only if small enough to be precise (less than 1/3 of 3.83 + maximum possible distance) */ 3.84 + if (s <= PGL_FADELIMIT) return s; 3.85 + /* calculate tunnel distance to antipodal point through scaled ellipsoid */ 3.86 + g = sqrt((x2+x1)*(x2+x1) + y2*y2 + (z2+z1)*(z2+z1)); 3.87 + /* convert tunnel distance to antipodal point through scaled ellipsoid to 3.88 + approximated surface distance to antipodal point on original ellipsoid */ 3.89 + if (g > 1.0) g = 1.0; 3.90 + t = PGL_DIAMETER * asin(g); 3.91 + /* surface distance between original points can now be approximated by 3.92 + substracting antipodal distance from maximum possible distance; 3.93 + return result only if small enough (less than 1/3 of maximum possible 3.94 + distance) */ 3.95 + if (t <= PGL_FADELIMIT) return PGL_MAXDIST-t; 3.96 + /* otherwise crossfade direct and antipodal result to ensure monotonicity */ 3.97 + return ( 3.98 + (s * (t-PGL_FADELIMIT) + (PGL_MAXDIST-t) * (s-PGL_FADELIMIT)) / 3.99 + (s + t - 2*PGL_FADELIMIT) 3.100 + ); 3.101 +} 3.102 + 3.103 +/* finite distance that can not be reached on earth */ 3.104 +#define PGL_ULTRA_DISTANCE (3 * PGL_MAXDIST) 3.105 + 3.106 + 3.107 +/*--------------------------------* 3.108 + * simple geographic data types * 3.109 + *--------------------------------*/ 3.110 + 3.111 +/* point on earth given by latitude and longitude in degrees */ 3.112 +/* (type "epoint" in SQL) */ 3.113 +typedef struct { 3.114 + double lat; /* between -90 and 90 (both inclusive) */ 3.115 + double lon; /* between -180 and 180 (both inclusive) */ 3.116 +} pgl_point; 3.117 + 3.118 +/* box delimited by two parallels and two meridians (all in degrees) */ 3.119 +/* (type "ebox" in SQL) */ 3.120 +typedef struct { 3.121 + double lat_min; /* between -90 and 90 (both inclusive) */ 3.122 + double lat_max; /* between -90 and 90 (both inclusive) */ 3.123 + double lon_min; /* between -180 and 180 (both inclusive) */ 3.124 + double lon_max; /* between -180 and 180 (both inclusive) */ 3.125 + /* if lat_min > lat_max, then box is empty */ 3.126 + /* if lon_min > lon_max, then 180th meridian is crossed */ 3.127 +} pgl_box; 3.128 + 3.129 +/* circle on earth surface (for radial searches with fixed radius) */ 3.130 +/* (type "ecircle" in SQL) */ 3.131 +typedef struct { 3.132 + pgl_point center; 3.133 + double radius; /* positive (including +0 but excluding -0), or -INFINITY */ 3.134 + /* A negative radius (i.e. -INFINITY) denotes nothing (i.e. no point), 3.135 + zero radius (0) denotes a single point, 3.136 + a finite radius (0 < radius < INFINITY) denotes a filled circle, and 3.137 + a radius of INFINITY is valid and means complete coverage of earth. */ 3.138 +} pgl_circle; 3.139 + 3.140 + 3.141 +/*----------------------------------* 3.142 + * geographic "cluster" data type * 3.143 + *----------------------------------*/ 3.144 + 3.145 +/* A cluster is a collection of points, paths, outlines, and polygons. If two 3.146 + polygons in a cluster overlap, the area covered by both polygons does not 3.147 + belong to the cluster. This way, a cluster can be used to describe complex 3.148 + shapes like polygons with holes. Outlines are non-filled polygons. Paths are 3.149 + open by default (i.e. the last point in the list is not connected with the 3.150 + first point in the list). Note that each outline or polygon in a cluster 3.151 + must cover a longitude range of less than 180 degrees to avoid ambiguities. 3.152 + Areas which are larger may be split into multiple polygons. */ 3.153 + 3.154 +/* maximum number of points in a cluster */ 3.155 +/* (limited to avoid integer overflows, e.g. when allocating memory) */ 3.156 +#define PGL_CLUSTER_MAXPOINTS 16777216 3.157 + 3.158 +/* types of cluster entries */ 3.159 +#define PGL_ENTRY_POINT 1 /* a point */ 3.160 +#define PGL_ENTRY_PATH 2 /* a path from first point to last point */ 3.161 +#define PGL_ENTRY_OUTLINE 3 /* a non-filled polygon with given vertices */ 3.162 +#define PGL_ENTRY_POLYGON 4 /* a filled polygon with given vertices */ 3.163 + 3.164 +/* Entries of a cluster are described by two different structs: pgl_newentry 3.165 + and pgl_entry. The first is used only during construction of a cluster, the 3.166 + second is used in all other cases (e.g. when reading clusters from the 3.167 + database, performing operations, etc). */ 3.168 + 3.169 +/* entry for new geographic cluster during construction of that cluster */ 3.170 +typedef struct { 3.171 + int32_t entrytype; 3.172 + int32_t npoints; 3.173 + pgl_point *points; /* pointer to an array of points (pgl_point) */ 3.174 +} pgl_newentry; 3.175 + 3.176 +/* entry of geographic cluster */ 3.177 +typedef struct { 3.178 + int32_t entrytype; /* type of entry: point, path, outline, polygon */ 3.179 + int32_t npoints; /* number of stored points (set to 1 for point entry) */ 3.180 + int32_t offset; /* offset of pgl_point array from cluster base address */ 3.181 + /* use macro PGL_ENTRY_POINTS to obtain a pointer to the array of points */ 3.182 +} pgl_entry; 3.183 + 3.184 +/* geographic cluster which is a collection of points, (open) paths, polygons, 3.185 + and outlines (non-filled polygons) */ 3.186 +typedef struct { 3.187 + char header[VARHDRSZ]; /* PostgreSQL header for variable size data types */ 3.188 + int32_t nentries; /* number of stored points */ 3.189 + pgl_circle bounding; /* bounding circle */ 3.190 + /* Note: bounding circle ensures alignment of pgl_cluster for points */ 3.191 + pgl_entry entries[FLEXIBLE_ARRAY_MEMBER]; /* var-length data */ 3.192 +} pgl_cluster; 3.193 + 3.194 +/* macro to determine memory alignment of points */ 3.195 +/* (needed to store pgl_point array after entries in pgl_cluster) */ 3.196 +typedef struct { char dummy; pgl_point aligned; } pgl_point_alignment; 3.197 +#define PGL_POINT_ALIGNMENT offsetof(pgl_point_alignment, aligned) 3.198 + 3.199 +/* macro to extract a pointer to the array of points of a cluster entry */ 3.200 +#define PGL_ENTRY_POINTS(cluster, idx) \ 3.201 + ((pgl_point *)(((intptr_t)cluster)+(cluster)->entries[idx].offset)) 3.202 + 3.203 +/* convert pgl_newentry array to pgl_cluster */ 3.204 +static pgl_cluster *pgl_new_cluster(int nentries, pgl_newentry *entries) { 3.205 + int i; /* index of current entry */ 3.206 + int npoints = 0; /* number of points in whole cluster */ 3.207 + int entry_npoints; /* number of points in current entry */ 3.208 + int points_offset = PGL_POINT_ALIGNMENT * ( 3.209 + ( offsetof(pgl_cluster, entries) + 3.210 + nentries * sizeof(pgl_entry) + 3.211 + PGL_POINT_ALIGNMENT - 1 3.212 + ) / PGL_POINT_ALIGNMENT 3.213 + ); /* offset of pgl_point array from base address (considering alignment) */ 3.214 + pgl_cluster *cluster; /* new cluster to be returned */ 3.215 + /* determine total number of points */ 3.216 + for (i=0; i<nentries; i++) npoints += entries[i].npoints; 3.217 + /* allocate memory for cluster (including entries and points) */ 3.218 + cluster = palloc(points_offset + npoints * sizeof(pgl_point)); 3.219 + /* re-count total number of points to determine offset for each entry */ 3.220 + npoints = 0; 3.221 + /* copy entries and points */ 3.222 + for (i=0; i<nentries; i++) { 3.223 + /* determine number of points in entry */ 3.224 + entry_npoints = entries[i].npoints; 3.225 + /* copy entry */ 3.226 + cluster->entries[i].entrytype = entries[i].entrytype; 3.227 + cluster->entries[i].npoints = entry_npoints; 3.228 + /* calculate offset (in bytes) of pgl_point array */ 3.229 + cluster->entries[i].offset = points_offset + npoints * sizeof(pgl_point); 3.230 + /* copy points */ 3.231 + memcpy( 3.232 + PGL_ENTRY_POINTS(cluster, i), 3.233 + entries[i].points, 3.234 + entry_npoints * sizeof(pgl_point) 3.235 + ); 3.236 + /* update total number of points processed */ 3.237 + npoints += entry_npoints; 3.238 + } 3.239 + /* set number of entries in cluster */ 3.240 + cluster->nentries = nentries; 3.241 + /* set PostgreSQL header for variable sized data */ 3.242 + SET_VARSIZE(cluster, points_offset + npoints * sizeof(pgl_point)); 3.243 + /* return newly created cluster */ 3.244 + return cluster; 3.245 +} 3.246 + 3.247 + 3.248 +/*----------------------------------------* 3.249 + * C functions on geographic data types * 3.250 + *----------------------------------------*/ 3.251 + 3.252 +/* round latitude or longitude to 12 digits after decimal point */ 3.253 +static inline double pgl_round(double val) { 3.254 + return round(val * 1e12) / 1e12; 3.255 +} 3.256 + 3.257 +/* compare two points */ 3.258 +/* (equality when same point on earth is described, otherwise an arbitrary 3.259 + linear order) */ 3.260 +static int pgl_point_cmp(pgl_point *point1, pgl_point *point2) { 3.261 + double lon1, lon2; /* modified longitudes for special cases */ 3.262 + /* use latitude as first ordering criterion */ 3.263 + if (point1->lat < point2->lat) return -1; 3.264 + if (point1->lat > point2->lat) return 1; 3.265 + /* determine modified longitudes (considering special case of poles and 3.266 + 180th meridian which can be described as W180 or E180) */ 3.267 + if (point1->lat == -90 || point1->lat == 90) lon1 = 0; 3.268 + else if (point1->lon == 180) lon1 = -180; 3.269 + else lon1 = point1->lon; 3.270 + if (point2->lat == -90 || point2->lat == 90) lon2 = 0; 3.271 + else if (point2->lon == 180) lon2 = -180; 3.272 + else lon2 = point2->lon; 3.273 + /* use (modified) longitude as secondary ordering criterion */ 3.274 + if (lon1 < lon2) return -1; 3.275 + if (lon1 > lon2) return 1; 3.276 + /* no difference found, points are equal */ 3.277 + return 0; 3.278 +} 3.279 + 3.280 +/* compare two boxes */ 3.281 +/* (equality when same box on earth is described, otherwise an arbitrary linear 3.282 + order) */ 3.283 +static int pgl_box_cmp(pgl_box *box1, pgl_box *box2) { 3.284 + /* two empty boxes are equal, and an empty box is always considered "less 3.285 + than" a non-empty box */ 3.286 + if (box1->lat_min> box1->lat_max && box2->lat_min<=box2->lat_max) return -1; 3.287 + if (box1->lat_min> box1->lat_max && box2->lat_min> box2->lat_max) return 0; 3.288 + if (box1->lat_min<=box1->lat_max && box2->lat_min> box2->lat_max) return 1; 3.289 + /* use southern border as first ordering criterion */ 3.290 + if (box1->lat_min < box2->lat_min) return -1; 3.291 + if (box1->lat_min > box2->lat_min) return 1; 3.292 + /* use northern border as second ordering criterion */ 3.293 + if (box1->lat_max < box2->lat_max) return -1; 3.294 + if (box1->lat_max > box2->lat_max) return 1; 3.295 + /* use western border as third ordering criterion */ 3.296 + if (box1->lon_min < box2->lon_min) return -1; 3.297 + if (box1->lon_min > box2->lon_min) return 1; 3.298 + /* use eastern border as fourth ordering criterion */ 3.299 + if (box1->lon_max < box2->lon_max) return -1; 3.300 + if (box1->lon_max > box2->lon_max) return 1; 3.301 + /* no difference found, boxes are equal */ 3.302 + return 0; 3.303 +} 3.304 + 3.305 +/* compare two circles */ 3.306 +/* (equality when same circle on earth is described, otherwise an arbitrary 3.307 + linear order) */ 3.308 +static int pgl_circle_cmp(pgl_circle *circle1, pgl_circle *circle2) { 3.309 + /* two circles with same infinite radius (positive or negative infinity) are 3.310 + considered equal independently of center point */ 3.311 + if ( 3.312 + !isfinite(circle1->radius) && !isfinite(circle2->radius) && 3.313 + circle1->radius == circle2->radius 3.314 + ) return 0; 3.315 + /* use radius as first ordering criterion */ 3.316 + if (circle1->radius < circle2->radius) return -1; 3.317 + if (circle1->radius > circle2->radius) return 1; 3.318 + /* use center point as secondary ordering criterion */ 3.319 + return pgl_point_cmp(&(circle1->center), &(circle2->center)); 3.320 +} 3.321 + 3.322 +/* set box to empty box*/ 3.323 +static void pgl_box_set_empty(pgl_box *box) { 3.324 + box->lat_min = INFINITY; 3.325 + box->lat_max = -INFINITY; 3.326 + box->lon_min = 0; 3.327 + box->lon_max = 0; 3.328 +} 3.329 + 3.330 +/* check if point is inside a box */ 3.331 +static bool pgl_point_in_box(pgl_point *point, pgl_box *box) { 3.332 + return ( 3.333 + point->lat >= box->lat_min && point->lat <= box->lat_max && ( 3.334 + (box->lon_min > box->lon_max) ? ( 3.335 + /* box crosses 180th meridian */ 3.336 + point->lon >= box->lon_min || point->lon <= box->lon_max 3.337 + ) : ( 3.338 + /* box does not cross the 180th meridian */ 3.339 + point->lon >= box->lon_min && point->lon <= box->lon_max 3.340 + ) 3.341 + ) 3.342 + ); 3.343 +} 3.344 + 3.345 +/* check if two boxes overlap */ 3.346 +static bool pgl_boxes_overlap(pgl_box *box1, pgl_box *box2) { 3.347 + return ( 3.348 + box2->lat_max >= box2->lat_min && /* ensure box2 is not empty */ 3.349 + ( box2->lat_min >= box1->lat_min || box2->lat_max >= box1->lat_min ) && 3.350 + ( box2->lat_min <= box1->lat_max || box2->lat_max <= box1->lat_max ) && ( 3.351 + ( 3.352 + /* check if one and only one box crosses the 180th meridian */ 3.353 + ((box1->lon_min > box1->lon_max) ? 1 : 0) ^ 3.354 + ((box2->lon_min > box2->lon_max) ? 1 : 0) 3.355 + ) ? ( 3.356 + /* exactly one box crosses the 180th meridian */ 3.357 + box2->lon_min >= box1->lon_min || box2->lon_max >= box1->lon_min || 3.358 + box2->lon_min <= box1->lon_max || box2->lon_max <= box1->lon_max 3.359 + ) : ( 3.360 + /* no box or both boxes cross the 180th meridian */ 3.361 + ( 3.362 + (box2->lon_min >= box1->lon_min || box2->lon_max >= box1->lon_min) && 3.363 + (box2->lon_min <= box1->lon_max || box2->lon_max <= box1->lon_max) 3.364 + ) || 3.365 + /* handle W180 == E180 */ 3.366 + ( box1->lon_min == -180 && box2->lon_max == 180 ) || 3.367 + ( box2->lon_min == -180 && box1->lon_max == 180 ) 3.368 + ) 3.369 + ) 3.370 + ); 3.371 +} 3.372 + 3.373 +/* check unambiguousness of east/west orientation of cluster entries and set 3.374 + bounding circle of cluster */ 3.375 +static bool pgl_finalize_cluster(pgl_cluster *cluster) { 3.376 + int i, j; /* i: index of entry, j: index of point in entry */ 3.377 + int npoints; /* number of points in entry */ 3.378 + int total_npoints = 0; /* total number of points in cluster */ 3.379 + pgl_point *points; /* points in entry */ 3.380 + int lon_dir; /* first point of entry west (-1) or east (+1) */ 3.381 + double lon_break = 0; /* antipodal longitude of first point in entry */ 3.382 + double lon_min, lon_max; /* covered longitude range of entry */ 3.383 + double value; /* temporary variable */ 3.384 + /* reset bounding circle center to empty circle at 0/0 coordinates */ 3.385 + cluster->bounding.center.lat = 0; 3.386 + cluster->bounding.center.lon = 0; 3.387 + cluster->bounding.radius = -INFINITY; 3.388 + /* if cluster is not empty */ 3.389 + if (cluster->nentries != 0) { 3.390 + /* iterate over all cluster entries and ensure they each cover a longitude 3.391 + range less than 180 degrees */ 3.392 + for (i=0; i<cluster->nentries; i++) { 3.393 + /* get properties of entry */ 3.394 + npoints = cluster->entries[i].npoints; 3.395 + points = PGL_ENTRY_POINTS(cluster, i); 3.396 + /* get longitude of first point of entry */ 3.397 + value = points[0].lon; 3.398 + /* initialize lon_min and lon_max with longitude of first point */ 3.399 + lon_min = value; 3.400 + lon_max = value; 3.401 + /* determine east/west orientation of first point and calculate antipodal 3.402 + longitude (Note: rounding required here) */ 3.403 + if (value < 0) { lon_dir = -1; lon_break = pgl_round(value + 180); } 3.404 + else if (value > 0) { lon_dir = 1; lon_break = pgl_round(value - 180); } 3.405 + else lon_dir = 0; 3.406 + /* iterate over all other points in entry */ 3.407 + for (j=1; j<npoints; j++) { 3.408 + /* consider longitude wrap-around */ 3.409 + value = points[j].lon; 3.410 + if (lon_dir<0 && value>lon_break) value = pgl_round(value - 360); 3.411 + else if (lon_dir>0 && value<lon_break) value = pgl_round(value + 360); 3.412 + /* update lon_min and lon_max */ 3.413 + if (value < lon_min) lon_min = value; 3.414 + else if (value > lon_max) lon_max = value; 3.415 + /* return false if 180 degrees or more are covered */ 3.416 + if (lon_max - lon_min >= 180) return false; 3.417 + } 3.418 + } 3.419 + /* iterate over all points of all entries and calculate arbitrary center 3.420 + point for bounding circle (best if center point minimizes the radius, 3.421 + but some error is allowed here) */ 3.422 + for (i=0; i<cluster->nentries; i++) { 3.423 + /* get properties of entry */ 3.424 + npoints = cluster->entries[i].npoints; 3.425 + points = PGL_ENTRY_POINTS(cluster, i); 3.426 + /* check if first entry */ 3.427 + if (i==0) { 3.428 + /* get longitude of first point of first entry in whole cluster */ 3.429 + value = points[0].lon; 3.430 + /* initialize lon_min and lon_max with longitude of first point of 3.431 + first entry in whole cluster (used to determine if whole cluster 3.432 + covers a longitude range of 180 degrees or more) */ 3.433 + lon_min = value; 3.434 + lon_max = value; 3.435 + /* determine east/west orientation of first point and calculate 3.436 + antipodal longitude (Note: rounding not necessary here) */ 3.437 + if (value < 0) { lon_dir = -1; lon_break = value + 180; } 3.438 + else if (value > 0) { lon_dir = 1; lon_break = value - 180; } 3.439 + else lon_dir = 0; 3.440 + } 3.441 + /* iterate over all points in entry */ 3.442 + for (j=0; j<npoints; j++) { 3.443 + /* longitude wrap-around (Note: rounding not necessary here) */ 3.444 + value = points[j].lon; 3.445 + if (lon_dir < 0 && value > lon_break) value -= 360; 3.446 + else if (lon_dir > 0 && value < lon_break) value += 360; 3.447 + if (value < lon_min) lon_min = value; 3.448 + else if (value > lon_max) lon_max = value; 3.449 + /* set bounding circle to cover whole earth if more than 180 degrees 3.450 + are covered */ 3.451 + if (lon_max - lon_min >= 180) { 3.452 + cluster->bounding.center.lat = 0; 3.453 + cluster->bounding.center.lon = 0; 3.454 + cluster->bounding.radius = INFINITY; 3.455 + return true; 3.456 + } 3.457 + /* add point to bounding circle center (for average calculation) */ 3.458 + cluster->bounding.center.lat += points[j].lat; 3.459 + cluster->bounding.center.lon += value; 3.460 + } 3.461 + /* count total number of points */ 3.462 + total_npoints += npoints; 3.463 + } 3.464 + /* determine average latitude and longitude of cluster */ 3.465 + cluster->bounding.center.lat /= total_npoints; 3.466 + cluster->bounding.center.lon /= total_npoints; 3.467 + /* normalize longitude of center of cluster bounding circle */ 3.468 + if (cluster->bounding.center.lon < -180) { 3.469 + cluster->bounding.center.lon += 360; 3.470 + } 3.471 + else if (cluster->bounding.center.lon > 180) { 3.472 + cluster->bounding.center.lon -= 360; 3.473 + } 3.474 + /* round bounding circle center (useful if it is used by other functions) */ 3.475 + cluster->bounding.center.lat = pgl_round(cluster->bounding.center.lat); 3.476 + cluster->bounding.center.lon = pgl_round(cluster->bounding.center.lon); 3.477 + /* calculate radius of bounding circle */ 3.478 + for (i=0; i<cluster->nentries; i++) { 3.479 + npoints = cluster->entries[i].npoints; 3.480 + points = PGL_ENTRY_POINTS(cluster, i); 3.481 + for (j=0; j<npoints; j++) { 3.482 + value = pgl_distance( 3.483 + cluster->bounding.center.lat, cluster->bounding.center.lon, 3.484 + points[j].lat, points[j].lon 3.485 + ); 3.486 + if (value > cluster->bounding.radius) cluster->bounding.radius = value; 3.487 + } 3.488 + } 3.489 + } 3.490 + /* return true (east/west orientation is unambiguous) */ 3.491 + return true; 3.492 +} 3.493 + 3.494 +/* check if point is inside cluster */ 3.495 +static bool pgl_point_in_cluster(pgl_point *point, pgl_cluster *cluster) { 3.496 + int i, j, k; /* i: entry, j: point in entry, k: next point in entry */ 3.497 + int entrytype; /* type of entry */ 3.498 + int npoints; /* number of points in entry */ 3.499 + pgl_point *points; /* array of points in entry */ 3.500 + int lon_dir = 0; /* first vertex west (-1) or east (+1) */ 3.501 + double lon_break = 0; /* antipodal longitude of first vertex */ 3.502 + double lat0 = point->lat; /* latitude of point */ 3.503 + double lon0; /* (adjusted) longitude of point */ 3.504 + double lat1, lon1; /* latitude and (adjusted) longitude of vertex */ 3.505 + double lat2, lon2; /* latitude and (adjusted) longitude of next vertex */ 3.506 + double lon; /* longitude of intersection */ 3.507 + int counter = 0; /* counter for intersections east of point */ 3.508 + /* points outside bounding circle are always assumed to be non-overlapping */ 3.509 + /* (necessary for consistent table and index scans) */ 3.510 + if ( 3.511 + pgl_distance( 3.512 + point->lat, point->lon, 3.513 + cluster->bounding.center.lat, cluster->bounding.center.lon 3.514 + ) > cluster->bounding.radius 3.515 + ) return false; 3.516 + /* iterate over all entries */ 3.517 + for (i=0; i<cluster->nentries; i++) { 3.518 + /* get properties of entry */ 3.519 + entrytype = cluster->entries[i].entrytype; 3.520 + npoints = cluster->entries[i].npoints; 3.521 + points = PGL_ENTRY_POINTS(cluster, i); 3.522 + /* determine east/west orientation of first point of entry and calculate 3.523 + antipodal longitude */ 3.524 + lon_break = points[0].lon; 3.525 + if (lon_break < 0) { lon_dir = -1; lon_break += 180; } 3.526 + else if (lon_break > 0) { lon_dir = 1; lon_break -= 180; } 3.527 + else lon_dir = 0; 3.528 + /* get longitude of point */ 3.529 + lon0 = point->lon; 3.530 + /* consider longitude wrap-around for point */ 3.531 + if (lon_dir < 0 && lon0 > lon_break) lon0 = pgl_round(lon0 - 360); 3.532 + else if (lon_dir > 0 && lon0 < lon_break) lon0 = pgl_round(lon0 + 360); 3.533 + /* iterate over all edges and vertices */ 3.534 + for (j=0; j<npoints; j++) { 3.535 + /* return true if point is on vertex of polygon */ 3.536 + if (pgl_point_cmp(point, &(points[j])) == 0) return true; 3.537 + /* calculate index of next vertex */ 3.538 + k = (j+1) % npoints; 3.539 + /* skip last edge unless entry is (closed) outline or polygon */ 3.540 + if ( 3.541 + k == 0 && 3.542 + entrytype != PGL_ENTRY_OUTLINE && 3.543 + entrytype != PGL_ENTRY_POLYGON 3.544 + ) continue; 3.545 + /* get latitude and longitude values of edge */ 3.546 + lat1 = points[j].lat; 3.547 + lat2 = points[k].lat; 3.548 + lon1 = points[j].lon; 3.549 + lon2 = points[k].lon; 3.550 + /* consider longitude wrap-around for edge */ 3.551 + if (lon_dir < 0 && lon1 > lon_break) lon1 = pgl_round(lon1 - 360); 3.552 + else if (lon_dir > 0 && lon1 < lon_break) lon1 = pgl_round(lon1 + 360); 3.553 + if (lon_dir < 0 && lon2 > lon_break) lon2 = pgl_round(lon2 - 360); 3.554 + else if (lon_dir > 0 && lon2 < lon_break) lon2 = pgl_round(lon2 + 360); 3.555 + /* return true if point is on horizontal (west to east) edge of polygon */ 3.556 + if ( 3.557 + lat0 == lat1 && lat0 == lat2 && 3.558 + ( (lon0 >= lon1 && lon0 <= lon2) || (lon0 >= lon2 && lon0 <= lon1) ) 3.559 + ) return true; 3.560 + /* check if edge crosses east/west line of point */ 3.561 + if ((lat1 < lat0 && lat2 >= lat0) || (lat2 < lat0 && lat1 >= lat0)) { 3.562 + /* calculate longitude of intersection */ 3.563 + lon = (lon1 * (lat2-lat0) + lon2 * (lat0-lat1)) / (lat2-lat1); 3.564 + /* return true if intersection goes (approximately) through point */ 3.565 + if (pgl_round(lon) == lon0) return true; 3.566 + /* count intersection if east of point and entry is polygon*/ 3.567 + if (entrytype == PGL_ENTRY_POLYGON && lon > lon0) counter++; 3.568 + } 3.569 + } 3.570 + } 3.571 + /* return true if number of intersections is odd */ 3.572 + return counter & 1; 3.573 +} 3.574 + 3.575 +/* calculate (approximate) distance between point and cluster */ 3.576 +static double pgl_point_cluster_distance(pgl_point *point, pgl_cluster *cluster) { 3.577 + int i, j, k; /* i: entry, j: point in entry, k: next point in entry */ 3.578 + int entrytype; /* type of entry */ 3.579 + int npoints; /* number of points in entry */ 3.580 + pgl_point *points; /* array of points in entry */ 3.581 + int lon_dir = 0; /* first vertex west (-1) or east (+1) */ 3.582 + double lon_break = 0; /* antipodal longitude of first vertex */ 3.583 + double lon_min = 0; /* minimum (adjusted) longitude of entry vertices */ 3.584 + double lon_max = 0; /* maximum (adjusted) longitude of entry vertices */ 3.585 + double lat0 = point->lat; /* latitude of point */ 3.586 + double lon0; /* (adjusted) longitude of point */ 3.587 + double lat1, lon1; /* latitude and (adjusted) longitude of vertex */ 3.588 + double lat2, lon2; /* latitude and (adjusted) longitude of next vertex */ 3.589 + double s; /* scalar for vector calculations */ 3.590 + double dist; /* distance calculated in one step */ 3.591 + double min_dist = INFINITY; /* minimum distance */ 3.592 + /* distance is zero if point is contained in cluster */ 3.593 + if (pgl_point_in_cluster(point, cluster)) return 0; 3.594 + /* iterate over all entries */ 3.595 + for (i=0; i<cluster->nentries; i++) { 3.596 + /* get properties of entry */ 3.597 + entrytype = cluster->entries[i].entrytype; 3.598 + npoints = cluster->entries[i].npoints; 3.599 + points = PGL_ENTRY_POINTS(cluster, i); 3.600 + /* determine east/west orientation of first point of entry and calculate 3.601 + antipodal longitude */ 3.602 + lon_break = points[0].lon; 3.603 + if (lon_break < 0) { lon_dir = -1; lon_break += 180; } 3.604 + else if (lon_break > 0) { lon_dir = 1; lon_break -= 180; } 3.605 + else lon_dir = 0; 3.606 + /* determine covered longitude range */ 3.607 + for (j=0; j<npoints; j++) { 3.608 + /* get longitude of vertex */ 3.609 + lon1 = points[j].lon; 3.610 + /* adjust longitude to fix potential wrap-around */ 3.611 + if (lon_dir < 0 && lon1 > lon_break) lon1 -= 360; 3.612 + else if (lon_dir > 0 && lon1 < lon_break) lon1 += 360; 3.613 + /* update minimum and maximum longitude of polygon */ 3.614 + if (j == 0 || lon1 < lon_min) lon_min = lon1; 3.615 + if (j == 0 || lon1 > lon_max) lon_max = lon1; 3.616 + } 3.617 + /* adjust longitude wrap-around according to full longitude range */ 3.618 + lon_break = (lon_max + lon_min) / 2; 3.619 + if (lon_break < 0) { lon_dir = -1; lon_break += 180; } 3.620 + else if (lon_break > 0) { lon_dir = 1; lon_break -= 180; } 3.621 + /* get longitude of point */ 3.622 + lon0 = point->lon; 3.623 + /* consider longitude wrap-around for point */ 3.624 + if (lon_dir < 0 && lon0 > lon_break) lon0 -= 360; 3.625 + else if (lon_dir > 0 && lon0 < lon_break) lon0 += 360; 3.626 + /* iterate over all edges and vertices */ 3.627 + for (j=0; j<npoints; j++) { 3.628 + /* get latitude and longitude values of current point */ 3.629 + lat1 = points[j].lat; 3.630 + lon1 = points[j].lon; 3.631 + /* consider longitude wrap-around for current point */ 3.632 + if (lon_dir < 0 && lon1 > lon_break) lon1 -= 360; 3.633 + else if (lon_dir > 0 && lon1 < lon_break) lon1 += 360; 3.634 + /* calculate distance to vertex */ 3.635 + dist = pgl_distance(lat0, lon0, lat1, lon1); 3.636 + /* store calculated distance if smallest */ 3.637 + if (dist < min_dist) min_dist = dist; 3.638 + /* calculate index of next vertex */ 3.639 + k = (j+1) % npoints; 3.640 + /* skip last edge unless entry is (closed) outline or polygon */ 3.641 + if ( 3.642 + k == 0 && 3.643 + entrytype != PGL_ENTRY_OUTLINE && 3.644 + entrytype != PGL_ENTRY_POLYGON 3.645 + ) continue; 3.646 + /* get latitude and longitude values of next point */ 3.647 + lat2 = points[k].lat; 3.648 + lon2 = points[k].lon; 3.649 + /* consider longitude wrap-around for next point */ 3.650 + if (lon_dir < 0 && lon2 > lon_break) lon2 -= 360; 3.651 + else if (lon_dir > 0 && lon2 < lon_break) lon2 += 360; 3.652 + /* go to next vertex and edge if edge is degenerated */ 3.653 + if (lat1 == lat2 && lon1 == lon2) continue; 3.654 + /* otherwise test if point can be projected onto edge of polygon */ 3.655 + s = ( 3.656 + ((lat0-lat1) * (lat2-lat1) + (lon0-lon1) * (lon2-lon1)) / 3.657 + ((lat2-lat1) * (lat2-lat1) + (lon2-lon1) * (lon2-lon1)) 3.658 + ); 3.659 + /* go to next vertex and edge if point cannot be projected */ 3.660 + if (!(s > 0 && s < 1)) continue; 3.661 + /* calculate distance from original point to projected point */ 3.662 + dist = pgl_distance( 3.663 + lat0, lon0, 3.664 + lat1 + s * (lat2-lat1), 3.665 + lon1 + s * (lon2-lon1) 3.666 + ); 3.667 + /* store calculated distance if smallest */ 3.668 + if (dist < min_dist) min_dist = dist; 3.669 + } 3.670 + } 3.671 + /* return minimum distance */ 3.672 + return min_dist; 3.673 +} 3.674 + 3.675 +/* estimator function for distance between box and point */ 3.676 +/* allowed to return smaller values than actually correct */ 3.677 +static double pgl_estimate_point_box_distance(pgl_point *point, pgl_box *box) { 3.678 + double dlon; /* longitude range of box (delta longitude) */ 3.679 + double h; /* half of distance along meridian */ 3.680 + double d; /* distance between both southern or both northern points */ 3.681 + double cur_dist; /* calculated distance */ 3.682 + double min_dist; /* minimum distance calculated */ 3.683 + /* return infinity if bounding box is empty */ 3.684 + if (box->lat_min > box->lat_max) return INFINITY; 3.685 + /* return zero if point is inside bounding box */ 3.686 + if (pgl_point_in_box(point, box)) return 0; 3.687 + /* calculate delta longitude */ 3.688 + dlon = box->lon_max - box->lon_min; 3.689 + if (dlon < 0) dlon += 360; /* 180th meridian crossed */ 3.690 + /* if delta longitude is greater than 180 degrees, perform safe fall-back */ 3.691 + if (dlon > 180) return 0; 3.692 + /* calculate half of distance along meridian */ 3.693 + h = pgl_distance(box->lat_min, 0, box->lat_max, 0) / 2; 3.694 + /* calculate full distance between southern points */ 3.695 + d = pgl_distance(box->lat_min, 0, box->lat_min, dlon); 3.696 + /* calculate maximum of full distance and half distance */ 3.697 + if (h > d) d = h; 3.698 + /* calculate distance from point to first southern vertex and substract 3.699 + maximum error */ 3.700 + min_dist = pgl_distance( 3.701 + point->lat, point->lon, box->lat_min, box->lon_min 3.702 + ) - d; 3.703 + /* return zero if estimated distance is smaller than zero */ 3.704 + if (min_dist <= 0) return 0; 3.705 + /* repeat procedure with second southern vertex */ 3.706 + cur_dist = pgl_distance( 3.707 + point->lat, point->lon, box->lat_min, box->lon_max 3.708 + ) - d; 3.709 + if (cur_dist <= 0) return 0; 3.710 + if (cur_dist < min_dist) min_dist = cur_dist; 3.711 + /* calculate full distance between northern points */ 3.712 + d = pgl_distance(box->lat_max, 0, box->lat_max, dlon); 3.713 + /* calculate maximum of full distance and half distance */ 3.714 + if (h > d) d = h; 3.715 + /* repeat procedure with northern vertices */ 3.716 + cur_dist = pgl_distance( 3.717 + point->lat, point->lon, box->lat_max, box->lon_max 3.718 + ) - d; 3.719 + if (cur_dist <= 0) return 0; 3.720 + if (cur_dist < min_dist) min_dist = cur_dist; 3.721 + cur_dist = pgl_distance( 3.722 + point->lat, point->lon, box->lat_max, box->lon_min 3.723 + ) - d; 3.724 + if (cur_dist <= 0) return 0; 3.725 + if (cur_dist < min_dist) min_dist = cur_dist; 3.726 + /* return smallest value (unless already returned zero) */ 3.727 + return min_dist; 3.728 +} 3.729 + 3.730 + 3.731 +/*----------------------------* 3.732 + * fractal geographic index * 3.733 + *----------------------------*/ 3.734 + 3.735 +/* number of bytes used for geographic (center) position in keys */ 3.736 +#define PGL_KEY_LATLON_BYTELEN 7 3.737 + 3.738 +/* maximum reference value for logarithmic size of geographic objects */ 3.739 +#define PGL_AREAKEY_REFOBJSIZE (PGL_DIAMETER/3.0) /* can be tweaked */ 3.740 + 3.741 +/* safety margin to avoid floating point errors in distance estimation */ 3.742 +#define PGL_FPE_SAFETY (1.0+1e-14) /* slightly greater than 1.0 */ 3.743 + 3.744 +/* pointer to index key (either pgl_pointkey or pgl_areakey) */ 3.745 +typedef unsigned char *pgl_keyptr; 3.746 + 3.747 +/* index key for points (objects with zero area) on the spheroid */ 3.748 +/* bit 0..55: interspersed bits of latitude and longitude, 3.749 + bit 56..57: always zero, 3.750 + bit 58..63: node depth in hypothetic (full) tree from 0 to 56 (incl.) */ 3.751 +typedef unsigned char pgl_pointkey[PGL_KEY_LATLON_BYTELEN+1]; 3.752 + 3.753 +/* index key for geographic objects on spheroid with area greater than zero */ 3.754 +/* bit 0..55: interspersed bits of latitude and longitude of center point, 3.755 + bit 56: always set to 1, 3.756 + bit 57..63: node depth in hypothetic (full) tree from 0 to (2*56)+1 (incl.), 3.757 + bit 64..71: logarithmic object size from 0 to 56+1 = 57 (incl.), but set to 3.758 + PGL_KEY_OBJSIZE_EMPTY (with interspersed bits = 0 and node depth 3.759 + = 113) for empty objects, and set to PGL_KEY_OBJSIZE_UNIVERSAL 3.760 + (with interspersed bits = 0 and node depth = 0) for keys which 3.761 + cover both empty and non-empty objects */ 3.762 + 3.763 +typedef unsigned char pgl_areakey[PGL_KEY_LATLON_BYTELEN+2]; 3.764 + 3.765 +/* helper macros for reading/writing index keys */ 3.766 +#define PGL_KEY_NODEDEPTH_OFFSET PGL_KEY_LATLON_BYTELEN 3.767 +#define PGL_KEY_OBJSIZE_OFFSET (PGL_KEY_NODEDEPTH_OFFSET+1) 3.768 +#define PGL_POINTKEY_MAXDEPTH (PGL_KEY_LATLON_BYTELEN*8) 3.769 +#define PGL_AREAKEY_MAXDEPTH (2*PGL_POINTKEY_MAXDEPTH+1) 3.770 +#define PGL_AREAKEY_MAXOBJSIZE (PGL_POINTKEY_MAXDEPTH+1) 3.771 +#define PGL_AREAKEY_TYPEMASK 0x80 3.772 +#define PGL_KEY_LATLONBIT(key, n) ((key)[(n)/8] & (0x80 >> ((n)%8))) 3.773 +#define PGL_KEY_LATLONBIT_DIFF(key1, key2, n) \ 3.774 + ( PGL_KEY_LATLONBIT(key1, n) ^ \ 3.775 + PGL_KEY_LATLONBIT(key2, n) ) 3.776 +#define PGL_KEY_IS_AREAKEY(key) ((key)[PGL_KEY_NODEDEPTH_OFFSET] & \ 3.777 + PGL_AREAKEY_TYPEMASK) 3.778 +#define PGL_KEY_NODEDEPTH(key) ((key)[PGL_KEY_NODEDEPTH_OFFSET] & \ 3.779 + (PGL_AREAKEY_TYPEMASK-1)) 3.780 +#define PGL_KEY_OBJSIZE(key) ((key)[PGL_KEY_OBJSIZE_OFFSET]) 3.781 +#define PGL_KEY_OBJSIZE_EMPTY 126 3.782 +#define PGL_KEY_OBJSIZE_UNIVERSAL 127 3.783 +#define PGL_KEY_IS_EMPTY(key) ( PGL_KEY_IS_AREAKEY(key) && \ 3.784 + (key)[PGL_KEY_OBJSIZE_OFFSET] == \ 3.785 + PGL_KEY_OBJSIZE_EMPTY ) 3.786 +#define PGL_KEY_IS_UNIVERSAL(key) ( PGL_KEY_IS_AREAKEY(key) && \ 3.787 + (key)[PGL_KEY_OBJSIZE_OFFSET] == \ 3.788 + PGL_KEY_OBJSIZE_UNIVERSAL ) 3.789 + 3.790 +/* set area key to match empty objects only */ 3.791 +static void pgl_key_set_empty(pgl_keyptr key) { 3.792 + memset(key, 0, sizeof(pgl_areakey)); 3.793 + /* Note: setting node depth to maximum is required for picksplit function */ 3.794 + key[PGL_KEY_NODEDEPTH_OFFSET] = PGL_AREAKEY_TYPEMASK | PGL_AREAKEY_MAXDEPTH; 3.795 + key[PGL_KEY_OBJSIZE_OFFSET] = PGL_KEY_OBJSIZE_EMPTY; 3.796 +} 3.797 + 3.798 +/* set area key to match any object (including empty objects) */ 3.799 +static void pgl_key_set_universal(pgl_keyptr key) { 3.800 + memset(key, 0, sizeof(pgl_areakey)); 3.801 + key[PGL_KEY_NODEDEPTH_OFFSET] = PGL_AREAKEY_TYPEMASK; 3.802 + key[PGL_KEY_OBJSIZE_OFFSET] = PGL_KEY_OBJSIZE_UNIVERSAL; 3.803 +} 3.804 + 3.805 +/* convert a point on earth into a max-depth key to be used in index */ 3.806 +static void pgl_point_to_key(pgl_point *point, pgl_keyptr key) { 3.807 + double lat = point->lat; 3.808 + double lon = point->lon; 3.809 + int i; 3.810 + /* clear latitude and longitude bits */ 3.811 + memset(key, 0, PGL_KEY_LATLON_BYTELEN); 3.812 + /* set node depth to maximum and type bit to zero */ 3.813 + key[PGL_KEY_NODEDEPTH_OFFSET] = PGL_POINTKEY_MAXDEPTH; 3.814 + /* iterate over all latitude/longitude bit pairs */ 3.815 + for (i=0; i<PGL_POINTKEY_MAXDEPTH/2; i++) { 3.816 + /* determine latitude bit */ 3.817 + if (lat >= 0) { 3.818 + key[i/4] |= 0x80 >> (2*(i%4)); 3.819 + lat *= 2; lat -= 90; 3.820 + } else { 3.821 + lat *= 2; lat += 90; 3.822 + } 3.823 + /* determine longitude bit */ 3.824 + if (lon >= 0) { 3.825 + key[i/4] |= 0x80 >> (2*(i%4)+1); 3.826 + lon *= 2; lon -= 180; 3.827 + } else { 3.828 + lon *= 2; lon += 180; 3.829 + } 3.830 + } 3.831 +} 3.832 + 3.833 +/* convert a circle on earth into a max-depth key to be used in an index */ 3.834 +static void pgl_circle_to_key(pgl_circle *circle, pgl_keyptr key) { 3.835 + /* handle special case of empty circle */ 3.836 + if (circle->radius < 0) { 3.837 + pgl_key_set_empty(key); 3.838 + return; 3.839 + } 3.840 + /* perform same action as for point keys */ 3.841 + pgl_point_to_key(&(circle->center), key); 3.842 + /* but overwrite type and node depth to fit area index key */ 3.843 + key[PGL_KEY_NODEDEPTH_OFFSET] = PGL_AREAKEY_TYPEMASK | PGL_AREAKEY_MAXDEPTH; 3.844 + /* check if radius is greater than (or equal to) reference size */ 3.845 + /* (treat equal values as greater values for numerical safety) */ 3.846 + if (circle->radius >= PGL_AREAKEY_REFOBJSIZE) { 3.847 + /* if yes, set logarithmic size to zero */ 3.848 + key[PGL_KEY_OBJSIZE_OFFSET] = 0; 3.849 + } else { 3.850 + /* otherwise, determine logarithmic size iteratively */ 3.851 + /* (one step is equivalent to a factor of sqrt(2)) */ 3.852 + double reference = PGL_AREAKEY_REFOBJSIZE / M_SQRT2; 3.853 + int objsize = 1; 3.854 + while (objsize < PGL_AREAKEY_MAXOBJSIZE) { 3.855 + /* stop when radius is greater than (or equal to) adjusted reference */ 3.856 + /* (treat equal values as greater values for numerical safety) */ 3.857 + if (circle->radius >= reference) break; 3.858 + reference /= M_SQRT2; 3.859 + objsize++; 3.860 + } 3.861 + /* set logarithmic size to determined value */ 3.862 + key[PGL_KEY_OBJSIZE_OFFSET] = objsize; 3.863 + } 3.864 +} 3.865 + 3.866 +/* check if one key is subkey of another key or vice versa */ 3.867 +static bool pgl_keys_overlap(pgl_keyptr key1, pgl_keyptr key2) { 3.868 + int i; /* key bit offset (includes both lat/lon and log. obj. size bits) */ 3.869 + /* determine smallest depth */ 3.870 + int depth1 = PGL_KEY_NODEDEPTH(key1); 3.871 + int depth2 = PGL_KEY_NODEDEPTH(key2); 3.872 + int depth = (depth1 < depth2) ? depth1 : depth2; 3.873 + /* check if keys are area keys (assuming that both keys have same type) */ 3.874 + if (PGL_KEY_IS_AREAKEY(key1)) { 3.875 + int j = 0; /* bit offset for logarithmic object size bits */ 3.876 + int k = 0; /* bit offset for latitude and longitude */ 3.877 + /* fetch logarithmic object size information */ 3.878 + int objsize1 = PGL_KEY_OBJSIZE(key1); 3.879 + int objsize2 = PGL_KEY_OBJSIZE(key2); 3.880 + /* handle special cases for empty objects (universal and empty keys) */ 3.881 + if ( 3.882 + objsize1 == PGL_KEY_OBJSIZE_UNIVERSAL || 3.883 + objsize2 == PGL_KEY_OBJSIZE_UNIVERSAL 3.884 + ) return true; 3.885 + if ( 3.886 + objsize1 == PGL_KEY_OBJSIZE_EMPTY || 3.887 + objsize2 == PGL_KEY_OBJSIZE_EMPTY 3.888 + ) return objsize1 == objsize2; 3.889 + /* iterate through key bits */ 3.890 + for (i=0; i<depth; i++) { 3.891 + /* every second bit is a bit describing the object size */ 3.892 + if (i%2 == 0) { 3.893 + /* check if object size bit is different in both keys (objsize1 and 3.894 + objsize2 describe the minimum index when object size bit is set) */ 3.895 + if ( 3.896 + (objsize1 <= j && objsize2 > j) || 3.897 + (objsize2 <= j && objsize1 > j) 3.898 + ) { 3.899 + /* bit differs, therefore keys are in separate branches */ 3.900 + return false; 3.901 + } 3.902 + /* increase bit counter for object size bits */ 3.903 + j++; 3.904 + } 3.905 + /* all other bits describe latitude and longitude */ 3.906 + else { 3.907 + /* check if bit differs in both keys */ 3.908 + if (PGL_KEY_LATLONBIT_DIFF(key1, key2, k)) { 3.909 + /* bit differs, therefore keys are in separate branches */ 3.910 + return false; 3.911 + } 3.912 + /* increase bit counter for latitude/longitude bits */ 3.913 + k++; 3.914 + } 3.915 + } 3.916 + } 3.917 + /* if not, keys are point keys */ 3.918 + else { 3.919 + /* iterate through key bits */ 3.920 + for (i=0; i<depth; i++) { 3.921 + /* check if bit differs in both keys */ 3.922 + if (PGL_KEY_LATLONBIT_DIFF(key1, key2, i)) { 3.923 + /* bit differs, therefore keys are in separate branches */ 3.924 + return false; 3.925 + } 3.926 + } 3.927 + } 3.928 + /* return true because keys are in the same branch */ 3.929 + return true; 3.930 +} 3.931 + 3.932 +/* combine two keys into new key which covers both original keys */ 3.933 +/* (result stored in first argument) */ 3.934 +static void pgl_unite_keys(pgl_keyptr dst, pgl_keyptr src) { 3.935 + int i; /* key bit offset (includes both lat/lon and log. obj. size bits) */ 3.936 + /* determine smallest depth */ 3.937 + int depth1 = PGL_KEY_NODEDEPTH(dst); 3.938 + int depth2 = PGL_KEY_NODEDEPTH(src); 3.939 + int depth = (depth1 < depth2) ? depth1 : depth2; 3.940 + /* check if keys are area keys (assuming that both keys have same type) */ 3.941 + if (PGL_KEY_IS_AREAKEY(dst)) { 3.942 + pgl_areakey dstbuf = { 0, }; /* destination buffer (cleared) */ 3.943 + int j = 0; /* bit offset for logarithmic object size bits */ 3.944 + int k = 0; /* bit offset for latitude and longitude */ 3.945 + /* fetch logarithmic object size information */ 3.946 + int objsize1 = PGL_KEY_OBJSIZE(dst); 3.947 + int objsize2 = PGL_KEY_OBJSIZE(src); 3.948 + /* handle special cases for empty objects (universal and empty keys) */ 3.949 + if ( 3.950 + objsize1 > PGL_AREAKEY_MAXOBJSIZE || 3.951 + objsize2 > PGL_AREAKEY_MAXOBJSIZE 3.952 + ) { 3.953 + if ( 3.954 + objsize1 == PGL_KEY_OBJSIZE_EMPTY && 3.955 + objsize2 == PGL_KEY_OBJSIZE_EMPTY 3.956 + ) pgl_key_set_empty(dst); 3.957 + else pgl_key_set_universal(dst); 3.958 + return; 3.959 + } 3.960 + /* iterate through key bits */ 3.961 + for (i=0; i<depth; i++) { 3.962 + /* every second bit is a bit describing the object size */ 3.963 + if (i%2 == 0) { 3.964 + /* increase bit counter for object size bits first */ 3.965 + /* (handy when setting objsize variable) */ 3.966 + j++; 3.967 + /* check if object size bit is set in neither key */ 3.968 + if (objsize1 >= j && objsize2 >= j) { 3.969 + /* set objsize in destination buffer to indicate that size bit is 3.970 + unset in destination buffer at the current bit position */ 3.971 + dstbuf[PGL_KEY_OBJSIZE_OFFSET] = j; 3.972 + } 3.973 + /* break if object size bit is set in one key only */ 3.974 + else if (objsize1 >= j || objsize2 >= j) break; 3.975 + } 3.976 + /* all other bits describe latitude and longitude */ 3.977 + else { 3.978 + /* break if bit differs in both keys */ 3.979 + if (PGL_KEY_LATLONBIT(dst, k)) { 3.980 + if (!PGL_KEY_LATLONBIT(src, k)) break; 3.981 + /* but set bit in destination buffer if bit is set in both keys */ 3.982 + dstbuf[k/8] |= 0x80 >> (k%8); 3.983 + } else if (PGL_KEY_LATLONBIT(src, k)) break; 3.984 + /* increase bit counter for latitude/longitude bits */ 3.985 + k++; 3.986 + } 3.987 + } 3.988 + /* set common node depth and type bit (type bit = 1) */ 3.989 + dstbuf[PGL_KEY_NODEDEPTH_OFFSET] = PGL_AREAKEY_TYPEMASK | i; 3.990 + /* copy contents of destination buffer to first key */ 3.991 + memcpy(dst, dstbuf, sizeof(pgl_areakey)); 3.992 + } 3.993 + /* if not, keys are point keys */ 3.994 + else { 3.995 + pgl_pointkey dstbuf = { 0, }; /* destination buffer (cleared) */ 3.996 + /* iterate through key bits */ 3.997 + for (i=0; i<depth; i++) { 3.998 + /* break if bit differs in both keys */ 3.999 + if (PGL_KEY_LATLONBIT(dst, i)) { 3.1000 + if (!PGL_KEY_LATLONBIT(src, i)) break; 3.1001 + /* but set bit in destination buffer if bit is set in both keys */ 3.1002 + dstbuf[i/8] |= 0x80 >> (i%8); 3.1003 + } else if (PGL_KEY_LATLONBIT(src, i)) break; 3.1004 + } 3.1005 + /* set common node depth (type bit = 0) */ 3.1006 + dstbuf[PGL_KEY_NODEDEPTH_OFFSET] = i; 3.1007 + /* copy contents of destination buffer to first key */ 3.1008 + memcpy(dst, dstbuf, sizeof(pgl_pointkey)); 3.1009 + } 3.1010 +} 3.1011 + 3.1012 +/* determine center(!) boundaries and radius estimation of index key */ 3.1013 +static double pgl_key_to_box(pgl_keyptr key, pgl_box *box) { 3.1014 + int i; 3.1015 + /* determine node depth */ 3.1016 + int depth = PGL_KEY_NODEDEPTH(key); 3.1017 + /* center point of possible result */ 3.1018 + double lat = 0; 3.1019 + double lon = 0; 3.1020 + /* maximum distance of real center point from key center */ 3.1021 + double dlat = 90; 3.1022 + double dlon = 180; 3.1023 + /* maximum radius of contained objects */ 3.1024 + double radius = 0; /* always return zero for point index keys */ 3.1025 + /* check if key is area key */ 3.1026 + if (PGL_KEY_IS_AREAKEY(key)) { 3.1027 + /* get logarithmic object size */ 3.1028 + int objsize = PGL_KEY_OBJSIZE(key); 3.1029 + /* handle special cases for empty objects (universal and empty keys) */ 3.1030 + if (objsize == PGL_KEY_OBJSIZE_EMPTY) { 3.1031 + pgl_box_set_empty(box); 3.1032 + return 0; 3.1033 + } else if (objsize == PGL_KEY_OBJSIZE_UNIVERSAL) { 3.1034 + box->lat_min = -90; 3.1035 + box->lat_max = 90; 3.1036 + box->lon_min = -180; 3.1037 + box->lon_max = 180; 3.1038 + return 0; /* any value >= 0 would do */ 3.1039 + } 3.1040 + /* calculate maximum possible radius of objects covered by the given key */ 3.1041 + if (objsize == 0) radius = INFINITY; 3.1042 + else { 3.1043 + radius = PGL_AREAKEY_REFOBJSIZE; 3.1044 + while (--objsize) radius /= M_SQRT2; 3.1045 + } 3.1046 + /* iterate over latitude and longitude bits in key */ 3.1047 + /* (every second bit is a latitude or longitude bit) */ 3.1048 + for (i=0; i<depth/2; i++) { 3.1049 + /* check if latitude bit */ 3.1050 + if (i%2 == 0) { 3.1051 + /* cut latitude dimension in half */ 3.1052 + dlat /= 2; 3.1053 + /* increase center latitude if bit is 1, otherwise decrease */ 3.1054 + if (PGL_KEY_LATLONBIT(key, i)) lat += dlat; 3.1055 + else lat -= dlat; 3.1056 + } 3.1057 + /* otherwise longitude bit */ 3.1058 + else { 3.1059 + /* cut longitude dimension in half */ 3.1060 + dlon /= 2; 3.1061 + /* increase center longitude if bit is 1, otherwise decrease */ 3.1062 + if (PGL_KEY_LATLONBIT(key, i)) lon += dlon; 3.1063 + else lon -= dlon; 3.1064 + } 3.1065 + } 3.1066 + } 3.1067 + /* if not, keys are point keys */ 3.1068 + else { 3.1069 + /* iterate over all bits in key */ 3.1070 + for (i=0; i<depth; i++) { 3.1071 + /* check if latitude bit */ 3.1072 + if (i%2 == 0) { 3.1073 + /* cut latitude dimension in half */ 3.1074 + dlat /= 2; 3.1075 + /* increase center latitude if bit is 1, otherwise decrease */ 3.1076 + if (PGL_KEY_LATLONBIT(key, i)) lat += dlat; 3.1077 + else lat -= dlat; 3.1078 + } 3.1079 + /* otherwise longitude bit */ 3.1080 + else { 3.1081 + /* cut longitude dimension in half */ 3.1082 + dlon /= 2; 3.1083 + /* increase center longitude if bit is 1, otherwise decrease */ 3.1084 + if (PGL_KEY_LATLONBIT(key, i)) lon += dlon; 3.1085 + else lon -= dlon; 3.1086 + } 3.1087 + } 3.1088 + } 3.1089 + /* calculate boundaries from center point and remaining dlat and dlon */ 3.1090 + /* (return values through pointer to box) */ 3.1091 + box->lat_min = lat - dlat; 3.1092 + box->lat_max = lat + dlat; 3.1093 + box->lon_min = lon - dlon; 3.1094 + box->lon_max = lon + dlon; 3.1095 + /* return radius (as a function return value) */ 3.1096 + return radius; 3.1097 +} 3.1098 + 3.1099 +/* estimator function for distance between point and index key */ 3.1100 +/* allowed to return smaller values than actually correct */ 3.1101 +static double pgl_estimate_key_distance(pgl_keyptr key, pgl_point *point) { 3.1102 + pgl_box box; /* center(!) bounding box of area index key */ 3.1103 + /* calculate center(!) bounding box and maximum radius of objects covered 3.1104 + by area index key (radius is zero for point index keys) */ 3.1105 + double distance = pgl_key_to_box(key, &box); 3.1106 + /* calculate estimated distance between bounding box of center point of 3.1107 + indexed object and point passed as second argument, then substract maximum 3.1108 + radius of objects covered by index key */ 3.1109 + /* (use PGL_FPE_SAFETY factor to cope with minor floating point errors) */ 3.1110 + distance = ( 3.1111 + pgl_estimate_point_box_distance(point, &box) / PGL_FPE_SAFETY - 3.1112 + distance * PGL_FPE_SAFETY 3.1113 + ); 3.1114 + /* truncate negative results to zero */ 3.1115 + if (distance <= 0) distance = 0; 3.1116 + /* return result */ 3.1117 + return distance; 3.1118 +} 3.1119 + 3.1120 + 3.1121 +/*---------------------------------* 3.1122 + * helper functions for text I/O * 3.1123 + *---------------------------------*/ 3.1124 + 3.1125 +#define PGL_NUMBUFLEN 64 /* buffer size for number to string conversion */ 3.1126 + 3.1127 +/* convert floating point number to string (round-trip safe) */ 3.1128 +static void pgl_print_float(char *buf, double flt) { 3.1129 + /* check if number is integral */ 3.1130 + if (trunc(flt) == flt) { 3.1131 + /* for integral floats use maximum precision */ 3.1132 + snprintf(buf, PGL_NUMBUFLEN, "%.17g", flt); 3.1133 + } else { 3.1134 + /* otherwise check if 15, 16, or 17 digits needed (round-trip safety) */ 3.1135 + snprintf(buf, PGL_NUMBUFLEN, "%.15g", flt); 3.1136 + if (strtod(buf, NULL) != flt) snprintf(buf, PGL_NUMBUFLEN, "%.16g", flt); 3.1137 + if (strtod(buf, NULL) != flt) snprintf(buf, PGL_NUMBUFLEN, "%.17g", flt); 3.1138 + } 3.1139 +} 3.1140 + 3.1141 +/* convert latitude floating point number (in degrees) to string */ 3.1142 +static void pgl_print_lat(char *buf, double lat) { 3.1143 + if (signbit(lat)) { 3.1144 + /* treat negative latitudes (including -0) as south */ 3.1145 + snprintf(buf, PGL_NUMBUFLEN, "S%015.12f", -lat); 3.1146 + } else { 3.1147 + /* treat positive latitudes (including +0) as north */ 3.1148 + snprintf(buf, PGL_NUMBUFLEN, "N%015.12f", lat); 3.1149 + } 3.1150 +} 3.1151 + 3.1152 +/* convert longitude floating point number (in degrees) to string */ 3.1153 +static void pgl_print_lon(char *buf, double lon) { 3.1154 + if (signbit(lon)) { 3.1155 + /* treat negative longitudes (including -0) as west */ 3.1156 + snprintf(buf, PGL_NUMBUFLEN, "W%016.12f", -lon); 3.1157 + } else { 3.1158 + /* treat positive longitudes (including +0) as east */ 3.1159 + snprintf(buf, PGL_NUMBUFLEN, "E%016.12f", lon); 3.1160 + } 3.1161 +} 3.1162 + 3.1163 +/* bit masks used as return value of pgl_scan() function */ 3.1164 +#define PGL_SCAN_NONE 0 /* no value has been parsed */ 3.1165 +#define PGL_SCAN_LAT (1<<0) /* latitude has been parsed */ 3.1166 +#define PGL_SCAN_LON (1<<1) /* longitude has been parsed */ 3.1167 +#define PGL_SCAN_LATLON (PGL_SCAN_LAT | PGL_SCAN_LON) /* bitwise OR of both */ 3.1168 + 3.1169 +/* parse a coordinate (can be latitude or longitude) */ 3.1170 +static int pgl_scan(char **str, double *lat, double *lon) { 3.1171 + double val; 3.1172 + int len; 3.1173 + if ( 3.1174 + sscanf(*str, " N %lf %n", &val, &len) || 3.1175 + sscanf(*str, " n %lf %n", &val, &len) 3.1176 + ) { 3.1177 + *str += len; *lat = val; return PGL_SCAN_LAT; 3.1178 + } 3.1179 + if ( 3.1180 + sscanf(*str, " S %lf %n", &val, &len) || 3.1181 + sscanf(*str, " s %lf %n", &val, &len) 3.1182 + ) { 3.1183 + *str += len; *lat = -val; return PGL_SCAN_LAT; 3.1184 + } 3.1185 + if ( 3.1186 + sscanf(*str, " E %lf %n", &val, &len) || 3.1187 + sscanf(*str, " e %lf %n", &val, &len) 3.1188 + ) { 3.1189 + *str += len; *lon = val; return PGL_SCAN_LON; 3.1190 + } 3.1191 + if ( 3.1192 + sscanf(*str, " W %lf %n", &val, &len) || 3.1193 + sscanf(*str, " w %lf %n", &val, &len) 3.1194 + ) { 3.1195 + *str += len; *lon = -val; return PGL_SCAN_LON; 3.1196 + } 3.1197 + return PGL_SCAN_NONE; 3.1198 +} 3.1199 + 3.1200 + 3.1201 +/*-----------------* 3.1202 + * SQL functions * 3.1203 + *-----------------*/ 3.1204 + 3.1205 +/* Note: These function names use "epoint", "ebox", etc. notation here instead 3.1206 + of "point", "box", etc. in order to distinguish them from any previously 3.1207 + defined functions. */ 3.1208 + 3.1209 +/* function needed for dummy types and/or not implemented features */ 3.1210 +PG_FUNCTION_INFO_V1(pgl_notimpl); 3.1211 +Datum pgl_notimpl(PG_FUNCTION_ARGS) { 3.1212 + ereport(ERROR, (errmsg("not implemented by pgLatLon"))); 3.1213 +} 3.1214 + 3.1215 +/* set point to latitude and longitude (including checks) */ 3.1216 +static void pgl_epoint_set_latlon(pgl_point *point, double lat, double lon) { 3.1217 + /* reject infinite or NaN values */ 3.1218 + if (!isfinite(lat) || !isfinite(lon)) { 3.1219 + ereport(ERROR, ( 3.1220 + errcode(ERRCODE_DATA_EXCEPTION), 3.1221 + errmsg("epoint requires finite coordinates") 3.1222 + )); 3.1223 + } 3.1224 + /* check latitude bounds */ 3.1225 + if (lat < -90) { 3.1226 + ereport(WARNING, (errmsg("latitude exceeds south pole"))); 3.1227 + lat = -90; 3.1228 + } else if (lat > 90) { 3.1229 + ereport(WARNING, (errmsg("latitude exceeds north pole"))); 3.1230 + lat = 90; 3.1231 + } 3.1232 + /* check longitude bounds */ 3.1233 + if (lon < -180) { 3.1234 + ereport(NOTICE, (errmsg("longitude west of 180th meridian normalized"))); 3.1235 + lon += 360 - trunc(lon / 360) * 360; 3.1236 + } else if (lon > 180) { 3.1237 + ereport(NOTICE, (errmsg("longitude east of 180th meridian normalized"))); 3.1238 + lon -= 360 + trunc(lon / 360) * 360; 3.1239 + } 3.1240 + /* store rounded latitude/longitude values for round-trip safety */ 3.1241 + point->lat = pgl_round(lat); 3.1242 + point->lon = pgl_round(lon); 3.1243 +} 3.1244 + 3.1245 +/* create point ("epoint" in SQL) from latitude and longitude */ 3.1246 +PG_FUNCTION_INFO_V1(pgl_create_epoint); 3.1247 +Datum pgl_create_epoint(PG_FUNCTION_ARGS) { 3.1248 + pgl_point *point = (pgl_point *)palloc(sizeof(pgl_point)); 3.1249 + pgl_epoint_set_latlon(point, PG_GETARG_FLOAT8(0), PG_GETARG_FLOAT8(1)); 3.1250 + PG_RETURN_POINTER(point); 3.1251 +} 3.1252 + 3.1253 +/* parse point ("epoint" in SQL) */ 3.1254 +/* format: '[NS]<float> [EW]<float>' */ 3.1255 +PG_FUNCTION_INFO_V1(pgl_epoint_in); 3.1256 +Datum pgl_epoint_in(PG_FUNCTION_ARGS) { 3.1257 + char *str = PG_GETARG_CSTRING(0); /* input string */ 3.1258 + char *strptr = str; /* current position within string */ 3.1259 + int done = 0; /* bit mask storing if latitude or longitude was read */ 3.1260 + double lat, lon; /* parsed values as double precision floats */ 3.1261 + pgl_point *point; /* return value (to be palloc'ed) */ 3.1262 + /* parse two floats (each latitude or longitude) separated by white-space */ 3.1263 + done |= pgl_scan(&strptr, &lat, &lon); 3.1264 + if (strptr != str && isspace(strptr[-1])) { 3.1265 + done |= pgl_scan(&strptr, &lat, &lon); 3.1266 + } 3.1267 + /* require end of string, and latitude and longitude parsed successfully */ 3.1268 + if (strptr[0] || done != PGL_SCAN_LATLON) { 3.1269 + ereport(ERROR, ( 3.1270 + errcode(ERRCODE_INVALID_TEXT_REPRESENTATION), 3.1271 + errmsg("invalid input syntax for type epoint: \"%s\"", str) 3.1272 + )); 3.1273 + } 3.1274 + /* allocate memory for result */ 3.1275 + point = (pgl_point *)palloc(sizeof(pgl_point)); 3.1276 + /* set latitude and longitude (and perform checks) */ 3.1277 + pgl_epoint_set_latlon(point, lat, lon); 3.1278 + /* return result */ 3.1279 + PG_RETURN_POINTER(point); 3.1280 +} 3.1281 + 3.1282 +/* create box ("ebox" in SQL) that is empty */ 3.1283 +PG_FUNCTION_INFO_V1(pgl_create_empty_ebox); 3.1284 +Datum pgl_create_empty_ebox(PG_FUNCTION_ARGS) { 3.1285 + pgl_box *box = (pgl_box *)palloc(sizeof(pgl_box)); 3.1286 + pgl_box_set_empty(box); 3.1287 + PG_RETURN_POINTER(box); 3.1288 +} 3.1289 + 3.1290 +/* set box to given boundaries (including checks) */ 3.1291 +static void pgl_ebox_set_boundaries( 3.1292 + pgl_box *box, 3.1293 + double lat_min, double lat_max, double lon_min, double lon_max 3.1294 +) { 3.1295 + /* if minimum latitude is greater than maximum latitude, return empty box */ 3.1296 + if (lat_min > lat_max) { 3.1297 + pgl_box_set_empty(box); 3.1298 + return; 3.1299 + } 3.1300 + /* otherwise reject infinite or NaN values */ 3.1301 + if ( 3.1302 + !isfinite(lat_min) || !isfinite(lat_max) || 3.1303 + !isfinite(lon_min) || !isfinite(lon_max) 3.1304 + ) { 3.1305 + ereport(ERROR, ( 3.1306 + errcode(ERRCODE_DATA_EXCEPTION), 3.1307 + errmsg("ebox requires finite coordinates") 3.1308 + )); 3.1309 + } 3.1310 + /* check latitude bounds */ 3.1311 + if (lat_max < -90) { 3.1312 + ereport(WARNING, (errmsg("northern latitude exceeds south pole"))); 3.1313 + lat_max = -90; 3.1314 + } else if (lat_max > 90) { 3.1315 + ereport(WARNING, (errmsg("northern latitude exceeds north pole"))); 3.1316 + lat_max = 90; 3.1317 + } 3.1318 + if (lat_min < -90) { 3.1319 + ereport(WARNING, (errmsg("southern latitude exceeds south pole"))); 3.1320 + lat_min = -90; 3.1321 + } else if (lat_min > 90) { 3.1322 + ereport(WARNING, (errmsg("southern latitude exceeds north pole"))); 3.1323 + lat_min = 90; 3.1324 + } 3.1325 + /* check if all longitudes are included */ 3.1326 + if (lon_max - lon_min >= 360) { 3.1327 + if (lon_max - lon_min > 360) ereport(WARNING, ( 3.1328 + errmsg("longitude coverage greater than 360 degrees") 3.1329 + )); 3.1330 + lon_min = -180; 3.1331 + lon_max = 180; 3.1332 + } else { 3.1333 + /* normalize longitude bounds */ 3.1334 + if (lon_min < -180) lon_min += 360 - trunc(lon_min / 360) * 360; 3.1335 + else if (lon_min > 180) lon_min -= 360 + trunc(lon_min / 360) * 360; 3.1336 + if (lon_max < -180) lon_max += 360 - trunc(lon_max / 360) * 360; 3.1337 + else if (lon_max > 180) lon_max -= 360 + trunc(lon_max / 360) * 360; 3.1338 + } 3.1339 + /* store rounded latitude/longitude values for round-trip safety */ 3.1340 + box->lat_min = pgl_round(lat_min); 3.1341 + box->lat_max = pgl_round(lat_max); 3.1342 + box->lon_min = pgl_round(lon_min); 3.1343 + box->lon_max = pgl_round(lon_max); 3.1344 + /* ensure that rounding does not change orientation */ 3.1345 + if (lon_min > lon_max && box->lon_min == box->lon_max) { 3.1346 + box->lon_min = -180; 3.1347 + box->lon_max = 180; 3.1348 + } 3.1349 +} 3.1350 + 3.1351 +/* create box ("ebox" in SQL) from min/max latitude and min/max longitude */ 3.1352 +PG_FUNCTION_INFO_V1(pgl_create_ebox); 3.1353 +Datum pgl_create_ebox(PG_FUNCTION_ARGS) { 3.1354 + pgl_box *box = (pgl_box *)palloc(sizeof(pgl_box)); 3.1355 + pgl_ebox_set_boundaries( 3.1356 + box, 3.1357 + PG_GETARG_FLOAT8(0), PG_GETARG_FLOAT8(1), 3.1358 + PG_GETARG_FLOAT8(2), PG_GETARG_FLOAT8(3) 3.1359 + ); 3.1360 + PG_RETURN_POINTER(box); 3.1361 +} 3.1362 + 3.1363 +/* create box ("ebox" in SQL) from two points ("epoint"s) */ 3.1364 +/* (can not be used to cover a longitude range of more than 120 degrees) */ 3.1365 +PG_FUNCTION_INFO_V1(pgl_create_ebox_from_epoints); 3.1366 +Datum pgl_create_ebox_from_epoints(PG_FUNCTION_ARGS) { 3.1367 + pgl_point *point1 = (pgl_point *)PG_GETARG_POINTER(0); 3.1368 + pgl_point *point2 = (pgl_point *)PG_GETARG_POINTER(1); 3.1369 + pgl_box *box = (pgl_box *)palloc(sizeof(pgl_box)); 3.1370 + double lat_min, lat_max, lon_min, lon_max; 3.1371 + double dlon; /* longitude range (delta longitude) */ 3.1372 + /* order latitude and longitude boundaries */ 3.1373 + if (point2->lat < point1->lat) { 3.1374 + lat_min = point2->lat; 3.1375 + lat_max = point1->lat; 3.1376 + } else { 3.1377 + lat_min = point1->lat; 3.1378 + lat_max = point2->lat; 3.1379 + } 3.1380 + if (point2->lon < point1->lon) { 3.1381 + lon_min = point2->lon; 3.1382 + lon_max = point1->lon; 3.1383 + } else { 3.1384 + lon_min = point1->lon; 3.1385 + lon_max = point2->lon; 3.1386 + } 3.1387 + /* calculate longitude range (round to avoid floating point errors) */ 3.1388 + dlon = pgl_round(lon_max - lon_min); 3.1389 + /* determine east-west direction */ 3.1390 + if (dlon >= 240) { 3.1391 + /* assume that 180th meridian is crossed and swap min/max longitude */ 3.1392 + double swap = lon_min; lon_min = lon_max; lon_max = swap; 3.1393 + } else if (dlon > 120) { 3.1394 + /* unclear orientation since delta longitude > 120 */ 3.1395 + ereport(ERROR, ( 3.1396 + errcode(ERRCODE_DATA_EXCEPTION), 3.1397 + errmsg("can not determine east/west orientation for ebox") 3.1398 + )); 3.1399 + } 3.1400 + /* use boundaries to setup box (and perform checks) */ 3.1401 + pgl_ebox_set_boundaries(box, lat_min, lat_max, lon_min, lon_max); 3.1402 + /* return result */ 3.1403 + PG_RETURN_POINTER(box); 3.1404 +} 3.1405 + 3.1406 +/* parse box ("ebox" in SQL) */ 3.1407 +/* format: '[NS]<float> [EW]<float> [NS]<float> [EW]<float>' 3.1408 + or: '[NS]<float> [NS]<float> [EW]<float> [EW]<float>' */ 3.1409 +PG_FUNCTION_INFO_V1(pgl_ebox_in); 3.1410 +Datum pgl_ebox_in(PG_FUNCTION_ARGS) { 3.1411 + char *str = PG_GETARG_CSTRING(0); /* input string */ 3.1412 + char *str_lower; /* lower case version of input string */ 3.1413 + char *strptr; /* current position within string */ 3.1414 + int valid; /* number of valid chars */ 3.1415 + int done; /* specifies if latitude or longitude was read */ 3.1416 + double val; /* temporary variable */ 3.1417 + int lat_count = 0; /* count of latitude values parsed */ 3.1418 + int lon_count = 0; /* count of longitufde values parsed */ 3.1419 + double lat_min, lat_max, lon_min, lon_max; /* see pgl_box struct */ 3.1420 + pgl_box *box; /* return value (to be palloc'ed) */ 3.1421 + /* lowercase input */ 3.1422 + str_lower = psprintf("%s", str); 3.1423 + for (strptr=str_lower; *strptr; strptr++) { 3.1424 + if (*strptr >= 'A' && *strptr <= 'Z') *strptr += 'a' - 'A'; 3.1425 + } 3.1426 + /* reset reading position to start of (lowercase) string */ 3.1427 + strptr = str_lower; 3.1428 + /* check if empty box */ 3.1429 + valid = 0; 3.1430 + sscanf(strptr, " empty %n", &valid); 3.1431 + if (valid && strptr[valid] == 0) { 3.1432 + /* allocate and return empty box */ 3.1433 + box = (pgl_box *)palloc(sizeof(pgl_box)); 3.1434 + pgl_box_set_empty(box); 3.1435 + PG_RETURN_POINTER(box); 3.1436 + } 3.1437 + /* demand four blocks separated by whitespace */ 3.1438 + valid = 0; 3.1439 + sscanf(strptr, " %*s %*s %*s %*s %n", &valid); 3.1440 + /* if four blocks separated by whitespace exist, parse those blocks */ 3.1441 + if (strptr[valid] == 0) while (strptr[0]) { 3.1442 + /* parse either latitude or longitude (whichever found in input string) */ 3.1443 + done = pgl_scan(&strptr, &val, &val); 3.1444 + /* store latitude or longitude in lat_min, lat_max, lon_min, or lon_max */ 3.1445 + if (done == PGL_SCAN_LAT) { 3.1446 + if (!lat_count) lat_min = val; else lat_max = val; 3.1447 + lat_count++; 3.1448 + } else if (done == PGL_SCAN_LON) { 3.1449 + if (!lon_count) lon_min = val; else lon_max = val; 3.1450 + lon_count++; 3.1451 + } else { 3.1452 + break; 3.1453 + } 3.1454 + } 3.1455 + /* require end of string, and two latitude and two longitude values */ 3.1456 + if (strptr[0] || lat_count != 2 || lon_count != 2) { 3.1457 + ereport(ERROR, ( 3.1458 + errcode(ERRCODE_INVALID_TEXT_REPRESENTATION), 3.1459 + errmsg("invalid input syntax for type ebox: \"%s\"", str) 3.1460 + )); 3.1461 + } 3.1462 + /* free lower case string */ 3.1463 + pfree(str_lower); 3.1464 + /* order boundaries (maximum greater than minimum) */ 3.1465 + if (lat_min > lat_max) { val = lat_min; lat_min = lat_max; lat_max = val; } 3.1466 + if (lon_min > lon_max) { val = lon_min; lon_min = lon_max; lon_max = val; } 3.1467 + /* allocate memory for result */ 3.1468 + box = (pgl_box *)palloc(sizeof(pgl_box)); 3.1469 + /* set boundaries (and perform checks) */ 3.1470 + pgl_ebox_set_boundaries(box, lat_min, lat_max, lon_min, lon_max); 3.1471 + /* return result */ 3.1472 + PG_RETURN_POINTER(box); 3.1473 +} 3.1474 + 3.1475 +/* set circle to given latitude, longitude, and radius (including checks) */ 3.1476 +static void pgl_ecircle_set_latlon_radius( 3.1477 + pgl_circle *circle, double lat, double lon, double radius 3.1478 +) { 3.1479 + /* set center point (including checks) */ 3.1480 + pgl_epoint_set_latlon(&(circle->center), lat, lon); 3.1481 + /* handle non-positive radius */ 3.1482 + if (isnan(radius)) { 3.1483 + ereport(ERROR, ( 3.1484 + errcode(ERRCODE_DATA_EXCEPTION), 3.1485 + errmsg("invalid radius for ecircle") 3.1486 + )); 3.1487 + } 3.1488 + if (radius == 0) radius = 0; /* avoids -0 */ 3.1489 + else if (radius < 0) { 3.1490 + if (isfinite(radius)) { 3.1491 + ereport(NOTICE, (errmsg("negative radius converted to minus infinity"))); 3.1492 + } 3.1493 + radius = -INFINITY; 3.1494 + } 3.1495 + /* store radius (round-trip safety is ensured by pgl_print_float) */ 3.1496 + circle->radius = radius; 3.1497 +} 3.1498 + 3.1499 +/* create circle ("ecircle" in SQL) from latitude, longitude, and radius */ 3.1500 +PG_FUNCTION_INFO_V1(pgl_create_ecircle); 3.1501 +Datum pgl_create_ecircle(PG_FUNCTION_ARGS) { 3.1502 + pgl_circle *circle = (pgl_circle *)palloc(sizeof(pgl_circle)); 3.1503 + pgl_ecircle_set_latlon_radius( 3.1504 + circle, PG_GETARG_FLOAT8(0), PG_GETARG_FLOAT8(1), PG_GETARG_FLOAT8(2) 3.1505 + ); 3.1506 + PG_RETURN_POINTER(circle); 3.1507 +} 3.1508 + 3.1509 +/* create circle ("ecircle" in SQL) from point ("epoint"), and radius */ 3.1510 +PG_FUNCTION_INFO_V1(pgl_create_ecircle_from_epoint); 3.1511 +Datum pgl_create_ecircle_from_epoint(PG_FUNCTION_ARGS) { 3.1512 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.1513 + double radius = PG_GETARG_FLOAT8(1); 3.1514 + pgl_circle *circle = (pgl_circle *)palloc(sizeof(pgl_circle)); 3.1515 + /* set latitude, longitude, radius (and perform checks) */ 3.1516 + pgl_ecircle_set_latlon_radius(circle, point->lat, point->lon, radius); 3.1517 + /* return result */ 3.1518 + PG_RETURN_POINTER(circle); 3.1519 +} 3.1520 + 3.1521 +/* parse circle ("ecircle" in SQL) */ 3.1522 +/* format: '[NS]<float> [EW]<float> <float>' */ 3.1523 +PG_FUNCTION_INFO_V1(pgl_ecircle_in); 3.1524 +Datum pgl_ecircle_in(PG_FUNCTION_ARGS) { 3.1525 + char *str = PG_GETARG_CSTRING(0); /* input string */ 3.1526 + char *strptr = str; /* current position within string */ 3.1527 + double lat, lon, radius; /* parsed values as double precision flaots */ 3.1528 + int valid = 0; /* number of valid chars */ 3.1529 + int done = 0; /* stores if latitude and/or longitude was read */ 3.1530 + pgl_circle *circle; /* return value (to be palloc'ed) */ 3.1531 + /* demand three blocks separated by whitespace */ 3.1532 + sscanf(strptr, " %*s %*s %*s %n", &valid); 3.1533 + /* if three blocks separated by whitespace exist, parse those blocks */ 3.1534 + if (strptr[valid] == 0) { 3.1535 + /* parse latitude and longitude */ 3.1536 + done |= pgl_scan(&strptr, &lat, &lon); 3.1537 + done |= pgl_scan(&strptr, &lat, &lon); 3.1538 + /* parse radius (while incrementing strptr by number of bytes parsed) */ 3.1539 + valid = 0; 3.1540 + if (sscanf(strptr, " %lf %n", &radius, &valid) == 1) strptr += valid; 3.1541 + } 3.1542 + /* require end of string and both latitude and longitude being parsed */ 3.1543 + if (strptr[0] || done != PGL_SCAN_LATLON) { 3.1544 + ereport(ERROR, ( 3.1545 + errcode(ERRCODE_INVALID_TEXT_REPRESENTATION), 3.1546 + errmsg("invalid input syntax for type ecircle: \"%s\"", str) 3.1547 + )); 3.1548 + } 3.1549 + /* allocate memory for result */ 3.1550 + circle = (pgl_circle *)palloc(sizeof(pgl_circle)); 3.1551 + /* set latitude, longitude, radius (and perform checks) */ 3.1552 + pgl_ecircle_set_latlon_radius(circle, lat, lon, radius); 3.1553 + /* return result */ 3.1554 + PG_RETURN_POINTER(circle); 3.1555 +} 3.1556 + 3.1557 +/* parse cluster ("ecluster" in SQL) */ 3.1558 +PG_FUNCTION_INFO_V1(pgl_ecluster_in); 3.1559 +Datum pgl_ecluster_in(PG_FUNCTION_ARGS) { 3.1560 + int i; 3.1561 + char *str = PG_GETARG_CSTRING(0); /* input string */ 3.1562 + char *str_lower; /* lower case version of input string */ 3.1563 + char *strptr; /* pointer to current reading position of input */ 3.1564 + int npoints_total = 0; /* total number of points in cluster */ 3.1565 + int nentries = 0; /* total number of entries */ 3.1566 + pgl_newentry *entries; /* array of pgl_newentry to create pgl_cluster */ 3.1567 + int entries_buflen = 4; /* maximum number of elements in entries array */ 3.1568 + int valid; /* number of valid chars processed */ 3.1569 + double lat, lon; /* latitude and longitude of parsed point */ 3.1570 + int entrytype; /* current entry type */ 3.1571 + int npoints; /* number of points in current entry */ 3.1572 + pgl_point *points; /* array of pgl_point for pgl_newentry */ 3.1573 + int points_buflen; /* maximum number of elements in points array */ 3.1574 + int done; /* return value of pgl_scan function */ 3.1575 + pgl_cluster *cluster; /* created cluster */ 3.1576 + /* lowercase input */ 3.1577 + str_lower = psprintf("%s", str); 3.1578 + for (strptr=str_lower; *strptr; strptr++) { 3.1579 + if (*strptr >= 'A' && *strptr <= 'Z') *strptr += 'a' - 'A'; 3.1580 + } 3.1581 + /* reset reading position to start of (lowercase) string */ 3.1582 + strptr = str_lower; 3.1583 + /* allocate initial buffer for entries */ 3.1584 + entries = palloc(entries_buflen * sizeof(pgl_newentry)); 3.1585 + /* parse until end of string */ 3.1586 + while (strptr[0]) { 3.1587 + /* require previous white-space or closing parenthesis before next token */ 3.1588 + if (strptr != str_lower && !isspace(strptr[-1]) && strptr[-1] != ')') { 3.1589 + goto pgl_ecluster_in_error; 3.1590 + } 3.1591 + /* ignore token "empty" */ 3.1592 + valid = 0; sscanf(strptr, " empty %n", &valid); 3.1593 + if (valid) { strptr += valid; continue; } 3.1594 + /* test for "point" token */ 3.1595 + valid = 0; sscanf(strptr, " point ( %n", &valid); 3.1596 + if (valid) { 3.1597 + strptr += valid; 3.1598 + entrytype = PGL_ENTRY_POINT; 3.1599 + goto pgl_ecluster_in_type_ok; 3.1600 + } 3.1601 + /* test for "path" token */ 3.1602 + valid = 0; sscanf(strptr, " path ( %n", &valid); 3.1603 + if (valid) { 3.1604 + strptr += valid; 3.1605 + entrytype = PGL_ENTRY_PATH; 3.1606 + goto pgl_ecluster_in_type_ok; 3.1607 + } 3.1608 + /* test for "outline" token */ 3.1609 + valid = 0; sscanf(strptr, " outline ( %n", &valid); 3.1610 + if (valid) { 3.1611 + strptr += valid; 3.1612 + entrytype = PGL_ENTRY_OUTLINE; 3.1613 + goto pgl_ecluster_in_type_ok; 3.1614 + } 3.1615 + /* test for "polygon" token */ 3.1616 + valid = 0; sscanf(strptr, " polygon ( %n", &valid); 3.1617 + if (valid) { 3.1618 + strptr += valid; 3.1619 + entrytype = PGL_ENTRY_POLYGON; 3.1620 + goto pgl_ecluster_in_type_ok; 3.1621 + } 3.1622 + /* error if no valid token found */ 3.1623 + goto pgl_ecluster_in_error; 3.1624 + pgl_ecluster_in_type_ok: 3.1625 + /* check if pgl_newentry array needs to grow */ 3.1626 + if (nentries == entries_buflen) { 3.1627 + pgl_newentry *newbuf; 3.1628 + entries_buflen *= 2; 3.1629 + newbuf = palloc(entries_buflen * sizeof(pgl_newentry)); 3.1630 + memcpy(newbuf, entries, nentries * sizeof(pgl_newentry)); 3.1631 + pfree(entries); 3.1632 + entries = newbuf; 3.1633 + } 3.1634 + /* reset number of points for current entry */ 3.1635 + npoints = 0; 3.1636 + /* allocate array for points */ 3.1637 + points_buflen = 4; 3.1638 + points = palloc(points_buflen * sizeof(pgl_point)); 3.1639 + /* parse until closing parenthesis */ 3.1640 + while (strptr[0] != ')') { 3.1641 + /* error on unexpected end of string */ 3.1642 + if (strptr[0] == 0) goto pgl_ecluster_in_error; 3.1643 + /* mark neither latitude nor longitude as read */ 3.1644 + done = PGL_SCAN_NONE; 3.1645 + /* require white-space before second, third, etc. point */ 3.1646 + if (npoints != 0 && !isspace(strptr[-1])) goto pgl_ecluster_in_error; 3.1647 + /* scan latitude (or longitude) */ 3.1648 + done |= pgl_scan(&strptr, &lat, &lon); 3.1649 + /* require white-space before second coordinate */ 3.1650 + if (strptr != str && !isspace(strptr[-1])) goto pgl_ecluster_in_error; 3.1651 + /* scan longitude (or latitude) */ 3.1652 + done |= pgl_scan(&strptr, &lat, &lon); 3.1653 + /* error unless both latitude and longitude were parsed */ 3.1654 + if (done != PGL_SCAN_LATLON) goto pgl_ecluster_in_error; 3.1655 + /* throw error if number of points is too high */ 3.1656 + if (npoints_total == PGL_CLUSTER_MAXPOINTS) { 3.1657 + ereport(ERROR, ( 3.1658 + errcode(ERRCODE_INVALID_TEXT_REPRESENTATION), 3.1659 + errmsg( 3.1660 + "too many points for ecluster entry (maximum %i)", 3.1661 + PGL_CLUSTER_MAXPOINTS 3.1662 + ) 3.1663 + )); 3.1664 + } 3.1665 + /* check if pgl_point array needs to grow */ 3.1666 + if (npoints == points_buflen) { 3.1667 + pgl_point *newbuf; 3.1668 + points_buflen *= 2; 3.1669 + newbuf = palloc(points_buflen * sizeof(pgl_point)); 3.1670 + memcpy(newbuf, points, npoints * sizeof(pgl_point)); 3.1671 + pfree(points); 3.1672 + points = newbuf; 3.1673 + } 3.1674 + /* append point to pgl_point array (includes checks) */ 3.1675 + pgl_epoint_set_latlon(&(points[npoints++]), lat, lon); 3.1676 + /* increase total number of points */ 3.1677 + npoints_total++; 3.1678 + } 3.1679 + /* error if entry has no points */ 3.1680 + if (!npoints) goto pgl_ecluster_in_error; 3.1681 + /* entries with one point are automatically of type "point" */ 3.1682 + if (npoints == 1) entrytype = PGL_ENTRY_POINT; 3.1683 + /* if entries have more than one point */ 3.1684 + else { 3.1685 + /* throw error if entry type is "point" */ 3.1686 + if (entrytype == PGL_ENTRY_POINT) { 3.1687 + ereport(ERROR, ( 3.1688 + errcode(ERRCODE_INVALID_TEXT_REPRESENTATION), 3.1689 + errmsg("invalid input syntax for type ecluster (point entry with more than one point)") 3.1690 + )); 3.1691 + } 3.1692 + /* coerce outlines and polygons with more than 2 points to be a path */ 3.1693 + if (npoints == 2) entrytype = PGL_ENTRY_PATH; 3.1694 + } 3.1695 + /* append entry to pgl_newentry array */ 3.1696 + entries[nentries].entrytype = entrytype; 3.1697 + entries[nentries].npoints = npoints; 3.1698 + entries[nentries].points = points; 3.1699 + nentries++; 3.1700 + /* consume closing parenthesis */ 3.1701 + strptr++; 3.1702 + /* consume white-space */ 3.1703 + while (isspace(strptr[0])) strptr++; 3.1704 + } 3.1705 + /* free lower case string */ 3.1706 + pfree(str_lower); 3.1707 + /* create cluster from pgl_newentry array */ 3.1708 + cluster = pgl_new_cluster(nentries, entries); 3.1709 + /* free pgl_newentry array */ 3.1710 + for (i=0; i<nentries; i++) pfree(entries[i].points); 3.1711 + pfree(entries); 3.1712 + /* set bounding circle of cluster and check east/west orientation */ 3.1713 + if (!pgl_finalize_cluster(cluster)) { 3.1714 + ereport(ERROR, ( 3.1715 + errcode(ERRCODE_DATA_EXCEPTION), 3.1716 + errmsg("can not determine east/west orientation for ecluster"), 3.1717 + errhint("Ensure that each entry has a longitude span of less than 180 degrees.") 3.1718 + )); 3.1719 + } 3.1720 + /* return cluster */ 3.1721 + PG_RETURN_POINTER(cluster); 3.1722 + /* code to throw error */ 3.1723 + pgl_ecluster_in_error: 3.1724 + ereport(ERROR, ( 3.1725 + errcode(ERRCODE_INVALID_TEXT_REPRESENTATION), 3.1726 + errmsg("invalid input syntax for type ecluster: \"%s\"", str) 3.1727 + )); 3.1728 +} 3.1729 + 3.1730 +/* convert point ("epoint") to string representation */ 3.1731 +PG_FUNCTION_INFO_V1(pgl_epoint_out); 3.1732 +Datum pgl_epoint_out(PG_FUNCTION_ARGS) { 3.1733 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.1734 + char latstr[PGL_NUMBUFLEN]; 3.1735 + char lonstr[PGL_NUMBUFLEN]; 3.1736 + pgl_print_lat(latstr, point->lat); 3.1737 + pgl_print_lon(lonstr, point->lon); 3.1738 + PG_RETURN_CSTRING(psprintf("%s %s", latstr, lonstr)); 3.1739 +} 3.1740 + 3.1741 +/* convert box ("ebox") to string representation */ 3.1742 +PG_FUNCTION_INFO_V1(pgl_ebox_out); 3.1743 +Datum pgl_ebox_out(PG_FUNCTION_ARGS) { 3.1744 + pgl_box *box = (pgl_box *)PG_GETARG_POINTER(0); 3.1745 + double lon_min = box->lon_min; 3.1746 + double lon_max = box->lon_max; 3.1747 + char lat_min_str[PGL_NUMBUFLEN]; 3.1748 + char lat_max_str[PGL_NUMBUFLEN]; 3.1749 + char lon_min_str[PGL_NUMBUFLEN]; 3.1750 + char lon_max_str[PGL_NUMBUFLEN]; 3.1751 + /* return string "empty" if box is set to be empty */ 3.1752 + if (box->lat_min > box->lat_max) PG_RETURN_CSTRING("empty"); 3.1753 + /* use boundaries exceeding W180 or E180 if 180th meridian is enclosed */ 3.1754 + /* (required since pgl_box_in orders the longitude boundaries) */ 3.1755 + if (lon_min > lon_max) { 3.1756 + if (lon_min + lon_max >= 0) lon_min -= 360; 3.1757 + else lon_max += 360; 3.1758 + } 3.1759 + /* format and return result */ 3.1760 + pgl_print_lat(lat_min_str, box->lat_min); 3.1761 + pgl_print_lat(lat_max_str, box->lat_max); 3.1762 + pgl_print_lon(lon_min_str, lon_min); 3.1763 + pgl_print_lon(lon_max_str, lon_max); 3.1764 + PG_RETURN_CSTRING(psprintf( 3.1765 + "%s %s %s %s", 3.1766 + lat_min_str, lon_min_str, lat_max_str, lon_max_str 3.1767 + )); 3.1768 +} 3.1769 + 3.1770 +/* convert circle ("ecircle") to string representation */ 3.1771 +PG_FUNCTION_INFO_V1(pgl_ecircle_out); 3.1772 +Datum pgl_ecircle_out(PG_FUNCTION_ARGS) { 3.1773 + pgl_circle *circle = (pgl_circle *)PG_GETARG_POINTER(0); 3.1774 + char latstr[PGL_NUMBUFLEN]; 3.1775 + char lonstr[PGL_NUMBUFLEN]; 3.1776 + char radstr[PGL_NUMBUFLEN]; 3.1777 + pgl_print_lat(latstr, circle->center.lat); 3.1778 + pgl_print_lon(lonstr, circle->center.lon); 3.1779 + pgl_print_float(radstr, circle->radius); 3.1780 + PG_RETURN_CSTRING(psprintf("%s %s %s", latstr, lonstr, radstr)); 3.1781 +} 3.1782 + 3.1783 +/* convert cluster ("ecluster") to string representation */ 3.1784 +PG_FUNCTION_INFO_V1(pgl_ecluster_out); 3.1785 +Datum pgl_ecluster_out(PG_FUNCTION_ARGS) { 3.1786 + pgl_cluster *cluster = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); 3.1787 + char latstr[PGL_NUMBUFLEN]; /* string buffer for latitude */ 3.1788 + char lonstr[PGL_NUMBUFLEN]; /* string buffer for longitude */ 3.1789 + char ***strings; /* array of array of strings */ 3.1790 + char *string; /* string of current token */ 3.1791 + char *res, *resptr; /* result and pointer to current write position */ 3.1792 + size_t reslen = 1; /* length of result (init with 1 for terminator) */ 3.1793 + int npoints; /* number of points of current entry */ 3.1794 + int i, j; /* i: entry, j: point in entry */ 3.1795 + /* handle empty clusters */ 3.1796 + if (cluster->nentries == 0) { 3.1797 + /* free detoasted cluster (if copy) */ 3.1798 + PG_FREE_IF_COPY(cluster, 0); 3.1799 + /* return static result */ 3.1800 + PG_RETURN_CSTRING("empty"); 3.1801 + } 3.1802 + /* allocate array of array of strings */ 3.1803 + strings = palloc(cluster->nentries * sizeof(char **)); 3.1804 + /* iterate over all entries in cluster */ 3.1805 + for (i=0; i<cluster->nentries; i++) { 3.1806 + /* get number of points in entry */ 3.1807 + npoints = cluster->entries[i].npoints; 3.1808 + /* allocate array of strings (one string for each point plus two extra) */ 3.1809 + strings[i] = palloc((2 + npoints) * sizeof(char *)); 3.1810 + /* determine opening string */ 3.1811 + switch (cluster->entries[i].entrytype) { 3.1812 + case PGL_ENTRY_POINT: string = (i==0)?"point (" :" point ("; break; 3.1813 + case PGL_ENTRY_PATH: string = (i==0)?"path (" :" path ("; break; 3.1814 + case PGL_ENTRY_OUTLINE: string = (i==0)?"outline (":" outline ("; break; 3.1815 + case PGL_ENTRY_POLYGON: string = (i==0)?"polygon (":" polygon ("; break; 3.1816 + default: string = (i==0)?"unknown" :" unknown"; 3.1817 + } 3.1818 + /* use opening string as first string in array */ 3.1819 + strings[i][0] = string; 3.1820 + /* update result length (for allocating result string later) */ 3.1821 + reslen += strlen(string); 3.1822 + /* iterate over all points */ 3.1823 + for (j=0; j<npoints; j++) { 3.1824 + /* create string representation of point */ 3.1825 + pgl_print_lat(latstr, PGL_ENTRY_POINTS(cluster, i)[j].lat); 3.1826 + pgl_print_lon(lonstr, PGL_ENTRY_POINTS(cluster, i)[j].lon); 3.1827 + string = psprintf((j == 0) ? "%s %s" : " %s %s", latstr, lonstr); 3.1828 + /* copy string pointer to string array */ 3.1829 + strings[i][j+1] = string; 3.1830 + /* update result length (for allocating result string later) */ 3.1831 + reslen += strlen(string); 3.1832 + } 3.1833 + /* use closing parenthesis as last string in array */ 3.1834 + strings[i][npoints+1] = ")"; 3.1835 + /* update result length (for allocating result string later) */ 3.1836 + reslen++; 3.1837 + } 3.1838 + /* allocate result string */ 3.1839 + res = palloc(reslen); 3.1840 + /* set write pointer to begin of result string */ 3.1841 + resptr = res; 3.1842 + /* copy strings into result string */ 3.1843 + for (i=0; i<cluster->nentries; i++) { 3.1844 + npoints = cluster->entries[i].npoints; 3.1845 + for (j=0; j<npoints+2; j++) { 3.1846 + string = strings[i][j]; 3.1847 + strcpy(resptr, string); 3.1848 + resptr += strlen(string); 3.1849 + /* free strings allocated by psprintf */ 3.1850 + if (j != 0 && j != npoints+1) pfree(string); 3.1851 + } 3.1852 + /* free array of strings */ 3.1853 + pfree(strings[i]); 3.1854 + } 3.1855 + /* free array of array of strings */ 3.1856 + pfree(strings); 3.1857 + /* free detoasted cluster (if copy) */ 3.1858 + PG_FREE_IF_COPY(cluster, 0); 3.1859 + /* return result */ 3.1860 + PG_RETURN_CSTRING(res); 3.1861 +} 3.1862 + 3.1863 +/* binary input function for point ("epoint") */ 3.1864 +PG_FUNCTION_INFO_V1(pgl_epoint_recv); 3.1865 +Datum pgl_epoint_recv(PG_FUNCTION_ARGS) { 3.1866 + StringInfo buf = (StringInfo)PG_GETARG_POINTER(0); 3.1867 + pgl_point *point = (pgl_point *)palloc(sizeof(pgl_point)); 3.1868 + point->lat = pq_getmsgfloat8(buf); 3.1869 + point->lon = pq_getmsgfloat8(buf); 3.1870 + PG_RETURN_POINTER(point); 3.1871 +} 3.1872 + 3.1873 +/* binary input function for box ("ebox") */ 3.1874 +PG_FUNCTION_INFO_V1(pgl_ebox_recv); 3.1875 +Datum pgl_ebox_recv(PG_FUNCTION_ARGS) { 3.1876 + StringInfo buf = (StringInfo)PG_GETARG_POINTER(0); 3.1877 + pgl_box *box = (pgl_box *)palloc(sizeof(pgl_box)); 3.1878 + box->lat_min = pq_getmsgfloat8(buf); 3.1879 + box->lat_max = pq_getmsgfloat8(buf); 3.1880 + box->lon_min = pq_getmsgfloat8(buf); 3.1881 + box->lon_max = pq_getmsgfloat8(buf); 3.1882 + PG_RETURN_POINTER(box); 3.1883 +} 3.1884 + 3.1885 +/* binary input function for circle ("ecircle") */ 3.1886 +PG_FUNCTION_INFO_V1(pgl_ecircle_recv); 3.1887 +Datum pgl_ecircle_recv(PG_FUNCTION_ARGS) { 3.1888 + StringInfo buf = (StringInfo)PG_GETARG_POINTER(0); 3.1889 + pgl_circle *circle = (pgl_circle *)palloc(sizeof(pgl_circle)); 3.1890 + circle->center.lat = pq_getmsgfloat8(buf); 3.1891 + circle->center.lon = pq_getmsgfloat8(buf); 3.1892 + circle->radius = pq_getmsgfloat8(buf); 3.1893 + PG_RETURN_POINTER(circle); 3.1894 +} 3.1895 + 3.1896 +/* TODO: binary receive function for cluster */ 3.1897 + 3.1898 +/* binary output function for point ("epoint") */ 3.1899 +PG_FUNCTION_INFO_V1(pgl_epoint_send); 3.1900 +Datum pgl_epoint_send(PG_FUNCTION_ARGS) { 3.1901 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.1902 + StringInfoData buf; 3.1903 + pq_begintypsend(&buf); 3.1904 + pq_sendfloat8(&buf, point->lat); 3.1905 + pq_sendfloat8(&buf, point->lon); 3.1906 + PG_RETURN_BYTEA_P(pq_endtypsend(&buf)); 3.1907 +} 3.1908 + 3.1909 +/* binary output function for box ("ebox") */ 3.1910 +PG_FUNCTION_INFO_V1(pgl_ebox_send); 3.1911 +Datum pgl_ebox_send(PG_FUNCTION_ARGS) { 3.1912 + pgl_box *box = (pgl_box *)PG_GETARG_POINTER(0); 3.1913 + StringInfoData buf; 3.1914 + pq_begintypsend(&buf); 3.1915 + pq_sendfloat8(&buf, box->lat_min); 3.1916 + pq_sendfloat8(&buf, box->lat_max); 3.1917 + pq_sendfloat8(&buf, box->lon_min); 3.1918 + pq_sendfloat8(&buf, box->lon_max); 3.1919 + PG_RETURN_BYTEA_P(pq_endtypsend(&buf)); 3.1920 +} 3.1921 + 3.1922 +/* binary output function for circle ("ecircle") */ 3.1923 +PG_FUNCTION_INFO_V1(pgl_ecircle_send); 3.1924 +Datum pgl_ecircle_send(PG_FUNCTION_ARGS) { 3.1925 + pgl_circle *circle = (pgl_circle *)PG_GETARG_POINTER(0); 3.1926 + StringInfoData buf; 3.1927 + pq_begintypsend(&buf); 3.1928 + pq_sendfloat8(&buf, circle->center.lat); 3.1929 + pq_sendfloat8(&buf, circle->center.lon); 3.1930 + pq_sendfloat8(&buf, circle->radius); 3.1931 + PG_RETURN_BYTEA_P(pq_endtypsend(&buf)); 3.1932 +} 3.1933 + 3.1934 +/* TODO: binary send functions for cluster */ 3.1935 + 3.1936 +/* cast point ("epoint") to box ("ebox") */ 3.1937 +PG_FUNCTION_INFO_V1(pgl_epoint_to_ebox); 3.1938 +Datum pgl_epoint_to_ebox(PG_FUNCTION_ARGS) { 3.1939 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.1940 + pgl_box *box = palloc(sizeof(pgl_box)); 3.1941 + box->lat_min = point->lat; 3.1942 + box->lat_max = point->lat; 3.1943 + box->lon_min = point->lon; 3.1944 + box->lon_max = point->lon; 3.1945 + PG_RETURN_POINTER(box); 3.1946 +} 3.1947 + 3.1948 +/* cast point ("epoint") to circle ("ecircle") */ 3.1949 +PG_FUNCTION_INFO_V1(pgl_epoint_to_ecircle); 3.1950 +Datum pgl_epoint_to_ecircle(PG_FUNCTION_ARGS) { 3.1951 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.1952 + pgl_circle *circle = palloc(sizeof(pgl_box)); 3.1953 + circle->center = *point; 3.1954 + circle->radius = 0; 3.1955 + PG_RETURN_POINTER(circle); 3.1956 +} 3.1957 + 3.1958 +/* cast point ("epoint") to cluster ("ecluster") */ 3.1959 +PG_FUNCTION_INFO_V1(pgl_epoint_to_ecluster); 3.1960 +Datum pgl_epoint_to_ecluster(PG_FUNCTION_ARGS) { 3.1961 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.1962 + pgl_newentry entry; 3.1963 + entry.entrytype = PGL_ENTRY_POINT; 3.1964 + entry.npoints = 1; 3.1965 + entry.points = point; 3.1966 + PG_RETURN_POINTER(pgl_new_cluster(1, &entry)); 3.1967 +} 3.1968 + 3.1969 +/* cast box ("ebox") to cluster ("ecluster") */ 3.1970 +#define pgl_ebox_to_ecluster_macro(i, a, b) \ 3.1971 + entries[i].entrytype = PGL_ENTRY_POLYGON; \ 3.1972 + entries[i].npoints = 4; \ 3.1973 + entries[i].points = points[i]; \ 3.1974 + points[i][0].lat = box->lat_min; \ 3.1975 + points[i][0].lon = (a); \ 3.1976 + points[i][1].lat = box->lat_min; \ 3.1977 + points[i][1].lon = (b); \ 3.1978 + points[i][2].lat = box->lat_max; \ 3.1979 + points[i][2].lon = (b); \ 3.1980 + points[i][3].lat = box->lat_max; \ 3.1981 + points[i][3].lon = (a); 3.1982 +PG_FUNCTION_INFO_V1(pgl_ebox_to_ecluster); 3.1983 +Datum pgl_ebox_to_ecluster(PG_FUNCTION_ARGS) { 3.1984 + pgl_box *box = (pgl_box *)PG_GETARG_POINTER(0); 3.1985 + double lon, dlon; 3.1986 + int nentries; 3.1987 + pgl_newentry entries[3]; 3.1988 + pgl_point points[3][4]; 3.1989 + if (box->lat_min > box->lat_max) { 3.1990 + nentries = 0; 3.1991 + } else if (box->lon_min > box->lon_max) { 3.1992 + if (box->lon_min < 0) { 3.1993 + lon = pgl_round((box->lon_min + 180) / 2.0); 3.1994 + nentries = 3; 3.1995 + pgl_ebox_to_ecluster_macro(0, box->lon_min, lon); 3.1996 + pgl_ebox_to_ecluster_macro(1, lon, 180); 3.1997 + pgl_ebox_to_ecluster_macro(2, -180, box->lon_max); 3.1998 + } else if (box->lon_max > 0) { 3.1999 + lon = pgl_round((box->lon_max - 180) / 2.0); 3.2000 + nentries = 3; 3.2001 + pgl_ebox_to_ecluster_macro(0, box->lon_min, 180); 3.2002 + pgl_ebox_to_ecluster_macro(1, -180, lon); 3.2003 + pgl_ebox_to_ecluster_macro(2, lon, box->lon_max); 3.2004 + } else { 3.2005 + nentries = 2; 3.2006 + pgl_ebox_to_ecluster_macro(0, box->lon_min, 180); 3.2007 + pgl_ebox_to_ecluster_macro(1, -180, box->lon_max); 3.2008 + } 3.2009 + } else { 3.2010 + dlon = pgl_round(box->lon_max - box->lon_min); 3.2011 + if (dlon < 180) { 3.2012 + nentries = 1; 3.2013 + pgl_ebox_to_ecluster_macro(0, box->lon_min, box->lon_max); 3.2014 + } else { 3.2015 + lon = pgl_round((box->lon_min + box->lon_max) / 2.0); 3.2016 + if ( 3.2017 + pgl_round(lon - box->lon_min) < 180 && 3.2018 + pgl_round(box->lon_max - lon) < 180 3.2019 + ) { 3.2020 + nentries = 2; 3.2021 + pgl_ebox_to_ecluster_macro(0, box->lon_min, lon); 3.2022 + pgl_ebox_to_ecluster_macro(1, lon, box->lon_max); 3.2023 + } else { 3.2024 + nentries = 3; 3.2025 + pgl_ebox_to_ecluster_macro(0, box->lon_min, -60); 3.2026 + pgl_ebox_to_ecluster_macro(1, -60, 60); 3.2027 + pgl_ebox_to_ecluster_macro(2, 60, box->lon_max); 3.2028 + } 3.2029 + } 3.2030 + } 3.2031 + PG_RETURN_POINTER(pgl_new_cluster(nentries, entries)); 3.2032 +} 3.2033 + 3.2034 +/* extract latitude from point ("epoint") */ 3.2035 +PG_FUNCTION_INFO_V1(pgl_epoint_lat); 3.2036 +Datum pgl_epoint_lat(PG_FUNCTION_ARGS) { 3.2037 + PG_RETURN_FLOAT8(((pgl_point *)PG_GETARG_POINTER(0))->lat); 3.2038 +} 3.2039 + 3.2040 +/* extract longitude from point ("epoint") */ 3.2041 +PG_FUNCTION_INFO_V1(pgl_epoint_lon); 3.2042 +Datum pgl_epoint_lon(PG_FUNCTION_ARGS) { 3.2043 + PG_RETURN_FLOAT8(((pgl_point *)PG_GETARG_POINTER(0))->lon); 3.2044 +} 3.2045 + 3.2046 +/* extract minimum latitude from box ("ebox") */ 3.2047 +PG_FUNCTION_INFO_V1(pgl_ebox_lat_min); 3.2048 +Datum pgl_ebox_lat_min(PG_FUNCTION_ARGS) { 3.2049 + PG_RETURN_FLOAT8(((pgl_box *)PG_GETARG_POINTER(0))->lat_min); 3.2050 +} 3.2051 + 3.2052 +/* extract maximum latitude from box ("ebox") */ 3.2053 +PG_FUNCTION_INFO_V1(pgl_ebox_lat_max); 3.2054 +Datum pgl_ebox_lat_max(PG_FUNCTION_ARGS) { 3.2055 + PG_RETURN_FLOAT8(((pgl_box *)PG_GETARG_POINTER(0))->lat_max); 3.2056 +} 3.2057 + 3.2058 +/* extract minimum longitude from box ("ebox") */ 3.2059 +PG_FUNCTION_INFO_V1(pgl_ebox_lon_min); 3.2060 +Datum pgl_ebox_lon_min(PG_FUNCTION_ARGS) { 3.2061 + PG_RETURN_FLOAT8(((pgl_box *)PG_GETARG_POINTER(0))->lon_min); 3.2062 +} 3.2063 + 3.2064 +/* extract maximum longitude from box ("ebox") */ 3.2065 +PG_FUNCTION_INFO_V1(pgl_ebox_lon_max); 3.2066 +Datum pgl_ebox_lon_max(PG_FUNCTION_ARGS) { 3.2067 + PG_RETURN_FLOAT8(((pgl_box *)PG_GETARG_POINTER(0))->lon_max); 3.2068 +} 3.2069 + 3.2070 +/* extract center point from circle ("ecircle") */ 3.2071 +PG_FUNCTION_INFO_V1(pgl_ecircle_center); 3.2072 +Datum pgl_ecircle_center(PG_FUNCTION_ARGS) { 3.2073 + PG_RETURN_POINTER(&(((pgl_circle *)PG_GETARG_POINTER(0))->center)); 3.2074 +} 3.2075 + 3.2076 +/* extract radius from circle ("ecircle") */ 3.2077 +PG_FUNCTION_INFO_V1(pgl_ecircle_radius); 3.2078 +Datum pgl_ecircle_radius(PG_FUNCTION_ARGS) { 3.2079 + PG_RETURN_FLOAT8(((pgl_circle *)PG_GETARG_POINTER(0))->radius); 3.2080 +} 3.2081 + 3.2082 +/* check if point is inside box (overlap operator "&&") in SQL */ 3.2083 +PG_FUNCTION_INFO_V1(pgl_epoint_ebox_overlap); 3.2084 +Datum pgl_epoint_ebox_overlap(PG_FUNCTION_ARGS) { 3.2085 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.2086 + pgl_box *box = (pgl_box *)PG_GETARG_POINTER(1); 3.2087 + PG_RETURN_BOOL(pgl_point_in_box(point, box)); 3.2088 +} 3.2089 + 3.2090 +/* check if point is inside circle (overlap operator "&&") in SQL */ 3.2091 +PG_FUNCTION_INFO_V1(pgl_epoint_ecircle_overlap); 3.2092 +Datum pgl_epoint_ecircle_overlap(PG_FUNCTION_ARGS) { 3.2093 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.2094 + pgl_circle *circle = (pgl_circle *)PG_GETARG_POINTER(1); 3.2095 + PG_RETURN_BOOL( 3.2096 + pgl_distance( 3.2097 + point->lat, point->lon, 3.2098 + circle->center.lat, circle->center.lon 3.2099 + ) <= circle->radius 3.2100 + ); 3.2101 +} 3.2102 + 3.2103 +/* check if point is inside cluster (overlap operator "&&") in SQL */ 3.2104 +PG_FUNCTION_INFO_V1(pgl_epoint_ecluster_overlap); 3.2105 +Datum pgl_epoint_ecluster_overlap(PG_FUNCTION_ARGS) { 3.2106 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.2107 + pgl_cluster *cluster = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); 3.2108 + bool retval = pgl_point_in_cluster(point, cluster); 3.2109 + PG_FREE_IF_COPY(cluster, 1); 3.2110 + PG_RETURN_BOOL(retval); 3.2111 +} 3.2112 + 3.2113 +/* check if point may be inside cluster (lossy overl. operator "&&+") in SQL */ 3.2114 +PG_FUNCTION_INFO_V1(pgl_epoint_ecluster_may_overlap); 3.2115 +Datum pgl_epoint_ecluster_may_overlap(PG_FUNCTION_ARGS) { 3.2116 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.2117 + pgl_cluster *cluster = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); 3.2118 + bool retval = pgl_distance( 3.2119 + point->lat, point->lon, 3.2120 + cluster->bounding.center.lat, cluster->bounding.center.lon 3.2121 + ) <= cluster->bounding.radius; 3.2122 + PG_FREE_IF_COPY(cluster, 1); 3.2123 + PG_RETURN_BOOL(retval); 3.2124 +} 3.2125 + 3.2126 +/* check if two boxes overlap (overlap operator "&&") in SQL */ 3.2127 +PG_FUNCTION_INFO_V1(pgl_ebox_overlap); 3.2128 +Datum pgl_ebox_overlap(PG_FUNCTION_ARGS) { 3.2129 + pgl_box *box1 = (pgl_box *)PG_GETARG_POINTER(0); 3.2130 + pgl_box *box2 = (pgl_box *)PG_GETARG_POINTER(1); 3.2131 + PG_RETURN_BOOL(pgl_boxes_overlap(box1, box2)); 3.2132 +} 3.2133 + 3.2134 +/* check if box and circle may overlap (lossy overl. operator "&&+") in SQL */ 3.2135 +PG_FUNCTION_INFO_V1(pgl_ebox_ecircle_may_overlap); 3.2136 +Datum pgl_ebox_ecircle_may_overlap(PG_FUNCTION_ARGS) { 3.2137 + pgl_box *box = (pgl_box *)PG_GETARG_POINTER(0); 3.2138 + pgl_circle *circle = (pgl_circle *)PG_GETARG_POINTER(1); 3.2139 + PG_RETURN_BOOL( 3.2140 + pgl_estimate_point_box_distance(&circle->center, box) <= circle->radius 3.2141 + ); 3.2142 +} 3.2143 + 3.2144 +/* check if box and cluster may overlap (lossy overl. operator "&&+") in SQL */ 3.2145 +PG_FUNCTION_INFO_V1(pgl_ebox_ecluster_may_overlap); 3.2146 +Datum pgl_ebox_ecluster_may_overlap(PG_FUNCTION_ARGS) { 3.2147 + pgl_box *box = (pgl_box *)PG_GETARG_POINTER(0); 3.2148 + pgl_cluster *cluster = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); 3.2149 + bool retval = pgl_estimate_point_box_distance( 3.2150 + &cluster->bounding.center, 3.2151 + box 3.2152 + ) <= cluster->bounding.radius; 3.2153 + PG_FREE_IF_COPY(cluster, 1); 3.2154 + PG_RETURN_BOOL(retval); 3.2155 +} 3.2156 + 3.2157 +/* check if two circles overlap (overlap operator "&&") in SQL */ 3.2158 +PG_FUNCTION_INFO_V1(pgl_ecircle_overlap); 3.2159 +Datum pgl_ecircle_overlap(PG_FUNCTION_ARGS) { 3.2160 + pgl_circle *circle1 = (pgl_circle *)PG_GETARG_POINTER(0); 3.2161 + pgl_circle *circle2 = (pgl_circle *)PG_GETARG_POINTER(1); 3.2162 + PG_RETURN_BOOL( 3.2163 + pgl_distance( 3.2164 + circle1->center.lat, circle1->center.lon, 3.2165 + circle2->center.lat, circle2->center.lon 3.2166 + ) <= circle1->radius + circle2->radius 3.2167 + ); 3.2168 +} 3.2169 + 3.2170 +/* check if circle and cluster overlap (overlap operator "&&") in SQL */ 3.2171 +PG_FUNCTION_INFO_V1(pgl_ecircle_ecluster_overlap); 3.2172 +Datum pgl_ecircle_ecluster_overlap(PG_FUNCTION_ARGS) { 3.2173 + pgl_circle *circle = (pgl_circle *)PG_GETARG_POINTER(0); 3.2174 + pgl_cluster *cluster = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); 3.2175 + bool retval = ( 3.2176 + pgl_point_cluster_distance(&(circle->center), cluster) <= circle->radius 3.2177 + ); 3.2178 + PG_FREE_IF_COPY(cluster, 1); 3.2179 + PG_RETURN_BOOL(retval); 3.2180 +} 3.2181 + 3.2182 +/* check if circle and cluster may be overlap (l. ov. operator "&&+") in SQL */ 3.2183 +PG_FUNCTION_INFO_V1(pgl_ecircle_ecluster_may_overlap); 3.2184 +Datum pgl_ecircle_ecluster_may_overlap(PG_FUNCTION_ARGS) { 3.2185 + pgl_circle *circle = (pgl_circle *)PG_GETARG_POINTER(0); 3.2186 + pgl_cluster *cluster = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); 3.2187 + bool retval = pgl_distance( 3.2188 + circle->center.lat, circle->center.lon, 3.2189 + cluster->bounding.center.lat, cluster->bounding.center.lon 3.2190 + ) <= circle->radius + cluster->bounding.radius; 3.2191 + PG_FREE_IF_COPY(cluster, 1); 3.2192 + PG_RETURN_BOOL(retval); 3.2193 +} 3.2194 + 3.2195 +/* check if two clusters may overlap (lossy overlap operator "&&+") in SQL */ 3.2196 +PG_FUNCTION_INFO_V1(pgl_ecluster_may_overlap); 3.2197 +Datum pgl_ecluster_may_overlap(PG_FUNCTION_ARGS) { 3.2198 + pgl_cluster *cluster1 = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); 3.2199 + pgl_cluster *cluster2 = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); 3.2200 + bool retval = pgl_distance( 3.2201 + cluster1->bounding.center.lat, cluster1->bounding.center.lon, 3.2202 + cluster2->bounding.center.lat, cluster2->bounding.center.lon 3.2203 + ) <= cluster1->bounding.radius + cluster2->bounding.radius; 3.2204 + PG_FREE_IF_COPY(cluster1, 0); 3.2205 + PG_FREE_IF_COPY(cluster2, 1); 3.2206 + PG_RETURN_BOOL(retval); 3.2207 +} 3.2208 + 3.2209 +/* calculate distance between two points ("<->" operator) in SQL */ 3.2210 +PG_FUNCTION_INFO_V1(pgl_epoint_distance); 3.2211 +Datum pgl_epoint_distance(PG_FUNCTION_ARGS) { 3.2212 + pgl_point *point1 = (pgl_point *)PG_GETARG_POINTER(0); 3.2213 + pgl_point *point2 = (pgl_point *)PG_GETARG_POINTER(1); 3.2214 + PG_RETURN_FLOAT8(pgl_distance( 3.2215 + point1->lat, point1->lon, point2->lat, point2->lon 3.2216 + )); 3.2217 +} 3.2218 + 3.2219 +/* calculate point to circle distance ("<->" operator) in SQL */ 3.2220 +PG_FUNCTION_INFO_V1(pgl_epoint_ecircle_distance); 3.2221 +Datum pgl_epoint_ecircle_distance(PG_FUNCTION_ARGS) { 3.2222 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.2223 + pgl_circle *circle = (pgl_circle *)PG_GETARG_POINTER(1); 3.2224 + double distance = pgl_distance( 3.2225 + point->lat, point->lon, circle->center.lat, circle->center.lon 3.2226 + ) - circle->radius; 3.2227 + PG_RETURN_FLOAT8((distance <= 0) ? 0 : distance); 3.2228 +} 3.2229 + 3.2230 +/* calculate point to cluster distance ("<->" operator) in SQL */ 3.2231 +PG_FUNCTION_INFO_V1(pgl_epoint_ecluster_distance); 3.2232 +Datum pgl_epoint_ecluster_distance(PG_FUNCTION_ARGS) { 3.2233 + pgl_point *point = (pgl_point *)PG_GETARG_POINTER(0); 3.2234 + pgl_cluster *cluster = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); 3.2235 + double distance = pgl_point_cluster_distance(point, cluster); 3.2236 + PG_FREE_IF_COPY(cluster, 1); 3.2237 + PG_RETURN_FLOAT8(distance); 3.2238 +} 3.2239 + 3.2240 +/* calculate distance between two circles ("<->" operator) in SQL */ 3.2241 +PG_FUNCTION_INFO_V1(pgl_ecircle_distance); 3.2242 +Datum pgl_ecircle_distance(PG_FUNCTION_ARGS) { 3.2243 + pgl_circle *circle1 = (pgl_circle *)PG_GETARG_POINTER(0); 3.2244 + pgl_circle *circle2 = (pgl_circle *)PG_GETARG_POINTER(1); 3.2245 + double distance = pgl_distance( 3.2246 + circle1->center.lat, circle1->center.lon, 3.2247 + circle2->center.lat, circle2->center.lon 3.2248 + ) - (circle1->radius + circle2->radius); 3.2249 + PG_RETURN_FLOAT8((distance <= 0) ? 0 : distance); 3.2250 +} 3.2251 + 3.2252 +/* calculate circle to cluster distance ("<->" operator) in SQL */ 3.2253 +PG_FUNCTION_INFO_V1(pgl_ecircle_ecluster_distance); 3.2254 +Datum pgl_ecircle_ecluster_distance(PG_FUNCTION_ARGS) { 3.2255 + pgl_circle *circle = (pgl_circle *)PG_GETARG_POINTER(0); 3.2256 + pgl_cluster *cluster = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); 3.2257 + double distance = ( 3.2258 + pgl_point_cluster_distance(&(circle->center), cluster) - circle->radius 3.2259 + ); 3.2260 + PG_FREE_IF_COPY(cluster, 1); 3.2261 + PG_RETURN_FLOAT8((distance <= 0) ? 0 : distance); 3.2262 +} 3.2263 + 3.2264 + 3.2265 +/*-----------------------------------------------------------* 3.2266 + * B-tree comparison operators and index support functions * 3.2267 + *-----------------------------------------------------------*/ 3.2268 + 3.2269 +/* macro for a B-tree operator (without detoasting) */ 3.2270 +#define PGL_BTREE_OPER(func, type, cmpfunc, oper) \ 3.2271 + PG_FUNCTION_INFO_V1(func); \ 3.2272 + Datum func(PG_FUNCTION_ARGS) { \ 3.2273 + type *a = (type *)PG_GETARG_POINTER(0); \ 3.2274 + type *b = (type *)PG_GETARG_POINTER(1); \ 3.2275 + PG_RETURN_BOOL(cmpfunc(a, b) oper 0); \ 3.2276 + } 3.2277 + 3.2278 +/* macro for a B-tree comparison function (without detoasting) */ 3.2279 +#define PGL_BTREE_CMP(func, type, cmpfunc) \ 3.2280 + PG_FUNCTION_INFO_V1(func); \ 3.2281 + Datum func(PG_FUNCTION_ARGS) { \ 3.2282 + type *a = (type *)PG_GETARG_POINTER(0); \ 3.2283 + type *b = (type *)PG_GETARG_POINTER(1); \ 3.2284 + PG_RETURN_INT32(cmpfunc(a, b)); \ 3.2285 + } 3.2286 + 3.2287 +/* macro for a B-tree operator (with detoasting) */ 3.2288 +#define PGL_BTREE_OPER_DETOAST(func, type, cmpfunc, oper) \ 3.2289 + PG_FUNCTION_INFO_V1(func); \ 3.2290 + Datum func(PG_FUNCTION_ARGS) { \ 3.2291 + bool res; \ 3.2292 + type *a = (type *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); \ 3.2293 + type *b = (type *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); \ 3.2294 + res = cmpfunc(a, b) oper 0; \ 3.2295 + PG_FREE_IF_COPY(a, 0); \ 3.2296 + PG_FREE_IF_COPY(b, 1); \ 3.2297 + PG_RETURN_BOOL(res); \ 3.2298 + } 3.2299 + 3.2300 +/* macro for a B-tree comparison function (with detoasting) */ 3.2301 +#define PGL_BTREE_CMP_DETOAST(func, type, cmpfunc) \ 3.2302 + PG_FUNCTION_INFO_V1(func); \ 3.2303 + Datum func(PG_FUNCTION_ARGS) { \ 3.2304 + int32_t res; \ 3.2305 + type *a = (type *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); \ 3.2306 + type *b = (type *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); \ 3.2307 + res = cmpfunc(a, b); \ 3.2308 + PG_FREE_IF_COPY(a, 0); \ 3.2309 + PG_FREE_IF_COPY(b, 1); \ 3.2310 + PG_RETURN_INT32(res); \ 3.2311 + } 3.2312 + 3.2313 +/* B-tree operators and comparison function for point */ 3.2314 +PGL_BTREE_OPER(pgl_btree_epoint_lt, pgl_point, pgl_point_cmp, <) 3.2315 +PGL_BTREE_OPER(pgl_btree_epoint_le, pgl_point, pgl_point_cmp, <=) 3.2316 +PGL_BTREE_OPER(pgl_btree_epoint_eq, pgl_point, pgl_point_cmp, ==) 3.2317 +PGL_BTREE_OPER(pgl_btree_epoint_ne, pgl_point, pgl_point_cmp, !=) 3.2318 +PGL_BTREE_OPER(pgl_btree_epoint_ge, pgl_point, pgl_point_cmp, >=) 3.2319 +PGL_BTREE_OPER(pgl_btree_epoint_gt, pgl_point, pgl_point_cmp, >) 3.2320 +PGL_BTREE_CMP(pgl_btree_epoint_cmp, pgl_point, pgl_point_cmp) 3.2321 + 3.2322 +/* B-tree operators and comparison function for box */ 3.2323 +PGL_BTREE_OPER(pgl_btree_ebox_lt, pgl_box, pgl_box_cmp, <) 3.2324 +PGL_BTREE_OPER(pgl_btree_ebox_le, pgl_box, pgl_box_cmp, <=) 3.2325 +PGL_BTREE_OPER(pgl_btree_ebox_eq, pgl_box, pgl_box_cmp, ==) 3.2326 +PGL_BTREE_OPER(pgl_btree_ebox_ne, pgl_box, pgl_box_cmp, !=) 3.2327 +PGL_BTREE_OPER(pgl_btree_ebox_ge, pgl_box, pgl_box_cmp, >=) 3.2328 +PGL_BTREE_OPER(pgl_btree_ebox_gt, pgl_box, pgl_box_cmp, >) 3.2329 +PGL_BTREE_CMP(pgl_btree_ebox_cmp, pgl_box, pgl_box_cmp) 3.2330 + 3.2331 +/* B-tree operators and comparison function for circle */ 3.2332 +PGL_BTREE_OPER(pgl_btree_ecircle_lt, pgl_circle, pgl_circle_cmp, <) 3.2333 +PGL_BTREE_OPER(pgl_btree_ecircle_le, pgl_circle, pgl_circle_cmp, <=) 3.2334 +PGL_BTREE_OPER(pgl_btree_ecircle_eq, pgl_circle, pgl_circle_cmp, ==) 3.2335 +PGL_BTREE_OPER(pgl_btree_ecircle_ne, pgl_circle, pgl_circle_cmp, !=) 3.2336 +PGL_BTREE_OPER(pgl_btree_ecircle_ge, pgl_circle, pgl_circle_cmp, >=) 3.2337 +PGL_BTREE_OPER(pgl_btree_ecircle_gt, pgl_circle, pgl_circle_cmp, >) 3.2338 +PGL_BTREE_CMP(pgl_btree_ecircle_cmp, pgl_circle, pgl_circle_cmp) 3.2339 + 3.2340 + 3.2341 +/*--------------------------------* 3.2342 + * GiST index support functions * 3.2343 + *--------------------------------*/ 3.2344 + 3.2345 +/* GiST "consistent" support function */ 3.2346 +PG_FUNCTION_INFO_V1(pgl_gist_consistent); 3.2347 +Datum pgl_gist_consistent(PG_FUNCTION_ARGS) { 3.2348 + GISTENTRY *entry = (GISTENTRY *) PG_GETARG_POINTER(0); 3.2349 + pgl_keyptr key = (pgl_keyptr)DatumGetPointer(entry->key); 3.2350 + StrategyNumber strategy = (StrategyNumber)PG_GETARG_UINT16(2); 3.2351 + bool *recheck = (bool *)PG_GETARG_POINTER(4); 3.2352 + /* demand recheck because index and query methods are lossy */ 3.2353 + *recheck = true; 3.2354 + /* strategy number aliases for different operators using the same strategy */ 3.2355 + strategy %= 100; 3.2356 + /* strategy number 11: equality of two points */ 3.2357 + if (strategy == 11) { 3.2358 + /* query datum is another point */ 3.2359 + pgl_point *query = (pgl_point *)PG_GETARG_POINTER(1); 3.2360 + /* convert other point to key */ 3.2361 + pgl_pointkey querykey; 3.2362 + pgl_point_to_key(query, querykey); 3.2363 + /* return true if both keys overlap */ 3.2364 + PG_RETURN_BOOL(pgl_keys_overlap(key, querykey)); 3.2365 + } 3.2366 + /* strategy number 13: equality of two circles */ 3.2367 + if (strategy == 13) { 3.2368 + /* query datum is another circle */ 3.2369 + pgl_circle *query = (pgl_circle *)PG_GETARG_POINTER(1); 3.2370 + /* convert other circle to key */ 3.2371 + pgl_areakey querykey; 3.2372 + pgl_circle_to_key(query, querykey); 3.2373 + /* return true if both keys overlap */ 3.2374 + PG_RETURN_BOOL(pgl_keys_overlap(key, querykey)); 3.2375 + } 3.2376 + /* for all remaining strategies, keys on empty objects produce no match */ 3.2377 + /* (check necessary because query radius may be infinite) */ 3.2378 + if (PGL_KEY_IS_EMPTY(key)) PG_RETURN_BOOL(false); 3.2379 + /* strategy number 21: overlapping with point */ 3.2380 + if (strategy == 21) { 3.2381 + /* query datum is a point */ 3.2382 + pgl_point *query = (pgl_point *)PG_GETARG_POINTER(1); 3.2383 + /* return true if estimated distance (allowed to be smaller than real 3.2384 + distance) between index key and point is zero */ 3.2385 + PG_RETURN_BOOL(pgl_estimate_key_distance(key, query) == 0); 3.2386 + } 3.2387 + /* strategy number 22: (point) overlapping with box */ 3.2388 + if (strategy == 22) { 3.2389 + /* query datum is a box */ 3.2390 + pgl_box *query = (pgl_box *)PG_GETARG_POINTER(1); 3.2391 + /* determine bounding box of indexed key */ 3.2392 + pgl_box keybox; 3.2393 + pgl_key_to_box(key, &keybox); 3.2394 + /* return true if query box overlaps with bounding box of indexed key */ 3.2395 + PG_RETURN_BOOL(pgl_boxes_overlap(query, &keybox)); 3.2396 + } 3.2397 + /* strategy number 23: overlapping with circle */ 3.2398 + if (strategy == 23) { 3.2399 + /* query datum is a circle */ 3.2400 + pgl_circle *query = (pgl_circle *)PG_GETARG_POINTER(1); 3.2401 + /* return true if estimated distance (allowed to be smaller than real 3.2402 + distance) between index key and circle center is smaller than radius */ 3.2403 + PG_RETURN_BOOL( 3.2404 + pgl_estimate_key_distance(key, &(query->center)) <= query->radius 3.2405 + ); 3.2406 + } 3.2407 + /* strategy number 24: overlapping with cluster */ 3.2408 + if (strategy == 24) { 3.2409 + bool retval; /* return value */ 3.2410 + /* query datum is a cluster */ 3.2411 + pgl_cluster *query = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); 3.2412 + /* return true if estimated distance (allowed to be smaller than real 3.2413 + distance) between index key and circle center is smaller than radius */ 3.2414 + retval = ( 3.2415 + pgl_estimate_key_distance(key, &(query->bounding.center)) <= 3.2416 + query->bounding.radius 3.2417 + ); 3.2418 + PG_FREE_IF_COPY(query, 1); /* free detoasted cluster (if copy) */ 3.2419 + PG_RETURN_BOOL(retval); 3.2420 + } 3.2421 + /* throw error for any unknown strategy number */ 3.2422 + elog(ERROR, "unrecognized strategy number: %d", strategy); 3.2423 +} 3.2424 + 3.2425 +/* GiST "union" support function */ 3.2426 +PG_FUNCTION_INFO_V1(pgl_gist_union); 3.2427 +Datum pgl_gist_union(PG_FUNCTION_ARGS) { 3.2428 + GistEntryVector *entryvec = (GistEntryVector *)PG_GETARG_POINTER(0); 3.2429 + pgl_keyptr out; /* return value (to be palloc'ed) */ 3.2430 + int i; 3.2431 + /* determine key size */ 3.2432 + size_t keysize = PGL_KEY_IS_AREAKEY( 3.2433 + (pgl_keyptr)DatumGetPointer(entryvec->vector[0].key) 3.2434 + ) ? sizeof (pgl_areakey) : sizeof(pgl_pointkey); 3.2435 + /* begin with first key as result */ 3.2436 + out = palloc(keysize); 3.2437 + memcpy(out, (pgl_keyptr)DatumGetPointer(entryvec->vector[0].key), keysize); 3.2438 + /* unite current result with second, third, etc. key */ 3.2439 + for (i=1; i<entryvec->n; i++) { 3.2440 + pgl_unite_keys(out, (pgl_keyptr)DatumGetPointer(entryvec->vector[i].key)); 3.2441 + } 3.2442 + /* return result */ 3.2443 + PG_RETURN_POINTER(out); 3.2444 +} 3.2445 + 3.2446 +/* GiST "compress" support function for indicis on points */ 3.2447 +PG_FUNCTION_INFO_V1(pgl_gist_compress_epoint); 3.2448 +Datum pgl_gist_compress_epoint(PG_FUNCTION_ARGS) { 3.2449 + GISTENTRY *entry = (GISTENTRY *) PG_GETARG_POINTER(0); 3.2450 + GISTENTRY *retval; /* return value (to be palloc'ed unless set to entry) */ 3.2451 + /* only transform new leaves */ 3.2452 + if (entry->leafkey) { 3.2453 + /* get point to be transformed */ 3.2454 + pgl_point *point = (pgl_point *)DatumGetPointer(entry->key); 3.2455 + /* allocate memory for key */ 3.2456 + pgl_keyptr key = palloc(sizeof(pgl_pointkey)); 3.2457 + /* transform point to key */ 3.2458 + pgl_point_to_key(point, key); 3.2459 + /* create new GISTENTRY structure as return value */ 3.2460 + retval = palloc(sizeof(GISTENTRY)); 3.2461 + gistentryinit( 3.2462 + *retval, PointerGetDatum(key), 3.2463 + entry->rel, entry->page, entry->offset, FALSE 3.2464 + ); 3.2465 + } else { 3.2466 + /* inner nodes have already been transformed */ 3.2467 + retval = entry; 3.2468 + } 3.2469 + /* return pointer to old or new GISTENTRY structure */ 3.2470 + PG_RETURN_POINTER(retval); 3.2471 +} 3.2472 + 3.2473 +/* GiST "compress" support function for indicis on circles */ 3.2474 +PG_FUNCTION_INFO_V1(pgl_gist_compress_ecircle); 3.2475 +Datum pgl_gist_compress_ecircle(PG_FUNCTION_ARGS) { 3.2476 + GISTENTRY *entry = (GISTENTRY *) PG_GETARG_POINTER(0); 3.2477 + GISTENTRY *retval; /* return value (to be palloc'ed unless set to entry) */ 3.2478 + /* only transform new leaves */ 3.2479 + if (entry->leafkey) { 3.2480 + /* get circle to be transformed */ 3.2481 + pgl_circle *circle = (pgl_circle *)DatumGetPointer(entry->key); 3.2482 + /* allocate memory for key */ 3.2483 + pgl_keyptr key = palloc(sizeof(pgl_areakey)); 3.2484 + /* transform circle to key */ 3.2485 + pgl_circle_to_key(circle, key); 3.2486 + /* create new GISTENTRY structure as return value */ 3.2487 + retval = palloc(sizeof(GISTENTRY)); 3.2488 + gistentryinit( 3.2489 + *retval, PointerGetDatum(key), 3.2490 + entry->rel, entry->page, entry->offset, FALSE 3.2491 + ); 3.2492 + } else { 3.2493 + /* inner nodes have already been transformed */ 3.2494 + retval = entry; 3.2495 + } 3.2496 + /* return pointer to old or new GISTENTRY structure */ 3.2497 + PG_RETURN_POINTER(retval); 3.2498 +} 3.2499 + 3.2500 +/* GiST "compress" support function for indices on clusters */ 3.2501 +PG_FUNCTION_INFO_V1(pgl_gist_compress_ecluster); 3.2502 +Datum pgl_gist_compress_ecluster(PG_FUNCTION_ARGS) { 3.2503 + GISTENTRY *entry = (GISTENTRY *) PG_GETARG_POINTER(0); 3.2504 + GISTENTRY *retval; /* return value (to be palloc'ed unless set to entry) */ 3.2505 + /* only transform new leaves */ 3.2506 + if (entry->leafkey) { 3.2507 + /* get cluster to be transformed (detoasting necessary!) */ 3.2508 + pgl_cluster *cluster = (pgl_cluster *)PG_DETOAST_DATUM(entry->key); 3.2509 + /* allocate memory for key */ 3.2510 + pgl_keyptr key = palloc(sizeof(pgl_areakey)); 3.2511 + /* transform cluster to key */ 3.2512 + pgl_circle_to_key(&(cluster->bounding), key); 3.2513 + /* create new GISTENTRY structure as return value */ 3.2514 + retval = palloc(sizeof(GISTENTRY)); 3.2515 + gistentryinit( 3.2516 + *retval, PointerGetDatum(key), 3.2517 + entry->rel, entry->page, entry->offset, FALSE 3.2518 + ); 3.2519 + /* free detoasted datum */ 3.2520 + if ((void *)cluster != (void *)DatumGetPointer(entry->key)) pfree(cluster); 3.2521 + } else { 3.2522 + /* inner nodes have already been transformed */ 3.2523 + retval = entry; 3.2524 + } 3.2525 + /* return pointer to old or new GISTENTRY structure */ 3.2526 + PG_RETURN_POINTER(retval); 3.2527 +} 3.2528 + 3.2529 +/* GiST "decompress" support function for indices */ 3.2530 +PG_FUNCTION_INFO_V1(pgl_gist_decompress); 3.2531 +Datum pgl_gist_decompress(PG_FUNCTION_ARGS) { 3.2532 + /* return passed pointer without transformation */ 3.2533 + PG_RETURN_POINTER(PG_GETARG_POINTER(0)); 3.2534 +} 3.2535 + 3.2536 +/* GiST "penalty" support function */ 3.2537 +PG_FUNCTION_INFO_V1(pgl_gist_penalty); 3.2538 +Datum pgl_gist_penalty(PG_FUNCTION_ARGS) { 3.2539 + GISTENTRY *origentry = (GISTENTRY *)PG_GETARG_POINTER(0); 3.2540 + GISTENTRY *newentry = (GISTENTRY *)PG_GETARG_POINTER(1); 3.2541 + float *penalty = (float *)PG_GETARG_POINTER(2); 3.2542 + /* get original key and key to insert */ 3.2543 + pgl_keyptr orig = (pgl_keyptr)DatumGetPointer(origentry->key); 3.2544 + pgl_keyptr new = (pgl_keyptr)DatumGetPointer(newentry->key); 3.2545 + /* copy original key */ 3.2546 + union { pgl_pointkey pointkey; pgl_areakey areakey; } union_key; 3.2547 + if (PGL_KEY_IS_AREAKEY(orig)) { 3.2548 + memcpy(union_key.areakey, orig, sizeof(union_key.areakey)); 3.2549 + } else { 3.2550 + memcpy(union_key.pointkey, orig, sizeof(union_key.pointkey)); 3.2551 + } 3.2552 + /* calculate union of both keys */ 3.2553 + pgl_unite_keys((pgl_keyptr)&union_key, new); 3.2554 + /* penalty equal to reduction of key length (logarithm of added area) */ 3.2555 + /* (return value by setting referenced value and returning pointer) */ 3.2556 + *penalty = ( 3.2557 + PGL_KEY_NODEDEPTH(orig) - PGL_KEY_NODEDEPTH((pgl_keyptr)&union_key) 3.2558 + ); 3.2559 + PG_RETURN_POINTER(penalty); 3.2560 +} 3.2561 + 3.2562 +/* GiST "picksplit" support function */ 3.2563 +PG_FUNCTION_INFO_V1(pgl_gist_picksplit); 3.2564 +Datum pgl_gist_picksplit(PG_FUNCTION_ARGS) { 3.2565 + GistEntryVector *entryvec = (GistEntryVector *)PG_GETARG_POINTER(0); 3.2566 + GIST_SPLITVEC *v = (GIST_SPLITVEC *)PG_GETARG_POINTER(1); 3.2567 + OffsetNumber i; /* between FirstOffsetNumber and entryvec->n (inclusive) */ 3.2568 + union { 3.2569 + pgl_pointkey pointkey; 3.2570 + pgl_areakey areakey; 3.2571 + } union_all; /* union of all keys (to be calculated from scratch) 3.2572 + (later cut in half) */ 3.2573 + int is_areakey = PGL_KEY_IS_AREAKEY( 3.2574 + (pgl_keyptr)DatumGetPointer(entryvec->vector[FirstOffsetNumber].key) 3.2575 + ); 3.2576 + int keysize = is_areakey ? sizeof(pgl_areakey) : sizeof(pgl_pointkey); 3.2577 + pgl_keyptr unionL = palloc(keysize); /* union of keys that go left */ 3.2578 + pgl_keyptr unionR = palloc(keysize); /* union of keys that go right */ 3.2579 + pgl_keyptr key; /* current key to be processed */ 3.2580 + /* allocate memory for array of left and right keys, set counts to zero */ 3.2581 + v->spl_left = (OffsetNumber *)palloc(entryvec->n * sizeof(OffsetNumber)); 3.2582 + v->spl_nleft = 0; 3.2583 + v->spl_right = (OffsetNumber *)palloc(entryvec->n * sizeof(OffsetNumber)); 3.2584 + v->spl_nright = 0; 3.2585 + /* calculate union of all keys from scratch */ 3.2586 + memcpy( 3.2587 + (pgl_keyptr)&union_all, 3.2588 + (pgl_keyptr)DatumGetPointer(entryvec->vector[FirstOffsetNumber].key), 3.2589 + keysize 3.2590 + ); 3.2591 + for (i=FirstOffsetNumber+1; i<entryvec->n; i=OffsetNumberNext(i)) { 3.2592 + pgl_unite_keys( 3.2593 + (pgl_keyptr)&union_all, 3.2594 + (pgl_keyptr)DatumGetPointer(entryvec->vector[i].key) 3.2595 + ); 3.2596 + } 3.2597 + /* check if trivial split is necessary due to exhausted key length */ 3.2598 + /* (Note: keys for empty objects must have node depth set to maximum) */ 3.2599 + if (PGL_KEY_NODEDEPTH((pgl_keyptr)&union_all) == ( 3.2600 + is_areakey ? PGL_AREAKEY_MAXDEPTH : PGL_POINTKEY_MAXDEPTH 3.2601 + )) { 3.2602 + /* half of all keys go left */ 3.2603 + for ( 3.2604 + i=FirstOffsetNumber; 3.2605 + i<FirstOffsetNumber+(entryvec->n - FirstOffsetNumber)/2; 3.2606 + i=OffsetNumberNext(i) 3.2607 + ) { 3.2608 + /* pointer to current key */ 3.2609 + key = (pgl_keyptr)DatumGetPointer(entryvec->vector[i].key); 3.2610 + /* update unionL */ 3.2611 + /* check if key is first key that goes left */ 3.2612 + if (!v->spl_nleft) { 3.2613 + /* first key that goes left is just copied to unionL */ 3.2614 + memcpy(unionL, key, keysize); 3.2615 + } else { 3.2616 + /* unite current value and next key */ 3.2617 + pgl_unite_keys(unionL, key); 3.2618 + } 3.2619 + /* append offset number to list of keys that go left */ 3.2620 + v->spl_left[v->spl_nleft++] = i; 3.2621 + } 3.2622 + /* other half goes right */ 3.2623 + for ( 3.2624 + i=FirstOffsetNumber+(entryvec->n - FirstOffsetNumber)/2; 3.2625 + i<entryvec->n; 3.2626 + i=OffsetNumberNext(i) 3.2627 + ) { 3.2628 + /* pointer to current key */ 3.2629 + key = (pgl_keyptr)DatumGetPointer(entryvec->vector[i].key); 3.2630 + /* update unionR */ 3.2631 + /* check if key is first key that goes right */ 3.2632 + if (!v->spl_nright) { 3.2633 + /* first key that goes right is just copied to unionR */ 3.2634 + memcpy(unionR, key, keysize); 3.2635 + } else { 3.2636 + /* unite current value and next key */ 3.2637 + pgl_unite_keys(unionR, key); 3.2638 + } 3.2639 + /* append offset number to list of keys that go right */ 3.2640 + v->spl_right[v->spl_nright++] = i; 3.2641 + } 3.2642 + } 3.2643 + /* otherwise, a non-trivial split is possible */ 3.2644 + else { 3.2645 + /* cut covered area in half */ 3.2646 + /* (union_all then refers to area of keys that go left) */ 3.2647 + /* check if union of all keys covers empty and non-empty objects */ 3.2648 + if (PGL_KEY_IS_UNIVERSAL((pgl_keyptr)&union_all)) { 3.2649 + /* if yes, split into empty and non-empty objects */ 3.2650 + pgl_key_set_empty((pgl_keyptr)&union_all); 3.2651 + } else { 3.2652 + /* otherwise split by next bit */ 3.2653 + ((pgl_keyptr)&union_all)[PGL_KEY_NODEDEPTH_OFFSET]++; 3.2654 + /* NOTE: type bit conserved */ 3.2655 + } 3.2656 + /* determine for each key if it goes left or right */ 3.2657 + for (i=FirstOffsetNumber; i<entryvec->n; i=OffsetNumberNext(i)) { 3.2658 + /* pointer to current key */ 3.2659 + key = (pgl_keyptr)DatumGetPointer(entryvec->vector[i].key); 3.2660 + /* keys within one half of the area go left */ 3.2661 + if (pgl_keys_overlap((pgl_keyptr)&union_all, key)) { 3.2662 + /* update unionL */ 3.2663 + /* check if key is first key that goes left */ 3.2664 + if (!v->spl_nleft) { 3.2665 + /* first key that goes left is just copied to unionL */ 3.2666 + memcpy(unionL, key, keysize); 3.2667 + } else { 3.2668 + /* unite current value of unionL and processed key */ 3.2669 + pgl_unite_keys(unionL, key); 3.2670 + } 3.2671 + /* append offset number to list of keys that go left */ 3.2672 + v->spl_left[v->spl_nleft++] = i; 3.2673 + } 3.2674 + /* the other keys go right */ 3.2675 + else { 3.2676 + /* update unionR */ 3.2677 + /* check if key is first key that goes right */ 3.2678 + if (!v->spl_nright) { 3.2679 + /* first key that goes right is just copied to unionR */ 3.2680 + memcpy(unionR, key, keysize); 3.2681 + } else { 3.2682 + /* unite current value of unionR and processed key */ 3.2683 + pgl_unite_keys(unionR, key); 3.2684 + } 3.2685 + /* append offset number to list of keys that go right */ 3.2686 + v->spl_right[v->spl_nright++] = i; 3.2687 + } 3.2688 + } 3.2689 + } 3.2690 + /* store unions in return value */ 3.2691 + v->spl_ldatum = PointerGetDatum(unionL); 3.2692 + v->spl_rdatum = PointerGetDatum(unionR); 3.2693 + /* return all results */ 3.2694 + PG_RETURN_POINTER(v); 3.2695 +} 3.2696 + 3.2697 +/* GiST "same"/"equal" support function */ 3.2698 +PG_FUNCTION_INFO_V1(pgl_gist_same); 3.2699 +Datum pgl_gist_same(PG_FUNCTION_ARGS) { 3.2700 + pgl_keyptr key1 = (pgl_keyptr)PG_GETARG_POINTER(0); 3.2701 + pgl_keyptr key2 = (pgl_keyptr)PG_GETARG_POINTER(1); 3.2702 + bool *boolptr = (bool *)PG_GETARG_POINTER(2); 3.2703 + /* two keys are equal if they are binary equal */ 3.2704 + /* (return result by setting referenced boolean and returning pointer) */ 3.2705 + *boolptr = !memcmp( 3.2706 + key1, 3.2707 + key2, 3.2708 + PGL_KEY_IS_AREAKEY(key1) ? sizeof(pgl_areakey) : sizeof(pgl_pointkey) 3.2709 + ); 3.2710 + PG_RETURN_POINTER(boolptr); 3.2711 +} 3.2712 + 3.2713 +/* GiST "distance" support function */ 3.2714 +PG_FUNCTION_INFO_V1(pgl_gist_distance); 3.2715 +Datum pgl_gist_distance(PG_FUNCTION_ARGS) { 3.2716 + GISTENTRY *entry = (GISTENTRY *)PG_GETARG_POINTER(0); 3.2717 + pgl_keyptr key = (pgl_keyptr)DatumGetPointer(entry->key); 3.2718 + StrategyNumber strategy = (StrategyNumber)PG_GETARG_UINT16(2); 3.2719 + bool *recheck = (bool *)PG_GETARG_POINTER(4); 3.2720 + double distance; /* return value */ 3.2721 + /* demand recheck because distance is just an estimation */ 3.2722 + /* (real distance may be bigger) */ 3.2723 + *recheck = true; 3.2724 + /* strategy number aliases for different operators using the same strategy */ 3.2725 + strategy %= 100; 3.2726 + /* strategy number 31: distance to point */ 3.2727 + if (strategy == 31) { 3.2728 + /* query datum is a point */ 3.2729 + pgl_point *query = (pgl_point *)PG_GETARG_POINTER(1); 3.2730 + /* use pgl_estimate_pointkey_distance() function to compute result */ 3.2731 + distance = pgl_estimate_key_distance(key, query); 3.2732 + /* avoid infinity (reserved!) */ 3.2733 + if (!isfinite(distance)) distance = PGL_ULTRA_DISTANCE; 3.2734 + /* return result */ 3.2735 + PG_RETURN_FLOAT8(distance); 3.2736 + } 3.2737 + /* strategy number 33: distance to circle */ 3.2738 + if (strategy == 33) { 3.2739 + /* query datum is a circle */ 3.2740 + pgl_circle *query = (pgl_circle *)PG_GETARG_POINTER(1); 3.2741 + /* estimate distance to circle center and substract circle radius */ 3.2742 + distance = ( 3.2743 + pgl_estimate_key_distance(key, &(query->center)) - query->radius 3.2744 + ); 3.2745 + /* convert non-positive values to zero and avoid infinity (reserved!) */ 3.2746 + if (distance <= 0) distance = 0; 3.2747 + else if (!isfinite(distance)) distance = PGL_ULTRA_DISTANCE; 3.2748 + /* return result */ 3.2749 + PG_RETURN_FLOAT8(distance); 3.2750 + } 3.2751 + /* strategy number 34: distance to cluster */ 3.2752 + if (strategy == 34) { 3.2753 + /* query datum is a cluster */ 3.2754 + pgl_cluster *query = (pgl_cluster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); 3.2755 + /* estimate distance to bounding center and substract bounding radius */ 3.2756 + distance = ( 3.2757 + pgl_estimate_key_distance(key, &(query->bounding.center)) - 3.2758 + query->bounding.radius 3.2759 + ); 3.2760 + /* convert non-positive values to zero and avoid infinity (reserved!) */ 3.2761 + if (distance <= 0) distance = 0; 3.2762 + else if (!isfinite(distance)) distance = PGL_ULTRA_DISTANCE; 3.2763 + /* free detoasted cluster (if copy) */ 3.2764 + PG_FREE_IF_COPY(query, 1); 3.2765 + /* return result */ 3.2766 + PG_RETURN_FLOAT8(distance); 3.2767 + } 3.2768 + /* throw error for any unknown strategy number */ 3.2769 + elog(ERROR, "unrecognized strategy number: %d", strategy); 3.2770 +} 3.2771 +