如何将具有异常维度的 netCDF 转换为标准 netCDF(ltime、lat、lon)(python)

How to convert netCDFs with unusual dimensions to a standard netCDF (ltime, lat, lon) (python)

我有多个最终想要合并的 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...

维度 'land' 是数字 1 到 67209 的计数,'tstep' 是 1 到 248 的计数。

变量'lat'和'lon'是经纬度值,形状为(67209,)

变量'time'是自月初以来的秒数(netcdf 是一个月)

接下来我将尺寸从 'tstep' 交换为 '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

从这里我想创建一个具有维度(时间、纬度、经度)和变量 Rainf 的 netCDF。

我已经尝试创建一个新的 netCDF(或改变这个)但是当我尝试传递 Rainf 变量时不起作用,因为它具有 (248, 67209) 的形状并且需要 (248, 67209) 的形状, 67209).尽管 'Rainf' 的当前 'land' 维度具有纬度和经度坐标。是否可以重塑此变量以具有时间、纬度和经度维度?

最后看来你想要的是将 "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()方法。但这不是你想要的,因为它只会回到原来的 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" 维度并提供了所需的输出。