气象编程 | cdo配合shape文件提取CMIP6指定区域数据

1、前言

前面已经给大家介绍过如何用cdo配合shape文件提取CMIP6指定区域数据，这个思路首先需要用gdal将shape文件转为netCDF格式的文件，然后再用cdo ifthen进行裁剪，今天要给大家介绍的是如何用gepandas+salem库直接用shape文件对CMIP6文件进行掩码，无需把shape文件转换为netCDF格式文件，最后再调用cdo的selbox操作进行裁剪。

2、数据处理

首先导入需要的库

import numpy as np
import xarray as xa
import geopandas as geo
import salem

然后，把中国中国区域的shape内容打印看一下：

China_shp  = geo.read_file('/mnt/e/shp/china.shp')
print(China_shp)

接下来，我们就可以用geopandas库直接对当前的CMIP6数据进行掩码。

dset = xa.open_dataset('/mnt/e/test/CMIP6.nc')
land  = dset.salem.roi(shape=China_shp)
land = land.sfcWind
land.to_netcdf(path = '/mnt/e/test/output.nc')

代码非常简单，第一行是读取原始CMIP6的数据，第二行代码就完成了转换工作，是不是非常简单？第三行，仅保留需要的风速数据，第四行，我们就把掩码后的数据另存为nc文件了。

注意：其实第三行的代码是可以去掉的，但是本人的shp文件可能存在一点问题，当没有提取sfcWind数据而直接保存时，第四行代码在执行时出现问题，故采取了折中的处理方式，大家在实际操作过程中，如何保存数据，需要依情况决定。

接下来通过对比掩码前后的图片检查效果。

掩码前是全球数据：

掩码后的图片，仅包括中国区域：

这里还需要注意一点的是，经过掩码之后，output.nc文件跟CMIP6.nc文件大小(基本上)还是一样大，至于原因，我想可能是这个操作只是根据shape文件，把shape文件之外的区域值全部设置为nan，如果想让output.nc文件进一步降低文件大小，还需要进行cdo的selbox操作。

导入cdo模块后，就可以直接调用cdo.sellonlatbox函数：

from cdo import Cdo
cdo = Cdo()
cdo.sellonlatbox( '70,140,0,60', input = '/mnt/e/test/output.nc', output = '/mnt/e/test/output2.nc')# '70,140,15,40'是东经70-140的范围，北纬0-60度的范围内

从图中可以看出CMIP6.nc和output.nc文件大小是一样大，但是经过selbox操作之后，output2.nc文件就很小了。

最后，我们再来看一下output2.nc的效果：

其实只看图片效果，可能还想不明白，为什么output2的文件大小会减少的如此之多，但结合nc文件的lat和lon数据，就能想明白了，这里的数据已经不是全球180x360的数据了。