Type: | Package |
Title: | Seismic Reflection and Scattering Coefficients |
Version: | 1.0-9 |
Date: | 2023-08-19 |
Author: | Jonathan M. Lees [aut, cre] |
Maintainer: | Jonathan M. Lees <jonathan.lees@unc.edu> |
Description: | Calculate and plot scattering matrix coefficients for plane waves at interface. |
License: | GPL-2 | GPL-3 [expanded from: GPL (≥ 2)] |
NeedsCompilation: | no |
Packaged: | 2023-08-19 10:48:49 UTC; lees |
Repository: | CRAN |
Date/Publication: | 2023-08-19 16:10:02 UTC |
Zoeppritz Equations
Description
Calculate and plot scattering matrix coefficients for plane waves at interface.
Author(s)
Jonathan M. Lees<jonathan.lees@unc.edu> Maintainer: Jonathan M. Lees<jonathan.lees@unc.edu>
References
Young,~G.B., Braile, L. W. 1976. A computer program for the application of Zoeppritz's amplitude equations and Knott's energy equations, Bulletin of the Seismological Society of America, vol.66, no.6,1881-1885.
K. Aki and P.G. Richards.Quantitative seismology. University Science Books, Sausalito, Calif., 2nd edition, 2002.
Examples
##### demo(ZOEP)
############# Incident wave in Low velocity layer
alpha1 = 4.98
beta1 = 2.9
rho1 = 2.667
alpha2 = 8.0
beta2 = 4.6
rho2 = 3.38
App = pzoeppritz( "Amplitude" , alpha1, alpha2, beta1, beta2, rho1 ,rho2, "P", "ALL");
App = pzoeppritz( "Amplitude" , alpha1, alpha2, beta1, beta2, rho1 ,rho2, "S", "ALL");
App = pzoeppritz( "Energy" , alpha1, alpha2, beta1, beta2, rho1 ,rho2, "P", "ALL");
App = pzoeppritz( "Potential" , alpha1, alpha2, beta1, beta2, rho1 ,rho2, "P", "ALL");
############# Incident wave in high velocity layer
alpha1 = 8.0
beta1 = 4.6
rho1 = 3.38
alpha2 = 4.98
beta2 = 2.9
rho2 = 2.667
App = pzoeppritz( "Amplitude" , alpha1, alpha2, beta1, beta2, rho1 ,rho2, "P", "ALL");
App = pzoeppritz( "Amplitude" , alpha1, alpha2, beta1, beta2, rho1 ,rho2, "S", "ALL");
Show Scattering Diagram Cartoon
Description
Adds a small diagram showing two layers and labeled scattered ray paths.
Usage
piczoeppritz(LL = list(x = c(0, 1), y = c(0, 1)), chincw = "P")
Arguments
LL |
Bounds of Box for plotting |
chincw |
character for incident wave |
Details
This code simply adds a small cartoon showing incoming and outgoing waves in scattering matrix.
Value
Graphical side effect.
Author(s)
Jonathan M. Lees<jonathan.lees@unc.edu>
See Also
pzoeppritz
Examples
plot(c(0,1), c(0,1), type='n')
piczoeppritz(LL = list(x = c(0.5, 0.75), y = c(0.5, 0.75)), chincw = "P")
Plot Scattering (Zoeppritz) Coefficients
Description
Plot Scattering (Zoeppritz) Coefficients
Usage
plotzoeppritz(A, zoepcols = c("red", "green", "blue", "purple"), zoeplty = c(1, 1, 1, 1))
Arguments
A |
list output of pzoeppritz or zoeppritz |
zoepcols |
vector of 4 colors |
zoeplty |
vector of 4 line types |
Details
Used to plot the matrix of scattering coefficients with different colros and/or line types.
Value
Graphical side effects.
Author(s)
Jonathan M. Lees<jonathan.lees@unc.edu>
See Also
zoeppritz, pzoeppritz, piczoeppritz
Examples
alpha1 = 4.98
beta1 = 2.9
rho1 = 2.667
alpha2 = 8.0
beta2 = 4.6
rho2 = 3.38
############### create the scattering matrix:
App = pzoeppritz( "Amplitude" , alpha1, alpha2, beta1, beta2, rho1 ,rho2, "P", "NONE");
######### plot
plotzoeppritz(App)
Plot Scattering (Zoeppritz) Coefficients
Description
Calculate and plot the P and S-wave scattering amplitudes for a plane wave at an interface.
Usage
pzoeppritz(chtype = "Amplitude", alpha1, alpha2,
beta1, beta2, rho1, rho2, chincw = "P",
choutkind = "ALL")
Arguments
chtype |
character, type of output, one of: Amplitude, Potential, Energy |
alpha1 |
P-wave Velocity of Upper Layer, km/s |
alpha2 |
P-wave Velocity of Lower Layer, km/s |
beta1 |
S-wave Velocity of Upper Layer, km/s |
beta2 |
S-wave Velocity of Lower Layer, km/s |
rho1 |
Density of Upper Layer, kg/m3 |
rho2 |
Density of Lower Layer, kg/m3 |
chincw |
Incident Wave: P, S |
choutkind |
character, type of out put one of: P, S, ALL, NONE |
Details
Front end for zoeppritz program.
Value
List output of zoeppritz call:
angle |
Incident angles, degrees |
rmat |
Matrix of 4 by n reflection coefficients for each angle |
rra |
Matrix of 4 by n real part of scattering matrix |
rra |
Matrix of 4 by n imaginary part of scattering matrix |
ang |
Matrix of 4 by n phase angle |
incw |
integer, from input parameter |
icoef |
integer, from input parameter |
alphacrit |
critical angle |
Note
This front end is easier to call because it is more verbose. Creates a plot of the coefficients versus incident angle. If coefficients are complex, they are replaced with NA and they are thus not plotted.
Author(s)
Jonathan M. Lees<jonathan.lees@unc.edu>
References
Young,~G.B., Braile, L. W. 1976. A computer program for the application of Zoeppritz's amplitude equations and Knott's energy equations, Bulletin of the Seismological Society of America, vol.66, no.6,1881-1885.
K. Aki and P. G. Richards. Quantitative seismology. University Science Books, Sausalito, Calif., 2nd edition, 2002.
See Also
zoeppritz, pzoeppritz, piczoeppritz
Examples
############# Incident wave in Low velocity layer
alpha1 = 4.98
beta1 = 2.9
rho1 = 2.667
alpha2 = 8.0
beta2 = 4.6
rho2 = 3.38
App = pzoeppritz( "Amplitude" , alpha1, alpha2,
beta1, beta2, rho1 ,rho2, "P", "ALL");
App = pzoeppritz( "Amplitude" , alpha1, alpha2,
beta1, beta2, rho1 ,rho2, "S", "ALL");
############# Incident wave in high velocity layer
alpha1 = 8.0
beta1 = 4.6
rho1 = 3.38
alpha2 = 4.98
beta2 = 2.9
rho2 = 2.667
App = pzoeppritz( "Amplitude" , alpha1, alpha2,
beta1, beta2, rho1 ,rho2, "P", "ALL");
App = pzoeppritz( "Amplitude" , alpha1, alpha2,
beta1, beta2, rho1 ,rho2, "S", "ALL");
Zoeppritz Equations
Description
Calculate the P and S-wave scattering amplitudes for a plane wave at an interface.
Usage
zoeppritz(icoef, vp1, vp2, vs1, vs2, rho1, rho2, incw)
Arguments
icoef |
type of out put Amplitude=1, Potential=2, Energy=3 |
vp1 |
P-wave Velocity of Upper Layer, km/s |
vp2 |
P-wave Velocity of Lower Layer, km/s |
vs1 |
S-wave Velocity of Upper Layer, km/s |
vs2 |
S-wave Velocity of Lower Layer, km/s |
rho1 |
Density of Upper Layer, kg/m3 |
rho2 |
Density of Lower Layer, kg/m3 |
incw |
integer,Incident Wave: P=1, S=2 |
Details
Coeficiants are calculated at angles from 0-90 degrees. Zero is returned where coefficients are imaginary.
Value
List:
angle |
Incident angles (degrees) |
rmat |
Matrix of 4 by n reflection coefficients for each angle |
rra |
Matrix of 4 by n real part of scattering matrix |
rra |
Matrix of 4 by n imaginary part of scattering matrix |
ang |
Matrix of 4 by n phase angle |
incw |
integer, from input parameter |
icoef |
integer, from input parameter |
Note
Based on the fortran algorithm in Young and Braile. Uses a linear approximation by Aki and Richards.
Author(s)
Jonathan M. Lees<jonathan.lees@unc.edu>
References
Young, G.B., Braile, L. W. 1976. A computer program for the application of Zoeppritz's amplitude equations and Knott's energy equations, Bulletin of the Seismological Society of America, vol.66, no.6,1881-1885.
K. Aki and P.G. Richards.Quantitative seismology. University Science Books, Sausalito, Calif., 2nd edition, 2002.
See Also
pzoeppritz, plotzoeppritz
Examples
######### set up 2-layer model
alpha1 = 4.98
beta1 = 2.9
rho1 = 2.667
alpha2 = 8.0
beta2 = 4.6
rho2 = 3.38
################### P-wave incident = 1
incw=1;
icoef=1
A = zoeppritz(icoef, alpha1, alpha2, beta1, beta2, rho1,rho2, incw)
plot(A$angle, A$rmat[,1], xlab="Incident Angle", ylab="Ratio of Amplitudes",
main="P-wave incident/P-wave Reflected" )
plot(A$angle, A$rmat[,2], xlab="Incident Angle", ylab="Ratio of Amplitudes",
main="P-wave incident/S-wave Reflected" )
plot(A$angle, A$rmat[,3], xlab="Incident Angle", ylab="Ratio of Amplitudes",
main="P-wave incident/P-wave Refracted" )
plot(A$angle, A$rmat[,4], xlab="Incident Angle", ylab="Ratio of Amplitudes",
main="P-wave incident/S-wave Refracted" )
################### S-wave incident = 2
incw=2
icoef=1
A = zoeppritz(icoef, alpha1, alpha2, beta1, beta2, rho1,rho2, incw)
plot(A$angle, A$rmat[,1], xlab="Incident Angle", ylab="Ratio of Amplitudes",
main="S-wave incident/P-wave Reflected" )
plot(A$angle, A$rmat[,2], xlab="Incident Angle", ylab="Ratio of Amplitudes",
main="S-wave incident/S-wave Reflected" )
plot(A$angle, A$rmat[,3], xlab="Incident Angle", ylab="Ratio of Amplitudes",
main="S-wave incident/P-wave Refracted" )
plot(A$angle, A$rmat[,4], xlab="Incident Angle", ylab="Ratio of Amplitudes",
main="S-wave incident/S-wave Refracted" )