??得益于 gdal 的强大功能,gma 继承了其对各类栅格/矢量(目前支持超过80种矢量文件格式)文件的读取支持,并在其基础上进行简化、优化等改造,以便其符合 gma 的整体建库逻辑。 ??本文基于 gma 定义的矢量读取逻辑和方式,重点展示如何结合 Matplotlib,实现矢量面文件的空间绘制。
绘制思路
开始
打开文件
读取图层/要素
读取折点坐标
生成并绘制多边形
结束
gma 矢量读取逻辑介绍
??在 gma 里,所有的打开的矢量数据可分为三级,分别为 数据资源—图层—要素。一般来说,普通的矢量数据均只含一个图层(例如:shp文件等)。一些特殊的文件可能会包含多个图层,这点与含有子数据集的多维栅格类似。可通过如下代码查看支持多图层的矢量格式:
import gma
gma.config.VectorFormatInfo().MultiLayersFormats
{‘AVCBIN’, ‘CAD’, ‘DGN’, ‘DXF’, ‘EDIGEO’, ‘GML’, ‘GPKG’, ‘GPX’, ‘KML’, ‘LIBKML’, ‘MSSQLSpatial’, ‘Memory’, ‘NGW’, ‘ODBC’, ‘ODS’, ‘OSM’, ‘OpenFileGDB’, ‘PDF’, ‘PGeo’, ‘PostgreSQL’, ‘S57’, ‘SQLite’, ‘VRT’, ‘XLS’, ‘XLSX’}
层级之间的关系及应用方式如下图:
注意:在 gdal/ogr 中,Feature 下还有 Geometry(几何)这一层(类),gma 舍弃了这一部分!!
环境与数据
系统: Window 10+ (X64) Python 版本: 3.8.8 + gma 版本: 1.1.2 + 其他 Python 第三方库: matplotlib、numpy
试验数据: https://gma.luosgeo.com/Open/China_Province_2022.7z
??gma:地理与气象分析库。安装和详细功能帮助见:地理与气象分析库。
绘制矢量的基础代码(详询:Luo_Suppe)
??基于 gma 最小矢量数据单元 Feature 的 PlotFeature (绘制要素)是绘制矢量面的最小算法,用以被 PlotLayer (绘制图层)调用。
注意:下述代码部分内容可忽略,因为其是为后续的标注和图例等功能准备的。 注意:未体现的导入在 from XXX import * 中已经引入。
import matplotlib.patches as ptc
from functools import partial
import numpy as np
from gma.algorithm.core.dataio import *
from gma.algorithm.core.gcreate import CreatePolygon
from gma.extend.mapplot.utils import *
class PlotFeature:
'''绘制矢量要素'''
def __init__(self, GMAFeature, Axes, **kwargs):
if isinstance(GMAFeature, Feature):
self.GMAFeature = GMAFeature
else:
raise
if issubclass(plt.Subplot, type(Axes)) is False:
raise
else:
self.Axes = Axes
plt.axis('equal')
self.ELECollection = self.ELECollections()
def ELECollections(self):
'''解包多要素集合'''
P = self.GMAFeature.Points
while isinstance(P[0], list):
if isinstance(P[0][0], list):
if isinstance(P[0][0][0], list) is False:
break
P = sum(P, [])
return P
def DrawPolygon(self, FaceColor = None, EdgeColor = None, Hatch = None, LineStyle = None, LineWidth = None, ADBoundary = True, **kwargs):
'''绘制多边形'''
for k in ['facecolor', 'edgecolor', 'hatch', 'linestyle', 'linewidth']:
if k in kwargs:
kwargs.pop(k)
XMIN, XMAX, YMIN, YMAX = [np.inf, -np.inf, np.inf, -np.inf]
FaceColor, EdgeColor, Hatch, LineStyle, LineWidth = PlotOptions(FaceColor, EdgeColor, Hatch, LineStyle, LineWidth)
for ELE in self.ELECollection:
ELE = np.array(ELE)[:, :2]
if len(ELE) < 3 or ELE[0] == []:
continue
XMIN, YMIN = np.concatenate([ELE, [[XMIN, YMIN]]], axis = 0).min(axis = 0)
XMAX, YMAX = np.concatenate([ELE, [[XMAX, YMAX]]], axis = 0).max(axis = 0)
Poly = ptc.Polygon(ELE,
facecolor = FaceColor,
edgecolor = EdgeColor,
hatch = Hatch,
linestyle = LineStyle,
linewidth = LineWidth,
**kwargs)
self.Axes.add_patch(Poly)
if ADBoundary == True:
self.Axes.set_xlim(XMIN - (XMAX - XMIN) * 0.05, XMAX + (XMAX - XMIN) * 0.05)
self.Axes.set_ylim(YMIN - (YMAX - YMIN) * 0.05, YMAX + (YMAX - YMIN) * 0.05)
class PlotLayer:
'''绘制矢量图层'''
def __init__(self, GMALayer, Axes = None, Boundary = None, **kwargs):
if isinstance(GMALayer, Layer):
self.GMALayer = GMALayer
else:
raise
self.BBox = self.GMALayer.Boundary
if Boundary is None:
self.Left = self.BBox[0] - (self.BBox[2] - self.BBox[0]) * 0.05
self.Right = self.BBox[2] + (self.BBox[2] - self.BBox[0]) * 0.05
self.Bottom = self.BBox[1] - (self.BBox[3] - self.BBox[1]) * 0.05
self.Top = self.BBox[3] + (self.BBox[3] - self.BBox[1]) * 0.05
else:
self.Left, self.Bottom, self.Right, self.Top = INITBoundary(Boundary)
self.BFeature = self._CreateBGeom()
if Axes is None:
X = self.Right - self.Left
Y = self.Top - self.Bottom
if X > Y:
Width = 10
Height = Y / X * Width
else:
Height = 10
Width = X / Y * Height
Fig = plt.figure(figsize = (Width, Height), dpi = 100)
Axes = plt.subplot(1, 1, 1, **kwargs)
elif issubclass(plt.Subplot, type(Axes)) is False:
raise
self.Axes = Axes
self.ColorOPT = []
self.DrawFeatureInfo = list(self._FeatureIterator())
def _CreateBGeom(self):
'''创建一个视图范围内的几何体'''
Points = [[self.Left, self.Top],
[self.Left, self.Bottom],
[self.Right, self.Bottom],
[self.Right, self.Top]]
self.BPolygon = CreatePolygon(Points)
def _FeatureIterator(self):
'''要素迭代器'''
i = -1
for FID in range(self.GMALayer.FeatureCount):
FU = self.GMALayer.GetFeature(FID)
if FU._Geometry.Intersects(self.BPolygon):
FEA = ogr.Feature(ogr.FeatureDefn())
FEA.SetGeometry(FU._Geometry.Intersection(self.BPolygon))
i = i + 1
yield i, FID, FU, Feature(FEA)
def DrawPolygon(self, FaceColor = None, EdgeColor = None, Hatch = None, LineStyle = None, LineWidth = None, **kwargs):
'''绘制图层'''
GetP = partial(GetPOptions, FaceColor = FaceColor, EdgeColor = EdgeColor,
Hatch = Hatch, LineStyle = LineStyle, LineWidth = LineWidth)
for i, FID, FU, FPlot in self.DrawFeatureInfo:
gmp = PlotFeature(FPlot, self.Axes)
FC, EC, H, LS, LW = GetP(i)
self.ColorOPT.append([FC, EC, H, LS, LW])
gmp.DrawPolygon(FaceColor = FC, EdgeColor = EC, Hatch = H, LineStyle = LS, LineWidth = LW, ADBoundary = False, **kwargs)
self.Axes.set_xlim(self.Left, self.Right)
self.Axes.set_ylim(self.Bottom, self.Top)
return self.Axes
应用及说明(详询:Luo_Suppe)
import gma
DS = gma.Open('China_Province_2022.shp')
LY = DS.GetLayer(0)
from gma.extend.mapplot import plot
Colors = plot.GetColorFromCMap('jet', N = LY.FeatureCount, Alpha=1)
PlotP = PlotLayer(LY, Boundary = None)
Axes = PlotP.DrawPolygon(FaceColor = Colors,
EdgeColor = 'gray',
Hatch = ['/', '\\', '|', '-', '+', 'x', 'o', 'O', '.', '*'] * 10,
LineStyle = '-',
LineWidth = 1)
PlotP = plot.PlotLayer(LY, Boundary = [109,31,118,38])
Axes = PlotP.DrawPolygon(FaceColor = 'lightblue',
EdgeColor = 'gray',
Hatch = None,
LineStyle = '-.',
LineWidth = None)
注意:左右边界处有缝隙,暂不清楚问题出在哪里。
参数说明
PlotLayer(GMALayer, Axes=None, Boundary=None, **kwargs)
GMALayer: 通过 gma.Open 打开并获取的矢量图层。
Axes: matplotlib 子图,不存在则自动生成。
Boundary: 绘图范围(左、下、右、上),默认是输入矢量范围外扩 0.05 倍。
**kwargs: 自动生成 matplotlib 子图时传递给 plt.subplot 的其他参数。
PlotP.DrawPolygon(FaceColor=None, EdgeColor=None, LineStyle=None, Hatch=None, LineWidth=None, **kwargs)
FaceColor: 填充颜色。可为每个多边形分配不同的颜色。详见 matplotlib.patches.Polygon。
EdgeColor: 边界线颜色。可为每个边分配不同的颜色。详见 matplotlib.patches.Polygon。
LineStyle: 边界线线形。可为每个边分配不同的线性。详见 matplotlib.patches.Polygon。
LineWidth: 边界线线宽。可为每个边分配不同的线宽。详见 matplotlib.patches.Polygon。
Hatch: 填充纹理。可为每个多边形分配不同的纹理。详见 matplotlib.patches.Polygon。
**kwargs: 绘制多边形时传递给 matplotlib.patches.Polygon 的其他参数。
|