<html xmlns:v="urn:schemas-microsoft-com:vml" xmlns:o="urn:schemas-microsoft-com:office:office" xmlns:w="urn:schemas-microsoft-com:office:word" xmlns:m="http://schemas.microsoft.com/office/2004/12/omml" xmlns="http://www.w3.org/TR/REC-html40"><head><meta http-equiv=Content-Type content="text/html; charset=utf-8"><meta name=Generator content="Microsoft Word 15 (filtered medium)"><style><!--
/* Font Definitions */
@font-face
        {font-family:"Cambria Math";
        panose-1:2 4 5 3 5 4 6 3 2 4;}
@font-face
        {font-family:Calibri;
        panose-1:2 15 5 2 2 2 4 3 2 4;}
@font-face
        {font-family:Consolas;
        panose-1:2 11 6 9 2 2 4 3 2 4;}
/* Style Definitions */
p.MsoNormal, li.MsoNormal, div.MsoNormal
        {margin:0in;
        margin-bottom:.0001pt;
        font-size:12.0pt;
        font-family:"Times New Roman",serif;}
a:link, span.MsoHyperlink
        {mso-style-priority:99;
        color:blue;
        text-decoration:underline;}
a:visited, span.MsoHyperlinkFollowed
        {mso-style-priority:99;
        color:purple;
        text-decoration:underline;}
pre
        {mso-style-priority:99;
        mso-style-link:"HTML Preformatted Char";
        margin:0in;
        margin-bottom:.0001pt;
        font-size:10.0pt;
        font-family:"Courier New";}
span.HTMLPreformattedChar
        {mso-style-name:"HTML Preformatted Char";
        mso-style-priority:99;
        mso-style-link:"HTML Preformatted";
        font-family:"Consolas",serif;}
span.EmailStyle19
        {mso-style-type:personal-reply;
        font-family:"Calibri",sans-serif;
        color:#1F497D;}
span.stringliteral
        {mso-style-name:stringliteral;}
span.keywordflow
        {mso-style-name:keywordflow;}
.MsoChpDefault
        {mso-style-type:export-only;
        font-size:10.0pt;}
@page WordSection1
        {size:8.5in 11.0in;
        margin:1.0in 1.0in 1.0in 1.0in;}
div.WordSection1
        {page:WordSection1;}
--></style><!--[if gte mso 9]><xml>
<o:shapedefaults v:ext="edit" spidmax="1026" />
</xml><![endif]--><!--[if gte mso 9]><xml>
<o:shapelayout v:ext="edit">
<o:idmap v:ext="edit" data="1" />
</o:shapelayout></xml><![endif]--></head><body lang=EN-US link=blue vlink=purple><div class=WordSection1><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>David,<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'><o:p> </o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>I took a cursory glance at the union code we have in place.<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>What it seems to do is two passes<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'><o:p> </o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>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)<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>2<sup>nd</sup> pass does an ST_Union using 'SUM'  <o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>And returns a new raster  where each pixel is  SUM/COUNT<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'><o:p> </o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>Basic code is here:<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'><o:p> </o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'><a href="http://postgis.net/docs/doxygen/2.4/da/dde/rtpg__mapalgebra_8c_a1d94065e6cef5d5d61417b82b2cf4fb6.html#a1d94065e6cef5d5d61417b82b2cf4fb6">http://postgis.net/docs/doxygen/2.4/da/dde/rtpg__mapalgebra_8c_a1d94065e6cef5d5d61417b82b2cf4fb6.html#a1d94065e6cef5d5d61417b82b2cf4fb6</a><o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>(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. )<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>I don't know why a count band would have a nodata value aside from when it's 0<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'><o:p> </o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>What does COUNT return for you?  <o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'><o:p> </o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>Thanks,<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'>Regina<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri",sans-serif;color:#1F497D'><o:p> </o:p></span></p><div><div style='border:none;border-top:solid #E1E1E1 1.0pt;padding:3.0pt 0in 0in 0in'><p class=MsoNormal style='margin-left:.5in'><b><span style='font-size:11.0pt;font-family:"Calibri",sans-serif'>From:</span></b><span style='font-size:11.0pt;font-family:"Calibri",sans-serif'> postgis-users [mailto:postgis-users-bounces@lists.osgeo.org] <b>On Behalf Of </b>David M. Kaplan<br><b>Sent:</b> Tuesday, March 06, 2018 10:29 AM<br><b>To:</b> postgis-users@lists.osgeo.org<br><b>Subject:</b> Re: [postgis-users] Aggregating rasters by adding and other confusions<o:p></o:p></span></p></div></div><p class=MsoNormal style='margin-left:.5in'><o:p> </o:p></p><div><p class=MsoNormal style='margin-left:.5in'>Hi,<o:p></o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'><o:p> </o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'>Yes, this does seem to be the solution. For completeness, the following did the same thing as my aggregate function:<o:p></o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'><o:p> </o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'>SELECT ST_Union(rast,1,'SUM') FROM blah;<o:p></o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'>This includes replacing 0's with NULL's - ST_Union had the same behavior at my aggregate function. <o:p></o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'><o:p> </o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'>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:<o:p></o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'><o:p> </o:p></p></div><pre style='margin-left:.5in'>db=# SELECT (ST_SummaryStats(ST_Union(rast,1,'SUM'))).sum FROM blah; sum  <o:p></o:p></pre><pre style='margin-left:.5in'>------<o:p></o:p></pre><pre style='margin-left:.5in'> 2792<o:p></o:p></pre><pre style='margin-left:.5in'>(1 row)<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>db=# SELECT (ST_SummaryStats(ST_Union(rast,1,'MEAN'))).sum FROM blah;<o:p></o:p></pre><pre style='margin-left:.5in'> sum  <o:p></o:p></pre><pre style='margin-left:.5in'>------<o:p></o:p></pre><pre style='margin-left:.5in'> 1346<o:p></o:p></pre><pre style='margin-left:.5in'>(1 row)<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><div><p class=MsoNormal style='margin-left:.5in'>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.<o:p></o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'><o:p> </o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'>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:<o:p></o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'><o:p> </o:p></p></div><pre style='margin-left:.5in'>fads=# SELECT (ST_SummaryStats(ST_Union(ST_SetBandNoDataValue(rast,1,0),1,'MEAN'))).sum FROM blah;<o:p></o:p></pre><pre style='margin-left:.5in'> sum  <o:p></o:p></pre><pre style='margin-left:.5in'>------<o:p></o:p></pre><pre style='margin-left:.5in'> 2025<o:p></o:p></pre><pre style='margin-left:.5in'>(1 row)<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>fads=# SELECT (ST_SummaryStats(ST_Union(ST_SetBandNoDataValue(rast,1,-999),1,'MEAN'))).sum FROM blah;<o:p></o:p></pre><pre style='margin-left:.5in'> sum  <o:p></o:p></pre><pre style='margin-left:.5in'>------<o:p></o:p></pre><pre style='margin-left:.5in'> 2025<o:p></o:p></pre><pre style='margin-left:.5in'>(1 row)<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><div><p class=MsoNormal style='margin-left:.5in'>Setting the no data value has no effect on the sum of the sum.<o:p></o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'><o:p> </o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'>Does any of this make sense?<o:p></o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'><o:p> </o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'>Thanks,<o:p></o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'>David<o:p></o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'><o:p> </o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'><o:p> </o:p></p></div><div><p class=MsoNormal style='margin-left:.5in'>On Mon, 2018-03-05 at 12:00 -0800, <a href="mailto:postgis-users-request@lists.osgeo.org">postgis-users-request@lists.osgeo.org</a> wrote:<o:p></o:p></p></div><blockquote style='border:none;border-left:solid #729FCF 1.5pt;padding:0in 0in 0in 6.0pt;margin-left:4.8pt;margin-right:0in'><pre style='margin-left:.5in'>Date: Mon, 5 Mar 2018 12:05:40 -0500<o:p></o:p></pre><pre style='margin-left:.5in'>From: "Regina Obe" <<a href="mailto:lr@pcorp.us">lr@pcorp.us</a>><o:p></o:p></pre><pre style='margin-left:.5in'>To: "'PostGIS Users Discussion'" <<a href="mailto:postgis-users@lists.osgeo.org">postgis-users@lists.osgeo.org</a>><o:p></o:p></pre><pre style='margin-left:.5in'>Subject: Re: [postgis-users] Aggregating rasters by adding and other<o:p></o:p></pre><pre style='margin-left:.5in'>        confusions<o:p></o:p></pre><pre style='margin-left:.5in'>Message-ID: <<a href="mailto:001f01d3b4a4$30d1fa20$9275ee60$@pcorp.us">001f01d3b4a4$30d1fa20$9275ee60$@pcorp.us</a>><o:p></o:p></pre><pre style='margin-left:.5in'>Content-Type: text/plain; charset="utf-8"<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>I think what you are looking for is ST_Union (.. SUM)  note this has union types – FIRST, MIN, MAX, COUNT, SUM, MEAN, RANGE<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'><a href="http://postgis.net/docs/manual-2.4/RT_ST_Union.html">http://postgis.net/docs/manual-2.4/RT_ST_Union.html</a><o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>From: postgis-users [<a href="mailto:postgis-users-bounces@lists.osgeo.org">mailto:postgis-users-bounces@lists.osgeo.org</a>] On Behalf Of David M. Kaplan<o:p></o:p></pre><pre style='margin-left:.5in'>Sent: Monday, March 05, 2018 10:03 AM<o:p></o:p></pre><pre style='margin-left:.5in'>To: <a href="mailto:postgis-users@lists.osgeo.org">postgis-users@lists.osgeo.org</a><o:p></o:p></pre><pre style='margin-left:.5in'>Subject: [postgis-users] Aggregating rasters by adding and other confusions<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>Hi,<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>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.<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>(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? <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>Given that I couldn't find one, I defined an aggregate function as follows:<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>CREATE OR REPLACE FUNCTION AddRasters(r1 raster, r2 raster)<o:p></o:p></pre><pre style='margin-left:.5in'>       RETURNS raster AS<o:p></o:p></pre><pre style='margin-left:.5in'>$BODY$<o:p></o:p></pre><pre style='margin-left:.5in'>SELECT ST_MapAlgebra($1,$2,'[rast1]+[rast2]');<o:p></o:p></pre><pre style='margin-left:.5in'>$BODY$ LANGUAGE 'sql' IMMUTABLE STRICT;<o:p></o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'>CREATE AGGREGATE sum (raster)<o:p></o:p></pre><pre style='margin-left:.5in'>(<o:p></o:p></pre><pre style='margin-left:.5in'>    sfunc = AddRasters,<o:p></o:p></pre><pre style='margin-left:.5in'>    stype = raster<o:p></o:p></pre><pre style='margin-left:.5in'>);<o:p></o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>(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):<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'># SELECT ST_BandMetadata(rast), ST_Metadata(rast), ST_SummaryStats(rast) FROM blah;<o:p></o:p></pre><pre style='margin-left:.5in'>-[ RECORD 1 ]---+-------------------------------------------------------<o:p></o:p></pre><pre style='margin-left:.5in'>st_bandmetadata | (16BUI,,f,)<o:p></o:p></pre><pre style='margin-left:.5in'>st_metadata     | (-180,90,360,180,1,-1,0,0,4326,1)<o:p></o:p></pre><pre style='margin-left:.5in'>st_summarystats | (64800,417,0.00643518518518519,0.223617719977485,0,46)<o:p></o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>I have not defined a nodataval for these layers and the original layers have no NULL values.<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>(3) Is there a postgis command to turn all the NULL values back into zeros?<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>(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?<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>(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?<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'> <o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>Thanks for the assistance,<o:p></o:p></pre><pre style='margin-left:.5in'><o:p> </o:p></pre><pre style='margin-left:.5in'>David<o:p></o:p></pre></blockquote></div></body></html>