[R-sig-Geo] ncdf 4-dim file to raster brick - level error??

Robert J. Hijmans r.hijmans at gmail.com
Wed Oct 14 04:53:13 CEST 2015


John
Thank you for reporting this. This has now been fixed in raster-devel
install.packages("raster", repos="http://R-Forge.R-project.org")  and
soon on CRAN.
The error (always using the first level on extraction of cell values)
occurred in ncdf files with "levels" as fourth dimension (rather than
time).
Robert

On Sun, Oct 4, 2015 at 12:48 PM, Gross, John <john_gross at nps.gov> wrote:
> I'm trying to read a 4-dimensional ncdf file (dims: lat lon time
> projection) into a raster brick and can't get it to read anything
> other than the first "level" (i.e. projection, the 4th dimension).
>
> Here's the code, output, and sessionInfo.  This is reminiscent of a
> problem a couple years ago with earlier versions of both ncdf4 and
> raster that was caused by the ncdf file not closing between calls  But
> I'm not sure how to determine it that's happening again.
>
> You can download the data file (10KB - it's only 3x3x59x4) from:
> https://www.dropbox.com/s/5jh5n8etpk1y1hu/Extraction_tasmin.nc?dl=0
>
> #===  Code ==========
> library(ncdf4)
> library(raster)
>
> nc <- nc_open("Extraction_tasmin.nc")    # this is a 3x3 text grid
> nc    # dims are lat long time and projection
>
> xy <- data.frame(x = min(nc$dim$lon$vals) + 1/32, y =
> max(nc$dim$lat$vals) - 1/32) #get top left point not on grid boundary
> coordinates(xy) <- ~x+y
>
> for(i in 1:5){
>   b1 <- brick("Extraction_tasmin.nc", lvar=4, level=i)
>   pts <- t(extract(b1, xy))
>   print(c("pts[1:5} -> ",pts[1:5])); flush.console()
> }
> #====  Output  =======
>
> File Extraction_tasmin.nc (NC_FORMAT_CLASSIC):
>
>      1 variables (excluding dimension variables):
>         float tasmin[longitude,latitude,time,projection]
>             typeConversion_op_ncl: double converted to float
>             _FillValue: 1.00000002004088e+20
>
>      4 dimensions:
>         latitude  Size:3
>             valid_max: 52.875
>             long_name: Latitude
>             valid_min: 25.125
>             units: degrees_north
>             axis: Y
>         longitude  Size:3
>             long_name: Longitude
>             axis: X
>             units: degrees_east
>             modulo: 360
>             topology: circular
>         time  Size:59
>             calendar: standard
>             units: days since 1950-01-01 00:00:00
>             standard_name: time
>             long_name: time
>             axis: T
>         projection  Size:4   *** is unlimited ***
>
>     2 global attributes:
>         NCO: 4.0.8
>         Projections: access1-0.1.rcp85, bcc-csm1-1.1.rcp85,
> canesm2.1.rcp85, ccsm4.1.rcp85,
>>
>
> Output:
> [1] "pts[1:5} -> "      "-5.08642387390137" "-3.23856282234192"
> "-5.46992635726929" "-3.83646535873413" "-3.00846982002258"
> [1] "pts[1:5} -> "      "-5.08642387390137" "-3.23856282234192"
> "-5.46992635726929" "-3.83646535873413" "-3.00846982002258"
> [1] "pts[1:5} -> "      "-5.08642387390137" "-3.23856282234192"
> "-5.46992635726929" "-3.83646535873413" "-3.00846982002258"
> [1] "pts[1:5} -> "      "-5.08642387390137" "-3.23856282234192"
> "-5.46992635726929" "-3.83646535873413" "-3.00846982002258"
> [1] "pts[1:5} -> "      "-5.08642387390137" "-3.23856282234192"
> "-5.46992635726929" "-3.83646535873413" "-3.00846982002258"
> Warning message:
> In .rasterObjectFromCDF(x, type = objecttype, band = band, ...) :
>   level set to: 4
> # ==== sessionInfo =========
> R version 3.2.2 (2015-08-14)
> Platform: x86_64-apple-darwin13.4.0 (64-bit)
> Running under: OS X 10.10.5 (Yosemite)
>
> locale:
> [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
>
> attached base packages:
> [1] stats     graphics  grDevices utils     datasets  methods   base
>
> other attached packages:
> [1] ncdf4_1.13    raster_2.4-18 sp_1.1-1
>
> loaded via a namespace (and not attached):
> [1] rgdal_1.0-4     tools_3.2.2     Rcpp_0.12.1     grid_3.2.2
> lattice_0.20-33
>
> --
> John Gross, PhD
> NPS, Ft Collins, CO
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo



More information about the R-sig-Geo mailing list