Как преобразовать netCDF с необычными размерами в стандартный netCDF (ltime, lat, lon) (python)

Question

У меня есть несколько файлов netCDF, которые в конечном итоге я хочу объединить. Пример netCDF выглядит следующим образом.

import xarray as xr
import numpy as np
import cftime

Rain_nc = xr.open_dataset('filepath.nc', decode_times=False)
print(Rain_nc)

<xarray.Dataset>
Dimensions: (land: 67209, tstep:248)
Dimensions without coordinates: land, tstep
Data variables:
    lon    (land) float32...
    lat    (land) float32...
    timestp(tstep) int32...
    time   (tstep) int32...
    Rainf  (tstep, land) float32...

измерение «земля» - это количество чисел от 1 до 67209, а «tstep» - это количество от 1 до 248.

переменная 'lat' и 'lon' - это значения широты и долготы в форме (67209,)

переменная 'time' - это время в секундах с начала месяца (netcdf - это месяц)

Затем я поменял местами измерения с tstep на time, преобразовал их для последующего объединения и установил координаты как lon, lat и time.

rain_nc = rain_nc.swap_dims({'tstep':'time'})
rain_nc = rain_nc.set_coords(['lon', 'lat', 'time'])

rain_nc['time'] = cftime.num2date(rain_nc['time'], units='seconds since 2016-01-01 00:00:00', calendar = 'standard')
rain_nc['time'] = cftime.date2num(rain_nc['time'], units='seconds since 1970-01-01 00:00:00', calendar = 'standard')

это оставило меня со следующим Набор данных:

print(rain_nc)

<xarray.Dataset>
Dimensions: (land: 67209, time: 248)
Coordinates:
    lon        (land)float32
    lat        (land)float32
  * time       (time)float64
Dimensions without coordinates: land
Data variables:
    timestp   (time)int32
    Rainf     (time, land)


print(rain_nc['land'])
<xarray.DataArray 'land' (land: 67209)>
array([    0,    1,    2,..., 67206, 67207, 67208])
Coordinates:
    lon     (land) float32 ...
    lat     (land) float32 ...
Dimensions without coordinates: land

интересующая меня переменная Rainf выглядит следующим образом:

<xarray.DataArray 'Rainf' (time: 248, land: 67209)>
[16667832 values with dtype=float32]
Coordinates:
    lon      (land) float32 -179.75 -179.75 -179.75 ... 179.75 179.75 
179.75
    lat      (land) float32 71.25 70.75 68.75 68.25 ... -16.25 -16.75 
-19.25
  * time     (time) float64 1.452e+09 1.452e+09 ... 1.454e+09 1.454e+09
Dimensions without coordinates: land
Attributes:
    title:       Rainf
    units:       kg/m2s
    long_name:   Mean rainfall rate over the \nprevious 3 hours
    actual_max:  0.008114143
    actual_min:  0.0
    Fill_value:  1e+20

Отсюда я хотел бы создать netCDF с размерами (время, широта, долгота) и переменная Rainf.

Я попытался создать новый netCDF (или изменить его), но когда я пытаюсь передать переменную Rainf, она не работает, так как она имеет форму (248, 67209) и нуждается в форма (248, 67209, 67209). Несмотря на то, что текущее «наземное» измерение «Rainf» имеет координаты широты и долготы. Можно ли изменить форму этой переменной, чтобы у нее были измерения времени, широты и долготы?

Answer 1

В конце концов, кажется, что вы хотите изменить размеры "land" на ("lat", "lon").

Итак, у вас есть DataArray, похожий на этот:

# Setting sizes and coordinates
lon_size, lat_size = 50, 80                                                                                                                                                                           
lon, lat = [arr.flatten() for arr in np.meshgrid(range(lon_size), range(lat_size))]                                                                                                                   
land_size = lon_size * lat_size                                                                                                                                                                       
time_size = 100 

da = xr.DataArray( 
    dims=("time", "land"), 
    data=np.random.randn(time_size, land_size), 
    coords=dict( 
        time=np.arange(time_size), 
        lon=("land", lon), 
        lat=("land", lat), 
    ) 
)  

, который выглядит так:

>>> da
<xarray.DataArray (time: 100, land: 4000)>
array([[...]])
Coordinates:
  * time     (time) int64 0 1 ... 98 99
    lon      (land) int64 0 1 ... 48 49
    lat      (land) int64 0 0 ... 79 79
Dimensions without coordinates: land

Сначала мы воспользуемся методом .set_index(), чтобы сообщить xarray, что индекс "land" должен быть представлен из "lon" и "lat" координаты:

>>> da.set_index(land=("lon", "lat"))                                                                                                                                                                    
<xarray.DataArray (time: 100, land: 4000)>
array([[...]])
Coordinates:
  * time     (time) int64 0 1 ... 98 99
  * land     (land) MultiIndex
  - lon      (land) int64 0 1 ... 48 49
  - lat      (land) int64 0 0 ... 79 79

Размеры все еще ("time", "land"), но теперь "land" это MultiIndex.

Обратите внимание, что если вы попытаетесь записать в NETCDF по адресу на этом этапе вы получите следующую ошибку:

>>> da.set_index(land=("lon", "lat")).to_netcdf("data.nc")   
NotImplementedError: variable 'land' is a MultiIndex, which cannot yet be serialized to netCDF files (https://github.com/pydata/xarray/issues/1077). Use reset_index() to convert MultiIndex levels into coordinate variables instead.

Он говорит вам использовать метод .reset_index(). Но это не то, что вы хотите здесь, потому что он просто go вернется в исходное состояние da.

Теперь вы хотите использовать метод .unstack():

>>> da.set_index(land=("lon", "lat")).unstack("land")                                                                                                                                                    
<xarray.DataArray (time: 100, lon: 50, lat: 80)>
array([[[...]]])
Coordinates:
  * time     (time) int64 0 1 ... 98 99
  * lon      (lon) int64 0 1 ... 48 49
  * lat      (lat) int64 0 1 ... 78 79

Эффективно уничтожает измерение "land" и дает желаемый результат.