[postgis-users] Need a method for "noding" a street network
Stephen Woodbridge
woodbri at swoodbridge.com
Thu May 9 11:53:29 PDT 2013
On 5/9/2013 1:27 PM, Nicolas Ribot wrote:
> Steve,
>
> I modified your script to include the final step gathering all needed
> input lines in the results table.
> I also corrected the linestring ends insertion queries:
>
> ...quote_ident(n_geom) || ') <> ''LINESTRING'' '; to become:
> ...quote_ident(n_geom) || ') = ''LINESTRING'' ';
>
> which inserted null values for geom column as st_startpoint(geom) is
> null for any geometry <> LINESTRING.
Duh, my bad on that and you even mentioned it. I'm pleading a nasty
sinus cold which has scramble my brains.
I would like to include this in pgRouting, if that is OK with you. I
would also like to credit you for your effort to develop this tool, so
if you would like to add a comment and copyright that would be great.
> It seems to work well on bdaways dataset and costa rica highways.
>
> I will add some snapping operations as suggested by Sandro and try to
> see if all segments are noded.
> (Topology still creating after 4270 ms)
The snapping would be nice. I'm wondering if Sandro could possible use a
variant of this to speed up the topology creation?
Thanks,
-Steve
> Nicolas
>
>
> On 9 May 2013 18:44, Nicolas Ribot <nicolas.ribot at gmail.com
> <mailto:nicolas.ribot at gmail.com>> wrote:
>
> Thanks Sandro, I will try to add these snapping steps in the process.
>
> The following code:
>
> select dropTopology('topo1');
> SELECT topology.CreateTopology('topo1',4326, 0.000001);
> select st_createTopoGeo('topo1', st_collect(geom))
> from ways;
>
> Still runs after 1790 sec.
>
> Nicolas
>
>
> On 9 May 2013 18:02, Sandro Santilli <strk at keybit.net
> <mailto:strk at keybit.net>> wrote:
>
> Note that ST_Intersect may sligthly move the inputs while attempting
> to catch robustness issues, so it's recommended to snap the final
> set of intersection points to a grid and remove duplicates, then
> finally snap each input segment to the so-found points to ensure
> the are found on the lines they initially were computed from.
>
> Out of curiosity: how much does ST_CreateTopoGeo takes when
> passed the input as a collection ?
>
> --strk;
>
> On Thu, May 09, 2013 at 05:52:54PM +0200, Nicolas Ribot wrote:
> > On 9 May 2013 15:43, Stephen Woodbridge
> <woodbri at swoodbridge.com <mailto:woodbri at swoodbridge.com>> wrote:
> >
> > > On 5/9/2013 7:20 AM, Nicolas Ribot wrote:
> > >
> > >> (I spammed this thread a bit with image attachment....)
> > >>
> > >> Hi Steve,
> > >>
> > >> We the given dataset, my approach is indeed slow compared
> to st_union
> > >> approach (though precision for the st_dwithin clause must
> be adapted to
> > >> the current dataset. I took the following precision: 0.000001)
> > >>
> > >> The st_union method generates 18322 segments in 7318 ms,
> though the
> > >> final association between original lines and new segment
> is not done
> > >> here.
> > >>
> > >> With the query I gave, the st_dwithin part takes 11.7 sec
> on a recent
> > >> laptop machine (1.8 Ghz Intel Core I7, 1024 mb of ram for
> shared_buffer,
> > >> 512 for work_mem)...
> > >>
> > >> The complete query returns 17292 segments in 17956 ms.
> > >>
> > >> As the lines are almost already noded, it generates a lot of
> > >> intersection points coincident with one line ends.
> > >>
> > >> As you noted, intermediate temp tables may help here:
> > >>
> > >> I decomposed the query into intermediate steps and the
> performance is
> > >> about the same as with st_union :
> > >>
> > >> -- First creates temp table with intersection points
> > >> drop table if exists intergeom;
> > >> create temp table intergeom as
> > >> select l1.id <http://l1.id> as l1id, l2.id <http://l2.id>
> as l2id,
> > >> st_intersection(l1.geom, l2.geom) as geom
> > >> from bdaways l1 join bdaways l2 on (st_dwithin(l1.geom,
> l2.geom,
> > >> 0.000001))
> > >> where l1.id <http://l1.id> <> l2.id <http://l2.id> ;
> > >>
> > >> -- keeps only true intersection points
> > >> -- must handle the case where lines intersects at a
> linestring...
> > >>
> > >
> > > Would it make sense to take all the geometryType(geom) <>
> 'LINESTRING' and
> > > just add the end points to the intergeom table. I think
> this would then add
> > > break points at the ends of overlapping segments insure
> that they get
> > > divided at those points and it would also add extraneous
> points at the
> > > other end, but these would get filter out later. So add
> these two lines
> > > here?
> > >
> >
> > You meant, geometryType(geom) = 'LINESTRING ?
> > Yes it would make sense. I will try.
> >
> >
> > > insert intergeom (l1id, l2id, geom)
> > > values (l1id, l2id, st_startpoint(geom))
> > > where geometryType(geom) <> 'LINESTRING';
> > >
> > > insert intergeom (l1id, l2id, geom)
> > > values (l1id, l2id, st_endpoint(geom))
> > > where geometryType(geom) <> 'LINESTRING';
> > >
> > >
> > >
> > > delete from intergeom where geometryType(geom) <> 'POINT';
> > >>
> > >
> > > I think some OSM data from a small central american country
> might make a
> > > good test case. I have not tried any of these:
> > >
> > >
> http://downloads.cloudmade.**com/americas/central_america<http://downloads.cloudmade.com/americas/central_america>
> > >
> > > but I notice that if you click into a country like Belize
> or Costa Rica,
> > > there is a shapefile for the country near the bottom of the
> list or you can
> > > pick just a city for a smaller data set.
> > >
> > > I'll start packaging this into a stored procedure. This is
> an awesome job.
> > >
> >
> > Thanks.
> > I tested the costa rica highways dataset and found the
> following results:
> >
> > The query runs in about 28 seconds to process 29187 ways that
> are not at
> > all connected between them.
> > It generates 48415 segments.
> >
> > As an extra steps, all original ways that did not need to be
> cut by
> > intersection points have to be added to the final results:
> these lines
> > are already clean but were dismissed during the cutting process:
> >
> > insert into res (gid, sub_id, geom, type)
> > select ways.gid, 1 as sub_id, ways.geom, geometryType(ways.geom)
> > from ways
> > where not exists (
> > select res.gid from res where res.gid = ways.gid
> > );
> >
> > The result looks good, as all ways seems to be segmented (cf.
> > http://sd-38177.dedibox.fr/img1.png and
> http://sd-38177.dedibox.fr/img2.png:
> > black: initial lines, blue: new segments. Lines ends are
> shown as circles).
> >
> > A problem remains with non-simple or dirty lines as they
> self-intersect or
> > have several ends and so are not noded properly.
> > Though it should be simple to process them by identifying the
> problem. here
> > is an example of such a line: http://sd-38177.dedibox.fr/img3.png
> >
> > I go looking at the procedure ;)
> >
> > Nicolas
> _______________________________________________
> postgis-users mailing list
> postgis-users at lists.osgeo.org <mailto:postgis-users at lists.osgeo.org>
> 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
>
More information about the postgis-users
mailing list