Extracting site-specific information from NetCDF f

2019-09-13 18:22发布

问题:

I got a NetCDF file from the German Meteorological Service concerning mean temperatures in Europe (CDC FDP SERVER). The only thing I want to extract is the daily mean temperature for Bornholm, which is an island in the central Baltic.

I know how to extract information for certain coordinates (see code sample below). The only problem is that the file specific coordinates are 'rotated' which is why the geographic coordinates for Bornholm (extracted from GoogleMaps) are kind of useless.

packages <- c("RNetCDF",
              "ncdf4",
              "raster")

lapply(packages, require, character.only = TRUE)

x <- mean(14.68,15.16)        #coordinates for a rectangle around 
y <- mean(54.987,55.299)      #Bornholm extracted from GoogleMaps

temp <- nc_open("tas_decreg_europe_v20140120_20030101_20030131.nc")
temp

var <- ncvar_get(temp, "tas")
point <- var[x,y,]
as.data.frame(point)

To cut it short - Google uses a close variant of the Mercator projection. So how can I either convert the NetCDF file or the coordinates from GoogleMaps, so that I can find what I need. I could have bet that there's a simple solution out there but unfortunately not - at least I couldn't find one.

For information about the file generated by print(temp) see below:

File tas_decreg_europe_v20140120_20030101_20030131.nc (NC_FORMAT_CLASSIC):

     2 variables (excluding dimension variables):
        char rotated_pole[]   
            grid_mapping_name: rotated_latitude_longitude
            grid_north_pole_latitude: 39.25
            grid_north_pole_longitude: -162
        float tas[lon,lat,time]   
            long_name: Near-Surface Air Temperature
            units: K
            grid_mapping: rotated_pole
            _FillValue: 1.00000002004088e+20
            missing_value: 1.00000002004088e+20

     3 dimensions:
        lon  Size:1056
            standard_name: grid_longitude
            long_name: longitude
            units: degrees_east
            axis: X
        lat  Size:1026
            standard_name: grid_latitude
            long_name: latitude
            units: degrees_north
            axis: Y
        time  Size:31   *** is unlimited ***
            standard_name: time
            units: days since 2003-01-01 00:00:00
            calendar: standard

Any help is appreciated. Thanks a lot...

回答1:

You load the raster package, but you do not use it. Have you tried something like the below?

library(raster)  
x <- mean(14.68,15.16)  
y <- mean(54.987,55.299)
temp <- brick("tas_decreg_europe_v20140120_20030101_20030131.nc", var='tas')
extract(temp, cbind(x,y))


标签: r netcdf netcdf4