From 5a2ad29783f6751ec418ec7b88e495474aebe9db Mon Sep 17 00:00:00 2001 From: Greg Hennessy Date: Tue, 17 Oct 2023 19:59:38 -0400 Subject: [PATCH] Added a q3c_sindist_bool function to allow for calculating a better Selectivity --- q3c.c | 80 +++++++++++++++++++++++++++++++++++++++++- scripts/q3c--2.0.0.sql | 28 +++++++++++++++ 2 files changed, 107 insertions(+), 1 deletion(-) diff --git a/q3c.c b/q3c.c index 6acb1a5..06c5236 100644 --- a/q3c.c +++ b/q3c.c @@ -33,6 +33,8 @@ #include "utils/array.h" #include "utils/geo_decls.h" #include "catalog/pg_type.h" +#include "nodes/supportnodes.h" +#include "optimizer/optimizer.h" #include "fmgr.h" #if PG_VERSION_NUM >= 90300 #include "access/tupmacs.h" @@ -60,6 +62,7 @@ Datum pgq3c_pixarea(PG_FUNCTION_ARGS); Datum pgq3c_dist(PG_FUNCTION_ARGS); Datum pgq3c_dist_pm(PG_FUNCTION_ARGS); Datum pgq3c_sindist(PG_FUNCTION_ARGS); +Datum pgq3c_sindist_bool(PG_FUNCTION_ARGS); Datum pgq3c_sindist_pm(PG_FUNCTION_ARGS); Datum q3c_strquery(PG_FUNCTION_ARGS); Datum pgq3c_nearby_it(PG_FUNCTION_ARGS); @@ -78,7 +81,7 @@ Datum pgq3c_get_version(PG_FUNCTION_ARGS); Datum pgq3c_sel(PG_FUNCTION_ARGS); Datum pgq3c_seljoin(PG_FUNCTION_ARGS); Datum pgq3c_seloper(PG_FUNCTION_ARGS); - +Datum pgq3c_join_selectivity(PG_FUNCTION_ARGS); /* Dummy function that implements the selectivity operator */ PG_FUNCTION_INFO_V1(pgq3c_seloper); @@ -131,6 +134,66 @@ Datum pgq3c_sel(PG_FUNCTION_ARGS) PG_RETURN_FLOAT8(ratio); } +PG_FUNCTION_INFO_V1(pgq3c_join_selectivity); +Datum pgq3c_join_selectivity(PG_FUNCTION_ARGS) +{ + Node *rawreq = (Node *) PG_GETARG_POINTER(0); + Node *ret = NULL; + Datum radDatum; + double rad; + double ratio; + Node *dummy = NULL; + + if (IsA(rawreq, SupportRequestSelectivity)) + { + SupportRequestSelectivity *req = (SupportRequestSelectivity *) rawreq; + /* we pass a constant as the 5th argument */ + dummy = (Node *) lfirst(list_nth_cell(req->args, 4)); + + /* some simple error checking */ + if (IsA(dummy,Const)) + { + Const *con = (Const *) dummy; + radDatum = PointerGetDatum(con->constvalue); + rad = DatumGetFloat8(radDatum); + ratio = 3.14 * rad * rad / 41252.; /* pi*r^2/whole_sky_area */ + /* clamp at 0, 1*/ + CLAMP_PROBABILITY(ratio); + //elog(WARNING, "HERE0.... %e", ratio); + req->selectivity = ratio; + ret = (Node *) req; + } + } + if (IsA(rawreq, SupportRequestSimplify)) + { + SupportRequestSimplify *req = (SupportRequestSimplify *) rawreq; + + /* not yet implemented, not sure what can usefully be done */ + + PG_RETURN_POINTER((Node *)0); + } + if (IsA(rawreq, SupportRequestCost)) + { + /* Provide some generic estimate */ + SupportRequestCost *req = (SupportRequestCost *) rawreq; + + /* in my experiments making these costs low doesn't change + the total cost very much at all */ + + req->startup = 1e-8; + req->per_tuple = 1e-8; + + ret = (Node *) req; + } + if (IsA(rawreq, SupportRequestRows)) + { + /* not yet implemented, not sure what can usefully be done */ + } + + PG_RETURN_POINTER(ret); +} + + PG_FUNCTION_INFO_V1(pgq3c_seljoin); Datum pgq3c_seljoin(PG_FUNCTION_ARGS) @@ -353,6 +416,21 @@ Datum pgq3c_sindist(PG_FUNCTION_ARGS) PG_RETURN_FLOAT8(res); } +PG_FUNCTION_INFO_V1(pgq3c_sindist_bool); +Datum pgq3c_sindist_bool(PG_FUNCTION_ARGS) +{ + q3c_coord_t ra1 = PG_GETARG_FLOAT8(0); + q3c_coord_t dec1 = PG_GETARG_FLOAT8(1); + q3c_coord_t ra2 = PG_GETARG_FLOAT8(2); + q3c_coord_t dec2 = PG_GETARG_FLOAT8(3); + double rad = PG_GETARG_FLOAT8(4); + q3c_coord_t res; + res = q3c_sindist(ra1, dec1, ra2, dec2); + rad = rad*Q3C_DEGRA/2; + rad = sin(rad)*sin(rad); + PG_RETURN_BOOL(res < rad); +} + PG_FUNCTION_INFO_V1(pgq3c_sindist_pm); Datum pgq3c_sindist_pm(PG_FUNCTION_ARGS) diff --git a/scripts/q3c--2.0.0.sql b/scripts/q3c--2.0.0.sql index 7b4ad27..a7f2ae1 100644 --- a/scripts/q3c--2.0.0.sql +++ b/scripts/q3c--2.0.0.sql @@ -24,6 +24,12 @@ CREATE OR REPLACE FUNCTION q3c_seljoin(internal, oid, internal, int2, internal) AS 'MODULE_PATHNAME', 'pgq3c_seljoin' LANGUAGE C IMMUTABLE STRICT ; +CREATE OR REPLACE FUNCTION q3c_join_selectivity(internal) + RETURNS internal + AS 'MODULE_PATHNAME', 'pgq3c_join_selectivity' + LANGUAGE C IMMUTABLE STRICT; + + -- distance operator with correct selectivity CREATE OPERATOR ==<<>>== ( @@ -96,6 +102,15 @@ COMMENT ON FUNCTION q3c_sindist(double precision, double precision, double precision, double precision) IS 'Function q3c_sindist(ra1, dec1, ra2, dec2) computing the sin(distance/2)^2 between points (ra1, dec1) and (ra2, dec2)'; +CREATE OR REPLACE FUNCTION q3c_sindist_bool(leftra double precision, leftde double precision, + rightra double precision, rightde double precision,rad double precision) + RETURNS bool + AS 'MODULE_PATHNAME', 'pgq3c_sindist_bool' + LANGUAGE C IMMUTABLE STRICT COST 1 SUPPORT q3c_join_selectivity; +COMMENT ON FUNCTION q3c_sindist_bool(double precision, double precision, + double precision, double precision,double precision) + IS 'Function q3c_sindist_bool(ra1, dec1, ra2, dec2, rad) is rad than the sin(distance/2)^2 between points (ra1, dec1) and (ra2, dec2)'; + CREATE OR REPLACE FUNCTION q3c_sindist_pm( ra1 double precision, dec1 double precision, pmra1 double precision, pmdec1 double precision, @@ -209,6 +224,19 @@ SELECT (((q3c_ang2ipix($3,$4)>=(q3c_nearby_it($1,$2,$5,0))) AND (q3c_ang2ipix($3 AND ($5::double precision ==<<>>== ($1,$2,$3,$4)::q3c_type) ' LANGUAGE SQL IMMUTABLE; +CREATE OR REPLACE FUNCTION q3c_join_bool(leftra double precision, leftdec double precision, + rightra double precision, rightdec double precision, + radius double precision) + RETURNS boolean AS +' +SELECT (((q3c_ang2ipix($3,$4)>=(q3c_nearby_it($1,$2,$5,0))) AND (q3c_ang2ipix($3,$4)<=(q3c_nearby_it($1,$2,$5,1)))) + OR ((q3c_ang2ipix($3,$4)>=(q3c_nearby_it($1,$2,$5,2))) AND (q3c_ang2ipix($3,$4)<=(q3c_nearby_it($1,$2,$5,3)))) + OR ((q3c_ang2ipix($3,$4)>=(q3c_nearby_it($1,$2,$5,4))) AND (q3c_ang2ipix($3,$4)<=(q3c_nearby_it($1,$2,$5,5)))) + OR ((q3c_ang2ipix($3,$4)>=(q3c_nearby_it($1,$2,$5,6))) AND (q3c_ang2ipix($3,$4)<=(q3c_nearby_it($1,$2,$5,7))))) + AND q3c_sindist_bool($1,$2,$3,$4,$5) + AND ($5::double precision ==<<>>== ($1,$2,$3,$4)::q3c_type) +' LANGUAGE SQL IMMUTABLE; + CREATE OR REPLACE FUNCTION q3c_join(leftra double precision, leftdec double precision, rightra real, rightdec real,