矢栅处理(4)：Python读取netCDF文件

发表于 分类于 阅读次数： Waline： 本文字数： 3.4k 阅读时长 ≈ 3 分钟

NC文件由于便捷存储多维数据（如大范围长时间序列），因而广泛用于气候，气象数据存储及交换。python文档netCDF4

安装：conda install -c conda-forge netcdf4

1
2
3
import numpy as np
import netCDF4 as nc
from netCDF4 import Dataset

NetCDF

NetCDF (Network Common Data Form)，网络通用数据格式

NC文件包括以下几个属性

  • 简单来说，netcdf是一个包含多自变量及其函数的文件。用公式来说就是f(x,y,z,…)=value，同时存储了xyz
  • 函数的自变量x,y,z等在netcdf中叫做维(dimension) 或坐标轴(axix),ds.dimensions可以看到维度的信息
  • 函数值与自变量均存储在(Variables)，ds.variables属性则存储了数据的详细信息及数组文件。

具体案例见下方示例：

Python读nc文件

netCDF文件有五种（NETCDF3_CLASSIC, NETCDF3_64BIT_OFFSET, NETCDF3_64BIT_DATA, NETCDF4_CLASSIC和NETCDF4）读取之后ds.data_model输出可以看到其实际存储格式。

1
2
3
4
ncfile = 'E:\\Datasets\\temp\\China_1km_maxtmp_2017.nc'
ds = Dataset(ncfile, "r", format="NETCDF4")
print(ds.data_model)
# ds.close()

常见属性举例：时间这里是从1800年开始按天算的，实际中都要自行判断并转换，这里的温度的单位是K

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
>>> ds.dimensions
{'lon': <class 'netCDF4._netCDF4.Dimension'>: name = 'lon', size = 7386,
 'lat': <class 'netCDF4._netCDF4.Dimension'>: name = 'lat', size = 4267,
 'time': <class 'netCDF4._netCDF4.Dimension'> (unlimited): name = 'time', size = 365}

>>> ds_dict = ds.variables
>>> dt_dict.keys()
dict_keys(['time', 'lon', 'lat', 'maxtmp'])

# 查看数据信息
>>> dt_dict['time']
<class 'netCDF4._netCDF4.Variable'>
float32 time(time)
    units: days since 1800-01-01 00:00:00
    standard_name: time
    calendar: proleptic_gregorian
    long_name: Time
    axis: T
unlimited dimensions: time
current shape = (365,)
filling on, default _FillValue of 9.969209968386869e+36 used

>>> dt_dict['lon']
<class 'netCDF4._netCDF4.Variable'>
float32 lon(lon)
    standard_name: longitude
    long_name: Longitude
    units: degrees_east
    axis: X
unlimited dimensions: 
current shape = (7386,)
filling on, default _FillValue of 9.969209968386869e+36 used

>>> dt_dict['maxtmp']
<class 'netCDF4._netCDF4.Variable'>
float32 maxtmp(time, lat, lon)
    _FillValue: nan
    scale_factor: 0.01
    add_offset: 273.15
    units: K
    long_name: Daily Maximum Near-Surface Air Temperature
    standard_name: air_temperature
    missing_value: nan
unlimited dimensions: time
current shape = (365, 4267, 7386)
filling on

数据矩阵访问，通常是掩码数组，这里是False，相当于无掩码，普通数组

1
2
3
4
5
6
7
8
9
10
>>> dt_dict['lon'][:]
masked_array(data=[ 73.44696,  73.45529,  73.46363, ..., 134.97195,
                   134.9803 , 134.98863],
             mask=False,
       fill_value=1e+20,
            dtype=float32)

>>> array_lon = dt_dict['lon'][:].data
array([ 73.44696,  73.45529,  73.46363, ..., 134.97195, 134.9803 ,
       134.98863], dtype=float32)

其他类似，后面的操作都是numy多维数组的操作了

  • 访问某天的全区域的数据[1,:,:]，也可以输入特定经纬度索引范围

  • 访问某点的时间序列数据[:,index_lon,index_lat]

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
>>> dt_maxtmp = dt_dict['maxtmp'][1,:,:]
>>> dt_maxtmp
masked_array(
  data=[[--, --, --, ..., --, --, --],
        [--, --, --, ..., --, --, --],
        [--, --, --, ..., --, --, --],
        ...,
        [--, --, --, ..., --, --, --],
        [--, --, --, ..., --, --, --],
        [--, --, --, ..., --, --, --]],
  mask=[[ True,  True,  True, ...,  True,  True,  True],
        [ True,  True,  True, ...,  True,  True,  True],
        [ True,  True,  True, ...,  True,  True,  True],
        ...,
        [ True,  True,  True, ...,  True,  True,  True],
        [ True,  True,  True, ...,  True,  True,  True],
        [ True,  True,  True, ...,  True,  True,  True]],
  fill_value=nan,
  dtype=float32)

>>> dt_maxtmp.data
array([[nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       ...,
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan]], dtype=float32)