<div dir="ltr"><div><div><div><div>As previously, <br></div>you need to separate the processing from the spatial filtering if you want to know where the slowness comes from (run the exact same querry, without the CASE).<br><br>What you numbers tall you is that your index is relatively useless : it reduce the data set, but there are still way too much work to do.<br>You may have few polygons with insane number of points,<br></div><div>that would take a lot of time.<br></div>(for instance you have a very detailed layer for sea, or forest etc etc).<br></div><div>If you do land cover, you may end up with big and detailed area.<br></div><br>It is usual to limit the number of points in polygon, for instance by cutting your big polygon into smaller pieces, using st_simplify,e tc.<br></div>There are many example on the web ("tiling" is one methods amongst others)<br><div><br>Another axis of improvement could be to analyze what you really need instead of <br><div><font face="Courier">(ST_Intersects(p.geom, n.the_geom) </font></div><div><font face="Courier"> AND NOT ST_Touches(p.geom, n.the_geom) );<br><br></font></div>Maybe you could do faster with custom st_relate (I doubt it).<br><br><br>Cheers,<br>Rémi-C<br><div><br></div></div></div><div class="gmail_extra"><br><div class="gmail_quote">2015-02-20 1:26 GMT+01:00 John Abraham <span dir="ltr"><<a href="mailto:jea@hbaspecto.com" target="_blank">jea@hbaspecto.com</a>></span>:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div style="word-wrap:break-word">Thanks for the hint, Rémi. I'll give it a try, to see if the st_intersection is slow (as compared to just st_intersects).<div><br></div><div>Here's my EXPLAIN ANALYZE for my 5 zone test:</div><div><br></div><div><div><font face="Courier">"Nested Loop (cost=1.56..80.91 rows=8 width=560) (actual time=3929258.738..19756123.479 rows=29 loops=1)"</font></div><span class=""><div><font face="Courier">" Join Filter: (_st_intersects(p.geom, n.the_geom) AND ((NOT (p.geom && n.the_geom)) OR (NOT _st_touches(p.geom, n.the_geom))))"</font></div></span><div><font face="Courier">" CTE sometaz"</font></div><div><font face="Courier">" -> Limit (cost=0.00..1.56 rows=5 width=9537) (actual time=0.433..0.576 rows=5 loops=1)"</font></div><div><font face="Courier">" -> Seq Scan on tazjan2 (cost=0.00..843.09 rows=2709 width=9537) (actual time=0.415..0.511 rows=5 loops=1)"</font></div><div><font face="Courier">" -> CTE Scan on sometaz n (cost=0.00..0.10 rows=5 width=80) (actual time=0.460..0.724 rows=5 loops=1)"</font></div><div><font face="Courier">" -> Index Scan using lancover_polygons_snap_geom_idx on lancover_polygons_snap p (cost=0.00..12.47 rows=5 width=480) (actual time=0.242..65.640 rows=16 loops=5)"</font></div><span class=""><div><font face="Courier">" Index Cond: (geom && n.the_geom)"</font></div></span><div><font face="Courier">"Total runtime: 19756143.935 ms"</font></div></div><div><br></div><div>Here it is explained in layspeak: <a href="http://explain.depesz.com/s/dbVi" target="_blank">http://explain.depesz.com/s/dbVi</a></div><div><br></div><div>So, am I reading this right, that it took 5.5 hours to st_intersects 5 (multi)polygons with 29 (multi)polygons, and find the 8 that intersected? That the index scan was fast (328ms), and returned only a few (16, or maybe 29) of the 998031 polygons, but that the actual intersection afterwards, between just a few polygons, was so very slow?</div><div><br></div><div>(The query was this:</div><div><br></div><div><div><font face="Courier">with sometaz as (select * from tazjan2 limit 5)</font></div><span class=""><div><font face="Courier">SELECT p.lc_class, n.taz</font></div><div><font face="Courier"> , CASE </font></div><div><font face="Courier"> WHEN ST_CoveredBy(p.geom, n.the_geom) </font></div><div><font face="Courier"> THEN p.geom </font></div><div><font face="Courier"> ELSE </font></div><div><font face="Courier"> ST_Multi(</font></div><div><font face="Courier"> ST_Intersection(p.geom,n.geom)</font></div><div><font face="Courier"> ) END AS geom </font></div><div><font face="Courier"> FROM lancover_polygons_snap AS p </font></div></span><div><font face="Courier"> INNER JOIN sometaz AS n </font></div><span class=""><div><font face="Courier"> ON (ST_Intersects(p.geom, n.the_geom) </font></div><div><font face="Courier"> AND NOT ST_Touches(p.geom, n.the_geom) );</font></div></span></div><div><br></div><div>)</div><div><br></div><div><div>
<span style="border-collapse:separate;border-spacing:0px"><div>--</div><div>John Abraham</div><div><a href="mailto:jea@hbaspecto.com" target="_blank">jea@hbaspecto.com</a></div><div>403-232-1060</div></span>
</div><div><div class="h5">
<br><div><div>On Feb 19, 2015, at 11:56 AM, Rémi Cura <<a href="mailto:remi.cura@gmail.com" target="_blank">remi.cura@gmail.com</a>> wrote:</div><br><blockquote type="cite"><div dir="ltr"><div><div><div><div><div>Hey,<br>you could try to not use CASE (so separate the spatial join from the processing, which is easy to parallelize (assuming you have more than one core )).<br></div><br>First generate the table with <br>----<br></div>CREATE TABLE psatial_mapping_between_lancover_and_taz AS <br><div><div><font face="Courier">SELECT row_number() over() as row_id, p.your_primary_key, n.your_primary_key<br></font></div><div><font face="Courier"> FROM lancover_polygons_snap AS p </font></div><div><font face="Courier"> INNER JOIN tazjan2 AS n </font></div><div><font face="Courier"> ON (ST_Intersects(p.geom, n.the_geom) </font></div><div><font face="Courier"> AND NOT ST_Touches(p.geom, n.the_geom) );</font></div>----<br><br></div><div>Then when you have your table, you can do your processing (put into a SQL function , it is cleaner),<br></div><div>and you parallelize on the<font face="Courier"> "row_id".<br><br></font></div>Basically you get min(row_id), ax(row_id), you spearate it into N parts, then execute the parts with your K processes.<br></div>(both bash and python have utilities for this).<br><br></div><div>Of course if it is actually the spatial join that is slow, you can also parallelise this.<br></div><div>(by cutting <font face="Courier">lancover</font>into N pieces for instance )<br></div>Cheers,<br></div>Rémi-C<br></div><div class="gmail_extra"><br><div class="gmail_quote">2015-02-19 19:45 GMT+01:00 John Abraham <span dir="ltr"><<a href="mailto:jea@hbaspecto.com" target="_blank">jea@hbaspecto.com</a>></span>:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div style="word-wrap:break-word">So I've was running this query for 866000 s (10 days) before I decided to kill it:<div><br></div><div><div><font face="Courier">create table taz_and_lancover_10_fast_2 as</font></div><div><font face="Courier">SELECT p.lc_class, n.taz</font></div><div><font face="Courier"> , CASE </font></div><div><font face="Courier"> WHEN ST_CoveredBy(p.geom, n.the_geom) </font></div><div><font face="Courier"> THEN p.geom </font></div><div><font face="Courier"> ELSE </font></div><div><font face="Courier"> ST_Multi(</font></div><div><font face="Courier"> ST_Intersection(p.geom,n.geom)</font></div><div><font face="Courier"> ) END AS geom </font></div><div><font face="Courier"> FROM lancover_polygons_snap AS p </font></div><div><font face="Courier"> INNER JOIN tazjan2 AS n </font></div><div><font face="Courier"> ON (ST_Intersects(p.geom, n.the_geom) </font></div><div><font face="Courier"> AND NOT ST_Touches(p.geom, n.the_geom) );</font></div><div><font face="Courier"><br></font></div><div>Explain shows this, it's using the spatial index:</div><div><br></div><div><div>"Nested Loop (cost=0.00..310492.35 rows=483105 width=9973)"</div><div>" Join Filter: (_st_intersects(p.geom, n.the_geom) AND ((NOT (p.geom && n.the_geom)) OR (NOT _st_touches(p.geom, n.the_geom))))"</div><div>" -> Seq Scan on tazjan2 n (cost=0.00..843.09 rows=2709 width=9493)"</div><div>" -> Index Scan using lancover_polygons_snap_geom_idx on lancover_polygons_snap p (cost=0.00..21.67 rows=5 width=480)"</div><div>" Index Cond: (geom && n.the_geom)"</div></div><div><br></div><div>There are 2709 rows in tazjan2 and 998031 rows in lancover_polygons_snap, so I appreciate that it's a bit of a large problem. But MapInfo was able to do it interactively in a few days and Geomedia was also able to do it in about a day. </div><div><br></div><div>Both MapInfo and Geomedia ran out of memory (8GB machines) until the problem was broken into two regions (North and South), but postgresql seems to be chewing on the problem using only 400MB. The interactive approach in MapInfo was to divide lancover_polygons_snap by lc_class, to further divide each region into about 10 subproblems. Perhaps subprobleming this is the way to go? Can't the query subproblem it based in the indices or would I have to do that manually?</div><div><br></div><div>One potential thing I've realized is that a few of the geometries in tazjan2 are multipolygons, not single polygons. But it's only a few. There are a few very large and complex polygons in lancover_polygons_snap, but again, most of the 998031 are simple small polygons, about half would be ST_CoveredBy the polygons in tazjan2 and most of the rest would only overlap two or three of the polygons in tazjan2.</div><div><br></div><div>I must be doing something wrong. Any hints?</div><div><br></div><div>I have</div><div> max_connections set to 100 (currently only 11 connections active)</div><div> work_mem was defaulting to 1MB, I just bumped it to 256MB (machine has 32GB but it's Postgresql9.1 (x86), i.e. 32bit version</div><div> shared_buffers was 32MB, I am trying 512MB, could go higher</div><div> effective_cache_size was defaulting to 128MB, I am trying 20GB.</div><div> random_page_cost was 4.0, this is VMWare virtual server. I am trying 2.0 </div><div><br></div><div>I'm trying an "explain analyze" with just 5 rows of tazjan2, a 5 x 998031 intersection problem with indices shouldn't take that long, should it? It's been running for a little over an hour. I'm wondering if the few complex polygons in lancover_polygons_snap are causing the slowness? Is there some fast way to divide complex polygons, e.g. apply a 1km grid over them? </div><div><br></div><div>PostGISfullversion(): "POSTGIS="2.0.3 r11128" GEOS="3.3.8-CAPI-1.7.8" PROJ="Rel. 4.8.0, 6 March 2012" GDAL="GDAL 1.9.2, released 2012/10/08" LIBXML="2.7.8" LIBJSON="UNKNOWN" (core procs from "2.0.3 r11132" need upgrade) TOPOLOGY (topology procs from "2.0.3 r11132" need upgrade) RASTER (raster procs from "2.0.3 r11132" need upgrade)"</div><div><br></div><div>
<span style="border-collapse:separate;border-spacing:0px"><div>--</div><div>John Abraham</div><div><br></div></span>
</div>
<br></div></div><br>_______________________________________________<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-bin/mailman/listinfo/postgis-users</a><br></blockquote></div><br></div>
_______________________________________________<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-bin/mailman/listinfo/postgis-users</a></blockquote></div><br></div></div></div></div><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/cgi-bin/mailman/listinfo/postgis-users" target="_blank">http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users</a><br></blockquote></div><br></div>