Автоматизируйте "ncvar_get"-читание различных имен переменных в файлах NetCDF

Question

У меня есть набор данных NetCDF с двумя климатическими сценариями (rcp и hist), каждый из которых содержит 25 файлов. Каждый файл содержит данные для переменной "pr", "tas", "tasmax" или "tasmin". Я написал цикл for для итеративного чтения файлов hist и rcp, чтения их с помощью nc_open, извлечения переменной с помощью ncvar_get и, наконец, вычисления в форме среднего (abs(hist - rcp) для получения среднего абсолютного расстояния между каждой парой. of листа и rcp. Проблема: поскольку ncvar_get требует точного имени переменной текущего файла, я написал блок if else (см. ниже), который должен найти имя переменной текущего файла и применить его для ncvar_get. Запуск кода, который я получаю следующая ошибка:

[1] "vobjtovarid4: error #F: I could not find the requsted var (or dimvar) in the file!"

[1] "var (or dimvar) name: tas"

[1] "file name: /data/historical/tasmax_ICHEC-EC-EARTH_DMI-HIRHAM5_r3i1p1.nc" Error in vobjtovarid4(nc, varid, verbose = verbose, allowdimvar = TRUE) : Variable not found



 #Extract of the files in the hist list. Same file names in the rcp list, but different directory

    > hist.files.cl <- list.files("/historical", full.names = TRUE)

    > hist.files.cl

 [1] "/historical/pr_CNRM-CERFACS-CNRM-CM5_ALADIN53_r1i1p1.nc"           
 [2] "/historical/pr_CNRM-CERFACS-CNRM-CM5_ALARO-0_r1i1p1.nc"            
 [3] "/historical/pr_ICHEC-EC-EARTH_HIRHAM5_r3i1p1.nc"                   
 [4] "/historical/pr_ICHEC-EC-EARTH_RACMO22E_r12i1p1.nc"                 
 [5] "/historical/pr_ICHEC-EC-EARTH_RCA4_r12i1p1.nc"                     
 [6] "/historical/pr_MPI-M-MPI-ESM-LR_RCA4_r1i1p1.nc"                    
 [7] "/historical/pr_MPI-M-MPI-ESM-LR_REMO2009_r1i1p1.nc"                
 [8] "/historical/pr_MPI-M-MPI-ESM-LR_REMO2009_r2i1p1.nc"                
 [9] "/historical/tas_CNRM-CERFACS-CNRM-CM5_CNRM-ALADIN53_r1i1p1.nc"     
[10] "/historical/tas_CNRM-CERFACS-CNRM-CM5_RMIB-UGent-ALARO-0_r1i1p1.nc"
[11] "/historical/tas_ICHEC-EC-EARTH_DMI-HIRHAM5_r3i1p1.nc"              
[12] "/historical/tas_ICHEC-EC-EARTH_KNMI-RACMO22E_r12i1p1.nc"           
[13] "/historical/tas_ICHEC-EC-EARTH_SMHI-RCA4_r12i1p1.nc"               
[14] "/historical/tas_MPI-M-MPI-ESM-LR_MPI-CSC-REMO2009_r1i1p1.nc"       
[15] "/historical/tas_MPI-M-MPI-ESM-LR_MPI-CSC-REMO2009_r2i1p1.nc"       
[16] "/historical/tasmax_ICHEC-EC-EARTH_DMI-HIRHAM5_r3i1p1.nc"           
[17] "/historical/tasmax_ICHEC-EC-EARTH_KNMI-RACMO22E_r12i1p1.nc"        
[18] "/historical/tasmax_ICHEC-EC-EARTH_SMHI-RCA4_r12i1p1.nc"        


euc.distance <- list()


for(i in 1:length(hist.files.cl)) {

#Open ith file in list of hist files as well as in list of rcp files

  hist.data <- nc_open(hist.files.cl[i])   
  rcp.data <- nc_open(rcp.files.cl[i])

  if(grepl("pr", hist.data$filename)){
    hist.var <- ncvar_get(hist.data, "pr")
    rcp.var <- ncvar_get(rcp.data, "pr")
    }else if (grepl("tas", hist.data$filename)){
    hist.var <- ncvar_get(hist.data, "tas")
    rcp.var <- ncvar_get(rcp.data, "tas")
    }else if (grepl("tasmax", hist.data$filename)){
    hist.var <- ncvar_get(hist.data, "tasmax")
    rcp.var <- ncvar_get(rcp.data, "tasmax")
    }else{
    hist.var <- ncvar_get(hist.data, "tasmin")
    rcp.var <- ncvar_get(rcp.data, "tasmin")
    }
 #Converting temperature variable from K to °C: 

  if(grepl("tas", hist.data$filename)){

    hist.var <- hist.var-273.15
    rcp.var <- rcp.var-273.15
  }

 #Find for the ith rcp file with dim=(1,1,360) in the ith hist file with dim=(385,373,360) the grid point with the best fitting distribution (each grid point consists of a distribution of 360 time steps).The calculation may contain errors...

  euc.distance[[i]] <- apply(hist.var, c(1,2), function(x) mean(abs(rcp.var - x)))
  min_values <-  which(rank(euc.distance[[i]], ties.method='min') <= 10)
}

Как отметил cath, возможная причина ошибки, но предложенный подход для извлечения интересующей части (= имя переменной) из имени файла не работает. Раньше я пытался автоматизировать извлечение имени переменной с помощью stringr("имя файла", начальная позиция, конечная позиция), пока не заметил, что в этом нет никакого смысла, потому что каждое имя переменной (pr, tas, tasmax, tasmin) имеет свою строку длина. Какие возможности вы видите для меня? Спасибо большое!

Answer 1

Чтобы завершить мой комментарий, если вам нужно работать с каждым файлом, вы можете сделать это сразу, поместив все в список.

Итак, сначала получите "ключевую часть" для каждого файла:

keyparts <- sub("^([a-z]+)_.+", "\\1", basename(hist.files.cl))
keyparts
# [1] "pr"     "pr"     "pr"     "pr"     "pr"     "pr"     "pr"     "pr"    
# [9] "tas"    "tas"    "tas"    "tas"    "tas"    "tas"    "tas"    "tasmax"
#[17] "tasmax" "tasmax"

Тогда вы можете использовать lapply сделать то, что вам нужно сделать для каждого файла одновременно:

my_res <- lapply(seq(keyparts), 
                 function(i){

         hist.data <- nc_open(hist.files.cl[i])   
         rcp.data <- nc_open(rcp.files.cl[i])

         hist.var <- ncvar_get(hist.data, keyparts[i])
         rcp.var <- ncvar_get(rcp.data, keyparts[i])

         if(keyparts[i]=="tas"){
           hist.var <- hist.var-273.15
           rcp.var <- rcp.var-273.15
         }


        euc.distance <- apply(hist.var, c(1,2), function(x) mean(abs(rcp.var - x)))
        min_values <-  which(rank(euc.distance[[i]], ties.method='min') <= 10)

        return(list(euc.distance=euc.distance, min.values=min.values))

                   })