使用 Python 获取 NetCDF 变量 min/max 的最快方法？答案

【问题标题】：fastest way to get NetCDF variable min/max using Python?使用 Python 获取 NetCDF 变量 min/max 的最快方法？
【发布时间】：2014-03-11 12:14:20
【问题描述】：

与scipy.io.netcdf相比，当切换到netCDF4 Python模块时，我从NetCDF文件中提取变量数据值的min/max的常用方法是慢一个数量级。

我正在处理相对较大的海洋模型输出文件（来自 ROMS），在给定的地图区域（夏威夷）上具有多个深度级别。当这些在 NetCDF-3 中时，我使用了scipy.io.netcdf。

现在这些文件在 NetCDF-4（“经典”）中，我不能再使用 scipy.io.netcdf，而是转而使用 netCDF4 Python 模块。但是，缓慢是一个问题，我想知道是否有更有效的方法来提取变量的数据范围（最小和最大数据值）？

这是我使用 scipy 的 NetCDF-3 方法：

import scipy.io.netcdf
netcdf = scipy.io.netcdf.netcdf_file(file)
var = netcdf.variables['sea_water_potential_temperature']
min = var.data.min()
max = var.data.max()

这是我使用 netCDF4 的 NetCDF-4 方法：

import netCDF4
netcdf = netCDF4.Dataset(file)
var = netcdf.variables['sea_water_potential_temperature']
var_array = var.data.flatten()
min = var_array.data.min()
max = var_array.data.max()

显着的区别是我必须先将netCDF4中的数据数组展平，而这个操作显然会减慢速度。

有没有更好/更快的方法？

【问题讨论】：

不确定真正的原因是什么，但是 scipy.io.netcdf 和 netcdf-python 的代码库是完全不同的。第一个基于 pupynere，如果我没记错的话，它有一些很好的优化。此外，netCDF 的规格与版本相比发生了很大变化。 3 -> 4. 我想你可以尝试类似np.array(var.data).max() 的方法来避免netCDF Variable 的扁平化。很难说，因为netCDF文件的结构是未知的。
为什么要展平阵列？ import numpy as np; np.max(var[:]) 有效吗？
谢谢。我想扁平化不是瓶颈。通过var[:] 访问数组是。 @abudis：为了工作，我不得不将你的命令修改为np.array(var[:].data).max()。 @SpencerHill：是的，这可行，但同样慢。这两个建议每个都花费与我上面的示例相同的时间：仍然很慢。我想@abudis 提到的 scipy 优化和 netCDF4 更改可能是罪魁祸首。
@JohnMaurer 仅供参考 var[:] 已经是 numpy.ndarray，因此无需在其上调用 np.array() 函数或访问其 data 属性。只需var[:].max() 就可以了。不过，这对计算速度没有帮助。不会立即想到更快的方法，但其他更专业的用户可能知道。
如何使用nco 运算符之一，称为子进程？比如：ncwa -y max ...

标签： python numpy scipy netcdf

【解决方案1】：

根据hpaulj 的建议，这里有一个使用subprocess 调用nco 命令ncwa 的函数。使用 OPeNDAP 地址时，它挂得很厉害，而且我手头上没有任何文件可以在本地对其进行测试。

您可以查看它是否适合您以及速度差异是什么。

这假设您已安装 nco 库。

def ncwa(path, fnames, var, op_type, times=None, lons=None, lats=None):
    '''Perform arithmetic operations on netCDF file or OPeNDAP data

    Args
    ----
    path: str
        prefix
    fnames: str or iterable
        Names of file(s) to perform operation on
    op_type: str
        ncwa arithmetic operation to perform. Available operations are:
        avg,mabs,mebs,mibs,min,max,ttl,sqravg,avgsqr,sqrt,rms,rmssdn
    times: tuple
        Minimum and maximum timestamps within which to perform the operation
    lons: tuple
        Minimum and maximum longitudes within which to perform the operation
    lats: tuple
        Minimum and maximum latitudes within which to perform the operation

    Returns
    -------
    result: float
        Result of the operation on the selected data

    Note
    ----
    Adapted from the OPeNDAP examples in the NCO documentation:
    http://nco.sourceforge.net/nco.html#OPeNDAP
    '''
    import os
    import netCDF4
    import numpy
    import subprocess

    output = 'tmp_output.nc'

    # Concatenate subprocess command
    cmd = ['ncwa']
    cmd.extend(['-y', '{}'.format(op_type)])
    if times:
        cmd.extend(['-d', 'time,{},{}'.format(times[0], times[1])])
    if lons:
        cmd.extend(['-d', 'lon,{},{}'.format(lons[0], lons[1])])
    if lats:
        cmd.extend(['-d', 'lat,{},{}'.format(lats[0], lats[1])])
    cmd.extend(['-p', path])
    cmd.extend(numpy.atleast_1d(fnames).tolist())
    cmd.append(output)

    # Run cmd and check for errors
    subprocess.run(cmd, stdout=subprocess.PIPE, check=True)

    # Load, read, close data and delete temp .nc file
    data = netCDF4.Dataset(output)
    result = float(data[var][:])
    data.close()
    os.remove(output)

    return result

path = 'https://ecowatch.ncddc.noaa.gov/thredds/dodsC/hycom/hycom_reg6_agg/'
fname = 'HYCOM_Region_6_Aggregation_best.ncd'

times = (0.0, 48.0)
lons = (201.5, 205.5)
lats = (18.5, 22.5)

smax = ncwa(path, fname, 'salinity', 'max', times, lons, lats)

【讨论】：

【解决方案2】：

如果您只是获取变量数组中的最小值/最大值，则可以使用 xarray。

%matplotlib inline
import xarray as xr

da = xr.open_dataset('infile/file.nc')
max = da.sea_water_potential_temperature.max()
min = da.sea_water_potential_temperature.min()

这应该分别给你一个最小值/最大值。您还可以在时间、经度、纬度等选定维度上获取变量的最小值/最大值。Xarray 非常适合处理多维数组，这就是为什么当您不使用其他操作工具时在 python 中处理 NetCDF 非常容易的原因比如 CDO 和 NCO。最后，xarray 还用于其他相关库中，这些库在 python (http://xarray.pydata.org/en/stable/related-projects.html) 中处理天气和气候数据。

【讨论】：

【解决方案3】：

Python 解决方案（使用 CDO 作为后端）是我的包 nctoolkit (https://pypi.org/project/nctoolkit/https://nctoolkit.readthedocs.io/en/latest/installing.html)。

它有许多用于计算不同类型的最小值/最大值的内置方法。

我们首先需要将文件作为数据集读取：

将 nctoolkit 导入为 nc

数据 = nc.open_data(文件)

如果您想要跨空间的最大值，对于每个时间步，您可以执行以下操作：

data.spatial_max()

每个网格单元和时间步长的所有深度的最大值将按如下方式计算：

data.vertical_max()

如果你想要跨时间的最大值，你会这样做：

data.max()

这些方法是可链接的，CDO 后端非常高效，因此应该是处理 ROMS 数据的理想选择。

【讨论】：