This case study reads information from the gstlearn repository. It is meant to show the possibilities on Polygons.
import os
import matplotlib.pyplot as plt
import gstlearn as gl
import gstlearn.plot as gp
import gstlearn.document as gdoc
import numpy as np
import scipy as sp
from scipy import interpolate
gdoc.setNoScroll()
We read the contents of the Polygon definition from a CSV file.
filename = gdoc.loadData("Alluvial", "Seine_Alluvial.csv")
csv = gl.CSVformat(False, 0, ";", ",", "9999")
poly = gl.Polygons.createFromCSV(filename, csv, False)
ax = poly.plot()
ax.decoration(title="Polygon")
The current polygon contains a large number of vertices.
poly.display(gl.AStringFormat(2))
Polygons -------- Number of Polygon Sets = 1 PolyElem #1 ........... Number of Vertices = 30301
We define a rotated grid which contains the whole alluvial plain.
grid = gl.DbGrid()
err = grid.reset([4000,3500],[100,100],[480000,6650000])
ax = grid.plot(alpha=0.3)
ax = poly.plot()
ax.decoration(title="Polygon in its alluvial plain")
Reduce the complexity of the Polygon
newpoly = poly.reduceComplexity(3000)
newpoly.display(gl.AStringFormat(2))
Polygons -------- Number of Polygon Sets = 1 PolyElem #1 ........... Number of Vertices = 1123
ax = grid.plot(alpha=0.3)
ax = newpoly.plot()
ax = poly.plot(edgecolor='red')
ax.decoration(title="(Simplified) Polygon in its alluvial plain")
Use the simplified polygon to perform a coarse selection. At the same time, dilate this selected area in order to recover the initial space as much as possible.
err = gl.db_polygon(grid,newpoly,
namconv = gl.NamingConvention("Coarse",True, True, True, gl.ELoc.SEL))
print("Number of active grid nodes =",grid.getSampleNumber(True))
Number of active grid nodes = 315868
grid.setLocator("Coarse",gl.ELoc.Z)
err = grid.morpho(gl.EMorpho.DILATION, option=1, radius = [30,30])
grid.setLocator("Morpho*",gl.ELoc.SEL)
print("Number of active grid nodes =",grid.getSampleNumber(True))
Number of active grid nodes = 1593571
ax = grid.plot(alpha=0.3)
ax = newpoly.plot()
ax = poly.plot(edgecolor='red')
ax.decoration(title="(Simplified) Polygon in its dilated selection")