1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
|
-- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
--
-- PostGIS - Spatial Types for PostgreSQL
-- http://postgis.refractions.net
--
-- Copyright (C) 2011 Andrea Peri <aperi2007@gmail.com>
--
-- This is free software; you can redistribute and/or modify it under
-- the terms of the GNU General Public Licence. See the COPYING file.
--
-- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
--{
--
-- Andrea Peri (27 Feb 2011) creation
--
-- GetFaceByPoint(atopology, point, tol)
--
-- Retrieve a Face ID given a POINT and a tolerance
-- tolerance = 0 mean exactly intersection
--
-- Returns return the integer ID if there is a face on the Point.
--
-- When the Point is even a Node it raise an exception.
-- This case is testable with the GetNodeByPoint(atopology, apoint, tol)
--
-- If there isn't any face in the Point, GetFaceByPoint return 0.
--
-- if near the point there are two or more faces it throw an exception.
--
CREATE OR REPLACE FUNCTION topology.GetFaceByPoint(atopology varchar, apoint geometry, tol1 float8)
RETURNS int
AS
$$
DECLARE
sql text;
idface int;
BEGIN
idface := -1;
--
-- Atopology and apoint are required
--
IF atopology IS NULL OR apoint IS NULL THEN
RAISE EXCEPTION 'Invalid null argument';
END IF;
--
-- Apoint must be a point
--
IF substring(geometrytype(apoint), 1, 5) != 'POINT'
THEN
RAISE EXCEPTION 'Node geometry must be a point';
END IF;
--
-- Tolerance must be >= 0
--
IF tol1 < 0
THEN
RAISE EXCEPTION 'Tolerance must be >=0';
END IF;
--
-- first test is to check if there is inside an mbr
--
if tol1 = 0 then
sql := 'SELECT a.face_id FROM '
|| quote_ident(atopology)
|| '.face as a WHERE '
|| '(a.mbr && ' || quote_literal(apoint::text)||'::geometry) '
|| 'LIMIT 1;';
else
sql := 'SELECT a.face_id FROM '
|| quote_ident(atopology)
|| '.face as a WHERE '
|| '(ST_DWithin(a.mbr,' || quote_literal(apoint::text)||'::geometry,' || tol1::text || ') ) '
|| 'LIMIT 1;';
end if;
BEGIN
EXECUTE sql INTO STRICT idface;
EXCEPTION
WHEN NO_DATA_FOUND THEN
idface = 0;
END;
if idface > 0 then
--
-- probably there is something so now check the exact test
--
if tol1 = 0 then
sql := 'SELECT e.face_id FROM ('
|| 'SELECT d.face_id,ST_BuildArea(ST_Union(geom)) as geom FROM ('
|| 'SELECT b.edge_id as edge_id,b.left_face as face_id,b.geom as geom FROM '
|| quote_ident(atopology) || '.edge_data as b,'
|| '(SELECT a.face_id FROM '
|| quote_ident(atopology) || '.face as a '
|| 'WHERE ST_Intersects(a.mbr,' || quote_literal(apoint::text)||'::geometry)=true'
|| ') as c '
|| 'WHERE (b.left_face = c.face_id) '
|| ' UNION ALL '
|| 'SELECT b.edge_id as edge_id, b.right_face as face_id, b.geom as geom FROM '
|| quote_ident(atopology) || '.edge_data as b,'
|| '(SELECT a.face_id FROM '
|| quote_ident(atopology) || '.face as a '
|| 'WHERE ST_Intersects(a.mbr,' || quote_literal(apoint::text)||'::geometry)=true'
|| ') as c '
|| 'WHERE (b.right_face = c.face_id) '
|| ') as d '
|| 'GROUP BY face_id '
|| ') as e '
|| 'WHERE ST_Intersects(e.geom, ' || quote_literal(apoint::text)||'::geometry)=true;';
else
sql := 'SELECT e.face_id FROM ('
|| 'SELECT d.face_id,ST_BuildArea(ST_Union(geom)) as geom FROM ('
|| 'SELECT b.edge_id as edge_id,b.left_face as face_id,b.geom as geom FROM '
|| quote_ident(atopology) || '.edge_data as b,'
|| '(SELECT a.face_id FROM '
|| quote_ident(atopology) || '.face as a '
|| 'WHERE ST_DWithin(a.mbr,' || quote_literal(apoint::text)||'::geometry,' || tol1::text || ')=true'
|| ') as c '
|| 'WHERE (b.left_face = c.face_id) '
|| ' UNION ALL '
|| 'SELECT b.edge_id as edge_id, b.right_face as face_id, b.geom as geom FROM '
|| quote_ident(atopology) || '.edge_data as b,'
|| '(SELECT a.face_id FROM '
|| quote_ident(atopology) || '.face as a '
|| 'WHERE ST_DWithin(a.mbr,' || quote_literal(apoint::text)||'::geometry,' || tol1::text || ')=true'
|| ') as c '
|| 'WHERE (b.right_face = c.face_id) '
|| ') as d '
|| 'GROUP BY face_id '
|| ') as e '
|| 'WHERE ST_DWithin(e.geom, ' || quote_literal(apoint::text)||'::geometry,' || tol1::text || ')=true;';
end if;
RAISE DEBUG ' ==> %',sql;
BEGIN
EXECUTE sql INTO STRICT idface;
EXCEPTION
WHEN NO_DATA_FOUND THEN
idface = 0;
WHEN TOO_MANY_ROWS THEN
RAISE EXCEPTION 'Two or more faces found';
END;
end if;
RETURN idface;
END
$$
LANGUAGE 'plpgsql' STABLE STRICT;
--} GetFaceByPoint
|
