1 /*** 2 * 3 * Copyright (C) 2011-2017 Regina Obe and Leo Hsu (Paragon Corporation) 4 **/ 5-- This function given a point try to determine the approximate street address (norm_addy form) 6-- and array of cross streets, as well as interpolated points along the streets 7-- Use case example an address at the intersection of 3 streets: SELECT pprint_addy(r.addy[1]) As st1, pprint_addy(r.addy[2]) As st2, pprint_addy(r.addy[3]) As st3, array_to_string(r.street, ',') FROM reverse_geocode(ST_GeomFromText('POINT(-71.057811 42.358274)',4269)) As r; 8--set search_path=tiger,public; 9 10CREATE OR REPLACE FUNCTION reverse_geocode(IN pt geometry, IN include_strnum_range boolean DEFAULT false, OUT intpt geometry[], OUT addy norm_addy[], OUT street character varying[]) 11 RETURNS record AS 12$BODY$ 13DECLARE 14 var_redge RECORD; 15 var_state text := NULL; 16 var_stusps text := NULL; 17 var_countyfp text := NULL; 18 var_addy NORM_ADDY; 19 var_addy_alt NORM_ADDY; 20 var_nstrnum numeric(10); 21 var_primary_line geometry := NULL; 22 var_primary_dist numeric(10,2) ; 23 var_pt geometry; 24 var_place varchar; 25 var_county varchar; 26 var_stmt text; 27 var_debug boolean = get_geocode_setting('debug_reverse_geocode')::boolean; 28 var_rating_highway integer = COALESCE(get_geocode_setting('reverse_geocode_numbered_roads')::integer,0);/**0 no preference, 1 prefer highway number, 2 prefer local name **/ 29 var_zip varchar := NULL; 30 var_primary_fullname varchar := ''; 31BEGIN 32 IF pt IS NULL THEN 33 RETURN; 34 ELSE 35 IF ST_SRID(pt) = 4269 THEN 36 var_pt := pt; 37 ELSIF ST_SRID(pt) > 0 THEN 38 var_pt := ST_Transform(pt, 4269); 39 ELSE --If srid is unknown, assume its 4269 40 var_pt := ST_SetSRID(pt, 4269); 41 END IF; 42 var_pt := ST_SnapToGrid(var_pt, 0.00005); /** Get rid of floating point junk that would prevent intersections **/ 43 END IF; 44 -- Determine state tables to check 45 -- this is needed to take advantage of constraint exclusion 46 IF var_debug THEN 47 RAISE NOTICE 'Get matching states start: %', clock_timestamp(); 48 END IF; 49 SELECT statefp, stusps INTO var_state, var_stusps FROM state WHERE ST_Intersects(the_geom, var_pt) LIMIT 1; 50 IF var_debug THEN 51 RAISE NOTICE 'Get matching states end: % - %', var_state, clock_timestamp(); 52 END IF; 53 IF var_state IS NULL THEN 54 -- We don't have any data for this state 55 RETURN; 56 END IF; 57 IF var_debug THEN 58 RAISE NOTICE 'Get matching counties start: %', clock_timestamp(); 59 END IF; 60 -- locate county 61 var_stmt := 'SELECT countyfp, name FROM county WHERE statefp = $1 AND ST_Intersects(the_geom, $2) LIMIT 1;'; 62 EXECUTE var_stmt INTO var_countyfp, var_county USING var_state, var_pt ; 63 64 --locate zip 65 var_stmt := 'SELECT zcta5ce FROM zcta5 WHERE statefp = $1 AND ST_Intersects(the_geom, $2) LIMIT 1;'; 66 EXECUTE var_stmt INTO var_zip USING var_state, var_pt; 67 -- locate city 68 IF var_zip > '' THEN 69 var_addy.zip := var_zip ; 70 END IF; 71 72 var_stmt := 'SELECT z.name FROM place As z WHERE z.statefp = $1 AND ST_Intersects(the_geom, $2) LIMIT 1;'; 73 EXECUTE var_stmt INTO var_place USING var_state, var_pt ; 74 IF var_place > '' THEN 75 var_addy.location := var_place; 76 ELSE 77 var_stmt := 'SELECT z.name FROM cousub As z WHERE z.statefp = $1 AND ST_Intersects(the_geom, $2) LIMIT 1;'; 78 EXECUTE var_stmt INTO var_place USING var_state, var_pt ; 79 IF var_place > '' THEN 80 var_addy.location := var_place; 81 -- ELSIF var_zip > '' THEN 82 -- SELECT z.city INTO var_place FROM zip_lookup_base As z WHERE z.statefp = var_state AND z.county = var_county AND z.zip = var_zip LIMIT 1; 83 -- var_addy.location := var_place; 84 END IF; 85 END IF; 86 87 IF var_debug THEN 88 RAISE NOTICE 'Get matching counties end: % - %',var_countyfp, clock_timestamp(); 89 END IF; 90 IF var_countyfp IS NULL THEN 91 -- We don't have any data for this county 92 RETURN; 93 END IF; 94 95 var_addy.stateAbbrev = var_stusps; 96 97 -- Find the street edges that this point is closest to with tolerance of 0.005 but only consider the edge if the point is contained in the right or left face 98 -- Then order addresses by proximity to road 99 IF var_debug THEN 100 RAISE NOTICE 'Get matching edges start: %', clock_timestamp(); 101 END IF; 102 103 var_stmt := ' 104 WITH ref AS ( 105 SELECT ' || quote_literal(var_pt::text) || '::geometry As ref_geom ) , 106 f AS 107 ( SELECT faces.* FROM faces CROSS JOIN ref 108 WHERE faces.statefp = ' || quote_literal(var_state) || ' AND faces.countyfp = ' || quote_literal(var_countyfp) || ' 109 AND ST_Intersects(faces.the_geom, ref_geom) 110 ), 111 e AS 112 ( SELECT edges.tlid , edges.statefp, edges.the_geom, CASE WHEN edges.tfidr = f.tfid THEN ''R'' WHEN edges.tfidl = f.tfid THEN ''L'' ELSE NULL END::varchar As eside, 113 ST_ClosestPoint(edges.the_geom,ref_geom) As center_pt, ref_geom 114 FROM edges INNER JOIN f ON (f.statefp = edges.statefp AND (edges.tfidr = f.tfid OR edges.tfidl = f.tfid)) 115 CROSS JOIN ref 116 WHERE edges.statefp = ' || quote_literal(var_state) || ' AND edges.countyfp = ' || quote_literal(var_countyfp) || ' 117 AND ST_DWithin(edges.the_geom, ref.ref_geom, 0.01) AND (edges.mtfcc LIKE ''S%'') --only consider streets and roads 118 ) , 119 ea AS 120 (SELECT e.statefp, e.tlid, a.fromhn, a.tohn, e.center_pt, ref_geom, a.zip, a.side, e.the_geom 121 FROM e LEFT JOIN addr As a ON (a.statefp = ' || quote_literal(var_state) || ' AND e.tlid = a.tlid and e.eside = a.side) 122 ) 123 SELECT * 124 FROM (SELECT DISTINCT ON(tlid,side) foo.fullname, foo.predirabrv, foo.streetname, foo.sufdirabrv, foo.streettypeabbrev, foo.zip, foo.center_pt, 125 side, to_number(CASE WHEN trim(fromhn) ~ ''^[0-9]+$'' THEN fromhn ELSE NULL END,''99999999'') As fromhn, to_number(CASE WHEN trim(tohn) ~ ''^[0-9]+$'' THEN tohn ELSE NULL END,''99999999'') As tohn, 126 ST_GeometryN(ST_Multi(line),1) As line, dist 127 FROM 128 (SELECT e.tlid, e.the_geom As line, n.fullname, COALESCE(n.prequalabr || '' '','''') || n.name AS streetname, n.predirabrv, COALESCE(suftypabrv, pretypabrv) As streettypeabbrev, 129 n.sufdirabrv, e.zip, e.side, e.fromhn, e.tohn , e.center_pt, 130 ST_DistanceSphere(ST_SetSRID(e.center_pt,4326),ST_SetSRID(ref_geom,4326)) As dist 131 FROM ea AS e 132 LEFT JOIN (SELECT featnames.* FROM featnames 133 WHERE featnames.statefp = ' || quote_literal(var_state) ||' ) AS n ON (n.statefp = e.statefp AND n.tlid = e.tlid) 134 ORDER BY dist LIMIT 50 ) As foo 135 ORDER BY foo.tlid, foo.side, '; 136 137 -- for numbered street/road use var_rating_highway to determine whether to prefer numbered or not (0 no pref, 1 prefer numbered, 2 prefer named) 138 var_stmt := var_stmt || ' CASE $1 WHEN 0 THEN 0 WHEN 1 THEN CASE WHEN foo.fullname ~ ''[0-9]+'' THEN 0 ELSE 1 END ELSE CASE WHEN foo.fullname > '''' AND NOT (foo.fullname ~ ''[0-9]+'') THEN 0 ELSE 1 END END '; 139 var_stmt := var_stmt || ', foo.fullname ASC NULLS LAST, dist LIMIT 50) As f ORDER BY f.dist, CASE WHEN fullname > '''' THEN 0 ELSE 1 END '; --don't bother penalizing for distance if less than 20 meters 140 141 IF var_debug = true THEN 142 RAISE NOTICE 'Statement 1: %', replace(var_stmt, '$1', var_rating_highway::text); 143 END IF; 144 145 FOR var_redge IN EXECUTE var_stmt USING var_rating_highway LOOP 146 IF var_debug THEN 147 RAISE NOTICE 'Start Get matching edges loop: %,%', var_primary_line, clock_timestamp(); 148 END IF; 149 IF var_primary_line IS NULL THEN --this is the first time in the loop and our primary guess 150 var_primary_line := var_redge.line; 151 var_primary_dist := var_redge.dist; 152 END IF; 153 154 IF var_redge.fullname IS NOT NULL AND COALESCE(var_primary_fullname,'') = '' THEN -- this is the first non-blank name we are hitting grab info 155 var_primary_fullname := var_redge.fullname; 156 var_addy.streetname = var_redge.streetname; 157 var_addy.streettypeabbrev := var_redge.streettypeabbrev; 158 var_addy.predirabbrev := var_redge.predirabrv; 159 var_addy.postDirAbbrev := var_redge.sufdirabrv; 160 END IF; 161 162 IF ST_Intersects(var_redge.line, var_primary_line) THEN 163 var_addy.streetname := var_redge.streetname; 164 165 var_addy.streettypeabbrev := var_redge.streettypeabbrev; 166 var_addy.address := var_nstrnum; 167 IF var_redge.fromhn IS NOT NULL THEN 168 --interpolate the number -- note that if fromhn > tohn we will be subtracting which is what we want 169 var_nstrnum := (var_redge.fromhn + ST_LineLocatePoint(var_redge.line, var_pt)*(var_redge.tohn - var_redge.fromhn))::numeric(10); 170 -- The odd even street number side of street rule 171 IF (var_nstrnum % 2) != (var_redge.tohn % 2) THEN 172 var_nstrnum := CASE WHEN var_nstrnum + 1 NOT BETWEEN var_redge.fromhn AND var_redge.tohn THEN var_nstrnum - 1 ELSE var_nstrnum + 1 END; 173 END IF; 174 var_addy.address := var_nstrnum; 175 END IF; 176 IF var_redge.zip > '' THEN 177 var_addy.zip := var_redge.zip; 178 ELSE 179 var_addy.zip := var_zip; 180 END IF; 181 -- IF var_redge.location > '' THEN 182 -- var_addy.location := var_redge.location; 183 -- ELSE 184 -- var_addy.location := var_place; 185 -- END IF; 186 187 -- This is a cross streets - only add if not the primary adress street 188 IF var_redge.fullname > '' AND var_redge.fullname <> var_primary_fullname THEN 189 street := array_append(street, (CASE WHEN include_strnum_range THEN COALESCE(var_redge.fromhn::varchar, '')::varchar || COALESCE(' - ' || var_redge.tohn::varchar,'')::varchar || ' '::varchar ELSE '' END::varchar || COALESCE(var_redge.fullname::varchar,''))::varchar); 190 END IF; 191 192 -- consider this a potential address 193 IF (var_redge.dist < var_primary_dist*1.1 OR var_redge.dist < 20) THEN 194 -- We only consider this a possible address if it is really close to our point 195 intpt := array_append(intpt,var_redge.center_pt); 196 -- note that ramps don't have names or addresses but they connect at the edge of a range 197 -- so for ramps the address of connecting is still useful 198 IF var_debug THEN 199 RAISE NOTICE 'Current addresses: %, last added, %, street: %, %', addy, var_addy, var_addy.streetname, clock_timestamp(); 200 END IF; 201 addy := array_append(addy, var_addy); 202 203 -- Use current values streetname for previous value if previous value has no streetname 204 IF var_addy.streetname > '' AND array_upper(addy,1) > 1 AND COALESCE(addy[array_upper(addy,1) - 1].streetname, '') = '' THEN 205 -- the match is probably an offshoot of some sort 206 -- replace prior entry with streetname of new if prior had no streetname 207 var_addy_alt := addy[array_upper(addy,1)- 1]; 208 IF var_debug THEN 209 RAISE NOTICE 'Replacing answer : %, %', addy[array_upper(addy,1) - 1], clock_timestamp(); 210 END IF; 211 var_addy_alt.streetname := var_addy.streetname; 212 var_addy_alt.streettypeabbrev := var_addy.streettypeabbrev; 213 var_addy_alt.predirabbrev := var_addy.predirabbrev; 214 var_addy_alt.postDirAbbrev := var_addy.postDirAbbrev; 215 addy[array_upper(addy,1) - 1 ] := var_addy_alt; 216 IF var_debug THEN 217 RAISE NOTICE 'Replaced with : %, %', var_addy_alt, clock_timestamp(); 218 END IF; 219 END IF; 220 221 IF var_debug THEN 222 RAISE NOTICE 'End Get matching edges loop: %', clock_timestamp(); 223 RAISE NOTICE 'Final addresses: %, %', addy, clock_timestamp(); 224 END IF; 225 226 END IF; 227 END IF; 228 229 END LOOP; 230 231 -- not matching roads or streets, just return basic info 232 IF NOT FOUND THEN 233 addy := array_append(addy,var_addy); 234 IF var_debug THEN 235 RAISE NOTICE 'No address found: adding: % street: %, %', var_addy, var_addy.streetname, clock_timestamp(); 236 END IF; 237 END IF; 238 IF var_debug THEN 239 RAISE NOTICE 'current array count : %, %', array_upper(addy,1), clock_timestamp(); 240 END IF; 241 242 RETURN; 243END; 244$BODY$ 245 LANGUAGE plpgsql STABLE 246 COST 1000; 247