Re: [postgis-users] postgis-users Digest, Vol 203, Issue 11
Hi, I tried adding a gist index to geom::geography and oddly enough it made the query a lot slower... Not really sure why, but the best option in my case still seems to be to first doing a bounding box reduction in lon-lat space and then apply ST_DWithin after casting to geography... Cheers, David On 15/01/2019 21:00, postgis-users-requ...@lists.osgeo.org wrote: From: Devdatta Tengshe If you are running a query on the Geography, you should create the index on Geography, using something like this: Create Index g_geog_idx on mySchema.myTable using GIST(geom::geography); Using this increased the speed of my query 1000 times. Source: https://twitter.com/postgis/status/675001071505383424 Regards, Devdatta On Tue, Jan 15, 2019 at 4:07 PM David M. Kaplan wrote: Hi, For a while it has not been clear in my head if and when GIST indexes are used when doing geography based calculations. If the data is stored in a table with a geometry column that has a GIST index on it, will that index be used if one does something like ST_DWithin(a.geom::geography,b.geom,1000)? If not, then if the data was stored in geography format instead of geometry, would this make the index more useful? To give a specific context, I have a table of point geometries with SRID 4326 and a GIST index, and I want to find the number of points in that table that are within a certain distance of each other point in that table, but I am not sure what is the most efficient way to do so. To do this, I have a query of the form: SELECT a.gid, a.geom, count(*) AS num_points FROM mytable a JOIN mytable b ON a.gid<>b.gid AND ST_DWithin(a.geom::geography,b.geom,1000,FALSE) GROUP BY a.gid,a.geom ; This is quite slow, so I presume the GIST index is not being used (and EXPLAIN didn't show anything that made it clear that it was being used though bounding box comparisons are included in the join filter). Also adding a lonlat bounding box comparison does speed the calculation up significantly: SELECT a.gid, a.geom, count(*) AS num_points FROM mytable a JOIN mytable b ON a.gid<>b.gid AND a.geom <#> b.geom < 0.02 -- in region where 0.02 degrees is > 1000 m AND ST_DWithin(a.geom::geography,b.geom,1000,FALSE) GROUP BY a.gid,a.geom ; Is this the best approach or am I missing something? Would using a LATERAL join improve things? Thanks, David -- ********** David M. Kaplan Charge de Recherche 1 Institut de Recherche pour le Developpement (IRD) UMR MARBEC (IRD/Ifremer/CNRS/UMII) av. Jean Monnet CS 30171 34203 Sete cedex France Email: david.kap...@ird.fr Phone: +33 (0)4 99 57 32 25 Fax: +33 (0)4 99 57 32 95 http://www.umr-marbec.fr/kaplan-david.html http://www.davidmkaplan.fr/ ** ___ postgis-users mailing list postgis-users@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/postgis-users
[postgis-users] confusion regarding use of indexes when doing geography calculations
Hi, For a while it has not been clear in my head if and when GIST indexes are used when doing geography based calculations. If the data is stored in a table with a geometry column that has a GIST index on it, will that index be used if one does something like ST_DWithin(a.geom::geography,b.geom,1000)? If not, then if the data was stored in geography format instead of geometry, would this make the index more useful? To give a specific context, I have a table of point geometries with SRID 4326 and a GIST index, and I want to find the number of points in that table that are within a certain distance of each other point in that table, but I am not sure what is the most efficient way to do so. To do this, I have a query of the form: SELECT a.gid, a.geom, count(*) AS num_points FROM mytable a JOIN mytable b ON a.gid<>b.gid AND ST_DWithin(a.geom::geography,b.geom,1000,FALSE) GROUP BY a.gid,a.geom ; This is quite slow, so I presume the GIST index is not being used (and EXPLAIN didn't show anything that made it clear that it was being used though bounding box comparisons are included in the join filter). Also adding a lonlat bounding box comparison does speed the calculation up significantly: SELECT a.gid, a.geom, count(*) AS num_points FROM mytable a JOIN mytable b ON a.gid<>b.gid AND a.geom <#> b.geom < 0.02 -- in region where 0.02 degrees is > 1000 m AND ST_DWithin(a.geom::geography,b.geom,1000,FALSE) GROUP BY a.gid,a.geom ; Is this the best approach or am I missing something? Would using a LATERAL join improve things? Thanks, David -- ********** David M. Kaplan Charge de Recherche 1 Institut de Recherche pour le Developpement (IRD) UMR MARBEC (IRD/Ifremer/CNRS/UMII) av. Jean Monnet CS 30171 34203 Sete cedex France Email: david.kap...@ird.fr Phone: +33 (0)4 99 57 32 25 Fax: +33 (0)4 99 57 32 95 http://www.umr-marbec.fr/kaplan-david.html http://www.davidmkaplan.fr/ ** ___ postgis-users mailing list postgis-users@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/postgis-users
Re: [postgis-users] very strange bug generating incorrect ST_Distance calculations
> >/Glad to hear that you seem to have found the issue. If there are any other tests I can do to try to diagnose the problem, let me know. / Can you open an issue in trac (https://trac.osgeo.org/postgis/) with the details, please? If you don't have an account I can do it myself. Done! Ticket #4290. ___ postgis-users mailing list postgis-users@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/postgis-users
Re: [postgis-users] very strange bug generating incorrect ST_Distance calculations
Glad to hear that you seem to have found the issue. If there are any other tests I can do to try to diagnose the problem, let me know. Cheers, David : You’ve found a bug, it’s in PostGIS (geodetic code is in postgis native), and the extra detail that you can get it to manifest with MULTIPOLYGON but the same rings in POLYGON don’t manifest is an excellent extra detail to help in fixing it. On 10/01/2019 12:26, David M. Kaplan wrote: Comparing the results of _ST_DistanceTree and _ST_DistanceUnCached is a good way to tease out whether there are differences between the implementations without having to fight the caching machinery along the way. P. Thanks for the responses. I tried with _ST_DistanceTree and _ST_DistanceUnCached and as predicted, tree gives 0 and uncached gives the correct answer. What is the appropriate work around to avoid this issue? Who should I report this problem to? The maintainers of libgeos? Also, how should I interpret that this error only occurs when one of the geometries is a multipolygon? In response to Raúl, I found that both queries give 0 for postgis 2.5.0 r16836, but only the multiline query gives zero for POSTGIS="2.4.3 r16312"... Cheers, David On Thu, Jan 10, 2019 at 7:00 AM Raúl Marín Rodríguez https://lists.osgeo.org/mailman/listinfo/postgis-users>> wrote: >//>/Hi, />//>/In my case I get dist_geography between gid 1 and gid 3 0 in both queries. />//>/This looks like a possible issue around the geos cache, something similar />/to https://trac.osgeo.org/postgis/ticket/4269. Can you run them with />/EXPLAIN ANALYZE? />//>/-- />/Raúl Marín Rodríguez />/carto.com />/___ />/postgis-users mailing list />/postgis-users at lists.osgeo.org <https://lists.osgeo.org/mailman/listinfo/postgis-users> />/https://lists.osgeo.org/mailman/listinfo/postgis-users / Hi, I have found what seems to be an extremely bizarre bug. I am not sure if it is a postgis issue or postgresql issue, but it leads ST_Distance(geography(Polygon,4326),geography(MultiPolygon,4326)) to return 0 (zero) when the true answer is not zero under certain peculiar situations. The table data for a set of polygons that generate the error can be downloaded here: http://www.davidmkaplan.fr/bad_geography_dist_calcs.csv Once you have that data, the code to generate the error is: CREATE TEMP TABLE tt (gid int PRIMARY KEY, geom geometry); \copy tt FROM 'bad_geography_dist_calcs.csv' WITH (FORMAT 'csv',HEADER TRUE) -- Incorrect dist_geography for t1.gid=1 and t2.gid=3 SELECT t1.gid AS gid1, t2.gid AS gid2, ST_Distance(t1.geom,t2.geom) AS dist_lonlat, ST_Distance(ST_Transform(t1.geom,26918), ST_Transform(t2.geom,26918)) AS dist_utm, ST_Distance(t1.geom::geography,t2.geom::geography) AS dist_geography, ST_Distance((ST_Dump(t1.geom)).geom::geography,(ST_Dump(t2.geom)).geom::geography) AS dist_geography_dump FROM tt t1 JOIN tt t2 ON t1.gid ST_Distance(t1.geom::geography,t2.geom::geography) AS dist_geography, ST_Distance((ST_Dump(t1.geom)).geom::geography,(ST_Dump(t2.geom)).geom::geography) AS dist_geography_dump FROM tt t1 JOIN tt t2 ON t1.gidHopefully, if I am not crazy, you will have a 0 value of dist_geography in the first query for t2.gid=3, but not in the second. Does anyone have any idea what could be causing this? To be honest, I am not really certain that I know where to start do diagnose the problem... The specifics of my installation are: POSTGIS="2.4.3 r16312" PGSQL="100" GEOS="3.7.0-CAPI-1.11.0 673b9939" PROJ="Rel. 4.9.3, 15 August 2016" GDAL="GDAL 2.3.2, released 2018/09/21" LIBXML="2.9.4" LIBJSON="0.12.1" LIBPROTOBUF="1.2.1" RASTER psql 10.6 (Ubuntu 10.6-0ubuntu0.18.04.1)) Ubuntu 18.04.1 LTS I have also tested this with postgis 2.5.0 r16836 and there the incorrect 0 is in both queries!?!?!? Thanks, David -- ** David M. Kaplan Charge de Recherche 1 Institut de Recherche pour le Developpement (IRD) UMR MARBEC (IRD/Ifremer/CNRS/UMII) av. Jean Monnet CS 30171 34203 Sete cedex France Email:david.kap...@ird.fr Phone: +33 (0)4 99 57 32 25 Fax: +33 (0)4 99 57 32 95 http://www.umr-marbec.fr/kaplan-david.html http://www.davidmkaplan.fr/ ** -- ** David M. Kaplan Charge de Recherche 1 Institut de Recherche pour le Developpement (IRD) UMR MARBEC (IRD/Ifremer/CNRS/UMII) av. Jean Monnet CS 30171 34203 Sete cedex France Email: david.kap...@ird.fr Phone: +33 (0)4 99 57 32 25 Fax: +33 (0)4 99 57 32 95 http://www.umr-marbec.fr/kaplan-david.html http://www.davidmkaplan.fr/ ** ___ postgis-users mailing list postgis-users@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/postgis-users
Re: [postgis-users] very strange bug generating incorrect ST_Distance calculations
On 10/01/2019 12:26, David M. Kaplan wrote: Comparing the results of _ST_DistanceTree and _ST_DistanceUnCached is a good way to tease out whether there are differences between the implementations without having to fight the caching machinery along the way. P. Thanks for the responses. I tried with _ST_DistanceTree and _ST_DistanceUnCached and as predicted, tree gives 0 and uncached gives the correct answer. What is the appropriate work around to avoid this issue? Who should I report this problem to? The maintainers of libgeos? Also, how should I interpret that this error only occurs when one of the geometries is a multipolygon? In response to Raúl, I found that both queries give 0 for postgis 2.5.0 r16836, but only the multiline query gives zero for POSTGIS="2.4.3 r16312"... Cheers, David On Thu, Jan 10, 2019 at 7:00 AM Raúl Marín Rodríguez https://lists.osgeo.org/mailman/listinfo/postgis-users>> wrote: >//>/Hi, />//>/In my case I get dist_geography between gid 1 and gid 3 0 in both queries. />//>/This looks like a possible issue around the geos cache, something similar />/to https://trac.osgeo.org/postgis/ticket/4269. Can you run them with />/EXPLAIN ANALYZE? />//>/-- />/Raúl Marín Rodríguez />/carto.com />/___ />/postgis-users mailing list />/postgis-users at lists.osgeo.org <https://lists.osgeo.org/mailman/listinfo/postgis-users> />/https://lists.osgeo.org/mailman/listinfo/postgis-users / Hi, I have found what seems to be an extremely bizarre bug. I am not sure if it is a postgis issue or postgresql issue, but it leads ST_Distance(geography(Polygon,4326),geography(MultiPolygon,4326)) to return 0 (zero) when the true answer is not zero under certain peculiar situations. The table data for a set of polygons that generate the error can be downloaded here: http://www.davidmkaplan.fr/bad_geography_dist_calcs.csv Once you have that data, the code to generate the error is: CREATE TEMP TABLE tt (gid int PRIMARY KEY, geom geometry); \copy tt FROM 'bad_geography_dist_calcs.csv' WITH (FORMAT 'csv',HEADER TRUE) -- Incorrect dist_geography for t1.gid=1 and t2.gid=3 SELECT t1.gid AS gid1, t2.gid AS gid2, ST_Distance(t1.geom,t2.geom) AS dist_lonlat, ST_Distance(ST_Transform(t1.geom,26918), ST_Transform(t2.geom,26918)) AS dist_utm, ST_Distance(t1.geom::geography,t2.geom::geography) AS dist_geography, ST_Distance((ST_Dump(t1.geom)).geom::geography,(ST_Dump(t2.geom)).geom::geography) AS dist_geography_dump FROM tt t1 JOIN tt t2 ON t1.gid ST_Distance(t1.geom::geography,t2.geom::geography) AS dist_geography, ST_Distance((ST_Dump(t1.geom)).geom::geography,(ST_Dump(t2.geom)).geom::geography) AS dist_geography_dump FROM tt t1 JOIN tt t2 ON t1.gidHopefully, if I am not crazy, you will have a 0 value of dist_geography in the first query for t2.gid=3, but not in the second. Does anyone have any idea what could be causing this? To be honest, I am not really certain that I know where to start do diagnose the problem... The specifics of my installation are: POSTGIS="2.4.3 r16312" PGSQL="100" GEOS="3.7.0-CAPI-1.11.0 673b9939" PROJ="Rel. 4.9.3, 15 August 2016" GDAL="GDAL 2.3.2, released 2018/09/21" LIBXML="2.9.4" LIBJSON="0.12.1" LIBPROTOBUF="1.2.1" RASTER psql 10.6 (Ubuntu 10.6-0ubuntu0.18.04.1)) Ubuntu 18.04.1 LTS I have also tested this with postgis 2.5.0 r16836 and there the incorrect 0 is in both queries!?!?!? Thanks, David -- ** David M. Kaplan Charge de Recherche 1 Institut de Recherche pour le Developpement (IRD) UMR MARBEC (IRD/Ifremer/CNRS/UMII) av. Jean Monnet CS 30171 34203 Sete cedex France Email: david.kap...@ird.fr Phone: +33 (0)4 99 57 32 25 Fax: +33 (0)4 99 57 32 95 http://www.umr-marbec.fr/kaplan-david.html http://www.davidmkaplan.fr/ ** ___ postgis-users mailing list postgis-users@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/postgis-users
[postgis-users] very strange bug generating incorrect ST_Distance calculations
Hi, I have found what seems to be an extremely bizarre bug. I am not sure if it is a postgis issue or postgresql issue, but it leads ST_Distance(geography(Polygon,4326),geography(MultiPolygon,4326)) to return 0 (zero) when the true answer is not zero under certain peculiar situations. The table data for a set of polygons that generate the error can be downloaded here: http://www.davidmkaplan.fr/bad_geography_dist_calcs.csv Once you have that data, the code to generate the error is: CREATE TEMP TABLE tt (gid int PRIMARY KEY, geom geometry); \copy tt FROM 'bad_geography_dist_calcs.csv' WITH (FORMAT 'csv',HEADER TRUE) -- Incorrect dist_geography for t1.gid=1 and t2.gid=3 SELECT t1.gid AS gid1, t2.gid AS gid2, ST_Distance(t1.geom,t2.geom) AS dist_lonlat, ST_Distance(ST_Transform(t1.geom,26918), ST_Transform(t2.geom,26918)) AS dist_utm, ST_Distance(t1.geom::geography,t2.geom::geography) AS dist_geography, ST_Distance((ST_Dump(t1.geom)).geom::geography,(ST_Dump(t2.geom)).geom::geography) AS dist_geography_dump FROM tt t1 JOIN tt t2 ON t1.gid ST_Distance(t1.geom::geography,t2.geom::geography) AS dist_geography, ST_Distance((ST_Dump(t1.geom)).geom::geography,(ST_Dump(t2.geom)).geom::geography) AS dist_geography_dump FROM tt t1 JOIN tt t2 ON t1.gidHopefully, if I am not crazy, you will have a 0 value of dist_geography in the first query for t2.gid=3, but not in the second. Does anyone have any idea what could be causing this? To be honest, I am not really certain that I know where to start do diagnose the problem... The specifics of my installation are: POSTGIS="2.4.3 r16312" PGSQL="100" GEOS="3.7.0-CAPI-1.11.0 673b9939" PROJ="Rel. 4.9.3, 15 August 2016" GDAL="GDAL 2.3.2, released 2018/09/21" LIBXML="2.9.4" LIBJSON="0.12.1" LIBPROTOBUF="1.2.1" RASTER psql 10.6 (Ubuntu 10.6-0ubuntu0.18.04.1)) Ubuntu 18.04.1 LTS I have also tested this with postgis 2.5.0 r16836 and there the incorrect 0 is in both queries!?!?!? Thanks, David -- ** David M. Kaplan Charge de Recherche 1 Institut de Recherche pour le Developpement (IRD) UMR MARBEC (IRD/Ifremer/CNRS/UMII) av. Jean Monnet CS 30171 34203 Sete cedex France Email: david.kap...@ird.fr Phone: +33 (0)4 99 57 32 25 Fax: +33 (0)4 99 57 32 95 http://www.umr-marbec.fr/kaplan-david.html http://www.davidmkaplan.fr/ ** ___ postgis-users mailing list postgis-users@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/postgis-users
Re: [postgis-users] Aggregating rasters by adding and other confusions
n place. > > What it seems to do is two passes > > > > 1 pass does an ST_Union using 'COUNT' (so in your case you'd get > numbers between 0 and 3, 0 being no rasters considered having any > data) > > 2nd pass does an ST_Union using 'SUM' > > And returns a new raster where each pixel is SUM/COUNT > > > > Basic code is here: > > > > http://postgis.net/docs/doxygen/2.4/da/dde/rtpg__mapalgebra_8c_a1d940 > 65e6cef5d5d61417b82b2cf4fb6.html#a1d94065e6cef5d5d61417b82b2cf4fb6 > > (note it only computes a value if the first band (which I presume to > be the count band) value > 0 and has no nodata value. ) > > I don't know why a count band would have a nodata value aside from > when it's 0 > > > > What does COUNT return for you? > > > > Thanks, > > Regina > > > > From: postgis-users [mailto:postgis-users-boun...@lists.osgeo.org] On > Behalf Of David M. Kaplan > Sent: Tuesday, March 06, 2018 10:29 AM > To: postgis-users@lists.osgeo.org > Subject: Re: [postgis-users] Aggregating rasters by adding and other > confusions > > > > Hi, > > > > Yes, this does seem to be the solution. For completeness, the > following did the same thing as my aggregate function: > > > > SELECT ST_Union(rast,1,'SUM') FROM blah; > > This includes replacing 0's with NULL's - ST_Union had the same > behavior at my aggregate function. > > > > Can you explain this? For my sum, this isn't a big deal, but > ST_Union(... 'MEAN') seems to be treating 0's at NULL values, thereby > removing them from the mean calculation: > > > > db=# SELECT (ST_SummaryStats(ST_Union(rast,1,'SUM'))).sum FROM blah; > sum > -- > 2792 > (1 row) > > db=# SELECT (ST_SummaryStats(ST_Union(rast,1,'MEAN'))).sum FROM blah; > sum > -- > 1346 > (1 row) > > > The table "blah" has 3 rasters with no null values, so I would expect > the sum of the mean to be 1/3 the sum of the sum, but clearly this is > not the case. > > > > I tried playing around with the no data value for my rasters to see > if this would change things. It does change things, but not in any > way I can explain: > > > > fads=# SELECT > (ST_SummaryStats(ST_Union(ST_SetBandNoDataValue(rast,1,0),1,'MEAN'))) > .sum FROM blah; > sum > -- > 2025 > (1 row) > > fads=# SELECT > (ST_SummaryStats(ST_Union(ST_SetBandNoDataValue(rast,1,- > 999),1,'MEAN'))).sum FROM blah; > sum > -- > 2025 > (1 row) > > > Setting the no data value has no effect on the sum of the sum. > > > > Does any of this make sense? > > > > Thanks, > > David > > > > > > On Mon, 2018-03-05 at 12:00 -0800, postgis-users-request@lists.osgeo. > org <mailto:postgis-users-requ...@lists.osgeo.org> wrote: > > Date: Mon, 5 Mar 2018 12:05:40 -0500 > From: "Regina Obe" <l...@pcorp.us <mailto:l...@pcorp.us> > > To: "'PostGIS Users Discussion'" <postgis-users@lists.osgeo.org to:postgis-users@lists.osgeo.org> > > Subject: Re: [postgis-users] Aggregating rasters by adding and other > confusions > Message-ID: <001f01d3b4a4$30d1fa20$9275ee60$@pcorp.us <mailto:001f01d > 3b4a4$30d1fa20$9275ee60$@pcorp.us> > > Content-Type: text/plain; charset="utf-8" > > I think what you are looking for is ST_Union (.. SUM) note this has > union types – FIRST, MIN, MAX, COUNT, SUM, MEAN, RANGE > > > > http://postgis.net/docs/manual-2.4/RT_ST_Union.html > > > > > > > > From: postgis-users [mailto:postgis-users-boun...@lists.osgeo.org] On > Behalf Of David M. Kaplan > Sent: Monday, March 05, 2018 10:03 AM > To: postgis-users@lists.osgeo.org <mailto:postgis-users@lists.osgeo.o > rg> > Subject: [postgis-users] Aggregating rasters by adding and other > confusions > > > > Hi, > > > > I have recently started working with the postgis raster > functionality. In general, I have found this really useful and have > been able to do some neat things fairly simply with this raster > functionality. Nevertheless, there are a few basic things that I am > confused about and I was hoping someone could give me a hand. > > > > (1) First of all, I have a table with a bunch of rasters that have > the same extent, alignment, scale, etc. and I want to aggregate them > together into a single raster using pixel-by-pixel addition. It seems > like there should be a function to do this, but I can't find one. Is > there an aggregate
Re: [postgis-users] Aggregating rasters by adding and other confusions
Hi, Yes, this does seem to be the solution. For completeness, the following did the same thing as my aggregate function: SELECT ST_Union(rast,1,'SUM') FROM blah; This includes replacing 0's with NULL's - ST_Union had the same behavior at my aggregate function. Can you explain this? For my sum, this isn't a big deal, but ST_Union(... 'MEAN') seems to be treating 0's at NULL values, thereby removing them from the mean calculation: db=# SELECT (ST_SummaryStats(ST_Union(rast,1,'SUM'))).sum FROM blah; sum -- 2792 (1 row) db=# SELECT (ST_SummaryStats(ST_Union(rast,1,'MEAN'))).sum FROM blah; sum -- 1346 (1 row) The table "blah" has 3 rasters with no null values, so I would expect the sum of the mean to be 1/3 the sum of the sum, but clearly this is not the case. I tried playing around with the no data value for my rasters to see if this would change things. It does change things, but not in any way I can explain: fads=# SELECT (ST_SummaryStats(ST_Union(ST_SetBandNoDataValue(rast,1,0),1,'MEAN'))).sum FROM blah; sum -- 2025 (1 row) fads=# SELECT (ST_SummaryStats(ST_Union(ST_SetBandNoDataValue(rast,1,-999),1,'MEAN'))).sum FROM blah; sum -- 2025 (1 row) Setting the no data value has no effect on the sum of the sum. Does any of this make sense? Thanks, David On Mon, 2018-03-05 at 12:00 -0800, postgis-users- requ...@lists.osgeo.org wrote: > Date: Mon, 5 Mar 2018 12:05:40 -0500 > From: "Regina Obe" <l...@pcorp.us> > To: "'PostGIS Users Discussion'" <postgis-users@lists.osgeo.org> > Subject: Re: [postgis-users] Aggregating rasters by adding and other > confusions > Message-ID: <001f01d3b4a4$30d1fa20$9275ee60$@pcorp.us> > Content-Type: text/plain; charset="utf-8" > > I think what you are looking for is ST_Union (.. SUM) note this has union > types – FIRST, MIN, MAX, COUNT, SUM, MEAN, RANGE > > > > http://postgis.net/docs/manual-2.4/RT_ST_Union.html > > > > > > > > From: postgis-users [mailto:postgis-users-boun...@lists.osgeo.org] On Behalf > Of David M. Kaplan > Sent: Monday, March 05, 2018 10:03 AM > To: postgis-users@lists.osgeo.org > Subject: [postgis-users] Aggregating rasters by adding and other confusions > > > > Hi, > > > > I have recently started working with the postgis raster functionality. In > general, I have found this really useful and have been able to do some neat > things fairly simply with this raster functionality. Nevertheless, there are > a few basic things that I am confused about and I was hoping someone could > give me a hand. > > > > (1) First of all, I have a table with a bunch of rasters that have the same > extent, alignment, scale, etc. and I want to aggregate them together into a > single raster using pixel-by-pixel addition. It seems like there should be a > function to do this, but I can't find one. Is there an aggregate > "ST_MapAlgebra" function? > > > > Given that I couldn't find one, I defined an aggregate function as follows: > > > > CREATE OR REPLACE FUNCTION AddRasters(r1 raster, r2 raster) >RETURNS raster AS > $BODY$ > SELECT ST_MapAlgebra($1,$2,'[rast1]+[rast2]'); > $BODY$ LANGUAGE 'sql' IMMUTABLE STRICT; > > CREATE AGGREGATE sum (raster) > ( > sfunc = AddRasters, > stype = raster > ); > > > (2) This seems to work, but it has the unexpected behavior that it replaces 0 > values with NULL. In my case, this is fine, but I am wondering why it does > this? I can't find anything that indicates that it should be replacing zeros > with NULL. Here is the metadata associated with one of my rasters (the others > are similar): > > > > # SELECT ST_BandMetadata(rast), ST_Metadata(rast), ST_SummaryStats(rast) FROM > blah; > -[ RECORD 1 ]---+--- > st_bandmetadata | (16BUI,,f,) > st_metadata | (-180,90,360,180,1,-1,0,0,4326,1) > st_summarystats | (64800,417,0.00643518518518519,0.223617719977485,0,46) > > > I have not defined a nodataval for these layers and the original layers have > no NULL values. > > > > (3) Is there a postgis command to turn all the NULL values back into zeros? > > > > (4) I was also considering just defining the '+' operator for raster + raster > to be pixel-by-pixel addition. Is there any reason that I wouldn't want to do > this? > > > > (5) Finally, I have been visualizing results with QGIS using the DB Manager. > However, I don't see how to select a row from a raster table and incorporate > just that row into the canvas. Is there a way to do this? > > > > Thanks for the assistance, > > David___ postgis-users mailing list postgis-users@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/postgis-users
[postgis-users] Aggregating rasters by adding and other confusions
Hi, I have recently started working with the postgis raster functionality. In general, I have found this really useful and have been able to do some neat things fairly simply with this raster functionality. Nevertheless, there are a few basic things that I am confused about and I was hoping someone could give me a hand. (1) First of all, I have a table with a bunch of rasters that have the same extent, alignment, scale, etc. and I want to aggregate them together into a single raster using pixel-by-pixel addition. It seems like there should be a function to do this, but I can't find one. Is there an aggregate "ST_MapAlgebra" function? Given that I couldn't find one, I defined an aggregate function as follows: CREATE OR REPLACE FUNCTION AddRasters(r1 raster, r2 raster) RETURNS raster AS$BODY$SELECT ST_MapAlgebra($1,$2,'[rast1]+[rast2]');$BODY$ LANGUAGE 'sql' IMMUTABLE STRICT; CREATE AGGREGATE sum (raster)(sfunc = AddRasters,stype = raster); (2) This seems to work, but it has the unexpected behavior that it replaces 0 values with NULL. In my case, this is fine, but I am wondering why it does this? I can't find anything that indicates that it should be replacing zeros with NULL. Here is the metadata associated with one of my rasters (the others are similar): # SELECT ST_BandMetadata(rast), ST_Metadata(rast), ST_SummaryStats(rast) FROM blah;-[ RECORD 1 ]---+ ---st_bandmetadata | (16BUI,,f,)st_metadata | (-180,90,360,180,1,- 1,0,0,4326,1)st_summarystats | (64800,417,0.00643518518518519,0.223617719977485,0,46) I have not defined a nodataval for these layers and the original layers have no NULL values. (3) Is there a postgis command to turn all the NULL values back into zeros? (4) I was also considering just defining the '+' operator for raster + raster to be pixel-by-pixel addition. Is there any reason that I wouldn't want to do this? (5) Finally, I have been visualizing results with QGIS using the DB Manager. However, I don't see how to select a row from a raster table and incorporate just that row into the canvas. Is there a way to do this? Thanks for the assistance, David -- **David M. KaplanCharge de Recherche 1 Institut de Recherche pour le Developpement (IRD)UMR MARBEC (IRD/Ifremer/CNRS/UMII)av. Jean MonnetCS 3017134203 Sete cedexFrance Email: david.kaplan@ird.frPhone: +33 (0)4 99 57 32 25Fax: +33 (0)4 99 57 32 95 http://www.umr-marbec.fr/kaplan-david.htmlhttp://www.davidmkaplan.fr/**___ postgis-users mailing list postgis-users@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/postgis-users
[postgis-users] error while loading dat for tiger geocoder
Hi, I am in the process of setting up the tiger geocoder for use with addresses from the state of virginia. I am using PostgreSQL 9.3.20, postgis 2.2.2 and postgis_tiger_geocoder 2.2.2. The state loader script (for linux shell) generates an SQL error towards the end of the script relating to the following line (line 210 in my script): ${PSQL} -c "ALTER TABLE tiger_staging.VA_tabblock RENAME geoid10 TO tabblock_id;" The error it generates is: ERROR: relation "tiger_staging.va_tabblock" does not exist I am under the impression that this error is harmless as at this stage the tiger_staging schema is empty because the function loader_load_staged_data() has already been run, but I want to make sure. If this is indeed a harmless error, then why is this line in the script? Thanks, David -- ********** David M. Kaplan Charge de Recherche 1 Institut de Recherche pour le Developpement (IRD) UMR MARBEC (IRD/Ifremer/CNRS/UMII) av. Jean Monnet CS 30171 34203 Sete cedex France Phone: +33 (0)4 99 57 32 25 Fax: +33 (0)4 99 57 32 95 http://www.umr-marbec.fr/kaplan-david.html ** ___ postgis-users mailing list postgis-users@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/postgis-users