11

I have a set of rasters (.vrt) with daily soil moisture data. I want to sum the pixels of all rasters in order to have a measure by month. However, the whole world is not covered each day, which results in nodata values at places where measures exist for the other days.

What I want to do, is sum the values of each raster. However, it seems that each time a nodata pixel is in the sum, the resulting pixel is directly categorized as nodata. I would like to have the opposite: ignoring all nodata values and summing the rest.

I thought of 2 ways of solving the problem:

  1. summing rasters ignoring nodata values
  2. converting nodata pixels to value 0, then sum all the rasters

enter image description here

Unfortunately, I can't find any tool to do this.

Can anyone help me?

|improve this question
  • I wonder how interpretable the result will be: after all, the sums will not include the values for the missing days, indicating they will be biased low by various amounts depending on the amount of missing data. If there is any appreciable amount of missing data, then you should consider instead imputing or predicting the values at the missing cells and then performing the sum. – whuber♦ Oct 16 '13 at 17:14
  • 1
    Thanks for this comment. I asked the managers of the data and they said they calculated the mean value by summing all available data, then divide by the number of days within the month when observations are available on that pixel. – Damien Oct 17 '13 at 13:09
  • 2
    That is identical to the local mean of a stack of rasters where the GIS is instructed just to ignore all NoData cells, which suggests another way for you to go. By the way, this procedure still potentially is biased when the missing observations are correlated with the values. For instance, when you are missing data due to cloud cover, it is plausible that on those days the soil moisture might be higher (on average) than usual. – whuber♦ Oct 17 '13 at 13:22
  • 1
    Alright, I just tested the cell statistics tool with the 'ignore nodata' checked and it does exactly what I wanted. I am aware of the shortcomings of such method. However, my work is at a rather large scale, both in space and time, so I assume that such biases will be less important. Anyway, I don't really have so much choice of data. – Damien Oct 17 '13 at 13:35
  • The scale will not be related to the amount of bias. The importance of the bias depends on its size and the sensitivity of your application to it. Although you might not be able to do anything about the data, you do have control over the procedures for analyzing them. There are plenty of ways to improve on this one. Among your options are interpolating over time and regressing against covariates that might be correlated with soil moisture and are available even on the missing dates. (There is an entire branch of statistics devoted to "imputing" missing values.) – whuber♦ Oct 17 '13 at 13:39
2

If you have access to ArcGIS then the Cell Statistics tool has the optional to ignore nodata which you sum\mean\min\max your rasters.

|improve this answer
  • If you modify this answer to replace "sum" by "average," it will be correct: see the comment thread to the question for the reason. – whuber♦ Oct 17 '13 at 13:40
  • Summing was only one part of my problem since I want to get the mean value. I just checked, and the Cell statitics tool is also the easiest way to sum rasters disregarding nodata pixels. – Damien Oct 17 '13 at 13:48
10

Firstly, you can use gdal_calc.py to change all -9999 to 0 and set the NoData value to 0.

For instance: 

gdal_calc.py -A input.tif --outfile=input_with_NoData.tif --calc="A+9999*(A==-9999)" --NoDataValue=0

Then, you can ignore NoData value using gdal_translate with the -a_nodata option followed by none.

-a_nodata value:

Assign a specified nodata value to output bands. Starting with GDAL 1.8.0, can be set to none to avoid setting a nodata value to the output file if one exists for the source file

Example:

gdal_translate -a_nodata none input_with_NoData.tif output_without_NoData.tif
|improve this answer
  • Thanks for the tip. I played a little with the -a_nodata parameter, and I manage to have the value of -9999 instead of NULL. However, I still need to change all -9999 pixels to 0 – Damien Oct 16 '13 at 7:43
  • 1
    You can use gdal_calc.py to change all -9999 to 0 before applying gdal_translate -a_nodata none .... For instance: gdal_calc.py -A input.tif --outfile=input_with_NoData.tif --type=Int16 --calc="A+9999*(A==-9999)" --NoDataValue=0 – Antonio Falciano Oct 16 '13 at 8:22
  • I tried several commands, but I really think I have a problem with gdal_calc.py. For instance, I just tried to add 1 to the whole raster (i.imgur.com/WiZG7MC.png) and I obtained this Maxrepeat error. I don't understand, I ran the command directly from the OSGeo shell, and this module seems to be correctly installed since it appears in the list when I open the shell (i.imgur.com/fgtMZQZ.png). Is my install broken? – Damien Oct 16 '13 at 8:57
  • It seems like a regex issue... Try changing your working directory using the cd command and retry in order to have only -A test.tif. – Antonio Falciano Oct 16 '13 at 9:28
  • 1
    It works for me, but I have installed all this stuff via the OSGeo4W installer. Instead, I see that you have installed QGIS via the standalone installer. So try to refresh/update your setup. – Antonio Falciano Oct 16 '13 at 10:54
5

In R:

library(raster)
s <- stack('file1', 'file2', ...)
ss <- sum(s, na.rm=TRUE)
|improve this answer
  • Could you add an explanation to this code (if you wrote it) or a link to where you found it, with a short summary from the site? – Paul Oct 15 '13 at 20:14
  • 1
    The above first loads the raster package, then creates a "stack" (a 3-dimensional array of raster files of identical extent and resolution, where each file becomes a slice of this stack) comprising the files indicated by the comma-separated file names. The final line performs cell-wise sums across all slices of stack s, with the argument na.rm=TRUE resulting in NA values being ignored. Output object ss is a raster object that can be exported with writeRaster (see ?writeRaster). (@RobertH is the creator of the R raster package.) – jbaums Jun 1 '14 at 0:08
1

The Conditional Tool in Spatial Analyst is an easy way to convert null values to zeros. Then taking the sums should be a snap.

|improve this answer
  • Thank you for your reply. I followed this post: support.esri.com/fr/knowledgebase/techarticles/detail/34932 and was able to convert all null pixels to 0. However, I would have liked to be able to do this operation for all the rasters within a directory. The 'Reclassify grid values' in QGIS seems to be able to do this, but all I don't know how it works. I choose "Simple table" method, and ask it to set null to 0, and don't change other values, but the resulting raster ranges from -0.99 to 0.08 while it originally spans from 428 to 3491 – Damien Oct 15 '13 at 14:33
  • Can you post an excerpt of your original raster somewhere, so we can have a look? – til_b Oct 15 '13 at 14:35
  • 2
    If you right click on the tool you can run a batch process on as many files as you'd like, i.e. an entire directory. – brock Oct 15 '13 at 14:49
  • Thanks for your answer. However, I can convert all NULLs to 0 for a single file using 'spatial analyst/Map algebra/Raster calculator'. If I right clic on this tool, I can batch a directory, but what will be the syntax in order to take into account each file automatically? (i.imgur.com/aYaUCzz.png) – Damien Oct 16 '13 at 15:44
  • I had a look at the 'spatial analyst/Reclass/Reclassify' tool and it does exactly what I want. I managed to use the batch window in order to process multiple files at once. I have a last problem: if I right-click/fill the 'output raster', it copies the same output files for all input rasters. How can I do to make it create a different output for each input? (i.imgur.com/jzTI2x9.png) – Damien Oct 16 '13 at 16:16
1

I had the same problem a while ago and I managed to solve it.

Just as a note about the display of nodata values in QGIS: nodata values are always shown as nodata values in the map window and when using the "object information" tool, no matter if they actually have the value -9999, 0, NULL etc. So after the reclassification of the nodata values to 0, if you go in the layer settings you will find in the Metadata that nodata have the value 0.

I finally found the perfect tool for the quest, the GRASS r.series. With the "aggregate operation" setting "maximum" and the "propagate NULLs" unchecked, the tool will create a Layer that fits your requirement.

|improve this answer
0

I remember having this problem a while back. As a solution I think I did +1 to all cells in the raster, added the rasters together, then did a final subtraction (of however many +1s you have added). It's pretty crude but if it gets the job done...!

Alternatively, use the 'Reclassify grid values' in the Processing toolbox in QGIS to convert your null values to zero.

|improve this answer
  • Hi, I tried many ways of using the 'reclassify grid values' tool with no success. I started on a raster where I converted all NULL pixels to -9999. I tried to reclassify using the method [0] single, and set 'old value = -9999', 'new value = 0.0', 'operator: [0]=', but all I get is a raster full of 'nan' pixels. What am I doing wrong? – Damien Oct 16 '13 at 15:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.