In [5]:
gp.plot(mydb, nameSize="z-1")
The module gstlearn.plot contains various plot functions for gstlearn objets: DbGrid, Db, Vario, Model, Polygons... These functions are also accessible as methods of each class. For example for a grid, we could use equivalently gp.grid(mygrid,...) or mygrid.plot(...), or for more specific functions: gp.point(mygrid,...) or mygrid.plot_point(...)
%%javascript
IPython.OutputArea.prototype._should_scroll = function(lines) {
return false;
}
import numpy as np
import pandas as pd
import sys
import os
import matplotlib
import matplotlib.pyplot as plt
import gstlearn as gl
import gstlearn.plot as gp
from scipy import ndimage, misc
We define the Space dimension and the defaulted figure dimension for Geographic plots. The command printDefault() provides the list of defaulted values.
gl.defineDefaultSpace(gl.ESpaceType.RN, 2)
gp.setDefault(dims=[6,6])
gp.printDefault()
Non geographical defaults: - Figure dimensions = [6, 6] - Limits along X (not defined) - Limits along Y (not defined) - Aspect = auto Geographical defaults: - Figure dimensions = [8, 8] - Limits along X (not defined) - Limits along Y (not defined) - Aspect = 1
The paragraph describes the deiferrent manners for calling the graphics contained in gstlearn.plot package. This is demonstrated on an example of a 2-D data set, organized as a set of isolated points (called mydb)
ndat = 30
mydb = gl.Db.createFillRandom(ndat=30, nvar=3)
mydb
Data Base Characteristics ========================= Data Base Summary ----------------- File is organized as a set of isolated points Space dimension = 2 Number of Columns = 6 Total number of samples = 30 Variables --------- Column = 0 - Name = rank - Locator = NA Column = 1 - Name = x-1 - Locator = x1 Column = 2 - Name = x-2 - Locator = x2 Column = 3 - Name = z-1 - Locator = z1 Column = 4 - Name = z-2 - Locator = z2 Column = 5 - Name = z-3 - Locator = z3
The graphic representation is performed in general, for each object of gstlearn, in many different ways:
This triggers the option attached specifically to each type of object.
gp.plot(mydb, nameSize="z-1")
You can call the generic function regrouping all the representation methods assigned to the target class (point).
This manner offers the possibility to check the list of options in the help.
help(gp.point)
Help on function point in module gstlearn.plot:
point(db, *args, **kwargs)
Construct a figure for plotting a point data base
ax: matplotlib.Axes
db: Db containing the variable to be plotted
nameColor: Name of the variable containing the color per sample
nameSize: Name of the variable containing the size per sample
nameLabel: Name of the variable containing the label per sample
nameCoorX: Name of the variable standing for X coordinate
nameCoorY: Name of the variable standing for Y coordinate
useSel : Boolean to indicate if the selection has to be considered
color: Constant color (used if 'nameColor' is not defined)
size: Constant size (used if 'nameSize' is not defined)
sizmin: Size corresponding to the smallest value (used if 'nameSize' is defined)
sizmax: Size corresponding to the largest value (used if 'nameSize' is defined)
flagAbsSize: Represent the Absolute value in Size representation
flagCst: When True, the size is kept constant (equel to 's')
cmap: Optional Color scale
flagGradient: Draw Gradient (if Gradients are defined)
colorGradient: Color attached to the Gradient representation
scaleGradient: Scale of the Gradient representation
flagTangent: Draw Tangent (if Tangents are defined)
colorTangent: Color attached to the Tangent representation
scaleTangent: Scale of the Tangent representation
flagLegendColor: Flag for representing the Color Legend (only if name-color is defined)
flagLegendSize: Flag for representing the Size legend (only if nameSize is defined)
flagLegendLabel: Flag for representing the Label Legend (only if nameLabel is defined)
legendNameColor: Title for the Color Legend (set to 'nameColor' if not defined)
legendNameSize: Title for the Size legend (set to 'nameSize' if not defined)
legendNameLabel: Title for the Label Legend (set to 'nameLabel' if not defined)
posX: rank of the first coordinate
posY: rank of the second coordinate
**kwargs : arguments passed to matplotllib.pyplot.scatter
For demonstration sake, we can use the generic function and ask for a combined display using symbols:
proportional to the variable z-2
with a color with relies on the variable z-1
The command point returns a pointer on the type of object (from matplotlib) which is generated internally. To get rid of this information (when considered as useless), you can simply store it in a dummy variable, as it is done in the next paragraph.
dum = gp.point(mydb, nameSize="z-2", nameColor="z-1")
You can also call the the same function, but in a object-based manner, where this representation is considered as a method of the object to be represented (here a Db).
dum = mydb.point(nameColor="z-1")
In addition, for this figure, we ask the display of a legend which combines the explanation for both representations (by size and by color).
Note that here, the display is initialized as a Geographic environment (the other one being called standard). This requires the system to use the default options registered for Geographic figures (dimensions, scale factor, ...).
fig, ax = gp.initGeographic()
ax.gstpoint(mydb, nameSize="z-2", nameColor="z-1", flagLegendSize=True, flagLegendColor=True)
plt.show()
In order to illustrate the use of this extensive version which allows overlay of several layers of information, we quickly create a polygon constituted as the (dilated) convex hull of the data (this option will be used later in this document).
mypoly = gl.Polygons.createFromDb(mydb, dilate=0.05)
The overlay can be handled as follows
fig, ax = gp.initGeographic()
ax.gstpoint(mydb, nameSize="z-2", nameColor="z-1", flagLegendColor=True, flagLegendSize=True)
ax.polygon(mypoly)
plt.show()
You can also be more specific by calling a function named explicitely. This is the case for displaying the literal values attached to each sample: the method is named literal.
Again the call to this function is covered using three different syntaxes.
dum = gp.literal(mydb, name="z-1")
The same function can be considered as a method of the object to be represented (Db)
dum = mydb.literal(name="z-2")
Fincally the same function can be specified as a method of the class Axis
fix, ax = gp.initGeographic()
ax.literal(mydb, name="z-1")
plt.show()
Creating a dummy Model used for simulating a random field and display it
ctxt = gl.CovContext() # use default space
mymodel = gl.Model(ctxt)
cova = gl.CovAniso(gl.ECov.CUBIC,10,1,1.5,mymodel.getContext())
mymodel.addCov(cova)
This first case gives the demonstration that we can combine gstlearn.plot functions with matplotlib functions. This is what is performed by representing a Model and adding some annotation. As an example, we wish to add the formula of the covariance represented, i.e.:
from IPython.display import display, Latex
display(Latex('$%s$'%cova.getFormula()))
fig, ax = gp.init() # gstlearn.plot entry function
ax.model(mymodel, flagLegend=True) # gstlearn.plot function
ax.decoration(title="Cubic model") # basic gstlean.plot function
ax.text(5,0.5,'$%s$'%cova.getFormula(),size='medium') # standard matplotlib function
plt.show() # standard matplotlib: flushing graphic
We create a rectangular non-rotated 2-D grid, and simulate random Gaussian field (using the Model previously defined). Two simulations are generated in order to emphasize the graphic posibilities in further parts of this note.
nx = [70,25]
dx = [1,2]
x0 = [-40, 20]
mygrid = gl.DbGrid.create(nx,dx,x0)
err = gl.simtub(None,mygrid,mymodel,None,2)
mygrid.display()
Data Base Grid Characteristics ============================== Data Base Summary ----------------- File is organized as a regular grid Space dimension = 2 Number of Columns = 5 Total number of samples = 1750 Grid characteristics: --------------------- Origin : -40.000 20.000 Mesh : 1.000 2.000 Number : 70 25 Variables --------- Column = 0 - Name = rank - Locator = NA Column = 1 - Name = x1 - Locator = x1 Column = 2 - Name = x2 - Locator = x2 Column = 3 - Name = Simu.1 - Locator = z1 Column = 4 - Name = Simu.2 - Locator = z2
Add a dummy selection to test visualization with Selection
dum = gp.symbol(mygrid, nameColor="Simu.1")
dum = mygrid.symbol(nameSize="Simu.1")
fig, ax = gp.initGeographic()
ax.symbol(mygrid, nameColor="Simu.1")
plt.show()
gp.plot(mygrid, "Simu.1")
mygrid["sel"] = 1. - (mygrid["x1"] > 0) * (mygrid["x1"] < 15) * (mygrid["x2"] > 40) * (mygrid["x2"] < 50)
mygrid.setLocator("sel",gl.ELoc.SEL)
gstlearn.plot gives the opportunity to store somme global parameters. Their list is given by the following function where we can clearly see those dedicated to geographical representations (where axes are dedicated to spatio / temporal coordinates). The other ones are called standard representations.
gp.printDefault()
Non geographical defaults: - Figure dimensions = [6, 6] - Limits along X (not defined) - Limits along Y (not defined) - Aspect = auto Geographical defaults: - Figure dimensions = [8, 8] - Limits along X (not defined) - Limits along Y (not defined) - Aspect = 1
Two functions (i.e. setDefault and setDefaultGeographic) allow changing the default values.
We simply represent the grid (using the defaulted color scale). As no variable nor representation is explicitely specified, the default variable is represented using the default representation (i.e. raster). The default variable is the first Z-locator variable (if any) or the variable created last in the file otherwise.
We also plot the legend on the right edge.
This is the opportunity to describe the mechanism for calling the functions of gstlearn.plot.
We have essentially 2 possible levels of calls:
fig, ax = gp.initGeographic()
res = gp.raster(mygrid, name="Simu.1", flagLegend=True)
plt.show()
fig, ax = gp.initGeographic()
ax.raster(mygrid, name="Simu.1", flagLegend=True)
ax.isoline(mygrid, name="Simu.2", levels=[-2,0,1,2,3,4,5,6,7,8], colors="yellow", flagLegend=False)
plt.show()
We can also choose to use a representation per grid node. As a matter of fact, we wish to represent the grid nodes as a set of points, whose size will correspond to the first simulation and whose color will correspond to the second simulation. In this last representation, we will also provide a title and some specific labels on axes.
fig, ax = gp.initGeographic()
ax.symbol(mygrid,nameSize="Simu.1",nameColor="Simu.2")
ax.decoration(title = "Representing Grid nodes as Points", xlabel="Easting", ylabel="Northing")
plt.show()
We can create a specific ColorScale, containing a limited number of colors, and sampling a given reference Color Scale.
For the next figure, we use the one defaulted by the system ('viridis') for sake of understanding. We simply reduce the number of colors
cmap = gp.getColorMap(5,'viridis')
fig, ax = gp.initGeographic()
ax.raster(mygrid, name="Simu.1", cmap=cmap)
ax.decoration(title="Non-conditional Simulation on Grid")
plt.show()
We can also change the reference color scale (using the one defaulted by the method getColorMap for example: 'gist rainbow') and increase the number of colors
cmap = gp.getColorMap(100)
fig, ax = gp.initGeographic()
ax.raster(mygrid,name="Simu.1",cmap=cmap)
ax.decoration(title="Non-conditional Simulation on Grid (other color scale)")
plt.show()
We have the same type of functions (as demonstrated for points and/or grids) for objects such as histograms of the values collected in the a data base (e.g. a simulation outcome over the grid), calling the relevant function in an object-based fashion
fig, ax = gp.init()
ax.histogram(mygrid,name="Simu.1", bins=30, color='yellow', edgecolor="blue")
plt.show()
Representing a scatter plot between two variables stored on the same Db.
fig, ax = gp.init()
ax.correlation(mygrid,"Simu.1","Simu.2", bins=50, cmin=1)
ax.decoration(title="Scatter Plot between two simulations")
plt.show()
A set of points is sampled from the previous Grid and stored in a new Point Db. The number of samples if fixed to 1% of the number of grid nodes.
mypoint = gl.Db()
mypoint.resetSamplingDb(mygrid,0.01)
mypoint.display()
Data Base Characteristics ========================= Data Base Summary ----------------- File is organized as a set of isolated points Space dimension = 2 Number of Columns = 7 Total number of samples = 17 Number of active samples = 17 Variables --------- Column = 0 - Name = rank - Locator = NA Column = 1 - Name = rank.1 - Locator = NA Column = 2 - Name = x1 - Locator = x1 Column = 3 - Name = x2 - Locator = x2 Column = 4 - Name = Simu.1 - Locator = z1 Column = 5 - Name = Simu.2 - Locator = z2 Column = 6 - Name = sel - Locator = sel
We create a polygon as the convex hull of the samples
mypoly = gl.Polygons.createFromDb(mypoint)
We now display the points (as no variable name is mentioned, the samples are posted using a constant size and color) and the polygon on top of the grid.
fig, ax = gp.initGeographic()
ax.raster(mygrid, name="Simu.1")
ax.polygon(mypoly,flagFace=False, edgecolor='yellow', linewidth=2)
ax.symbol(mypoint, c="black")
plt.show()
We create the same grid as before but with a rotation of 20 degrees.
mygrid = gl.DbGrid.create(nx,dx,x0,[20,0])
err = gl.simtub(None,mygrid,mymodel,nbsimu=2)
fig, ax = gp.initGeographic()
ax.raster(mygrid,name="Simu.1")
plt.show()
A new set of Points is sampled from the rotated Grid. As the same seed is used, the ranks of the selected samples within the grid are the same. Furthermore, we generate the Polygon as the convex hull of the newly created Point db.
mypoint = gl.Db()
mypoint.resetSamplingDb(mygrid,0.01)
mypoly = gl.Polygons.createFromDb(mypoint)
We represent again the three components (grid, points and polygon) on the same view
fig, ax = gp.initGeographic()
ax.raster(mygrid, name="Simu.1")
ax.polygon(mypoly,flagFace=False,edgecolor='yellow')
ax.symbol(mypoint,c="black")
plt.show()
Let us now add a selection in order to restrict the previous representation to the only non-masked samples
tab = mygrid.getColumn("x1")
sel = (np.asarray(tab) < 0).astype(float)
mygrid.addSelection(tuple(sel),'sel')
fig, ax = gp.initGeographic()
ax.raster(mygrid,name="Simu.1",useSel=True)
plt.show()
This paragraph is meant to present the possibility of splitting a figure in two scenes, to represent a grid in each scene (for example) and share the (same) color scale for the two scenes.
fig = plt.figure(figsize=(20,6))
vmin = -4
vmax = +4
ax1 = fig.add_subplot(1,2,1)
ax1.raster(mygrid,name="Simu.1", useSel=False, flagLegend = False, vmin=vmin, vmax=vmax)
ax2 = fig.add_subplot(1,2,2)
ax2.raster(mygrid,name="Simu.2", useSel=False, flagLegend = False, vmin=vmin, vmax=vmax)
fig.subplots_adjust(right=0.7)
cbar_ax = fig.add_axes([0.75, 0.1, 0.02, 0.75])
im = ax.collections[0] # get mappable described by the colorbar
err = fig.colorbar(im, cax = cbar_ax)
In this paragraph, we wish to display sample points with given colors.
The values at sample locations (coordinates and variable values) are provided sample per sample in array called 'tab'.
tab = [1., 1., 1., 2., 2., 3., 3., 3., 5.]
dat1 = gl.Db.createFromSamples(3, gl.ELoadBy.SAMPLE, tab, names=["x","y","z"], locatorNames=["x1","x2","z"])
dbfmt = gl.DbStringFormat()
dbfmt.setFlags(flag_resume=True, flag_array=True)
dat1.display(dbfmt)
Data Base Characteristics
=========================
Data Base Summary
-----------------
File is organized as a set of isolated points
Space dimension = 2
Number of Columns = 4
Total number of samples = 3
Data Base Contents
------------------
rank x y z
[ 0,] 1.000 1.000 1.000 1.000
[ 1,] 2.000 2.000 2.000 3.000
[ 2,] 3.000 3.000 3.000 5.000
Variables
---------
Column = 0 - Name = rank - Locator = NA
Column = 1 - Name = x - Locator = x1
Column = 2 - Name = y - Locator = x2
Column = 3 - Name = z - Locator = z1
We represent the samples without using any pre-specified color map. The system uses the default color scale and assigns the lowest value to be represented to the first color and the largest value to the last color. For all subsequent graphics, the dimensions of the non-geographical plots is fixed.
gp.setDefaultGeographic(dims=[5,4])
fig, ax = gp.initGeographic()
ax.symbol(dat1,nameColor="z",s=200)
plt.show()
Representing using a given color map (based on few colors [5]). The color scale is now discrete but the system still assigns the lowest value (i.e. 1) to the first color and the largest value (i.e. 5) to the last color
ncol = 5
cmap = gp.getColorMap(ncol)
fig, ax = gp.initGeographic()
ax.symbol(dat1,nameColor="z",s=200,cmap=cmap)
plt.show()
We use a new Db where the values at first sample has been modified (from 1 to 4) while the other have been left unchanged. We use the same color scale as before. Again the lowest value (i.e. 3) is assigned to the first color and the largest value (i.e. 5) to the last color.
tab = [1., 1., 4., 2., 2., 3., 3., 3., 5.]
dat2 = gl.Db.createFromSamples(3, gl.ELoadBy.SAMPLE, tab, ["x","y","z"], ["x1","x2","z"])
fig, ax = gp.initGeographic()
ax.symbol(dat2,nameColor="z",s=200,cmap=cmap)
plt.show()
Default dimensions for Geographical plots are set back to larger dimensions for subsequent graphics
gp.setDefaultGeographic(dims=[8,8])
In this section, we demonstrate the possibilities offered by the graphics for working with multiple figures and overlaying graphics. This is described through the use of variograms and models. For this reason we consider the two non-conditional simulations created earllier on the existing grid. We calculate the simple and cross variograms along the two main axes of the grid and fit a model automatically.
varioparam = gl.VarioParam.createMultipleFromGrid(mygrid, npas=10)
vario = gl.Vario(varioparam, mygrid)
err = vario.compute(gl.ECalcVario.VARIOGRAM)
model = gl.Model()
err = model.fit(vario,[gl.ECov.CUBIC])
In the next graphic, we produce the simple variogram of the first variable calculated in the first direction. Note that in all subsequent plots, we will use the same dimension for the elementary plot (which is set as a default)
gp.setDefault(dims=[6,5])
fig, ax = gp.init()
ax.variogram(vario, ivar=0, jvar=0, idir=0)
ax.decoration(title="First (simple) Variable - First Direction")
plt.show()
In the next graphic, we produce a single figure where the variograms of the first variable calculated in the first direction (black) and the second direction (red) are overlaid. The overlay is performed manually.
fig, ax = gp.init()
ax.variogram(vario,idir=0, ivar=0, jvar=0, color="black")
ax.variogram(vario,idir=1, ivar=0, jvar=0, color='red')
ax.decoration(title="First (simple) Variable - Two directions")
plt.show()
In the next graphic, we produce a single graphic where the cross-variograms between first and second variables are displayed for all directions. The colors are extracted from the Color Map provided as argument.
fig, ax = gp.init()
ax.variogram(vario,ivar=1,jvar=0,idir=-1)
ax.decoration(title="Cross-variogram - Two directions")
ax.legend()
plt.show()
In the next figure, we draw the first direction and overlay the second direction (on purpose using two statements). Moreover we also change the dimensions of the plot.
fig, ax = gp.init()
ax.variogram(vario,idir=0,ivar=0,jvar=0,cmap=cmap)
ax.variogram(vario,idir=1,ivar=0,jvar=0)
ax.decoration(title="Overlay test")
plt.show()
All graphic representations that were produced up to now were essentially performed with either a single view (or Axes). The exercise performed earlier with two view in the same figure was handled by the user, calling subplots facility explicaitly.
In the next examples, we will directly consider a figure consituted of many views. This is the case of the multivariate variograms and/or models. In other words, we wish to call the graphic representation, providing a single multivariate variogram and let the facility create all the necessary subplots. Moreover, we want to be able to consider this set of subplots (containing the experimental variograms for example) and overlay the model.
This requires using the other argument (fig) returned by the graphic initialization function (gp.init) which is called passing the expected number of rows and columns to be present in the figure. In our case, as the variogram is calculated for 2 variables, the figure should contain 2 rows and 2 columns.
Representing all simple and cross variograms for all directions.
nv = vario.getVariableNumber()
fig, axs = gp.init(nv,nv)
fig.variogram(vario,ivar=-1,jvar=-1,idir=-1,cmap=cmap)
fig.decoration(title="Simple and Cross Variograms in all directions")
fig.geometry(dims=[6,5])
plt.show()
Represent the Model calculated for the second variable. If the Model is not isotropic, the plot should differ per direction: as direction has not been mentionned, the first direction (of the geographic sysrtem) is used by default.
fig, ax = gp.init()
ax.model(model,ivar=1,jvar=1)
ax.decoration(title="Model for the Second Variable in First Direction")
plt.show()
Representing all simple and cross variograms together with the fitted model for the all directions. This is directly provided by the function varmod.
nv = vario.getVariableNumber()
fig, axs = gp.init(nv,nv)
fig.varmod(vario=vario, model=model, cmap=cmap)
fig.decoration(title="All variograms and Model")
plt.show()
The next figure is meant to demonstrate the overlay possibilities. We first represent the experimental variograms for all variable (in the first direction only to be legible). Then we overlay the model ... only over the experimental simple variogram of the second variable (in dashed blue).
nv = vario.getVariableNumber()
fig, axs = gp.init(nv,nv)
fig.variogram(vario=vario,idir=0,ivar=-1,jvar=-1,cmap=cmap)
axs[1,1].model(model,ivar=1,jvar=1,codir=vario.getCodirs(0),hmax = vario.getHmax(),
linestyle = 'dashed', color='blue')
fig.decoration(title="Overlay test")
plt.show()
