How to convert a rotated NetCDF back to normal lat

2019-05-01 04:40发布

I have a NetCDF file with rotated coordinates. I need to convert it to normal lat/long coordinates (-180 to 180 for long and -90 to 90 for lat).

library(ncdf4)
nc_open('dat.nf')

For the dimensions, it shows:

[1] "     5 variables (excluding dimension variables):"
[1] "        double time_bnds[bnds,time]   "
[1] "        double lon[rlon,rlat]   "
[1] "            long_name: longitude"
[1] "            units: degrees_east"
[1] "        double lat[rlon,rlat]   "
[1] "            long_name: latitude"
[1] "            units: degrees_north"
[1] "        char rotated_pole[]   "
[1] "            grid_mapping_name: rotated_latitude_longitude"
[1] "            grid_north_pole_longitude: 83"
[1] "            grid_north_pole_latitude: 42.5"
[1] "        float tasmax[rlon,rlat,time]   "
[1] "            long_name: Daily Maximum Near-Surface Air Temperature"
[1] "            standard_name: air_temperature"
[1] "            units: K"
[1] "            cell_methods: time:maximum within days time:mean over days"
[1] "            coordinates: lon lat"
[1] "            grid_mapping: rotated_pole"
[1] "            _FillValue: 1.00000002004088e+20"

[1] "     4 dimensions:"
[1] "        rlon  Size:310"
[1] "            long_name: longitude in rotated pole grid"
[1] "            units: degrees"
[1] "            axis: X"
[1] "            standard_name: grid_longitude"
[1] "        rlat  Size:260"
[1] "            long_name: latitude in rotated pole grid"
[1] "            units: degrees"
[1] "            axis: Y"
[1] "            standard_name: grid_latitude"
[1] "        bnds  Size:2"

Could anyone show me how to convert the rotated coordinates back to normal lat/long? Thanks.

4条回答
Ridiculous、
2楼-- · 2019-05-01 05:16

I went through the CDO link as suggested by @kakk11, but somehow that could not work for me. Afte much research, I found a way

First, convert the rotated grid to curvilinear grid

cdo setgridtype,curvilinear Sin.nc out.nc

Next transform to your desired grid e.g. for global 1X1 degree

cdo remapbil,global_1 out.nc out2.nc

or for a grid like below

gridtype = lonlat

xsize = 320 # replace by your value

ysize = 180 # replace by your value

xfirst = 1 # replace by your value

xinc = 0.0625 # replace by your value

yfirst = 43 # replace by your value

yinc = 0.0625 # replace by your value

save this info as target_grid.txt and then run

cdo remapbil,target_grid.txt out.nc out2.nc

In my case there was additional issue that my variables did not have the grid information. so CDO assumed it to be regular lat-long grid. So before all the above-mentioned steps, I had to add grid information attribute to all the variables (in my cases all the variables ended with _ave) using nco

ncatted -a coordinates,'_ave$',c,c,'lon lat' in.nc
ncatted -a grid_mapping,'_ave$',c,c,'rotated_pole' in.nc

Please note that your should have a variable called rotated_pole in your nc file with the lat long information of rotated pole.

查看更多
Luminary・发光体
3楼-- · 2019-05-01 05:18

NCO's ncks can probably do this in two commands using MSA

ncks -O -H --msa -d Lon,0.,180. -d Lon,-180.,-1.0 in.nc out.nc ncap2 -O -s 'where(Lon < 0) Lon=Lon+360' out.nc out.nc

查看更多
Ridiculous、
4楼-- · 2019-05-01 05:18

There is also the possibility to do that in R (as the User is referring to it in the question). Of course, NCO and CDO are more efficient (way faster). Please, look also at this answer.

library(ncdf4)
library(raster)

nsat<- stack (air_temperature.nc)

##check the extent
extent(nsat)
## this will be in the form 0-360 degrees

#change the coordinates
nsat1<-rotate(nstat)

#check result:
extent(nsat1)
##this should be in the format you are looking for: -180/180

Hope this helps.

查看更多
男人必须洒脱
5楼-- · 2019-05-01 05:28

I would use cdo for this purpose https://code.zmaw.de/boards/2/topics/102

Another option is just create a mapping between rotated and geographic coordinates and use the original data without interpolation. I can find the equations if necessary.

查看更多
登录 后发表回答