Python interface

From GeopsyWiki
Revision as of 21:21, 11 January 2024 by Marc (talk | contribs) (→‎Example)
Jump to navigation Jump to search

Some functionalities are accessible from a Python script. Currently there are two modules:

  • GeopsyPyCoreWave: surface wave computation
  • GeopsyPyScifigs: graphical tools to plot data

Other modules will be probably proposed in the future (e.g. Neighborhood inversion).

Installation

The modules are compiled with Geopsy package if Python include files are in the include paths. At configure step, add Python include path (among other options not listed here, see general instructions specific to your platform):

 ./configure -I /usr/include/python3.9

Currently, modules are tested on Linux only. Under Windows and Mac OS, it is probably not working.

Upon completion of the build process, the modules are dynamic libraries located in the lib directory of the install directory. Their names start with "lib" which prevents the recognition by Python. We suggest to create a new directory at the same level as lib and to modify the library names as follow:

 mkdir python
 cp lib/libGeopsyPyCoreWave.so python/GeopsyPyCoreWave.so
 cp lib/libGeopsyPySciFigs.so python/GeopsyPySciFigs.so

Finally, you have to add the path to Python

 export PYTHONPATH=/path/where/geopsy/is/installed/python:$PYTHONPATH

Interface documentation

The list of available functions in each module is provided with the help function. Their required arguments are also shortly described.

 import GeopsyPyCoreWave as gp
 help(gp)

Example

In this section, we describe a small example (SciFigs-example-01.py) to compute dispersion and ellipticity curves, and to plot them on the screen, in a .page file file and in a PDF file.

First import of the required modules:

 import GeopsyPySciFigs as sf
 import numpy as np
 import GeopsyPyCoreWave as gp

Construction of the NumPY arrays which define the 1D layered model and the frequency vector where to compute the curves.

 nmodes=7
 freq=np.logspace(np.log10(1), np.log10(50), 500)
 omega=freq.copy()
 omega*=2*np.pi
 h=np.array([10,20])
 vp=np.array([500,1000,3000])
 vs=np.array([200,600,1500])
 rho=np.array([2000,2000,2500])

Computation of dispersion and ellipticity curves.

 slowRayleigh=gp.rayleighDispersionCurve(nmodes, 0, h, vp, vs, rho, omega)
 slowRayleighGroup=gp.rayleighDispersionCurve(nmodes, 1, h, vp, vs, rho, omega)
 ellRayleigh=gp.rayleighEllipticityCurve(nmodes, h, vp, vs, rho, omega)
 ellRayleigh=np.arctan(ellRayleigh)*(180/np.pi)
 slowLove=gp.loveDispersionCurve(nmodes, 0, h, vs, rho, omega)

Definition of useful functions to set plot attributes.

 def layoutAttributes(p, x, y, w, h):
   p.setAnchor("TopRight")
   p.xAxis().setSizeInfo(w)
   p.yAxis().setSizeInfo(h)
   p.setPrintX(x)
   p.setPrintY(y)
 def frequencyAttributes(p):
   p.xAxis().setMin(1)
   p.xAxis().setMax(50)
   p.xAxis().setScaleType("Log")
   p.xAxis().setTitle("Frequency (Hz)")
 def dispersionAttributes(p):
   frequencyAttributes(p)
   p.yAxis().setScaleType("InverseLog")
   p.yAxis().setMin(1/1500)
   p.yAxis().setMax(1/150)

Creation of a blank GraphicSheet.

 s=sf.newSheet()

Add a dispersion curve plot with Rayleigh phase velocity.

 prp=sf.newPlot(s)
 sf.addCurves(prp, freq, slowRayleigh)
 layoutAttributes(prp, 12, 1, 11 ,6)
 dispersionAttributes(prp)
 prp.yAxis().setTitle("Rayleigh phase velocity (m/s)")

Add a dispersion curve plot with Rayleigh group velocity.

 prg=sf.newPlot(s)
 sf.addCurves(prg, freq, slowRayleighGroup)
 layoutAttributes(prg, 23, 1, 11 ,6)
 dispersionAttributes(prg)
 prg.yAxis().setTitle("Rayleigh group velocity (m/s)")

Add an ellipticity curve plot.

 per=sf.newPlot(s)
 sf.addCurves(per, freq, ellRayleigh)
 layoutAttributes(per, 12, 7, 11 ,6)
 frequencyAttributes(per)
 per.yAxis().setTitle("Angular ellipticity (deg)")
 per.yAxis().setMin(-90)
 per.yAxis().setMax(90)
 per.graphContents().layer(0).setSignThreshold(45)

Add a dispersion curve with Love phase velocity

 pl=sf.newPlot(s)
 sf.addCurves(pl, freq, slowLove)
 layoutAttributes(pl, 23, 7, 11 ,6)
 dispersionAttributes(pl)
 pl.yAxis().setTitle("Love phase velocity (m/s)")

Output results/

 s.setPaperOrientation("Landscape")
 s.fileSave_2("/tmp/dc.page")
 s.exportImage_2("/tmp/dc.pdf")