0

I have a .txt file containing a number of stations accompanied with their coordinates. So, the .txt file has three columns: One with the station id, one with lat and one with lon, as shown:

station , lat , lon 
ABTR2100 ,39.13,34.52
GRMR0100 ,20.18,49.00
DDDD0100 ,23.22,46.81
SLPT0100 ,26.91,32.23
NDRT0100 ,29.55,48.97

Also, I have an .nc file containing hourly temperature data for March of 2011. The structure of the .nc file is the following:

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 var11[rlon,rlat,height,time]   
        table: 2
        grid_mapping: rotated_pole
 4 dimensions:
    rlon  Size:848
        standard_name: grid_longitude
        long_name: longitude in rotated pole grid
        units: degrees
        axis: X
    rlat  Size:824
        standard_name: grid_latitude
        long_name: latitude in rotated pole grid
        units: degrees
        axis: Y
    height  Size:1
        standard_name: height
        long_name: height
        units: m
        positive: up
        axis: Z
    time  Size:744   *** is unlimited ***
        standard_name: time
        units: hours since 2011-02-28 18:00:00
        calendar: proleptic_gregorian

What I wish to do is extract the time series for each one of the stations contained in the .txt file for each time step, and store the time series to a separate .txt file that will have the name of the station id. So, this means that I will have a .txt file named "ABTR2100_temp" that will contain the time series for that specific location. The same will happen with all stations. Ηοw could I do this in R?

Tatus
  • 141
  • 1
  • 7
  • 1
    You could try `ncdf` package for reading in the *.nc file properly. – statespace Dec 02 '16 at 14:21
  • Paul Hiemstra, I have tried to apply the extraction using the raster package in R, reading the .nc file as a brick, and then use the "extract" command to exract my time series for my points, but this solution doesn't give what I want. I want to store the time series to separate files for each location and name the file with the station's id. – Tatus Dec 02 '16 at 14:22
  • A.Val, I have no problem reading the .nc file correctly, but to extract the time series and saving them to separate files, assigning to the file's name the station's id name! – Tatus Dec 02 '16 at 14:32
  • If you could provide a small sample (using the dput function) of the time series data, it would be easier to understand your problem. Initially I would look at using the dplyr package and performing a "join" between the 2 datasets. Then use a `for` loop a to subset and save. – Dave2e Dec 02 '16 at 16:45

1 Answers1

1

Read your station coordinates as a data frame. Suppress station names conversion to factor disable stringsAsFactor argument.

yourstations <- read.csv("yourstations.txt", stringsAsFactor = FALSE)

Load your data as a connection into your workspace with the help of ncdf4 package. Cheque everything is OK?

library(ncdf4)
nc.filename <- "yourdata.nc"
yourconn.nc <- nc_open(nc.filename)
yourconn.nc

Get the value of grid latitude and longitude.

lon <- ncvar_get(yourconn.nc,"rlon")
lat <- ncvar_get(yourconn.nc,"rlat")

Use a for loop for every row of data frame of stations. With the for loop find the closest latitude and longitude and based on the row-column number get all the values for one grid point (sta.temp). Write with the station name.

for(ttnum in 1:nrow(yourstations)){
    ## Find the closes lat    
    lat.diff <- abs(lat - yourstations[ttnum, "lat"])
    lat.nr <- which(lat.diff == min(lat.diff))[1]
    ## Find the closes lon
    lon.diff <- abs(lon - yourstations[ttnum, "lon"])
    lon.nr <- which(lon.diff == min(lon.diff))[1]
    ## Copy data of the gridpoint
    sta.temp <- ncvar_get(yourconn.nc, "var11",  c(lon.nr,lat.nr,1,1),c(1,1,1,744))
    ## Write the vector
    write(sta.temp, file = paste0(yourstations[ttnum, "station"],"_temp.txt"), ncol = 1)
}

I am not sure that you can work smoothly with your rotated coordinates (see the question Extracting site-specific information from NetCDF file in R where netcdf4 package used through raster).

Community
  • 1
  • 1
kaliczp
  • 457
  • 1
  • 12
  • 18