1

我正在使用 siphon 从 Unidata Thredds 服务器下载 GFS 数据，因此我可以使用 MetPy 绘制它。我写了一个脚本来做到这一点，昨天它工作得很好：

#Get data using siphon
best_gfs = TDSCatalog('http://thredds.ucar.edu/thredds/catalog/grib/NCEP/GFS/Global_0p25deg/catalog.xml?dataset=grib/NCEP/GFS/Global_0p25deg/Best')
best_ds = best_gfs.datasets[0]
ncss = best_ds.subset()
query = ncss.query()
query.lonlat_box(north=55, south=20, east=-60, west=-120).time(datetime.utcnow())
query.accept('netcdf4')
query.variables('Geopotential_height_isobaric')

data = ncss.get_data(query)

#Parse data using MetPy
ds = xr.open_dataset(NetCDF4DataStore(data))
data = ds.metpy.parse_cf()

time_of_run = data['reftime'][0].dt.strftime('%Y%m%d_%H%MZ').values
print(time_of_run)

当我在美国东部时间下午 2 点左右运行它时，这段代码输出了 2020-03-29 12:00Z，一切都很好。

今天早上运行它时，我得到了一个错误：

Traceback (most recent call last):
  File "C:\Users\jacks\Anaconda3\envs\metpy_test\lib\site-packages\xarray\core\dataset.py", line 1155, in _construct_dataarray
    variable = self._variables[name]
KeyError: 'reftime'

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "h5_rh_wind_gph_temp.py", line 51, in <module>
    time_of_run = data['reftime'][0].dt.strftime('%Y%m%d_%H%MZ').values
  File "C:\Users\jacks\Anaconda3\envs\metpy_test\lib\site-packages\xarray\core\dataset.py", line 1245, in __getitem__
    return self._construct_dataarray(key)
  File "C:\Users\jacks\Anaconda3\envs\metpy_test\lib\site-packages\xarray\core\dataset.py", line 1158, in _construct_dataarray
    self._variables, name, self._level_coords, self.dims
  File "C:\Users\jacks\Anaconda3\envs\metpy_test\lib\site-packages\xarray\core\dataset.py", line 165, in _get_virtual_variable
    ref_var = variables[ref_name]
KeyError: 'reftime'

这表明对“reftime”键的引用是无效的。为了调查，然后我打印了“数据”xarray：

<xarray.Dataset>
Dimensions:                       (isobaric6: 34, lat: 141, lon: 241, time1: 1)
Coordinates:
    reftime1                      (time1) datetime64[ns] ...
  * time1                         (time1) datetime64[ns] 2020-03-30T12:00:00
  * isobaric6                     (isobaric6) float32 40.0 100.0 ... 100000.0
  * lat                           (lat) float32 55.0 54.75 54.5 ... 20.25 20.0
  * lon                           (lon) float32 240.0 240.25 ... 299.75 300.0
    crs                           object Projection: latitude_longitude
Data variables:
    Geopotential_height_isobaric  (time1, isobaric6, lat, lon) float32 ...
    LatLon_Projection             int32 ...
Attributes:
    Originating_or_generating_Center:                                        ...
    Originating_or_generating_Subcenter:                                     ...
    GRIB_table_version:                                                      ...
    Type_of_generating_process:                                              ...
    Analysis_or_forecast_generating_process_identifier_defined_by_originating...
    Conventions:                                                             ...
    history:                                                                 ...
    featureType:                                                             ...
    History:                                                                 ...
    geospatial_lat_min:                                                      ...
    geospatial_lat_max:                                                      ...
    geospatial_lon_min:                                                      ...
    geospatial_lon_max:                                                      ...

这表明我想要的信息（模型的运行时间）现在存储为“reftime1”。为什么这个 1 突然出现在 reftime 键的末尾？它的出现有规律吗？我希望最终将此脚本作为 cron 作业运行以自动生成绘图，因此最好找出一种方法来预测名称的这种变化或完全规避键名。

4

1 回答 1

1

reftime从到的变化来自reftime1于 THREDDS 和 netCDF-java 如何处理生成底层 GRIB 数据的 netCDF 表示。GRIB 基本上是作为单独的 2D 数据切片到达的。为了创建时间、参考时间和各种垂直维度，netCDF-java 正在查看可用于给定字段（例如Geopotential_height_isobaric）的一组 GRIB 消息。如果字段具有不同的时间/垂直维度集，则使用唯一名称创建单独的维度，例如reftime、reftime1、reftime2。哪个字段以哪个维度名称结束取决于 netCDF-java 在集合中遇到特定 GRIB 消息的顺序。

解决这个问题的方法是避免依赖于名称，而是使用元数据来找出您要查找的内容。MetPy 可以通过为各种维度/坐标提供适当的别名来完成其中的一些工作：

# Will point to appropriate time1, time2, etc.
time = data.Geopotential_height_isobaric.metpy.time

这适用于特定变量的坐标。在 的情况下reftime，由于这不是变量的坐标，您还可以通过查找气候和预报 (CF) 元数据标准名称来查找forecast_reference_time：

filtered_ds = data_.filter_by_attrs(standard_name='forecast_reference_time')

这仍然留下了一个 xarray Dataset，您需要找到某种方法来提取内部唯一的变量——我不确定从那里向下走的最佳方法是什么。

于 2020-03-30T20:09:25.480 回答