Extract a series of point from NETCDF files and generate output in CSV format

This is a function in R helping retrieve a series of points in a NETCDF files and generate a well-structured CSV format on the output side.

rm(list=ls())

#==Load packages======

library("dplyr")
library("tibble")
library('ncdf4')
library("ncdf4.helpers")

getwd()

## [1] "C:/Users/De-RoOt/Documents/Blog/N1/CRU_data"

For this exercise, I am going to use the High-resolution gridded datasets - Climatic Research Unit and retrieve some stations values in different climatic zones.

The latest version of CRU (Ts 4.03) provide gridded data from Jan. 1901 to Dec. 2018 which can be downloaded from this link.

#=====Get data path
(data<-list.files(pattern="*.nc"))

## [1] "CRU_WA_pre.nc"   "CRU_WA_tmax.nc"  "CRU_WA_tmean.nc"

The stations data are arranged as follows

#===stations informations =======
zone_points<- as_tibble(read.table("Climatic_Zone_Points.csv", header=TRUE,sep=","))

zone_points$Stations<-as.character(zone_points$Stations)#Format colunm to character
zone_points$Climatic.Zones<-as.character(zone_points$Climatic.Zones)#Format colunm to character

#print
zone_points

## # A tibble: 3 x 4
##   Stations Climatic.Zones  Longitude Latitude
##   <chr>    <chr>               <dbl>    <dbl>
## 1 Cotonou  Guinean              2.31     6.40
## 2 SAVE     Soudano-Guinean      2.66     8.18
## 3 Kandi    Soudanien            3.22    11.0

Extraction Function

#Function of extraction

Extraction <-function(data_path,grid_info){
  
  nfiles<-length(data_path) #number of files in the path
   
  # Repeat operation for each files in the data path
  for (nfil in 1:nfiles){
    
    modExt1<-NULL #final ouptut variable
    print(data_path[nfil])
    
    nc_brick<-NULL
    fileName<-basename(data_path[nfil])# File name without the path
    
    #Variable name
    spl<-unlist(strsplit(fileName, "[.]"))#split the filename in different compoenents
    varlabel<-paste(spl[1])
    
    #Netcdf  multi-layer raster object 
    nc_brick<-brick(data_path[nfil]) 
    
    dataget<-nc_open(data_path[nfil],write=FALSE)
    
    #time axis of dataget
    ts<-nc.get.time.series(dataget)
    
    #Get a list of names of data variables
    var<-nc.get.variable.list(dataget)
    
    # Retrieve, day, month and year from  ts
    day<-format(ts, "%d")
    year<-format(ts, "%Y")
    month<-format(ts,"%m")
    date<-format(ts,"%Y-%m-%d")
    
    #New dataset with   date,year, month, day colunms
    times<-as.data.frame(cbind(date,year,month,day))
    
    #length of the time series
    tslen<-length(ts)
    
    
    #=Extract variable values for each station and given lon and lat===========  
    
    points<-length(unique(grid_info$Stations)) # Number of grid point to extract
    
    for (p in 1:points){
      
      print(grid_info$Stations[p])
      consol_tmp<-NULL
      grid_val<-NULL
      
      #Extract values
      grid_val<-extract(nc_brick,matrix(c(grid_info$Longitude[p],grid_info$Latitude[p]), ncol = 2))
      grid_val<-as.data.frame(t(grid_val), row.names = F)
      
      #descriptive colunms (Stations, Climate.Zones,Longitude,Longitude)
      info_val<-cbind(rep(grid_info$Stations[p],tslen),rep(grid_info$Climatic.Zones[p],tslen),
                      rep(grid_info$Longitude[p],tslen),rep(grid_info$Latitude[p],tslen))
    
      #Consolidated file output
      
      consol_tmp<-cbind(times,info_val,grid_val)
      
      #Model output file
      modExt1<-rbind(modExt1,consol_tmp)
      #rm(grid_val,consol_tmp)
    }
    
    #Rename colnames
    colnames(modExt1)<-c('Date','Year','Month', 'Day', 'Stations','Climatic_Zones','Longitude','Latitude',varlabel)
    
    modExt1<-as_tibble(modExt1)
    
    #Output file name
    filelabel<-paste(spl[1],'_daily','.csv',sep='')
    
    #Output file in csv file
    write.csv(modExt1,filelabel)
    
    print (paste('Output file ',nfil,filelabel))
    
    print ("======Content printing===========================")
    print (modExt1[,-1])
    print ("======End printing===============================")
  }

}

Operate the function Extraction on the netcdf files in the directory CRU Data

Extraction(data,zone_points)

## [1] "CRU_WA_pre.nc"
## [1] "Cotonou"
## [1] "SAVE"
## [1] "Kandi"
## [1] "Output file  1 CRU_WA_pre_daily.csv"
## [1] "======Content printing==========================="
## # A tibble: 4,176 x 8
##    Year  Month Day   Stations Climatic_Zones Longitude Latitude CRU_WA_pre
##    <fct> <fct> <fct> <fct>    <fct>          <fct>     <fct>         <dbl>
##  1 1901  01    16    Cotonou  Guinean        2.30763   6.39958        14.2
##  2 1901  02    15    Cotonou  Guinean        2.30763   6.39958        18.4
##  3 1901  03    16    Cotonou  Guinean        2.30763   6.39958       155. 
##  4 1901  04    16    Cotonou  Guinean        2.30763   6.39958        73.6
##  5 1901  05    16    Cotonou  Guinean        2.30763   6.39958       185. 
##  6 1901  06    16    Cotonou  Guinean        2.30763   6.39958       180. 
##  7 1901  07    16    Cotonou  Guinean        2.30763   6.39958       243. 
##  8 1901  08    16    Cotonou  Guinean        2.30763   6.39958       132. 
##  9 1901  09    16    Cotonou  Guinean        2.30763   6.39958       368  
## 10 1901  10    16    Cotonou  Guinean        2.30763   6.39958       251. 
## # ... with 4,166 more rows
## [1] "======End printing==============================="
## [1] "CRU_WA_tmax.nc"
## [1] "Cotonou"
## [1] "SAVE"
## [1] "Kandi"
## [1] "Output file  2 CRU_WA_tmax_daily.csv"
## [1] "======Content printing==========================="
## # A tibble: 4,176 x 8
##    Year  Month Day   Stations Climatic_Zones Longitude Latitude CRU_WA_tmax
##    <fct> <fct> <fct> <fct>    <fct>          <fct>     <fct>          <dbl>
##  1 1901  01    16    Cotonou  Guinean        2.30763   6.39958         32.2
##  2 1901  02    15    Cotonou  Guinean        2.30763   6.39958         33.5
##  3 1901  03    16    Cotonou  Guinean        2.30763   6.39958         32.4
##  4 1901  04    16    Cotonou  Guinean        2.30763   6.39958         32.3
##  5 1901  05    16    Cotonou  Guinean        2.30763   6.39958         30.7
##  6 1901  06    16    Cotonou  Guinean        2.30763   6.39958         29.4
##  7 1901  07    16    Cotonou  Guinean        2.30763   6.39958         28.5
##  8 1901  08    16    Cotonou  Guinean        2.30763   6.39958         27.9
##  9 1901  09    16    Cotonou  Guinean        2.30763   6.39958         28.8
## 10 1901  10    16    Cotonou  Guinean        2.30763   6.39958         30.5
## # ... with 4,166 more rows
## [1] "======End printing==============================="
## [1] "CRU_WA_tmean.nc"
## [1] "Cotonou"
## [1] "SAVE"
## [1] "Kandi"
## [1] "Output file  3 CRU_WA_tmean_daily.csv"
## [1] "======Content printing==========================="
## # A tibble: 4,176 x 8
##    Year  Month Day   Stations Climatic_Zones Longitude Latitude
##    <fct> <fct> <fct> <fct>    <fct>          <fct>     <fct>   
##  1 1901  01    16    Cotonou  Guinean        2.30763   6.39958 
##  2 1901  02    15    Cotonou  Guinean        2.30763   6.39958 
##  3 1901  03    16    Cotonou  Guinean        2.30763   6.39958 
##  4 1901  04    16    Cotonou  Guinean        2.30763   6.39958 
##  5 1901  05    16    Cotonou  Guinean        2.30763   6.39958 
##  6 1901  06    16    Cotonou  Guinean        2.30763   6.39958 
##  7 1901  07    16    Cotonou  Guinean        2.30763   6.39958 
##  8 1901  08    16    Cotonou  Guinean        2.30763   6.39958 
##  9 1901  09    16    Cotonou  Guinean        2.30763   6.39958 
## 10 1901  10    16    Cotonou  Guinean        2.30763   6.39958 
## # ... with 4,166 more rows, and 1 more variable: CRU_WA_tmean <dbl>
## [1] "======End printing==============================="

Let me know your comments below!

if you like this post, Don’t forget to share it!