Среднее для другого слоя и нескольких файлов netcdf с R - PullRequest
0 голосов
/ 02 октября 2018

У меня есть 15 файлов netCDF (.nc) для каждого года с 2000 по 2014 год. В одном файле nc у меня есть почасовые данные одной переменной в 8760 слоях.3 измерения: время (размер 8760), широта (размер 90) и долгота (размер 180) (разрешение 2 °).

Я хочу вычислить среднее значение моей переменной в период с 8:00 до 19:00 с апреля посентябрь и за период 2000-2014 гг.

Для одного файла .nc это соответствует среднему значению между

  • временем слоя с 2169 (т. е. 01/04/2000 8 утра) до2180 (т.е. 01/04/2000 7 вечера) (для i = 2169 до i + 11),
  • , затем с 2193 (т.е. 02/04/2000 8:00) до 2204 (т.е. 02/04/2000 19:00)(i + 22, i + 33)
  • и т. д. *
  • ... и с 6537 (т. е. 30/09/2000 8 утра) до 6548 (т. е. 30/09/20007 вечера)
  • А потом среднее по всем нс.файлы.

Результат должен быть представлен в одном файле .nc с 3 измерениями: - время (только одно значение как среднее), - широта (размер 90) и - долгота (размер 180) (2° разрешение)

тогда я могу нарисовать карту переменной, усредненной за 2000-2014 годы (с апреля по сентябрь, с 8:00 до 19:00).Я могу читать каждый файл nc, составлять карту для каждого часа каждого файла nc, но я знаю, как сделать среднее значение необходимым.Если кто-нибудь может мне помочь, это было бы здорово.

имя моей переменной: dname <- "sfvmro3" </p>

Вот мой код для чтения:

ncin <- nc_open("sfvmro3_hourly_2000.nc")
print(ncin)

lon <- ncvar_get(ncin, "lon")
lon[lon > 180] <- lon[lon > 180] - 360
nlon <- dim(lon)
head(lon)

lat <- ncvar_get(ncin, "lat", verbose = F)
nlat <- dim(lat)
head(lat)

print(c(nlon, nlat))

t <- ncvar_get(ncin, "time")
tunits <- ncatt_get(ncin, "time", "units")
nt <- dim(t)

dname <- "sfvmro3"
var.array <- ncvar_get(ncin, dname)*10^9  # from mol.mol-1 to ppb
dlname <- ncatt_get(ncin, dname, "long_name")
dunits <- ncatt_get(ncin, dname, "units")
fillvalue <- ncatt_get(ncin, dname, "_FillValue")
var.array[var.array == fillvalue$value] <- NA
dim(var.array)

tustr <- strsplit(tunits$value, " ")
tdstr <- strsplit(unlist(tustr)[3], "-")
tyear = as.integer(unlist(tdstr)[1])
tmonth = as.integer(unlist(tdstr)[2])
tday = as.integer(unlist(tdstr)[3])
chron = chron(t, origin = c(tmonth, tday, tyear))

Вот подробности на один годовой файл. Nc:

 4 variables (excluding dimension variables):
    double time_bnds[bnds,time]   
    double lat_bnds[bnds,lat]   
    double lon_bnds[bnds,lon]   
    float sfvmro3[lon,lat,time]   
        standard_name: mole_fraction_of_ozone_in_air
        long_name: Ozone Volume Mixing Ratio in the Lowest Model Layer
        units: mole mole-1
        original_name: O_x
        original_units: 1
        history: 2016-04-22T05:20:31Z altered by CMOR: Converted units from '1' to 'mole mole-1'.
        cell_methods: time: point (interval: 30 minutes)
        cell_measures: area: areacella
        missing_value: 1.00000002004088e+20
        _FillValue: 1.00000002004088e+20
        associated_files: ...

 4 dimensions:
    time  Size:8760   *** is unlimited ***
        bounds: time_bnds
        units: days since 1850-01-01
        calendar: noleap
        axis: T
        long_name: time
        standard_name: time
    lat  Size:90
        bounds: lat_bnds
        units: degrees_north
        axis: Y
        long_name: latitude
        standard_name: latitude
    lon  Size:180
        bounds: lon_bnds
        units: degrees_east
        axis: X
        long_name: longitude
        standard_name: longitude
    bnds  Size:2

26 global attributes:
    institution: aaaa
    institute_id: aaaa
    experiment_id: aaaa
    source: aaaa
    model_id: aaaa
    forcing: HG, SA, S
    parent_experiment_id: N/A
    parent_experiment_rip: N/A
    branch_time: 0
    contact: aaa
    history: aaa
    initialization_method: 1
    physics_version: 1
    tracking_id: aaa
    product: output
    experiment: aaa
    frequency: hr
    creation_date: 2016-04-22T05:20:31Z
    Conventions: aaa
    project_id: aaa
    table_id:aaa
    title: aaaa
    parent_experiment: N/A
    modeling_realm: aaa
    realization: 1
    cmor_version: 2.7.1

1 Ответ

0 голосов
/ 16 ноября 2018

Я знаю два разных возможных решения вашей проблемы.Один из них основан на получении среднего значения для каждого файла .nc, а затем - среднем весе этого значения, другой - на получении действительно большого массива и среднего значения с использованием этого массива.

  • Первое возможное решение

Каждый .nc, который вы прочитаете, даст вам и массив, массив1, массив2 и так далее.Также для каждого массива у вас будет временной ряд, связанный с одним измерением массива.Это означает, что time_serie1 имеет все разное время в формате POSIXct для array1.Итак, сначала вы должны встроить этот вектор.Один у вас есть, что вы можете получить векторный индекс времени, которое вы хотите использовать для среднего.Для этого я бы использовал пакет lubridate, но это не обязательно.

index1 <- month(time_serie1) < 10 & month(time_serie1) > 3 # this make an index from april to septembre
index1 <- index1 & hour(time_serie1) <= 19 & hour(time_serie1) >= 8 # then you add the hour restriction
mean1 <- apply(array1[,,index1],1:2,mean)

Этот код даст вам двумерный массив со средним значением за первый год, вы можете поместить свои массивы и time_series в список и зациклить его,Тогда вы будете иметь для каждого года двумерный массив среднего значения за этот год, и вы можете усреднить эти массивы.Часть среднего веса, которую я сказал, состоит в том, что если вы сделаете это, и в среднем вы включите февраль, то ваши средства будут потрачены на разное количество дней, для вашего примера это не обязательно, но если вы используете февраль, то вынужно взвесить количество данных, используемых для каждого среднего значения.

  • Второе возможное решение

Для этого решения почти то же самое, что и для другого, но мне оно нравитсяБольше.Вы можете объединить все свои массивы в большой массив, выполняя это по порядку, чтобы индекс времени находился в возрастающем порядке, я назову этот массив BigArray.Затем объедините временные ряды, связанные с каждым массивом, я назову его BigTime.И ищите индексы, которые вы хотите усреднить, и все готово.Большим преимуществом является то, что вам не нужно делать цикл с данными в списке, и вам не нужно заботиться об изменении размера февраля.

Index <- month(BigTime) < 10 & month(BigTime) > 3 # this make an index from april to septembre
Index <- Index & hour(BigTime) <= 19 & hour(BigTime) >= 8 # then you add the hour restriction
Mean <- apply(BigArray[,,Index],1:2,mean) 

И тогда это делаетсясреднее значение для ваших значений.

В обоих возможных вариантах строится двухмерный массив, если вам нужен трехмерный массив с одним измерением (временем), имеющий только одно значение, то добавьте это измерение.И если вы хотите получить дополнительную информацию, то для получения информации о среднем значении времени обычно называйте составную технику в науке о метеорологии.

Я надеюсь, что это решит вашу проблему.

...