NetCDF数据的标准化处理（二）_netcdf转规约-CSDN博客

本文链接：https://blog.csdn.net/liuniu1101/article/details/105515280

NetCDF数据的标准化处理（二）

在上篇“NetCDF数据的标准化处理（一）”中我们介绍了简单的标准格网的NetCDF文件的标准化方法，那么本篇我们将介绍复杂结构的NetCDF文件的标准化处理方法。

一、复杂NetCDF数据结构

NetCDF的结构如下图所示，包含了多个group组，我们需要空间化的变量为“nitrogendioxide_tropospheric_column”，这个NetCDF文件是不能直接入镶嵌数据集的：
在这里插入图片描述
我们可以看到“nitrogendioxide_tropospheric_column”的空间维度变量并不是经纬度，而是自定义的格网，同时经纬度“longitude”，“latitude”也是具有同样格网维度信息的多维变量

![在这里插入图片描述](https://img-blog.csdnimg.cn/2020041416404125.jpg?x-oss-process=image/watermark,type_ZmFuZ3poZW5naGVpdGk,shadow_10,text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2xpdW5pdTExMDE=,size_16,color_FFFFFF,t_70#pic_center
经过简单计算相邻格网之间的纬度差，发现这个NetCDF数据不是一个规则的方格网数据，那么我们应该怎么处理呢？
在这里插入图片描述

二、复杂结构的NetCDF空间化处理

对于这种不规则格网的NetCDF文件，我们可以采用下边的思路进行处理。
在这里插入图片描述
NetCDF转换为空间点要素的流程：

第一步：读取原始的NetCDF相关信息

oldnc=Dataset(inncfile,'r',format='NETCDF4')
#该nc数据有group属性
productgroup = oldnc.groups["PRODUCT"]
#读取原始NC文件的维度信息
len_scanline=len(productgroup.dimensions["scanline"])
len_groundpixel=len(productgroup.dimensions["ground_pixel"])
len_time=len(productgroup.dimensions["time"])

第二步：读取原始NC文件的变量值，并变换维度信息

lonarr=productgroup.variables["longitude"][:].reshape(len_scanline*len_groundpixel,len_time)
latarr=productgroup.variables["latitude"][:].reshape(len_scanline*len_groundpixel,len_time)
timearr=productgroup.variables["time"][:]
troposphericarr=productgroup.variables["nitrogendioxide_tropospheric_column"][:].reshape(len_scanline*len_groundpixel,len_time)

第三步：通过numpy变换，获取对应行列索引处的时间、精度、纬度、变量的值，并转换成空间点要素

#定义输出点图层的投影
sr = arcpy.SpatialReference(4326)
#以单个时刻的数据为例
for timeindex in range(len_time):
	#获取当前时刻的各个变量的numpy数组
	itemtimearr = numpy.array([timeindex+1]*len_scanline*len_groundpixel).reshape(len_scanline*len_groundpixel,1)
	newtimearr = itemtimearr[:,0].reshape(len_scanline*len_groundpixel,1)
	newlonarr = lonarr[:,0].reshape(len_scanline*len_groundpixel,1)
	newlatarr = latarr[:,0].reshape(len_scanline*len_groundpixel,1)
	newtroposphericarr = troposphericarr[:,0].reshape(len_scanline*len_groundpixel,1).data
	#将各个变量的数组水平拼接成    time,lon,lat,变量
	temparray = numpy.hstack((newtimearr,newlonarr,newlatarr,newtroposphericarr))
	datatype = numpy.dtype({'names':['timeindex','X','Y','tropospheric'],'formats':['i','f','f','f']})
	#剔除数组中原始NONE的填充值9.9692e+36
	array = temparray[numpy.where(temparray[:,3]<9.9692e+36)].tolist()
	#numpy数组需要转换成tuple的数组
	resultarray = numpy.array([tuple(x) for x in array],dtype=datatype)
	#将numpy转换成点要素图层
	arcpy.da.NumPyArrayToFeatureClass(resultarray, outfc, ('X','Y'), sr)
	break

生成的空间点要素图层结构如下图所示，我们可以看到这是一个不规则的格网的数据：
在这里插入图片描述
对点要素进行插值处理
上面得到了NetCDF转换生成的点要素图层，我们就可以采用插值的方法进行栅格化处理，下面是采用idw插值生成的栅格，如图：

对插值生成的栅格进行裁剪
大家有没有发现熟悉的问题，我们需要一个上面点要素图层实际范围的掩膜，那么面状范围制作请参见“环境镶嵌数据集的Boundary定制”中介绍的方法，下面就是裁剪之后的效果了。
在这里插入图片描述
上述示例的完整代码已经上传至网盘链接：链接：https://pan.baidu.com/s/13NRAs4pK56k78PxzqnmhGQ
提取码：onle，欢迎大家留言，互相交流学习。