I am trying to create a dataset that shows future climatic condition in Europe aggregates to NUT2 administration level in R. To do so, i work with cordex data (domain EU 11), that is coming withrotated coordinates.
I used the data hurs_EUR-11_ICHEC-EC-EARTH_rcp85_r1i1p1_SMHI-RCA4_v1_day_20310101-20351231.nc (data here ) which is the relative humidity Its looks as follows :
6 variables (excluding dimension variables): double height axis: Z long_name: height positive: up standard_name: height units: m double time_bnds[bnds,time] char rotated_pole grid_mapping_name: rotated_latitude_longitude grid_north_pole_latitude: 39.25 grid_north_pole_longitude: -162 float hurs[rlon,rlat,time] grid_mapping: rotated_pole _FillValue: 1.00000002004088e+20 missing_value: 1.00000002004088e+20 standard_name: relative_humidity long_name: Near-Surface Relative Humidity units: % coordinates: lon lat height cell_methods: time: mean double lon[rlon,rlat] standard_name: longitude long_name: longitude units: degrees_east double lat[rlon,rlat] standard_name: latitude long_name: latitude units: degrees_north
4 dimensions: time Size:1826 *** is unlimited *** standard_name: time units: days since 1949-12-01 00:00:00 calendar: standard long_name: time bounds: time_bnds axis: T bnds Size:2
2 "vobjtovarid4: **** WARNING **** I was asked to get a varid for dimension named bnds BUT this dimension HAS NO DIMVAR! Code will probably fail at this point" rlon Size:424 standard_name: grid_longitude long_name: longitude in rotated pole grid units: degrees axis: X rlat Size:412 standard_name: grid_latitude long_name: latitude in rotated pole grid units: degrees axis: Y
First, i created an array
hurs.array <- ncvar_get(data2, "hurs")
Then i create a grid with the projection as suggested in this post Choosing CRS / PROJ4 string for EURO-Cordex rotated pole projection? Note that the issue raised about this post in about an simple feature and not about raster.
r_brick <- brick(hurs.array, xmn=min(lat), xmx=max(lat), ymn=min(lon), ymx=max(lon), crs=CRS("+proj=ob_tran +o_proj=longlat +o_lon_p=-162.0 +o_lat_p=39.25 +lon_0=180 +lat_0=0"))
Now my question is how to transform this data into the EU projection? with
r_brickP<- projectRaster(r_brick, crs="+init=epsg:3035")
not finiteno non-missing arguments to min;
What is going on?
Is this the right approach?
Is there a better approach to unrotate the data?
If yes which one?
What is my most efficient way to then perform a zonal statics with as layer with epsg:3035 projection?
I also tried the function rotate, but i got a line suggesting that this is not this method was not appropriate.