利用geometry函数和shp文件进行数据的精确筛选

jl2587t

利用geometry函数和shp文件进行数据的精确筛选

本帖主要介绍的是利用shp文件和geometry函数实现数据的精确筛选，意思就是加入作者只需要一个省的数据，手上确是全国或者全球的数据，怎么去清洗数据，只得到符合省界经纬度圈内的数据。
先上结果图（以雷电数据为例）

图中是全国及周边的FY4闪电数据，但是作者只需要湖南省的，就需要数据清洗


说明：
1.这里只介绍重点部分，对于画图将不再介绍，源代码里面可以看懂
2.geometry函数其实就在shapely库中
3.该方法经作者实测应该不是最优解，因为作者i7十一代的CUP依然运行了很久，如果进行逐个点的判断，会比较慢，需要先进行处理，减少数据量。
4.需要进行初步的数据清洗，利用经纬度先做大致的裁剪，然后再进行精细的判断。


具体实现：
#本期重点，利用shp文件判断数据是否再内部
sf = shapefile.Reader('Hunan_Province')#读取shp文件
sf = sf.shapeRecords()[0].shape#获得经纬度的边缘，这里简写了，实际步骤为读取shapeRecords()，再读取第0个，再读取shape

flat_points = np.column_stack((lon, lat))#将一维的经纬度整合成二维
in_shape_points = []#建立一个list存储所需数据
for pt in flat_points:
    if geometry.Point(pt).within(geometry.shape(sf)):#判断点是否再shp形状内部
        in_shape_points.append(pt)#存入数组
selected_lon = [elem[0] for elem in in_shape_points]#重新展开为一维
selected_lat = [elem[1] for elem in in_shape_points]#重新展开为一维


代码解释：
1.lon, lat经纬度已经进行了初步的筛选，得到了一小区域的所有点
2.sf代表shp文件的shape属性，代码中进行了缩写，可以分开一步一步写
3.代码中用的shp文件是单省的，有的人可能是全国的文件，需要做循环进行判断
例如：
for city_rcd in sz_shp.shapeRecords():  # 遍历每一条shaperecord
         if city_rcd.record[6] == 'Shenzhen':  # 遍历时，record字段是地区的信息（由字符串表示）
             sz_shp = city_rcd.shape  # 遍历时，shape字段是shape——形状（由点组成）



4.geometry.Point(pt).within(geometry.shape(sf))就是精细的判断数据是否再shp这个圈里
5.将得到的数据再次展开成一维的，用于画散点图




具体代码里面有详细的解释，但是因为版权的关系，没有附上测试数据，往见谅

Masterpiece

既然有shp边界文件，那直接用matplotlib建立一个path类，用path的类函数contains_points去判断复数个点在不在边界内就行，算上import库的代码，这程序不会超过10行

jl2587t

Masterpiece 发表于 2021-04-02 11:00
既然有shp边界文件，那直接用matplotlib建立一个path类，用path的类函数contains_points去判断复数个点在不在边界内就行，算上import库的代码，这程序不会超过10行

数了一下，7行

来自安卓APP客户端

kingjames

A2021011103078357900001文件方便分享一下吗？

jl2587t

kingjames 发表于 2021-04-04 02:00
A2021011103078357900001文件方便分享一下吗？

这个不方便，网友分享给我的，让我读一下，不然分享。可以看读取fy4a闪电数据那一贴，里面有说哪里能下载

来自安卓APP客户端

kingjames

jl2587t 发表于 2021-4-4 07:14
这个不方便，网友分享给我的，让我读一下，不然分享。可以看读取fy4a闪电数据那一贴，里面有说哪里能下载

ok  理解，谢谢答复