【可视化】中国地级市PM2.5浓度分区统计

晨昏三变

已于 2024-08-15 17:06:48 修改

阅读量209

点赞数 8

分类专栏：可视化文章标签：信息可视化 python

于 2024-08-15 16:50:49 首次发布

本文链接：https://blog.csdn.net/m0_38065162/article/details/141226377

版权

可视化专栏收录该内容

5 篇文章 0 订阅

订阅专栏

利用 GeoPandas 和xarray(rioxarray)库对中国PM2.5 浓度数据开展城市尺度的分区统计，得到我国所有地级市的平均PM2.5浓度水平。

代码实现

导入相关库

# -*- coding: UTF-8 -*-
# 导入相关库
import pandas as pd
import geopandas as gpd
import numpy as np
import matplotlib.pyplot as plt
from shapely.geometry import Point
import xarray as xr
from shapely.geometry import box
from shapely.geometry import mapping
import rioxarray

加载数据

长时序无缝高分辨率空气污染物浓度数据集（简称LGHAP）对环境管理和地球系统科学分析具有重要意义。在当前发布的 LGHAP 数据集（LGHAP v2）中，提供了长达 22 年的无缝气溶胶光学深度（AOD），近地表 PM2.5 和PM10浓度格网数据，空间分辨率为1km，覆盖中国陆地。详细数据介绍和获取方式见-zenodo。

读取中国区域的PM2.5浓度数据、裁减区域：

#读取全球PM25数据，其以nc文件存储
pm25_data = xr.open_dataset('LGHAP.Global_PM25.Y001.A2021.nc')
# 提取 PM2.5 浓度值和坐标
pm25 = pm25_data['PM25']  
lat = pm25_data['lat']  
lon = pm25_data['lon']
#读取中国城市区划shp
china_gdf=gpd.read_file('市.shp')
# 使用 unary_union 将中国的多个行政区划合并为一个
china_shape = china_gdf.geometry.union_all()
#根据行政区划对象创建裁剪范围
min_lon, min_lat, max_lon, max_lat = china_shape.bounds
bbox = box(min_lon, min_lat, max_lon, max_lat)
# 裁剪 PM2.5 数据
china_pm25_data = pm25_data.sel(lon=slice(min_lon, max_lon), lat=slice(min_lat,max_lat))
china_pm25_values = china_pm25_data['PM25']
lon2d, lat2d = np.meshgrid(china_pm25_data['lon'], china_pm25_data['lat'])

PM2.5浓度分区统计

#编写自定义的分区统计函数
def stas_PM25(geometry,nc):
    #创建城市几何对象的外接矩形
    min_lon, min_lat, max_lon, max_lat = geometry.bounds
    bbox = box(min_lon, min_lat, max_lon, max_lat)
    #根据城市的几何边界对原始nc文件进行规则裁剪
    clip_nc=nc.sel(lon=slice(min_lon, max_lon), lat=slice(min_lat,max_lat))
    ## 向NC文件中写入WGS84坐标系信息
    clip_nc.rio.write_crs("epsg:4326", inplace=True)
    # 指定NC文件的xy识别名称 
    clip_nc.rio.set_spatial_dims(x_dim="lon", y_dim="lat", inplace=True) 
    #把Geoseries对象转换为GeoDataframe对象
    city_shp=gpd.GeoDataFrame(geometry=[geometry])
    #对nc进行不规则掩膜
    ds=clip_nc.rio.clip(city_shp.geometry.apply(mapping),city_shp.crs,drop=False)
    #返回当前城市掩膜下的PM2.5均值，使用nanmean函数来忽略空值nan
    return np.nanmean(ds['PM25'].values)

全国所有地级市PM2.5浓度分区统计

#定义用于统计的PM2.5格点数据
nc=china_pm25_data
#使用apply函数对各城市对象进行分区统计，返回对应的PM2.5浓度
china_gdf['PM25']=china_gdf.apply(lambda city:stas_PM25(city['geometry'],nc),axis=1)