使用地理信息系统光栅的工具
georasters的Python项目详细描述
^ {tt1}$包是一个提供快速灵活的Python模块。 用于处理gis光栅文件的工具。它提供了georaster类，这使得使用raster非常透明和容易。 在某种程度上，它试图像geopandas对几何学那样处理光栅。
它包括用于
- 合并光栅
- 绘图光栅
- 从光栅中提取信息
- 给定一个点（纬度，经度）在光栅中找到它的位置
- 聚合光栅以降低分辨率
- 将两个不同大小的光栅与公共区域和大小对齐
- 获取光栅的所有地理信息
- 轻松创建geotiff文件
- 将geotiff文件加载为屏蔽的numpy raster
- 使用几何图形剪裁光栅
- 使用几何图形获取分区统计信息
- 空间分析工具
安装
pipinstallgit+git://github.com/ozak/georasters.gitpipinstallgeorasters
示例用法：georasters
importgeorastersasgr# Load dataraster='./data/slope.tif'data=gr.from_file(raster)# Plot datadata.plot()# Get some statsdata.mean()data.sum()data.std()# Convert to Pandas DataFramedf=data.to_pandas()# Save transformed data to GeoTiffdata2=data**2data2.to_tiff('./data2')# Algebra with rastersdata3=np.sin(data.raster)/data2data3.plot()# Notice that by using the data.raster object,# you can do any mathematical operation that handles# Numpy Masked Arrays# Find value at point (x,y) or at vectors (X,Y)value=data.map_pixel(x,y)Value=data.map_pixel(X,Y)
合并georaster示例：
importosimportgeorastersasgrimportmatplotlib.pyplotaspltDATA="/path/to/tiff/files"# Import rasterraster=os.path.join(DATA,'pre1500.tif')data=gr.from_file(raster)(xmin,xsize,x,ymax,y,ysize)=data.geot# Split raster in twodata1=gr.GeoRaster(data.raster[:data.shape[0]/2,:],data.geot,nodata_value=data.nodata_value,projection=data.projection,datatype=data.datatype,)data2=gr.GeoRaster(data.raster[data.shape[0]/2:,:],(xmin,xsize,x,ymax+ysize*data.shape[0]/2,y,ysize),nodata_value=data.nodata_value,projection=data.projection,datatype=data.datatype,)# Plot both parts and save themplt.figure(figsize=(12,8))data1.plot()plt.savefig(os.path.join(DATA,'data1.png'),bbox_inches='tight')
plt.figure(figsize=(12,8))data2.plot()plt.savefig(os.path.join(DATA,'data2.png'),bbox_inches='tight')
# Generate merged rasterdata3=data1.union(data2)# Plot it and save the figureplt.figure(figsize=(12,8))data3.plot()plt.savefig(os.path.join(DATA,'data3.png'),bbox_inches='tight')
另一个合并：
示例用法：其他函数
importgeorastersasgr# Get info on rasterNDV,xsize,ysize,GeoT,Projection,DataType=gr.get_geo_info(raster)# Load rasterdata=load_tiff(raster)# Find location of point (x,y) on raster, e.g. to extract info at that locationcol,row=gr.map_pixel(x,y,GeoT[1],GeoT[-1],GeoT[0],GeoT[3])value=data[row,col]# Agregate raster by summing over cells in order to increase pixel size by e.g. 10gr.aggregate(data,NDV,(10,10))# Align two rastersdata2=load_tiff(raster2)(alignedraster_o,alignedraster_a,GeoT_a)=gr.align_rasters(raster,raster2,how=np.mean)# Create GeoRasterA=gr.GeoRaster(data,GeoT,nodata_value=NDV)# Load another rasterNDV,xsize,ysize,GeoT,Projection,DataType=gr.get_geo_info(raster2)data=load_tiff(raster2)B=gr.GeoRaster(data2,GeoT,nodata_value=NDV)# Plot RasterA.plot()# Merge both rasters and plotC=B.merge(A)C.plot()
问题
找到虫子了吗？通过github问题报告
- 下载再现错误所需的最小可能光栅和矢量数据集的链接
- 再现错误的python代码或命令
- 有关您的环境的信息：python、gdal和numpy的版本以及系统内存
欢迎加入QQ群-->： 979659372
