[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