2

I'm trying to aggregate data up from 0.05 to a coarser 0.25 resolution and am having trouble getting gdalwarp to work correctly. From comparing the original input data with the aggregated data and manually calculating the averages, I can see that the result I'm getting is not what I would expect. Many (not all) aggregated pixel values are just wrong (FYI, I have specifically been looking where the water meets land on the horn of Africa). Also, it's my understanding that gdal is supposed to ignore NoData pixels in the average calculation, but in some areas, the 25 pixels underneath the aggregated pixels were all NoData, yet the aggregated pixel reported a number, when it should have been NoData.

The command I used below:

gdalwarp -tr .25 .25 -r average -srcnodata -15000 -dstnodata -15000 /path/to/temp_CHIRPS_201512.tif /path/to/CHIRPS_201512_ave.tif

I can post the entire script if need be.

I'm using gdal 2.0.1 on linux.

Here is a link to the input file I used, in addition to the result I got by calling the above command: https://www.dropbox.com/sh/dk0c41638jg8fir/AAD9z9IXeEfJ2PyEUVnjG4aya?dl=0

I tried running the same command on a different dataset (I am doing the same process for several datasets) and it appears gdalwarp did in fact work.

Can anyone help me understand why it is working for one dataset and not the other?

  • 1
    Interesting, do you have a representative datasource to reconstruct the behavior? – huckfinn Feb 24 '16 at 3:34
  • @huckfinn I do. How can I get it to you? – user20408 Feb 24 '16 at 15:30
  • Download, Dropbox, Github – huckfinn Feb 24 '16 at 15:34
  • @huckfinn here is a link to the input file I used and the output I got with the above command: dropbox.com/sh/dk0c41638jg8fir/AAD9z9IXeEfJ2PyEUVnjG4aya?dl=0 – user20408 Feb 24 '16 at 15:38
  • Hi @huckfinn were you able to receive my data from drop box? Thanks so much! – user20408 Feb 25 '16 at 17:48
2

I think the index numbers for the row of the averaging window function is shifted and that seems to be an error. To get this opinion, I made the following test using QGIS and a subsample of the both datasets you left in the dropbox.

In first step I cut out a windows in the region of your interest ..some where to address an area in the northern part of the Red Sea:

gdal_translate -projwin 0 50 50 0 temp_CHIRPS_201512.tif hires.tif
gdal_translate -projwin 0 50 50 0 CHIRPS_201512_ave.tif lowres.tif

subset maps

I call the full resolution image hires.tif and the averaged one lowres.tif. To see what the gdalwarp tool produces I built xyz files to get frames and numbers corresponding to the cells

gdal_translate hires.tif hires.xyz -of xyz
gdal_translate lowres.tif lowres.xyz -of xyz

and import the stuff into QGIS via Layer -> Add Layer -> Add Delimited Text Layer with the following settings:

import hires.xyz

The same procedure was applied to the lowres stuff. I zoom to the paninsula of Sinai and set the symbology of the raster to a squared marker without background and a size of 0.05 map units for the hires positions and 0.25 map units for the lowres points (your produced raster sizes).

raster and frames

I choose colored symbology for the hires raster:

hires raster settings

and a gray setup for the lowres (averaged) raster to bring out the artefacts.

lowres raster settings

At 35° East and 27° North we find a nice artefact that should be a -15000 cell (NA in GNU-R speech).

artefact region

To see whats happend I zoom in and label the cells:

wrong cell

It seems, that we have a cell shift and the average is calculated from a wrong set of samples,

shift window

because if I use to north shifted window, I will calculate:

> cell <- c(
+   NA,  6, 27, 25, NA, 
+   NA, NA, NA, NA, NA,
+   NA, NA, NA, NA, NA,
+   NA, NA, NA, NA, NA,
+   NA, NA, NA, NA, NA
+ );
> mean(cell,na.rm=TRUE)
[1] 19.33333

and get the magic 19 at a wrong place in the raster. I think you should send a bug report to the gdal community.

Sorry forgot to use python...

In [1]: import numpy

In [2]: cell = numpy.array ([
   ...:   None,  6, 27, 25, None, 
   ...:   None, None, None, None, None,
   ...:   None, None, None, None, None,
   ...:   None, None, None, None, None,
   ...:   None, None, None, None, None
   ...: ])

In [3]: numpy.mean(filter(None, cell))
Out[3]: 19.333333333333332

In this dropbox you will find the test environment with:

  • avg-error.qgs the qgis project
  • temp_CHIRPS_201512.tif your hires file
  • CHIRPS_201512_ave.tif your lowres file
  • hires.tif the 0 50 50 0 cut version of the hires file
  • hires.xyz the xyz export of the hires cut version
  • lowres.tif the cut version of the lowres file
  • lowres.xyz the xyz export of the hires cut version
  • hires-point.qml style to draw the cells hires.xyz
  • lowres-point.qml style to draw the cells lowres.xyz
  • hires-raster.qml style to draw the hires raster
  • lowres-raster.qml style to draw the lowres raster
  • howto code snippets

enjoy...

  • WOW. thank you. I'll review this in more detail later. Do you have any idea why it works for other inputs (different dataset) but not this one? – user20408 Feb 26 '16 at 17:04
  • I just looked at my results and unfortunately I don't see the shift. Have you taken a look at the other raster I left on dropbox (the result I got)? Maybe we are getting different results. Thanks again. – user20408 Feb 26 '16 at 17:35
  • Hm... strange. I use the both datasets from the dropbox and choose exactly the cell samples. The position of the cell with the value 6 is (35,05:27,8) (35,1:27,75) in the file temp_CRISP_201512.tif. – huckfinn Feb 27 '16 at 5:42
  • 1
    I put the test environment in the dropbox dropbox.com/sh/jrsxye0ss24woz6/AAAoIf7pw1Fx6DpQtCJbkMnWa?dl=0. Look at the comments I added at the end of my answer... – huckfinn Feb 27 '16 at 6:22

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.