[prev in list] [next in list] [prev in thread] [next in thread]
List: postgis-users
Subject: Re: [postgis-users] split line at polygon edge
From: Rémi_Cura <remi.cura () gmail ! com>
Date: 2016-03-18 8:04:56
Message-ID: CAJvUf_tuzNHHGG2sn0LEcGbynxw5aJsNq1=rH7dTbh8rZp+qaA () mail ! gmail ! com
[Download RAW message or body]
[Attachment #2 (multipart/alternative)]
Hey, in theory st_intersection already uses index and st_intersects under
the hood,
I personally prefer to explicitly add it, as I found it more easy to
understand.
What you didn't take into account is that you may get multiline,
for instance if your polygon is a U, and the line cross it left to right,
the intersection result would be two lines (each in a vertical part of the
U).
If you don't mind having multilinestring :
--------------------------------------
SELECT
ST_Multi(ST_CollectionExtract(ST_Intersection(t.geom,d.geom),2))::geometry(multilinestring,4326)
AS geom
, t.trails_id, d.id AS trail_system
FROM temp_trails AS t, divisions AS d
WHERE ST_Intersects(t.geom,d.geom) = TRUE;
-------------------------------------
if you want only line, you have to break multilines into simple lines
-------------------------------------
SELECT row_number() over() as qgis_id, dmp.path AS line_id
, dmp.geom::geometry(linestring,4326) AS geom
, t.trails_id, d.id AS trail_system
FROM temp_trails AS t, divisions AS d
, ST_Dump(ST_CollectionExtract(ST_Intersection(t.geom,d.geom),2)) as dmp
WHERE ST_Intersects(t.geom,d.geom) = TRUE;
--------------------------------------
Cheers,
Rémi-C
2016-03-17 22:17 GMT+01:00 François Hugues <hugues.francois@irstea.fr>:
> Hi,
>
> Empty geometries are returned when there is no intersection and I think we
> forgot something obvious. When you want to intersect geometries you need to
> add WHERE ST_Intersects (a.geom,b.geom).
>
> Things should work better and faster.
>
> HugThanks Remi-C and Hugues for your suggestions, they got me what I
> needed!
>
> I first tried Remi-C's example, since I was curious about how it would turn
> out. It gave me an error mentioning that it could not convert
> GeometryCollection to LineString. This error brought me back to what Hugues
> mentioned. So I used ST_Summary() to verify the GeometryCollections, which
> appeared to be empty (0 elements), and mixed in I noticed the LineStrings,
> MultiLineStrings. Since the Collections seemed to be empty I opted to
> separate out the linestrings I using a function Hugues mentioned
> ST_GeometryType()
>
> Specifically I used:
> ST_GeometryType(geom) like '%Line%'
>
> to get both linestring and multilinestrings.
>
> In the end it took 2 statements, even though I knew someone much more
> proficient then myself could do it in one.
>
> My final statements where:
> #create table public.temp_trail_div1 as select st_intersection(t.geom,
> d.geom) as geom,t.trails_id,
> d.id as trail_system from public.temp_trails as t, public.divisions as d;
>
> #create table public.temp_trail_div_sep as select * from
> public.temp_trail_div1 where ST_GeometryType(geom) like '%Line%';
>
> This seems to have done the trick, for now. Could someone enlighten me on
> how that might be done in one statement?
>
> Thanks again,
> Garret
>
>
> On Thu, Mar 17, 2016 at 5:45 AM, Rémi Cura <remi.cura@gmail.com> wrote:
>
> > Hey,
> > two things :
> > recent version of QGIS are boringly strict about geometry type,
> > so if you want to be able to add the corresponding postgis layer to qgis,
> > you may have to explicitely cast the result.
> > QGIS also require a unique identifier per row,
> > which you can fabricate with row_number() for instance
> >
> > ----------------------------------------------------------------
> > CREATE TABLE my_table AS
> > SELECT row_number() over() AS qgis_unique_id,
> > st_intersection(t.geom, d.geom)::geometry(linestring,4326) AS geom
> > ,t.trails_id, d.id
> > FROM public.temp_trails as t, public.polys as d;
> > ----------------------------------------------------------------
> >
> > Cheers,
> > Rémi-C
> >
> > 2016-03-17 8:15 GMT+01:00 François Hugues <hugues.francois@irstea.fr>:
> >
> >> Hello,
> >>
> >>
> >>
> >> Dis you take a look at the query result ? I think you should first try
> >> to see what is the type of geometry returned using ST_GeometryType().
> You
> >> may have some geometrycollections and I'm not sure QGis can handle it.
> In
> >> this case you could extract lines using ST_CollectionExtract().
> >>
> >>
> >>
> >> To achieve what you want to do, you'll be able to compare your original
> >> lines table with the result of your query using ST_Difference().
> >>
> >>
> >>
> >> Regards,
> >>
> >>
> >>
> >> Hugues.
> >>
> >>
> >>
> >> *De :* postgis-users [mailto:postgis-users-bounces@lists.osgeo.org] *De
> >> la part de* Garret W
> >> *Envoyé :* jeudi 17 mars 2016 04:11
> >> *À :* postgis-users@lists.osgeo.org
> >> *Objet :* [postgis-users] split line at polygon edge
> >>
> >>
> >>
> >> Hi Ive been looking for a way to take several hundred lines and split
> >> them where they intersect a polygon while also giving them the ID of the
> >> polygon they fall in. Ive seen many posts on splitting polygons. But its
> >> been difficult for me to adapt those examples.
> >>
> >> Ive been able to get an output from this:
> >>
> >> select st_intersection(t.geom, d.geom),t.trails_id, d.id
> >> from public.temp_trails as t, public.polys as d;
> >>
> >> Its giving me the line and IDs that I wanted but the geom is unreadable
> >> for some reason by QGIS.
> >>
> >> 99.9% of the lines fall within a polygon. Id like to still hang on to
> >> those few lines that arent contained in a polygon. They should just be
> >> split with no ID added
> >>
> >> Im using; postgis 2.2, postgresql 9.5
> >>
> >> Thank you
> >> Garret
> >>
> >> _______________________________________________
> >> postgis-users mailing list
> >> postgis-users@lists.osgeo.org
> >> http://lists.osgeo.org/mailman/listinfo/postgis-users
> >>
> >
> >
> > _______________________________________________
> > postgis-users mailing list
> > postgis-users@lists.osgeo.org
> > http://lists.osgeo.org/mailman/listinfo/postgis-users
> >
> _______________________________________________
> postgis-users mailing list
> postgis-users@lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/postgis-users
>
[Attachment #5 (text/html)]
<div dir="ltr"><div class="gmail_default" \
style="font-family:monospace,monospace"><div \
style="font-family:arial,sans-serif;font-size:12.8px">Hey, in theory st_intersection \
already uses index and st_intersects under the hood, </div><div \
style="font-family:arial,sans-serif;font-size:12.8px">I personally prefer to \
explicitly add it, as I found it more easy to understand.</div><div \
style="font-family:arial,sans-serif;font-size:12.8px"><br></div><div \
style="font-family:arial,sans-serif;font-size:12.8px">What you didn't take into \
account is that you may get multiline,</div><div \
style="font-family:arial,sans-serif;font-size:12.8px">for instance if your polygon is \
a U, and the line cross it left to right, </div><div \
style="font-family:arial,sans-serif;font-size:12.8px">the intersection result would \
be two lines (each in a vertical part of the U).</div><div \
style="font-family:arial,sans-serif;font-size:12.8px">If you don't mind having \
multilinestring :</div><div \
style="font-family:arial,sans-serif;font-size:12.8px"><br></div><div \
style="font-family:arial,sans-serif;font-size:12.8px">--------------------------------------</div><div \
style="font-family:arial,sans-serif;font-size:12.8px">SELECT \
ST_Multi(ST_CollectionExtract(ST_Intersection(t.geom,d.geom),2))::geometry(multilinestring,4326) \
AS geom</div><div style="font-family:arial,sans-serif;font-size:12.8px"> , \
t.trails_id, <a href="http://d.id">d.id</a> AS trail_system</div><div \
style="font-family:arial,sans-serif;font-size:12.8px">FROM temp_trails AS t, \
divisions AS d</div><div style="font-family:arial,sans-serif;font-size:12.8px">WHERE \
ST_Intersects(<span style="font-size:12.8px">t.geom,d.geom</span><span \
style="font-size:12.8px">) = TRUE; </span></div><div \
style="font-family:arial,sans-serif;font-size:12.8px">-------------------------------------</div><div \
style="font-family:arial,sans-serif;font-size:12.8px">if you want only line, you have \
to break multilines into simple lines</div><div \
style="font-family:arial,sans-serif;font-size:12.8px">-------------------------------------</div><div \
style="font-family:arial,sans-serif;font-size:12.8px"><div \
style="font-size:12.8px">SELECT row_number() over() as qgis_id, dmp.path AS \
line_id</div><div style="font-size:12.8px"> , \
dmp.geom::geometry(linestring,4326) AS geom</div><div style="font-size:12.8px"> \
, t.trails_id, <a href="http://d.id">d.id</a> AS trail_system</div><div \
style="font-size:12.8px">FROM temp_trails AS t, divisions AS d</div><div \
style="font-size:12.8px"> , <span \
style="font-size:12.8px">ST_Dump(ST_CollectionExtract(ST_Intersection(t.geom,d.geom),2)) \
as dmp</span></div><div style="font-size:12.8px"><span style="font-size:12.8px">WHERE \
ST_Intersects(</span><span style="font-size:12.8px">t.geom,d.geom</span><span \
style="font-size:12.8px">) = TRUE</span><span style="font-size:12.8px">; \
</span></div></div><div \
style="font-family:arial,sans-serif;font-size:12.8px">--------------------------------------</div><div \
style="font-family:arial,sans-serif;font-size:12.8px"><br></div><div \
style="font-family:arial,sans-serif;font-size:12.8px">Cheers,</div><div \
style="font-family:arial,sans-serif;font-size:12.8px">Rémi-C</div><div \
style="font-family:arial,sans-serif;font-size:12.8px"><br></div><div \
style="font-family:arial,sans-serif;font-size:12.8px"><br></div></div></div><div \
class="gmail_extra"><br><div class="gmail_quote">2016-03-17 22:17 GMT+01:00 François \
Hugues <span dir="ltr"><<a href="mailto:hugues.francois@irstea.fr" \
target="_blank">hugues.francois@irstea.fr</a>></span>:<br><blockquote \
class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc \
solid;padding-left:1ex">Hi,<br> <br>
Empty geometries are returned when there is no intersection and I think we forgot \
something obvious. When you want to intersect geometries you need to add WHERE \
ST_Intersects (a.geom,b.geom).<br> <br>
Things should work better and faster.<br>
<br>
HugThanks Remi-C and Hugues for your suggestions, they got me what I needed!<br>
<div><div class="h5"><br>
I first tried Remi-C's example, since I was curious about how it would turn<br>
out. It gave me an error mentioning that it could not convert<br>
GeometryCollection to LineString. This error brought me back to what Hugues<br>
mentioned. So I used ST_Summary() to verify the GeometryCollections, which<br>
appeared to be empty (0 elements), and mixed in I noticed the LineStrings,<br>
MultiLineStrings. Since the Collections seemed to be empty I opted to<br>
separate out the linestrings I using a function Hugues mentioned<br>
ST_GeometryType()<br>
<br>
Specifically I used:<br>
ST_GeometryType(geom) like '%Line%'<br>
<br>
to get both linestring and multilinestrings.<br>
<br>
In the end it took 2 statements, even though I knew someone much more<br>
proficient then myself could do it in one.<br>
<br>
My final statements where:<br>
#create table public.temp_trail_div1 as select st_intersection(t.geom,<br>
d.geom) as geom,t.trails_id,<br>
<a href="http://d.id" rel="noreferrer" target="_blank">d.id</a> as trail_system from \
public.temp_trails as t, public.divisions as d;<br> <br>
#create table public.temp_trail_div_sep as select * from<br>
public.temp_trail_div1 where ST_GeometryType(geom) like '%Line%';<br>
<br>
This seems to have done the trick, for now. Could someone enlighten me on<br>
how that might be done in one statement?<br>
<br>
Thanks again,<br>
Garret<br>
<br>
<br>
On Thu, Mar 17, 2016 at 5:45 AM, Rémi Cura <<a \
href="mailto:remi.cura@gmail.com">remi.cura@gmail.com</a>> wrote:<br> <br>
> Hey,<br>
> two things :<br>
> recent version of QGIS are boringly strict about geometry type,<br>
> so if you want to be able to add the corresponding postgis layer to qgis,<br>
> you may have to explicitely cast the result.<br>
> QGIS also require a unique identifier per row,<br>
> which you can fabricate with row_number() for instance<br>
> <br>
> ----------------------------------------------------------------<br>
> CREATE TABLE my_table AS<br>
> SELECT row_number() over() AS qgis_unique_id,<br>
> st_intersection(t.geom, d.geom)::geometry(linestring,4326) AS geom<br>
> ,t.trails_id, <a href="http://d.id" rel="noreferrer" \
target="_blank">d.id</a><br> > FROM public.temp_trails as t, public.polys as \
d;<br> > ----------------------------------------------------------------<br>
><br>
> Cheers,<br>
> Rémi-C<br>
><br>
> 2016-03-17 8:15 GMT+01:00 François Hugues <<a \
href="mailto:hugues.francois@irstea.fr">hugues.francois@irstea.fr</a>>:<br> \
><br> >> Hello,<br>
>><br>
>><br>
>><br>
>> Dis you take a look at the query result ? I think you should first try<br>
>> to see what is the type of geometry returned using ST_GeometryType(). \
You<br> >> may have some geometrycollections and I'm not sure QGis can handle \
it. In<br> >> this case you could extract lines using \
ST_CollectionExtract().<br> >><br>
>><br>
>><br>
>> To achieve what you want to do, you'll be able to compare your original<br>
>> lines table with the result of your query using ST_Difference().<br>
>><br>
>><br>
>><br>
>> Regards,<br>
>><br>
>><br>
>><br>
>> Hugues.<br>
>><br>
>><br>
>><br>
</div></div>>> *De :* postgis-users [mailto:<a \
href="mailto:postgis-users-bounces@lists.osgeo.org">postgis-users-bounces@lists.osgeo.org</a>] \
*De<br> >> la part de* Garret W<br>
>> *Envoyé :* jeudi 17 mars 2016 04:11<br>
>> *À :* <a href="mailto:postgis-users@lists.osgeo.org">postgis-users@lists.osgeo.org</a><br>
>> *Objet :* [postgis-users] split line at polygon edge<br>
<div class="HOEnZb"><div class="h5">>><br>
>><br>
>><br>
>> Hi Ive been looking for a way to take several hundred lines and split<br>
>> them where they intersect a polygon while also giving them the ID of the<br>
>> polygon they fall in. Ive seen many posts on splitting polygons. But its<br>
>> been difficult for me to adapt those examples.<br>
>><br>
>> Ive been able to get an output from this:<br>
>><br>
>> select st_intersection(t.geom, d.geom),t.trails_id, <a href="http://d.id" \
rel="noreferrer" target="_blank">d.id</a><br> >> from public.temp_trails as \
t, public.polys as d;<br> >><br>
>> Its giving me the line and IDs that I wanted but the geom is unreadable<br>
>> for some reason by QGIS.<br>
>><br>
>> 99.9% of the lines fall within a polygon. Id like to still hang on to<br>
>> those few lines that arent contained in a polygon. They should just be<br>
>> split with no ID added<br>
>><br>
>> Im using; postgis 2.2, postgresql 9.5<br>
>><br>
>> Thank you<br>
>> Garret<br>
>><br>
>> _______________________________________________<br>
>> postgis-users mailing list<br>
>> <a href="mailto:postgis-users@lists.osgeo.org">postgis-users@lists.osgeo.org</a><br>
>> <a href="http://lists.osgeo.org/mailman/listinfo/postgis-users" \
rel="noreferrer" target="_blank">http://lists.osgeo.org/mailman/listinfo/postgis-users</a><br>
>><br>
><br>
><br>
> _______________________________________________<br>
> postgis-users mailing list<br>
> <a href="mailto:postgis-users@lists.osgeo.org">postgis-users@lists.osgeo.org</a><br>
> <a href="http://lists.osgeo.org/mailman/listinfo/postgis-users" \
rel="noreferrer" target="_blank">http://lists.osgeo.org/mailman/listinfo/postgis-users</a><br>
><br>
_______________________________________________<br>
postgis-users mailing list<br>
<a href="mailto:postgis-users@lists.osgeo.org">postgis-users@lists.osgeo.org</a><br>
<a href="http://lists.osgeo.org/mailman/listinfo/postgis-users" rel="noreferrer" \
target="_blank">http://lists.osgeo.org/mailman/listinfo/postgis-users</a></div></div></blockquote></div><br></div>
[Attachment #6 (text/plain)]
_______________________________________________
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/postgis-users
[prev in list] [next in list] [prev in thread] [next in thread]
Configure |
About |
News |
Add a list |
Sponsored by KoreLogic