[postgis-users] Need a method for "noding" a street network

Nicolas Ribot nicolas.ribot at gmail.com
Wed May 8 13:05:32 PDT 2013


Is it possible for you to share the dataset or a subset of it, to give it a
try ?

Nicolas


On 8 May 2013 22:04, Nicolas Ribot <nicolas.ribot at gmail.com> wrote:

>
>
>
> On 8 May 2013 21:50, Stephen Woodbridge <woodbri at swoodbridge.com> wrote:
>
>> I tried to change my first with to just query and existing table, like:
>>
>> with lines as (
>>         select id as gid, * from bdaways
>>         /*
>>
>>         select 1 as gid, 'MULTILINESTRING((0 1,2 1))'::geometry as geom
>>         union all
>>         select 2 as gid, 'MULTILINESTRING((1 0,1 2))'::geometry as geom
>>         union all
>>         select 3 as gid, 'LINESTRING(1 1.5,2 2)'::geometry as geom
>>         */
>> ),
>>
>> but when I run it I get:
>>
>> ERROR:  line_locate_point: 2st arg isnt a point
>>
>> ********** Error **********
>>
>> ERROR: line_locate_point: 2st arg isnt a point
>> SQL state: XX000
>>
>> This looks like st_intersection(l1.geom, l2.geom) is not returning a
>> point.
>
>
> Yes it may return an empty geometryCollection. It could be filtered in the
> cut_locations cte by adding a where clause geometryType = 'POINT'.
>
>
>>
>> On 5/8/2013 3:32 PM, Nicolas Ribot wrote:
>>
>>> Hi,
>>>
>>> Glad it helps.
>>>
>>> CTE are supported since 8.4 according to the doc:
>>> http://www.postgresql.org/**docs/8.4/static/queries-with.**html<http://www.postgresql.org/docs/8.4/static/queries-with.html>
>>>
>>
>> ok, cool, will give that a try also.
>>
>>
>>  Otherwise, the CTE's can be replaced by subqueries.
>>> The stored procedure may also help by creating temp tables instead of
>>> chaining subqueries, though I don't know if it will run faster.
>>>
>>
>> One thing I noticed is that CTE's can not be indexed, so I might get
>> better performance is I create tables and index them based on the needs of
>> successive queries in a stored procedure. I'll need to play with this a bit
>> to figure out what works best.
>>
>>
>>  Using topology should be pretty simple in your case: build a new
>>> topology based on the lines table, then query the topology.edge table,
>>> keeping initial line gid. It may be worth trying it.
>>>
>>
>> I need to find the how to for working with topologies. I have seem some
>> in the past and just need to give it a try now that I have a real use case
>> for it.
>>
>> -Steve
>>
>>  Nicolas
>>>
>>>
>>> On 8 May 2013 21:01, Stephen Woodbridge <woodbri at swoodbridge.com
>>> <mailto:woodbri at swoodbridge.**com <woodbri at swoodbridge.com>>> wrote:
>>>
>>>     Hi Nicolas,
>>>
>>>     Wow! thank you for an excellent example. This is very help. Since I
>>>     want this to work on pg 8.4+, I'll convert this into a stored
>>>     procedure since I can't use CTE subqueries.
>>>
>>>     Now I have some work cut out to do on this. :)
>>>
>>>     Thanks again,
>>>        -Steve
>>>
>>>
>>>     On 5/8/2013 2:39 PM, Nicolas Ribot wrote:
>>>
>>>         Hi Stephen,
>>>
>>>         Building a topology would definitively help in this situation,
>>>         though it
>>>         may take some time on very large dataset I guess.
>>>         If you plan to use some topological functions on the dataset in
>>>         addition
>>>         with pgRouting functions, it may be worth the effort.
>>>
>>>         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 ?
>>>
>>>         Looking at this subject recently (cutting lines by points, cf.
>>>         http://trac.osgeo.org/postgis/**__wiki/__**
>>> UsersWikiSplitPolygonWithPoint**__s<http://trac.osgeo.org/postgis/__wiki/__UsersWikiSplitPolygonWithPoint__s>
>>>         <http://trac.osgeo.org/**postgis/wiki/**
>>> UsersWikiSplitPolygonWithPoint**s<http://trac.osgeo.org/postgis/wiki/UsersWikiSplitPolygonWithPoints>
>>> >)
>>>
>>>         I
>>>         found that linear referencing functions can help in such a case.
>>>
>>>         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.
>>>         It appears to be very efficient, using st_dwithin to trigger
>>> spatial
>>>         index, then joining on the lines primary keys, which should be
>>> fast.
>>>
>>>         In your usecase, intersection nodes between lines have to be
>>>         identified
>>>         before their locations can be computed.
>>>
>>>         Concerning the tolerance, I'm pretty sure snapping the input
>>>         dataset to
>>>         a grid would help to run a precise st_intersection between lines.
>>>
>>>         Based on the linestring sample data, here is  the query using
>>> linear
>>>         referencing. It uses CTE subqueries to identify each step:
>>>
>>>         with lines as (
>>>                   select 1 as gid, 'MULTILINESTRING((0 1,2
>>>         1))'::geometry as geom
>>>                   union all
>>>                   select 2 as gid, 'MULTILINESTRING((1 0,1
>>>         2))'::geometry as geom
>>>                   union all
>>>                   select 3 as gid, 'LINESTRING(1 1.5,2 2)'::geometry as
>>> geom
>>>         ),
>>>         -- multilinestrings are dumped into simple objects
>>>         -- if multilinestrings have several parts, one should generate a
>>>         unique
>>>         id based
>>>         -- on their gid and path into the collection.
>>>         dumped_lines as (
>>>         select gid, (st_dump(l.geom)).geom
>>>         from lines l
>>>         ),
>>>         -- This query computes the locations, for each input line, of the
>>>         intersection points with other lines.
>>>         -- this will be used to cut lines based on these locations.
>>>         -- to be able to cut lines from their beginning to their end, we
>>>         generate the 0 and 1 location index
>>>         cut_locations as (
>>>         select l1.gid as lgid, st_line_locate_point(l1.geom,
>>>         st_intersection(l1.geom, l2.geom)) as locus
>>>         from dumped_lines l1 join dumped_lines l2 on (st_dwithin(l1.geom,
>>>         l2.geom, 0.01))
>>>         where l1.gid <> l2.gid
>>>         -- then generates start and end locus for each line, to be able
>>>         to cut them
>>>         UNION ALL
>>>         select l.gid as lgid, 0 as locus
>>>         from dumped_lines l
>>>         UNION ALL
>>>         select l.gid as lgid, 1 as locus
>>>         from dumped_lines l
>>>         order by lgid, locus
>>>         ),
>>>         -- This query generates a row_number index column for each input
>>>         line
>>>         and intersection point.
>>>         -- This index will be used to self-join the table to cut a line
>>>         between
>>>         two consecutive locations
>>>         -- (idx, idx+1) pairs.
>>>         -- window function is used to generate the index inside each
>>>         line partition
>>>         loc_with_idx as (
>>>         select lgid, locus, row_number() over (partition by lgid order
>>>         by locus)
>>>         as idx
>>>         from cut_locations
>>>         )
>>>         -- finally, each original line is cut with consecutive locations
>>>         using
>>>         linear referencing function.
>>>         -- a filtering is done to eliminate points produced when lines
>>>         connect
>>>         at their ends
>>>         select l.gid, loc1.idx as sub_id, st_line_substring(l.geom,
>>>         loc1.locus,
>>>         loc2.locus) as geom ,
>>>         st_geometryType(st_line___**substring(l.geom, loc1.locus,
>>>
>>>         loc2.locus)) as type
>>>         from loc_with_idx loc1 join loc_with_idx loc2 using (lgid) join
>>>         dumped_lines l on (l.gid = loc1.lgid)
>>>         where loc2.idx = loc1.idx+1
>>>         -- filter out point geometries occuring if intersection point is
>>> at
>>>         line's start or end point.
>>>         -- there must be a faster way to filter out theses geometries.
>>>         and st_geometryType(st_line___**substring(l.geom, loc1.locus,
>>>
>>>         loc2.locus))
>>>         <> 'ST_Point';
>>>
>>>
>>>         A new unique ID key can be computed based on line gid and subgid
>>>         generated by the query.
>>>         Initial line attributes can be moved to the new segments using
>>>         the line
>>>         gid key.
>>>
>>>         Nicolas
>>>
>>>
>>>         On 8 May 2013 16:27, Stephen Woodbridge <woodbri at swoodbridge.com
>>>         <mailto:woodbri at swoodbridge.**com <woodbri at swoodbridge.com>>
>>>         <mailto:woodbri at swoodbridge.__**com
>>>
>>>         <mailto:woodbri at swoodbridge.**com <woodbri at swoodbridge.com>>>>
>>> wrote:
>>>
>>>              Hi all,
>>>
>>>              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?"
>>>
>>>              It seems to me that this would be a good function to create
>>> in
>>>              either postgis or pgrouting.
>>>
>>>              THE PROBLEM:
>>>
>>>              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.
>>>
>>>              THE RESULTS:
>>>
>>>              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.
>>>
>>>              POTENTIAL SOLUTIONS:
>>>
>>>              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.
>>>
>>>              2. Another way seems to be using a query like:
>>>
>>>              select (st_dump(bar.the_geom)).* from (
>>>                   select st_union(foo.the_geom) as the_geom from mytable
>>> foo
>>>              ) as bar;
>>>
>>>              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.
>>>
>>>              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.
>>>
>>>              3. Is there a good recipe for doing this somewhere that I
>>>         have not
>>>              found? or other better approaches to this problem?
>>>
>>>              What would be the best way to add tolerance to the problem?
>>>         using
>>>              snap to grid?
>>>
>>>              Thoughts on how to do this efficiently?
>>>
>>>              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.
>>>
>>>              Thanks,
>>>                 -Steve
>>>
>>>
>>>              -- Example to demonstrate st_union above
>>>              select st_astext((st_dump(bar.the____**_geom)).geom) from (
>>>
>>>
>>>                   select st_union(foo.the_geom) as the_geom from (
>>>                       select 'MULTILINESTRING((0 1,2 1))'::geometry as
>>>         the_geom
>>>                       union all
>>>                       select 'MULTILINESTRING((1 0,1 2))'::geometry as
>>>         the_geom
>>>                       union all
>>>                       select 'LINESTRING(1 1.5,2 2)'::geometry as
>>> the_geom
>>>                       ) as foo
>>>                   ) as bar;
>>>
>>>              "LINESTRING(1 1.5,2 2)"
>>>              "LINESTRING(1 0,1 1)"
>>>              "LINESTRING(1 1,1 1.5)"
>>>              "LINESTRING(1 1.5,1 2)"
>>>              "LINESTRING(0 1,1 1)"
>>>              "LINESTRING(1 1,2 1)"
>>>              ______________________________**_____________________
>>>
>>>              postgis-users mailing list
>>>         postgis-users at lists.osgeo.org
>>>         <mailto:postgis-users at lists.**osgeo.org<postgis-users at lists.osgeo.org>
>>> >
>>>         <mailto:postgis-users at lists.__**osgeo.org <http://osgeo.org>
>>>         <mailto:postgis-users at lists.**osgeo.org<postgis-users at lists.osgeo.org>
>>> >>
>>>         http://lists.osgeo.org/cgi-___**_bin/mailman/listinfo/postgis-**
>>> ____users<http://lists.osgeo.org/cgi-____bin/mailman/listinfo/postgis-____users>
>>>         <http://lists.osgeo.org/cgi-__**bin/mailman/listinfo/postgis-_**
>>> _users<http://lists.osgeo.org/cgi-__bin/mailman/listinfo/postgis-__users>
>>> >
>>>
>>>         <http://lists.osgeo.org/cgi-__**bin/mailman/listinfo/postgis-_**
>>> _users<http://lists.osgeo.org/cgi-__bin/mailman/listinfo/postgis-__users>
>>>         <http://lists.osgeo.org/cgi-**bin/mailman/listinfo/postgis-**
>>> users <http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>         ______________________________**___________________
>>>         postgis-users mailing list
>>>         postgis-users at lists.osgeo.org <mailto:postgis-users at lists.**
>>> osgeo.org <postgis-users at lists.osgeo.org>>
>>>         http://lists.osgeo.org/cgi-__**bin/mailman/listinfo/postgis-_**
>>> _users<http://lists.osgeo.org/cgi-__bin/mailman/listinfo/postgis-__users><
>>> http://lists.osgeo.org/cgi-**bin/mailman/listinfo/postgis-**users<http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users>
>>> >
>>>
>>>
>>>     ______________________________**___________________
>>>     postgis-users mailing list
>>>     postgis-users at lists.osgeo.org <mailto:postgis-users at lists.**
>>> osgeo.org <postgis-users at lists.osgeo.org>>
>>>     http://lists.osgeo.org/cgi-__**bin/mailman/listinfo/postgis-_**
>>> _users<http://lists.osgeo.org/cgi-__bin/mailman/listinfo/postgis-__users>
>>>     <http://lists.osgeo.org/cgi-**bin/mailman/listinfo/postgis-**users<http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users>
>>> >
>>>
>>>
>>>
>>>
>>> ______________________________**_________________
>>> postgis-users mailing list
>>> postgis-users at lists.osgeo.org
>>> http://lists.osgeo.org/cgi-**bin/mailman/listinfo/postgis-**users<http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users>
>>>
>>>
>> ______________________________**_________________
>> postgis-users mailing list
>> postgis-users at lists.osgeo.org
>> http://lists.osgeo.org/cgi-**bin/mailman/listinfo/postgis-**users<http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users>
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/postgis-users/attachments/20130508/730cdde9/attachment.html>


More information about the postgis-users mailing list