用Python进行栅格数据的分区统计和批量提取

有时候我们会有这样的想法,就是针对某个区域的栅格数据,要提取它的平均值或者其他统计指标,比如在一个省内提取多年的降雨数据,最后分区域地计算一些统计值,或者从多个栅格数据中提取某个区域的数值形成一个序列。为了方便,画一个示意图看看,比如就像提取这个区域中的某一个市的区域,然后形成一个序列数据,这就可以使用rasterstats库了,此外的分区统计也可以用这个库

这个实验使用的数据格式分别是栅格(*.tif)和矢量(.shp),之后的分区统计操作和栅格数据的提取都是源于这两类数据。为了能使用上这个rasterstats库,选择了在google colab平台运行脚本,因为安装库实在是太方便了,在win上老是安装不上的,在google notebook立马就搞定了,而且可以把数据存储到谷歌云盘,直接在notebook中就是可以链接使用的

那么现在就开始做测试,使用的数据就是左侧的栅格和矢量数据集
导入相关的模块

import geopandas as gpd
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import rasterio
import rasterstats
from rasterio.plot import show
# show()方法用来展示栅格图形
from rasterio.plot import show_hist
# 用来展示直方图
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter

使用geopandas和rasterio分别读取矢量和栅格数据

# 使用geopandas读取矢量数据
districts = gpd.read_file('/content/drive/MyDrive/Datashpraster/Data/Districts/districts.shp')

# 使用rasterio读取栅格数据,栅格数据和矢量数据的坐标投影需要一致
raster = rasterio.open('/content/drive/MyDrive/Datashpraster/Data/Rainfall Data Rasters/2020-4-1.tif')
# 把矢量数据和栅格数据绘制到一个axis上,这个axis不是坐标轴,而是图形
plt.rcParams['font.family'] = 'Times New Roman'
plt.rcParams['font.size'] = 20

fig, (ax1,ax2) = plt.subplots(1,2,figsize=(15,6))

show(raster, ax=ax1,title='Rainfall')
# 读取进来的矢量数据可以直接调用gpd的plot()方法绘制
districts.plot(ax=ax1, facecolor='None', edgecolor='red')
show_hist(raster,ax=ax2,title='hist')

plt.show()

先绘制一下结果看看

读取栅格数据:

# 提取雨量栅格值到numpy数组
# 遵循GDAL规则从第一波段读取
rainfall_data = raster.read(1)
rainfall_data

开始分区统计:

# 设置坐标变换信息
affine = raster.transform

# 准备开始进行空间分区计算
# 第一个参数是矢量分区,第二个是栅格,第三个是坐标变换信息,第四个是统计均值
avg_rallrain = rasterstats.zonal_stats(districts,rainfall_data,affine=affine,stats=['mean'],geojson_out=True)
# avg_rallrain

# 除了统计平均值之外,还有最大最小值那些

绘制一下,只是一个简单的图形而已

当然第二部分更有意思,就是从多个分散的栅格数据中提取数据形成一个序列

,就是这些tif数据

loop这些栅格数据集:

获得提取到的结果,没错,就是这么一个序列数据,然后就是绘图了

转换数据格式

# 将Date列转为时间型
data['Date'] = pd.to_datetime(data['Date'], infer_datetime_format=True)

# print(data)

data['Date'] = data['Date'].dt.date
print(data)

绘图结果就是简单的图形而已

# 准备绘制图形
fig,(ax1,ax2)= plt.subplots(2,1,figsize=(18,6))
plt.rcParams['font.size'] = 15

data.plot(x='Date', y='Average_RF_Porto', ax=ax1, kind='bar', title='Avg_Rail_Porto')
data.plot(x='Date', y='Average_RF_Faro', ax=ax2, kind='bar', title='Avg_Rail_Faro',color='red')

#自动调整图形的分布
plt.tight_layout()
plt.show()

结果就这样一个序列图,目的就是从栅格提取指定的研究区,然后提取栅格的值,再来绘图

虽然感觉不是那么花里胡哨的图,但这个应该还是比较实用的,特别是大批量提取栅格值的时候。由于在google colab里面操作的步骤比较多,中间可能有省略的地方,但重要的应该都在文中了,当然也可以迁移运用到其他地方,也可以查看一下这个第三方库的教程,比如read(1)是什么意思,官网的docs就写得有,实在是很方便的

以上就是用Python进行栅格数据的分区统计和批量提取的详细内容,更多关于Python 栅格数据的分区统计和批量提取 的资料请关注我们其它相关文章!

(0)

相关推荐

  • python 矢量数据转栅格数据代码实例

    这篇文章主要介绍了python 矢量数据转栅格数据代码实例,文中通过示例代码介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值,需要的朋友可以参考下 投影包osr与proj4的使用 osr投影转换示例 from osgeo import osr,ogr #定义投影 #wgs84 source=osr.SpatialReference() source.ImportFromEPSG(4326) #google target=osr.SpatialReference() target.Imp

  • Python叠加两幅栅格图像的实现方法

    目的 现有两幅栅格图像,一个是某地区道路栅格图,一个是某地区土地利用类型图,需要将道路叠加到土地利用类型图中,即叠加后,重合的像元值以道路图为准,其余的像元值仍是土地利用类型图原有的像元值. 图1 道路信息图 图2 土地利用类型图 图3 结果图 具体实现 from gdalconst import * from osgeo import gdal import osr import sys import copy #叠加两个栅格图像(一个道路栅格图,一个土地利用类型图),两幅图像重叠的像元值都是

  • 在Python中用GDAL实现矢量对栅格的切割实例

    概述: 本文讲述如何在Python中用GDAL实现根据输入矢量边界对栅格数据的裁剪. 效果: 裁剪前 矢量边界 裁剪后 实现代码: # -*- coding: utf-8 -*- """ @author lzugis @date 2017-06-02 @brief 利用shp裁剪影像 """ from osgeo import gdal, gdalnumeric, ogr from PIL import Image, ImageDraw impo

  • python获取栅格点和面值的实现

    1.获取指定位置的点值: OutputFile = open(statisticResultTXT, 'w') cellvalue=arcpy.GetCellValue_management(inputfilepath+filenname+".tif",staionXY, "1") OutputFile.write(stationID+"_"+filenname+""+str(cellvalue)+'\n') OutputFi

  • 关于Python 的简单栅格图像边界提取方法

    在GIS中,栅格属性里有关于栅格自身的信息,背景(nodata value)对于识别一张图像的边界像元尤为重要,我们目的只要把每行每列中的第一次出现不是nodata的像元和最后一次出现nodata的前一个像元就可以了. 对于栅格,可以用ArcPy中的RasterToNumpyArray函数将将栅格转成numpy数组,然后就可以按照所想读取出每行列中首尾像元. 以下是部分代码提取边界像元的核心算法,其实是很简单的一个思路(假设0是nodata value). a=[[0 for col in ra

  • 用Python进行栅格数据的分区统计和批量提取

    有时候我们会有这样的想法,就是针对某个区域的栅格数据,要提取它的平均值或者其他统计指标,比如在一个省内提取多年的降雨数据,最后分区域地计算一些统计值,或者从多个栅格数据中提取某个区域的数值形成一个序列.为了方便,画一个示意图看看,比如就像提取这个区域中的某一个市的区域,然后形成一个序列数据,这就可以使用rasterstats库了,此外的分区统计也可以用这个库 这个实验使用的数据格式分别是栅格(*.tif)和矢量(.shp),之后的分区统计操作和栅格数据的提取都是源于这两类数据.为了能使用上这个r

  • python获取文件后缀名及批量更新目录下文件后缀名的方法

    本文实例讲述了python获取文件后缀名及批量更新目录下文件后缀名的方法.分享给大家供大家参考.具体实现方法如下: 1. 获取文件后缀名: 复制代码 代码如下: #!/usr/bin/python import os dict = {} for d, fd, fl in os.walk('/home/ahda/Program/'):         for f in fl:                 sufix = os.path.splitext(f)[1][1:]           

  • Python爬虫框架Scrapy实战之批量抓取招聘信息

    网络爬虫抓取特定网站网页的html数据,但是一个网站有上千上万条数据,我们不可能知道网站网页的url地址,所以,要有个技巧去抓取网站的所有html页面.Scrapy是纯Python实现的爬虫框架,用户只需要定制开发几个模块就可以轻松的实现一个爬虫,用来抓取网页内容以及各种图片,非常之方便- Scrapy 使用wisted这个异步网络库来处理网络通讯,架构清晰,并且包含了各种中间件接口,可以灵活的完成各种需求.整体架构如下图所示: 绿线是数据流向,首先从初始URL 开始,Scheduler 会将其

  • 利用Python批量提取Win10锁屏壁纸实战教程

    前言 相信使用Win10的朋友会发现,每次开机锁屏界面都会有不一样的漂亮图片,这些图片通常选自优秀的摄影作品,十分精美. 但是由于系统会自动更换这些图片,所以就算再好看的图片,也许下次开机之后就被替换掉了. 借助Python,我们可以用简单的几行代码,批量提取这些精美的锁屏图片.把喜欢的图片设置成桌面背景,就不用担心被替换掉啦. 下面话不多说了,来一起看看详细的介绍吧. 提取原理 Win10系统会自动下载最新的锁屏壁纸,并将他们保存在一个系统文件夹中,路径是C:\Users\[用户名]\AppD

  • Python批量提取PDF文件中文本的脚本

    本文实例为大家分享了Python批量提取PDF文件中文本的具体代码,供大家参考,具体内容如下 首先需要执行命令pip install pdfminer3k来安装处理PDF文件的扩展库. import os import sys import time pdfs = (pdfs for pdfs in os.listdir('.') if pdfs.endswith('.pdf')) for pdf1 in pdfs: pdf = pdf1.replace(' ', '_').replace('-

  • Python从数据库读取大量数据批量写入文件的方法

    使用机器学习训练数据时,如果数据量较大可能我们不能够一次性将数据加载进内存,这时我们需要将数据进行预处理,分批次加载进内存. 下面是代码作用是将数据从数据库读取出来分批次写入txt文本文件,方便我们做数据的预处理和训练机器学习模型. #%% import pymssql as MySQLdb #这里是python3 如果你是python2.x的话,import MySQLdb #数据库连接属性 hst = '188.10.34.18' usr = 'sa' passwd = 'p@ssw0rd'

  • Python将list中的string批量转化成int/float的方法

    最近在处理词向量这块,因为平时习惯把处理的词向量保存成文件,但是txt文件读取出来的都是string格式的数字,有必要转成float型 上网查了一下教程,在这记录一下: data = ['1','3.2','2'] data = map(eval, data) print data 不知道map函数怎么实现的,没看官方文档,反正实现了就好. 输出:[1, 3.2, 2] 原有string格式的数字是整形就输出整形,是浮点就输出浮点. 以上这篇Python将list中的string批量转化成int

  • Python如何使用Gitlab API实现批量的合并分支

    这篇文章主要介绍了Python如何使用Gitlab API实现批量的合并分支,文中通过示例代码介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值,需要的朋友可以参考下 1.需求:每次大批量上线完成后,都会进行将hotfix合并到Master,合并到test/uat等等重复操作(上线发布后自动合并master已完成). 2.现实:在完成发布后自动合并master后,可能还有的项目人员忘记合并到其他分支的情况,so #!/usr/bin/python3 #coding=utf-8 # 自动合

  • Python使用扩展库pywin32实现批量文档打印实例

    本文代码需要正确安装Python扩展库pywin32,建议下载whl文件进行离线安装.然后调用win32api的ShellExecute()函数来实现文档打印,系统会根据文档类型自动选择不同的软件进行打开并自动打印,如果要打印的是图片的话,需要手工确认一下. 关于ShellExecute()函数的参数含义请查阅Windows API或pywin32帮助文档. import win32print import win32api for fn in ['1.txt', '2.txt', '3.txt

  • Python用5行代码实现批量抠图的示例代码

    前言 对于会PhotoShop的人来说,抠图是非常简单的操作了,有时候几秒钟就能扣好一张图.不过一些比较复杂的图,有时候还是要画点时间的,今天就给大家带了一个非常快速简单的办法,用Python来批量抠取人像. 效果展示 开始吧,我也不看好什么自动抠图,总觉得不够精确,抠不出满意的图.下面我就直接展示一下效果图吧.我们先看看原图 这张图片背景未纯色,我们平时用PhotoShop抠起来也比较简单,对我们计算机来说也不是什么难题,下面是效果图: 因为本身是PNG图片,而且原图是白色背景,所以看不出什么

随机推荐