Scientific Code

agilegeo.wavelet

agilegeo.wavelet.ricker(duration, dt, f)

Also known as the mexican hat wavelet, models the function: A = (1-2 pi^2 f^2 t^2) e^{-pi^2 f^2 t^2}

Parameters:
  • duration – The length in seconds of the wavelet.
  • dt – is the sample interval in seconds (usually 0.001, 0.002, 0.004)
Params f:

Center frequency of the wavelet (in Hz). If a list or tuple is passed, the first element will be used.

Returns:

ricker wavelets with center frequency f sampled at t.

agilegeo.wavelet.ormsby(duration, dt, f)

The Ormsby wavelet requires four frequencies: f1 = low-cut frequency f2 = low-pass frequency f3 = high-pass frequency f4 = hi-cut frequency Together, the frequencies define a trapezoid shape in the spectrum. The Ormsby wavelet has several sidelobes, unlike Ricker wavelets which only have two, one either side.

Parameters:
  • duration – The length in seconds of the wavelet.
  • dt – is the sample interval in seconds (usually 0.001, 0.002, 0.004)
Params f:

Tuple of form (f1,f2,f3,f4), or a similar list.

Returns:

A vector containing the ormsby wavelet

agilegeo.wavelet.sweep(duration, dt, f, method='linear', phi=0, vertex_zero=True, autocorrelate=True)

Generates a linear frequency modulated wavelet (sweep) Does a wrapping of scipy.signal.chirp

Parameters:
  • duration – The length in seconds of the wavelet.
  • dt – is the sample interval in seconds (usually 0.001, 0.002, 0.004)
  • f – Tuple of (f1, f2), or a similar list. A list of lists will create a wavelet bank.
  • method – {‘linear’,’quadratic’,’logarithmic’}, optional
  • phi – float, phase offset in degrees
  • vertex_zero – bool, optional This parameter is only used when method is ‘quadratic’. It determines whether the vertex of the parabola that is the graph of the frequency is at t=0 or t=t1.
Returns:

An LFM waveform.

agilegeo.wavelet.rotate_phase(w, phi)

Performs a phase rotation of wavelet using: A = w(t)Cos(phi) + h(t)Sin(phi) Where w(t) is the wavelet and h(t) is it’s hilbert transform.

Params w:The wavelet vector.
Params phi:The phase rotation angle (in Radians) to apply.
Returns:The phase rotated signal.

agilegeo.avo

agilegeo.avo.zoeppritz(vp1, vs1, rho1, vp0, vs0, rho0, theta1)

Full Zoeppritz solution, considered the definitive solution. Calculates the angle dependent p-wave reflectivity of an interface between two mediums.

Parameters:
  • vp1 – The p-wave velocity of the upper medium.
  • vs1 – The s-wave velocity of the upper medium.
  • rho1 – The density of the upper medium.
  • vp0 – The p-wave velocity of the lower medium.
  • vs0 – The s-wave velocity of the lower medium.
  • rho0 – The density of the lower medium.
  • theta1 – An array of incident angles to use for reflectivity calculation [degrees].
Returns:

a vector of len(theta1) containing the reflectivity value corresponding to each angle.

agilegeo.avo.akirichards(vp1, vs1, rho1, vp2, vs2, rho2, theta1)

This is the formulation from Avseth et al., Quantitative seismic interpretation, Cambridge University Press, 2006. Adapted for a 4-term formula. See http://subsurfwiki.org/wiki/Aki-Richards_equation

Parameters:
  • vp1 – The p-wave velocity of the upper medium.
  • vs1 – The s-wave velocity of the upper medium.
  • rho1 – The density of the upper medium.
  • vp2 – The p-wave velocity of the lower medium.
  • vs2 – The s-wave velocity of the lower medium.
  • rho2 – The density of the lower medium.
  • theta1 – An array of incident angles to use for reflectivity calculation [degrees].
Returns:

a vector of len(theta1) containing the reflectivity value corresponding to each angle.

agilegeo.avo.akirichards_alt(vp1, vs1, rho1, vp2, vs2, rho2, theta1)

This is another formulation of the Aki-Richards solution. See http://subsurfwiki.org/wiki/Aki-Richards_equation

Parameters:
  • vp1 – The p-wave velocity of the upper medium.
  • vs1 – The s-wave velocity of the upper medium.
  • rho1 – The density of the upper medium.
  • vp2 – The p-wave velocity of the lower medium.
  • vs2 – The s-wave velocity of the lower medium.
  • rho2 – The density of the lower medium.
  • theta1 – An array of incident angles to use for reflectivity calculation [degrees].
Returns:

a vector of len(theta1) containing the reflectivity value corresponding to each angle.

agilegeo.avo.fatti(vp1, vs1, rho1, vp2, vs2, rho2, theta1)

Compute reflectivities with Fatti’s formulation of the Aki-Richards equation, which does not account for the critical angle. Fatti et al (1994), Geophysics 59 (9).

Parameters:
  • vp1 – The p-wave velocity of the upper medium.
  • vs1 – The s-wave velocity of the upper medium.
  • rho1 – The density of the upper medium.
  • vp2 – The p-wave velocity of the lower medium.
  • vs2 – The s-wave velocity of the lower medium.
  • rho2 – The density of the lower medium.
  • theta1 – An array of incident angles to use for reflectivity calculation [degrees].
Returns:

a vector of len(theta1) containing the reflectivity value corresponding to each angle.

agilegeo.avo.shuey2(vp1, vs1, rho1, vp2, vs2, rho2, theta1)

Compute Shuey approximation with 2 terms. http://subsurfwiki.org/wiki/Shuey_equation

Parameters:
  • vp1 – The p-wave velocity of the upper medium.
  • vs1 – The s-wave velocity of the upper medium.
  • rho1 – The density of the upper medium.
  • vp2 – The p-wave velocity of the lower medium.
  • vs2 – The s-wave velocity of the lower medium.
  • rho2 – The density of the lower medium.
  • theta1 – An array of incident angles to use for reflectivity calculation [degrees].
Returns:

a vector of len(theta1) containing the reflectivity value corresponding to each angle.

agilegeo.avo.shuey3(vp1, vs1, rho1, vp2, vs2, rho2, theta1)

Compute Shuey approximation with 3 terms. http://subsurfwiki.org/wiki/Shuey_equation

Parameters:
  • vp1 – The p-wave velocity of the upper medium.
  • vs1 – The s-wave velocity of the upper medium.
  • rho1 – The density of the upper medium.
  • vp2 – The p-wave velocity of the lower medium.
  • vs2 – The s-wave velocity of the lower medium.
  • rho2 – The density of the lower medium.
  • theta1 – An array of incident angles to use for reflectivity calculation [degrees].
Returns:

a vector of len(theta1) containing the reflectivity value corresponding to each angle.

agilegeo.avo.bortfeld2(vp1, vs1, rho1, vp2, vs2, rho2, theta1)

The 2-term Bortfeld approximation for ava analysis.

Parameters:
  • vp1 – The p-wave velocity of the upper medium.
  • vs1 – The s-wave velocity of the upper medium.
  • rho1 – The density of the upper medium.
  • vp2 – The p-wave velocity of the lower medium.
  • vs2 – The s-wave velocity of the lower medium.
  • rho2 – The density of the lower medium.
  • theta1 – An array of incident angles to use for reflectivity calculation [degrees].
Returns:

a vector of len(theta1) containing the reflectivity value corresponding to each angle.

agilegeo.avo.bortfeld3(vp1, vs1, rho1, vp2, vs2, rho2, theta1)

Compute Bortfeld approximation with three terms. http://sepwww.stanford.edu/public/docs/sep111/marie2/paper_html/node2.html

Parameters:
  • vp1 – The p-wave velocity of the upper medium.
  • vs1 – The s-wave velocity of the upper medium.
  • rho1 – The density of the upper medium.
  • vp2 – The p-wave velocity of the lower medium.
  • vs2 – The s-wave velocity of the lower medium.
  • rho2 – The density of the lower medium.
  • theta1 – An array of incident angles to use for reflectivity calculation [degrees].
Returns:

a vector of len(theta1) containing the reflectivity value corresponding to each angle.

agilegeo.avo.time_to_depth(data, vmodel, dt, dz)

Converts data from the time domain to the depth domain given a velocity model.

Parameters:
  • data – The data to convert, will work with a 1 or 2D numpy numpy array. array(samples,traces).
  • vmodel – P-wave velocity model that corresponds to the data. Must be the same shape as data.
  • dt – The sample interval of the input data [s].
  • dz – The sample interval of the output data [m].
Returns:

The data resampled in the depth domain.

agilegeo.avo.depth_to_time(data, vmodel, dz, dt)

Converts data from the time domain to the depth domain given a velocity model.

Parameters:
  • data – The data to convert, will work with a 1 or 2D numpy numpy array. array(samples,traces).
  • vmodel – P-wave velocity model that corresponds to the data. Must be the same shape as data.
  • dz – The sample interval of the input data [m].
  • dt – The sample interval of the output data [s].
Returns:

The data resampled in the time domain.

modelr.rock_properties

Container for physical rock properties

class modelr.rock_properties.RockProperties(vp, vs=None, rho=None, vp_sig=0, vs_sig=0, rho_sig=0)[source]

Class to store rock properties.

Parameters:
  • vp – pressure wave velocity
  • vs – shear wave velocity
  • rho – bulk density

Table Of Contents

Previous topic

Models and Calculations

Next topic

Plotting Code

This Page