Как извлечь подмножество из файла netCDF, используя границы широты / долготы в R

У меня есть файл netCDF, из которого я хочу извлечь подмножество, определяемое границами широты / долготы (т. Е. Поле широты / долготы), используя пакет 'ncdf' в R.

Резюме моего файла netCDF ниже. Он имеет два измерения (широта и долгота) и 1 переменная (10U_GDS4_SFC). По сути, это сетка широта / долгота, содержащая значения ветра:

[1] "file example.nc has 2 dimensions:"
[1] "lat_0   Size: 1280"
[1] "lon_1   Size: 2560"
[1] "------------------------"
[1] "file example.nc has 1 variables:"
[1] "float 10U_GDS4_SFC[lon_1,lat_0]  Longname:10 metre U wind component Missval:1e+30"

Переменная широты варьируется от +90 до -90, а переменная долготы от 0 до 360.

Я хочу извлечь подмножество общей сетки, используя следующие географические угловые границы:

нижний левый угол: широта: 34,5˚, длина: 355˚, левый верхний угол: широта: 44,5˚, длина: 355˚, правый верхний угол: широта: 44,5˚, длина: 12˚, правый нижний угол: широта: 34,5˚, Длинный: 12˚

Я знаю, что части переменной могут быть извлечены с помощью get.var.ncdf() команда (пример ниже):

z1 = get.var.ncdf(example.nc, "10U_GDS4_SFC", start=c(11,26), count=c(5,5))

Тем не менее, я не могу понять, как можно включить широту / долготу, чтобы в итоге я получил заданную пространственную сетку, содержащую значения переменных. Я новичок в работе со значениями netCDF в R, и любые советы будут с благодарностью. Большое спасибо!

2 ответа

Решение

В принципе вы 2/3 пути туда. Конечно, вы можете создать начальные индексы, используя что-то вроде этого:

require(ncdf4)

ncFile <- nc_open( MyNetCDF )
LonStartIdx <- which( ncFile$dim$lon$vals == 355)
LatStartIdx <- which( ncFile$dim$lat$vals == 34.5)

Сделайте то же самое для подсчета. Затем прочитайте переменную, которую вы хотите

MyVariable <- ncvar_get( ncFile, varName, start=c( LonStartIdx, LatStartIdx), count=...)

Однако в твоем случае тебе не повезло, насколько я знаю. Процедуры чтения / записи netcdf делают свое дело последовательно. Ваша сетка оборачивается, так как у вас есть координаты, которые идут от 0 до 360 по долготе, и вы заинтересованы в поле, которое содержит нулевой меридиан.

Для вас (при условии, что у вас не слишком много данных) было бы более целесообразно прочитать всю таблицу в R, а затем использовать либо subset или создать индексы, используя which и вырезать свою "коробку" в R.

ncFile <- nc_open( MyNetCDF )
LonIdx <- which( ncFile$dim$lon$vals > 355 | ncFile$dim$lon$vals < 10)
LatIdx <- which( ncFile$dim$lat$vals > 34.5 & ncFile$dim$lat$vals < 44.5)
MyVariable <- ncvar_get( ncFile, varName)[ LonIdx, LatIdx]
nc_close(ncFile)

Замечание: я предпочитаю ncdf4Я нахожу синтаксис немного легче запомнить (и было еще одно преимущество перед старым R-пакетом netcdf, которое я забыл...)

Хорошо. Комментарии не могут быть до тех пор, пока они мне понадобятся, поэтому я обновил ответ. Не беспокойтесь. Давайте рассмотрим вопросы шаг за шагом.

  • which Функциональный способ будет работать. Я использую это сам.
  • Данные будут в том же формате, что и в файле netcf, но я не слишком уверен, есть ли какая-то проблема с меридианом 0 (думаю, да). Возможно, вам придется поменять местами две половинки, выполнив что-то вроде этого (замените соответствующую строку во втором примере)

    LonIdx <- c(which( ncFile$dim$lon$vals > 355) , which( ncFile$dim$lon$vals < 10) )
    

    Это меняет порядок координатных индексов, так что сначала идет западная часть, а затем восточная.

  • Переформатирование всего в кадр данных 2x3 возможно. Возьмите данные, которые возвращает мой 2-й пример кода (будет матрица, [lon x lat]. Также получите значения координат из

    lon <- ncFile$dim$lon$val[LonIdx]
    

    (или как называется долгота в вашем примере, то же самое для lat). Затем собрать матрицу, используя

    cbind( rep(lat, each=length(lon)), rep(lon,length(lat)), c(myVariable) )
    
  • Координаты, конечно, будут такими же, как в файле netcdf...

Вам необходимо проверить последний cbind, так как я уверен только на 98%, что я не перепутал координаты. В сценариях R, которые я нашел на своем рабочем столе, я использую циклы, которые являются... злыми... Это должно быть (немного?) Быстрее, а также более разумно.

Вы также можете использовать CDO, чтобы сначала извлечь область из командной строки bash, а затем прочитать файл в R:

cdo sellonlatbox,-5,12,34.5,44.5 in.nc out.nc 

Я отмечаю в вышеупомянутой дискуссии, что была проблема, касающаяся порядка широт. Вы также можете использовать команду CDO "invertlat", чтобы разобраться в этом за вас.

Если вы используете Linux, этого легко добиться с помощью nctoolkit ( https://nctoolkit.readthedocs.io/en/latest/):

import nctoolkit as nc
data = nc.open_data("example.nc")    
data.clip(lon = [-12, -5], lat = [35.4, 44.5])
Другие вопросы по тегам