bogus interpretation by raster of north pole chl-a grid

[mdsumner]

Four variables, chla, longitude, latitude, (the 2D map with a bio-physical and a lon and a lat value for every cell), and stereographic (a placeholder for metadata coments in the file).

Details

https://stat.ethz.ch/pipermail/r-sig-geo/2017-January/025256.html

library(ncdf4)
ncvar_get(nc_open("subset.nc"), "stereographic")
## [1] -2147483647  ## a classic boilerplate value, the largest 32-bit integer:

.Machine$integer.max

fourcornervalues <- function(x) x[cbind(c(1, 1, nrow(x), nrow(x)), c(1, ncol(x), ncol(x), 1))]

fourcornervalues(ncvar_get(nc_open("subset.nc"), "longitude"))
[1] -50.000645 -91.041626  39.924549  -1.468801
 fourcornervalues(ncvar_get(nc_open("subset.nc"), "latitude"))
[1] 55.05590 85.63445 56.70048 43.22274

r <- raster("subset.nc")  ## assumes varname = chla

print(r)

extent(r)
class       : Extent 
xmin        : -3.5625 
xmax        : 38.0625 
ymin        : -40.0625 
ymax        : -3.3125 


##Can you see the connection?   Me neither. 

## get the raw longitude latitude values
xy <- cbind(values(raster("subset.nc", varname = "longitude")), 
                 values(raster("subset.nc", varname = "latitude")))


## from the grid_mapping I'd guess (note that datum is not specified in the file): 
prj <- "+proj=stere +lat_0=90 +lon_0=-45 +ellps=WGS84 +no_defs"

## looks close, right?

plot(project(xy, prj), pch = ".")

## guess again
prj <- "+proj=stere +lat_0=90 +lon_0=-45 +ellps=sphere +no_defs"

• where is the file from, can you share it?

• what R package/s are you using for NetCDF? Don't know? See the next item

• enter sessionInfo() and post it in these back ticks to provide details on your R session

## enter your sessionInfo() print out here

[niknikos]

Hi Mike,
For the record I am not a github expert I just have some minor working knowledge :)

• The file is from http://marine.copernicus.eu/ where you need to create an account to be able to download their products. However to see the documentation you shouldn't have any problem following this link
  http://marine.copernicus.eu/services-portfolio/access-to-products/?option=com_csw&view=details&product_id=ARCTIC_ANALYSIS_FORECAST_BIO_002_004
  I have downloaded several .nc files for the last 5 years and each one consists of bands representing days. I can share the file(s) but each one is approx. 90Mb that's why I provided only a subset (i.e. chla for one date)... Here is the link to one of them in case it comes handy
  https://dl.dropboxusercontent.com/u/2156982/BIOGEO2012.nc

• ncdf4 is the package I used to access the .nc just to explore it the first time, but then I utilized the raster package to be able to work with brick in the same way you used it above 'r <- raster("subset.nc") ## assumes varname = chla'

• The sessionInfo output in case you still need it

R version 3.3.1 (2016-06-21)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows >= 8 x64 (build 9200)

locale:
[1] LC_COLLATE=English_United States.1252 LC_CTYPE=English_United States.1252 LC_MONETARY=English_United States.1252
[4] LC_NUMERIC=C LC_TIME=English_United States.1252

attached base packages:
[1] stats graphics grDevices utils datasets methods base

other attached packages:
[1] proj4_1.0-8 rgdal_1.1-10 raster_2.5-8 sp_1.2-3 RevoUtilsMath_8.0.3

loaded via a namespace (and not attached):
[1] RevoUtils_10.0.1 tools_3.3.1 Rcpp_0.12.7 ncdf4_1.15 grid_3.3.1 lattice_0.20-34

[mdsumner]

thanks this is helpful, I'll try to get to it soonish-ish :)