<div dir="ltr">Hi Stephen,<div><br></div><div style>Building a topology would definitively help in this situation, though it may take some time on very large dataset I guess.</div><div style>If you plan to use some topological functions on the dataset in addition with pgRouting functions, it may be worth the effort. </div>
<div style><br></div><div style>Concerning st_union and its magic "segmentize" feature, would it be possible to divide the initial set of lines into smaller areas and process these subsets to avoid filling up the memory ?</div>
<div style><br></div><div style>Looking at this subject recently (cutting lines by points, cf. <a href="http://trac.osgeo.org/postgis/wiki/UsersWikiSplitPolygonWithPoints">http://trac.osgeo.org/postgis/wiki/UsersWikiSplitPolygonWithPoints</a>) I found that linear referencing functions can help in such a case.</div>
<div style><br></div><div style>The principle is to get the location index of intersection points for each line, and then to cut this line by its locations, using st_line_substring.</div><div style>It appears to be very efficient, using st_dwithin to trigger spatial index, then joining on the lines primary keys, which should be fast.</div>
<div style><br></div><div style>In your usecase, intersection nodes between lines have to be identified before their locations can be computed.</div><div style><br></div><div style>Concerning the tolerance, I'm pretty sure snapping the input dataset to a grid would help to run a precise st_intersection between lines.</div>
<div style><br></div><div style>Based on the linestring sample data, here is the query using linear referencing. It uses CTE subqueries to identify each step:</div><div style><br></div><div style><div><font face="courier new, monospace">with lines as (</font></div>
<div><font face="courier new, monospace"> select 1 as gid, 'MULTILINESTRING((0 1,2 1))'::geometry as geom</font></div><div><font face="courier new, monospace"> union all</font></div><div><font face="courier new, monospace"> select 2 as gid, 'MULTILINESTRING((1 0,1 2))'::geometry as geom</font></div>
<div><font face="courier new, monospace"> union all</font></div><div><font face="courier new, monospace"> select 3 as gid, 'LINESTRING(1 1.5,2 2)'::geometry as geom</font></div><div><font face="courier new, monospace">), </font></div>
<div><font face="courier new, monospace">-- multilinestrings are dumped into simple objects</font></div><div style><font face="courier new, monospace">-- if multilinestrings have several parts, one should generate a unique id based</font></div>
<div style><font face="courier new, monospace">-- on their gid and path into the collection.</font></div><div><font face="courier new, monospace">dumped_lines as (</font></div><div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>select gid, (st_dump(l.geom)).geom</font></div>
<div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>from lines l</font></div><div><font face="courier new, monospace">), </font></div><div style><font face="courier new, monospace">-- This query computes th</font><span style="font-family:'courier new',monospace">e locations, for each input line, of the intersection points with other lines.</span></div>
<div><font face="courier new, monospace">-- this will be used to cut lines based on these locations.</font></div><div style><font face="courier new, monospace">-- to be able to cut lines from their beginning to their end, we generate the 0 and 1 location index</font></div>
<div><font face="courier new, monospace">cut_locations as (</font></div><div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>select l1.gid as lgid, st_line_locate_point(l1.geom, st_intersection(l1.geom, l2.geom)) as locus </font></div>
<div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>from dumped_lines l1 join dumped_lines l2 on (st_dwithin(l1.geom, l2.geom, 0.01))</font></div><div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>where l1.gid <> l2.gid</font></div>
<div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>-- then generates start and end locus for each line, to be able to cut them</font></div><div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>UNION ALL</font></div>
<div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>select l.gid as lgid, 0 as locus</font></div><div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>from dumped_lines l</font></div>
<div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>UNION ALL</font></div><div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>select l.gid as lgid, 1 as locus</font></div>
<div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>from dumped_lines l</font></div><div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>order by lgid, locus</font></div>
<div><font face="courier new, monospace">), </font></div><div><font face="courier new, monospace">-- This query generates a row_number index column for each input line and intersection point. </font></div><div><font face="courier new, monospace">-- This index will be used to self-join the table to cut a line between two consecutive locations </font></div>
<div style><font face="courier new, monospace">-- (idx, idx+1) pairs.</font></div><div style><font face="courier new, monospace">-- window function is used to generate the index inside each line partition</font></div><div>
<font face="courier new, monospace">loc_with_idx as (</font></div><div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>select lgid, locus, row_number() over (partition by lgid order by locus) as idx</font></div>
<div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>from cut_locations</font></div><div><font face="courier new, monospace">) </font></div><div><font face="courier new, monospace">-- finally, each original line is cut with consecutive locations using linear referencing function.</font></div>
<div style><font face="courier new, monospace">-- a filtering is done to eliminate points produced when lines connect at their ends</font></div><div><font face="courier new, monospace">select l.gid, loc1.idx as sub_id, st_line_substring(l.geom, loc1.locus, loc2.locus) as geom ,</font></div>
<div><font face="courier new, monospace"><span class="" style="white-space:pre"> </span>st_geometryType(st_line_substring(l.geom, loc1.locus, loc2.locus)) as type</font></div><div><font face="courier new, monospace">from loc_with_idx loc1 join loc_with_idx loc2 using (lgid) join dumped_lines l on (l.gid = loc1.lgid)</font></div>
<div><font face="courier new, monospace">where loc2.idx = loc1.idx+1</font></div><div><font face="courier new, monospace">-- filter out point geometries occuring if intersection point is at line's start or end point.</font></div>
<div><font face="courier new, monospace">-- there must be a faster way to filter out theses geometries.</font></div><div><font face="courier new, monospace">and st_geometryType(st_line_substring(l.geom, loc1.locus, loc2.locus)) <> 'ST_Point';</font></div>
<div><br></div><div><br></div><div style>A new unique ID key can be computed based on line gid and subgid generated by the query.</div><div style>Initial line attributes can be moved to the new segments using the line gid key.</div>
<div style><br></div><div style>Nicolas</div></div></div><div class="gmail_extra"><br><br><div class="gmail_quote">On 8 May 2013 16:27, Stephen Woodbridge <span dir="ltr"><<a href="mailto:woodbri@swoodbridge.com" target="_blank">woodbri@swoodbridge.com</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Hi all,<br>
<br>
This question comes up reasonably often on the pgRouting list and has been posted he on occasion under titles like "How to break streets at intersections?"<br>
<br>
It seems to me that this would be a good function to create in either postgis or pgrouting.<br>
<br>
THE PROBLEM:<br>
<br>
I have a table of 10's of thousands of street segments to 10's of millions of street segments. These street segments are LINSTRING or MULTILINESTRING geometries with some arbitrary number of attribute columns. The geometries may cross one another and are not noded correctly for use with pgRouting.<br>
<br>
THE RESULTS:<br>
<br>
We want to process the table and create a new table with the same structure (see comment about primary key below), and in the new table all the geometries are broken at intersections and all the new pieces of the original segment that have been broken have the original attributes propagated to them. So if the original segment has column foo='abc' and was split into 3 new segments, each of the three new segments would also have foo='abc'. The exception to this might be that the new table needs a new primary column as the old primary key will now be duplicated for the multiple parts.<br>
<br>
POTENTIAL SOLUTIONS:<br>
<br>
1. I think one way to do this would be to create a topology and load the table into it, then extra a new table from the topology. Although I'm not sure of the specifics for doing this or the efficency of doing it this way.<br>
<br>
2. Another way seems to be using a query like:<br>
<br>
select (st_dump(bar.the_geom)).* from (<br>
select st_union(foo.the_geom) as the_geom from mytable foo<br>
) as bar;<br>
<br>
And then taking each of the dump.geom objects and using st_contains to find which original segment it belonged to so we can move the attributes to the new segment. This method also loose any association to the original record and forces the use of st_contains to re-associate the new segments to the original segments.<br>
<br>
My concern with this is that the st_union has to load the whole table which may be 10's of millions of street segments and this will likely be a memory problem. Also running the st_contains() does not seems to me to be optimal.<br>
<br>
3. Is there a good recipe for doing this somewhere that I have not found? or other better approaches to this problem?<br>
<br>
What would be the best way to add tolerance to the problem? using snap to grid?<br>
<br>
Thoughts on how to do this efficiently?<br>
<br>
Since I'm working on the pgRouting 2.0 release I thought this might be a nice function to add to that if we can come up with a generic way to do this.<br>
<br>
Thanks,<br>
-Steve<br>
<br>
<br>
-- Example to demonstrate st_union above<br>
select st_astext((st_dump(bar.the_<u></u>geom)).geom) from (<br>
select st_union(foo.the_geom) as the_geom from (<br>
select 'MULTILINESTRING((0 1,2 1))'::geometry as the_geom<br>
union all<br>
select 'MULTILINESTRING((1 0,1 2))'::geometry as the_geom<br>
union all<br>
select 'LINESTRING(1 1.5,2 2)'::geometry as the_geom<br>
) as foo<br>
) as bar;<br>
<br>
"LINESTRING(1 1.5,2 2)"<br>
"LINESTRING(1 0,1 1)"<br>
"LINESTRING(1 1,1 1.5)"<br>
"LINESTRING(1 1.5,1 2)"<br>
"LINESTRING(0 1,1 1)"<br>
"LINESTRING(1 1,2 1)"<br>
______________________________<u></u>_________________<br>
postgis-users mailing list<br>
<a href="mailto:postgis-users@lists.osgeo.org" target="_blank">postgis-users@lists.osgeo.org</a><br>
<a href="http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users" target="_blank">http://lists.osgeo.org/cgi-<u></u>bin/mailman/listinfo/postgis-<u></u>users</a><br>
</blockquote></div><br></div>