[R-sig-Geo] raster pkg: overwriting the dimensions detected by netcdf4 (PERSIANN-CDR)

Chris Reudenbach reudenbach at uni-marburg.de
Wed Feb 3 17:01:21 CET 2016


If you look at the file with ncview and analyse the header with ncdump - 
h you will find that en extent is missing.
The standard tool ncview shows also a "flipped version"

Try to use gdal more directly like this snippet

library(raster)
library(gdalUtils)
library(mapview)

# necessary because of nc file see configuration options 
http://www.gdal.org/frmt_netcdf.html
Sys.setenv(GDAL_NETCDF_BOTTOMUP="NO")
extentNC <- extent(0, 360,-60, 60)
proj4str <- '+proj=longlat +datum=WGS84 +no_defs'
ncfilename <- "PERSIANN-CDR_v01r01_20090101_c20140523.nc"
dataset<- "precipitation"
band <-1

x<-gdal_translate(paste0("NETCDF:",ncfilename,":",dataset), 
paste0(ncfilename,".tif"),
                b=band,
                of="GTiff",
                output_Raster=TRUE,
                verbose=TRUE,
                overwrite=TRUE)
# due to GDAL_NETCDF_BOTTOMUP="NO" we have to flip the raster
x <- flip(x, direction='y')
extent(x) <- extentNC
projection(x) <- proj4str

mapview(x)
plot(x)


cheers Chris


Am 03.02.2016 um 16:20 schrieb MAURICIO ZAMBRANO BIGIARINI:
> On 3 February 2016 at 11:34, Michael Sumner <mdsumner at gmail.com> wrote:
>>
>> On Thu, 4 Feb 2016 at 01:26 MAURICIO ZAMBRANO BIGIARINI
>> <mauricio.zambrano at ufrontera.cl> wrote:
>>> Dear spatial community,
>>>
>>> While reading some netCDFfiles with the raster package , I got a
>>> "rotated" file, where the columns where read as rows and vice versa:
>>>
>>> ------------- START ----------------
>>> x <- raster("PERSIANN-CDR_v01r01_20090101_c20140523.nc")
>>> Loading required namespace: ncdf4
>>> plot(x)
>>> x
>>> class       : RasterLayer
>>> dimensions  : 1440, 480, 691200  (nrow, ncol, ncell)
>>> resolution  : 0.25, 0.25  (x, y)
>>> extent      : -60, 60, 0, 360  (xmin, xmax, ymin, ymax)
>>> coord. ref. : NA
>>> data source : /home/hzambran/PERSIANN-CDR_v01r01_20090101_c20140523.nc
>>> names       : NOAA.Climate.Data.Record.of.PERSIANN.CDR.daily.precipitation
>>> z-value     : 2009-01-01
>>> zvar        : precipitation
>>> ------------- END ----------------
>>>
>> Transpose might be sufficient?
>>
>> tx <- t(x)
>>
>> plot(tx)
>> tx
>
> Thanks Michale,
>
> In fact, that is what I'm doing for reading the files, and the spatial
> extent is correct and the precipitation values also seems to be right.
>
> However, my question was more "conceptual", in the sense that the
> provider says that the original file is NOT rotated, while the only
> way  I have to get the correct files is rotating the original file
> with the 't' command you mentioned.
>
> In addition, while reading the file with the latest version of GRASS
> GIS (7.0.3), and I got the following error message:
>
> "
> (Tue Feb  2 19:41:29 2016)
> r.in.gdal input=/home/hzambran/PERSIANN-CDR_v01r01_20090101_c20140523.nc
> output=PERSIANN_CDR_v01r01_20090101_c20140523 -e
> Warning 1: dimension #2 (lat) is not a Longitude/X
> dimension.
> Warning 1: dimension #1 (lon) is not a Latitude/Y dimension.
> ERROR: Input raster map is flipped or rotated - cannot import. You may
> use 'gdalwarp' to transform the map to North-up.
> "
>
>
> So, I'm quite sure that the original file IS rotated, but I didn't
> have more technical arguments to give to the data provider to
> demonstrate that the file is rotated, because they insist that the
> rotation is because R is not reading the file in the proper way....
>
> mzb
>
> =====================================
> "In the end, it's not the years in your life that count.
>   It's the life in your years". (Abraham Lincoln)
> =====================================
> Linux user #454569 -- Linux Mint user
>
>>> 60W-60E and 0-360N, while the correct extent should be: 60S, 60N, 0E,
>>> 360E.
>>>
>>> The file can be downloaded from:
>>>
>>> ftp://data.ncdc.noaa.gov/cdr/persiann/files/2009/PERSIANN-CDR_v01r01_20090101_c20140523.nc
>>>
>>>
>>> After contacting the providers of the file, they mentioned that when
>>> reading the file I have to provide the dimension information, i.e.,
>>> 480 rows and 1440 columns, while the raster package detects the other
>>> way around.
>>>
>>> I would highly appreciate if you could tell me:
>>>
>>> 1) is the previous situation a bug of my netCDF driver (netCDF 4.1.3,
>>> GDAL 1.11.2, released 2015/02/10) or a problem of the dimensions
>>> actually recorded in the netCDF file ?
>>>
>>> 2) is there any way to  overwrite the dimensions detected by netcdf4
>>> when reading the file with the 'raster' command ?
>>>
>>>
>>> My sessionInfo():
>>> R version 3.2.3 (2015-12-10)
>>> Platform: x86_64-pc-linux-gnu (64-bit)
>>> Running under: Ubuntu 14.04.2 LTS
>>>
>>> locale:
>>>   [1] LC_CTYPE=en_GB.UTF-8       LC_NUMERIC=C
>>>   [3] LC_TIME=en_GB.UTF-8        LC_COLLATE=en_GB.UTF-8
>>>   [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_GB.UTF-8
>>>   [7] LC_PAPER=en_US.UTF-8       LC_NAME=C
>>>   [9] LC_ADDRESS=C               LC_TELEPHONE=C
>>> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
>>>
>>> attached base packages:
>>> [1] stats     graphics  grDevices utils     datasets  methods   base
>>>
>>> other attached packages:
>>> [1] raster_2.5-2 sp_1.2-1
>>>
>>> loaded via a namespace (and not attached):
>>> [1] tools_3.2.3     Rcpp_0.12.3     ncdf4_1.15      grid_3.2.3
>>> [5] lattice_0.20-33
>>>
>>>
>>> Thanks in advance and kind regards
>>>
>>> Mauricio Zambrano-Bigiarini, PhD
>>>
>>> =====================================
>>> Dept. of Civil Engineering
>>> Faculty of Engineering and Sciences
>>> Universidad de La Frontera
>>> PO Box 54-D, Temuco, Chile
>>> http://ingenieriacivil.ufro.cl/
>>> =====================================
>>> mailto     : mauricio.zambrano at ufrontera.cl
>>> work-phone : +56 45 259 2812
>>> =====================================
>>> "In the end, it's not the years in your life that count.
>>>   It's the life in your years". (Abraham Lincoln)
>>> =====================================
>>> Linux user #454569 -- Linux Mint user
>>>
>>> --
>>> La información contenida en este correo electrónico y cualquier anexo o
>>> respuesta relacionada, puede contener datos e información confidencial y
>>> no
>>> puede ser usada o difundida por personas distintas a su(s)
>>> destinatario(s).
>>> Si usted no es el destinatario de esta comunicación, le informamos que
>>> cualquier divulgación, distribución o copia de esta información constituye
>>> un delito conforme a la ley chilena. Si lo ha recibido por error, por
>>> favor
>>> borre el mensaje y todos sus anexos y notifique al remitente.
>>>
>>> _______________________________________________
>>> R-sig-Geo mailing list
>>> R-sig-Geo at r-project.org
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>> --
>> Dr. Michael Sumner
>> Software and Database Engineer
>> Australian Antarctic Division
>> 203 Channel Highway
>> Kingston Tasmania 7050 Australia
>>



More information about the R-sig-Geo mailing list