density_functions module

class density_functions.Convert(profile, profile_type=None, region_name=None, seismic_method_name=None, geotherm=None)[source]

Bases: object

Convert seismic velocity profiles to various parameters.

This class provides methods to read seismic velocity profiles, convert them to other parameters using different approaches, and write the converted data to output files.

Parameters:

  • profile (str) – File path to the seismic profile data.

  • profile_type (str) – Type of the seismic profile, “Vp” or “Vs.”

  • region_name (str, optional) – Geographic location of the profile. Default is None.

  • seismic_method_name (str) – Method used to acquire the velocity profile. If set to None, the read_data method will pick up the argument from the file string. Note if set to None the strict file convention must be set (see README.md) and tutorial_1.ipynb.

  • geotherm (class instance) – Instance of the Geotherm class containing information about the temperature profile at the site of the seismic profile.

data

Dictionary containing parsed seismic profile data.

Type:

dict

profile_type

Type of profile. “Vp” or “Vs”.

Type:

str

moho

Moho depth parsed from the profile data.

Type:

float

geotherm

Instance of the Geotherm class (if using temperature-dependent conversion. Default None)

Type:

class instance

region_name

Geographic location of the file (parsed from file string or given as argument)

Type:

str

method

Method used to collect the velocity profile. Parsed from file string or given as argument (e.g. ‘REFRACTION’).

Type:

str

rho, av_rho, rho_hi_res

Attributes generated by the Vs_to_Vp_brocher, Vp_to_density_brocher, and V_to_density_stephenson methods.

Type:

various types

vp_calc

Calculated vp profile generated by Vs_to_Vp_brocher method.

Type:

np.ndarray

V_to_density_stephenson(parameters, dz=0.1, constant_depth=None, constant_density=None, T_dependence=True, plot=False)[source]

Convert seismic velocity profiles to density using the Stephenson method.

This method takes seismic velocity profiles and converts them to density profiles using the Stephenson density conversion approach. It provides options for specifying conversion parameters, depth increment, constant density values, temperature dependence, and plotting. The resulting density fields are added to the ‘self.data’ dictionary.

Parameters:

  • parameters (Constants or TemperatureDependentConstants class instance) – Parameters for density conversion.

  • dz (float, optional) – Depth increment for density calculation (default is 0.1 km).

  • constant_depth (float, optional) – Depth for constant density (km).

  • constant_density (float, optional) – Constant density value.

  • T_dependence (bool, optional) – Include temperature dependence (default is True).

  • plot (bool, optional) – Whether to plot the density and pressure profiles (default is False).

Return type:

None

Vp_to_density_brocher()[source]

Convert Vp profile to density using Brocher (2005) method.

This method converts a seismic Vp profile to density using the Brocher (2005) method. It checks whether a Vp or calculated Vp profile exists and performs the conversion accordingly. The resulting density fields are added to the self.data dictionary.

Raises:

KeyError – If no Vp or calculated Vp profile is found when converting a Vs profile. In such cases, it reminds the user to convert Vs to Vp first.

Return type:

None

Vs_to_Vp_brocher()[source]

Convert Vs profile to Vp profile using Brocher’s (2005) approach.

This method converts a seismic Vs profile to a Vp profile using Brocher’s (2005) approach. It updates the self.data dictionary to include the Vp_calc profile field, average Vp, and Vp/Vs ratio.

Raises:

  • Exception – If the profile type is not “Vs,” indicating that you must be using a Vs profile for conversion. In such cases, it advises the user to re-initiate the class.

  • NameError – If the data dictionary has not been created yet, it reminds the user to run read_data first to create the data dictionary.

Return type:

None

read_data()[source]

Read in data file and parse it into a data dictionary.

This method reads the seismic profile data file specified by self.profile, extracts relevant information, and organizes it into a dictionary stored in self.data.

The method handles the following tasks:

  1. Extracts method type from the file path (e.g., refraction,

    reflection, RF).

  2. Parses the record header data, including station name, location

    coordinates, and Moho depth.

  3. Converts depth-velocity data into a 2D array and reorders columns.

  4. Ensures the velocity profile extends to the Moho depth.

  5. Checks for any velocity step at the Moho depth and removes it

    so that we do not sample mantle velocities.

  6. Corrects for potential data capture errors by ensuring depth is

    monotonically increasing.

  7. Calculates the average velocity for the profile.

  8. Constructs a data dictionary based on the profile type (Vs or Vp)

    and stores it in self.data.

Returns:

The parsed data is stored in self.data for further use.

Return type:

None

Notes

This method assumes that self.profile_type has been set to “Vs” or “Vp” to indicate the type of seismic profile being read.

write_data(path, file_structure=None, approach='stephenson', T_dependence=False)[source]

Write seismic profile data to appropriate file locations.

This method takes the seismic profile data stored in the class instance and writes it to the correct file locations based on the specified density conversion approach and temperature dependence settings.

Parameters:

  • path (str) – The path to the directory containing all velocity data. See README for directory structure details.

  • file_structure (str, optional) – If set to None (Default) then we will construct the path from the information in the metadata (i.e. following the default file structure). Otherwise a manual outpath needs to be used that leads to the output location. We will then append relevant method information to the output filename to bookmark the output profiles.

  • approach (str, optional) – The density (and Vp-Vs if used) conversion approach, which is needed for the file path (default is “stephenson”).

  • T_dependence (bool, optional) – Specifies whether temperature dependence is included (default is False).

Return type:

None

class density_functions.MultiConversion(path, which_location='ALL', write_data=False, approach='stephenson', parameters=None, master_geotherm=None, constant_depth=None, constant_density=None, T_dependence=False)[source]

Bases: object

Wrapper class to extract multiple files from path to send to Convert class using specified density conversion approach. Check that all required arguments have been provided.

Parameters:

  • path (str) – The master directory where all data are stored in directories named after their location.

  • which_location (str or list, optional) – Determines which locations the user wants to convert. Defaults to “ALL,” indicating that all locations will be converted. If specific locations are desired, provide a list of location names.

  • write_data (bool, optional) – Specifies whether to write the converted data to files. Defaults to False.

  • approach (str, optional) – The density conversion approach to use. Options are “stephenson” or “brocher.” Defaults to “stephenson.”

  • parameters (class instance, optional) – Class instance of the Constants class. Must be provided if using approach ‘stephenson’. Must contain a material_constants attribute if T_dependence is True.

  • master_geotherm (instance of Geotherm class, optional) – Used as a reference or template for other operations. When the master attribute of master_geotherm is True, deep copies are made and parameters are updated for all individual profiles.

  • constant_depth (float, optional) – The depth (from the surface) over which to use a constant density value (in kilometers). Defaults to None.

  • constant_density (float, optional) – The value of the constant density to use for the uppermost few kilometers (in Mg/m3), if constant_depth is set. Defaults to None.

  • T_dependence (bool, optional) – Determines whether to include temperature dependence of velocity to density conversion, including thermal expansion and compressibility. Defaults to False.

path

The master directory where all data are stored.

Type:

str

which_location

The locations to convert.

Type:

str or list

write_data

Whether to write the converted data to files.

Type:

bool

approach

The density conversion approach to use.

Type:

str

parameters

Instance of the Constants class.

Type:

class instance

master_geotherm

Used as a reference or template for other operations.

Type:

instance of Geotherm class

constant_depth

The depth over which to use a constant density value.

Type:

float

constant_density

The value of the constant density to use for the uppermost few kilometers.

Type:

float

T_dependence

Whether to include temperature dependence of velocity to density conversion.

Type:

bool

convert_metadata

Metadata for the conversion process. This is a list of dictionaries that is appended to in the process_file_list method. It is used to store information about the conversion process for each profile.

Type:

list

assemble_file_lists()[source]

Assemble file lists and necessary parameters for density conversion.

This method prepares the necessary data and parameters for density conversion by assembling file lists, checking provided arguments, and setting up parameters based on the chosen conversion approach. It also initializes the convert_metadata attribute for storing conversion metadata.

Parameters:

None –

Returns:

None

profiles_to_dataframe()[source]

Convert station profiles to a pandas DataFrame for bulk information.

This method iterates over the station profiles stored in the station_profiles attribute of the class instance. For each profile, it extracts the ‘station’, ‘location’, ‘av_rho’, ‘moho’, ‘region’, ‘av_Vp’, ‘av_Vs’, and ‘av_Vp_calc’ attributes (if they exist), and stores them in separate lists. These lists are then converted to numpy arrays.

The method then creates a pandas DataFrame from these arrays, with each array forming a column in the DataFrame. The DataFrame is stored in the station_profiles_df attribute of the class instance, and is also returned by the method.

Returns:

A pandas DataFrame where each row represents a station profile, and the columns represent the ‘station’, ‘lon’, ‘lat’, ‘moho’, ‘av_vp’, ‘av_vs’, ‘av_rho’, and ‘region’ attributes of the profile. ‘lon’ and ‘lat’ are derived from the ‘location’ attribute of the profile, which is assumed to be a two-element array or list. If ‘av_Vp’, ‘av_Vs’, or ‘av_Vp_calc’ do not exist in a profile, their values in the DataFrame are set to NaN.

Return type:

pd.DataFrame

send_to_conversion_function(parallel=False)[source]

Send data to the velocity conversion function and obtain station profiles.

This method initiates the velocity conversion process for each profile and assembles a list of data dictionaries containing velocity, Moho, density information, etc. The method provides an option for parallel processing but currently has limitations due to pickling issues. It works with both parallel and single-processor modes.

Parameters:

parallel (bool, optional) – Whether to use parallel processing (default is False). Note that this functionality is not currently available owing to issues with the class-based architecture for handling constants and some non-pickleable objects. This bug will be fixed in a later release. Setting parallel=True will trigger an exception.

Returns:

An array containing station profiles for further use.

Return type:

np.ndarray

class density_functions.V2RhoStephenson(data, parameters, profile='Vp', constant_depth=None, constant_density=None, T_dependence=True, geotherm=None)[source]

Bases: object

Class for density conversion of a given seismic profile using the Stephenson method. Only handles single profiles.

If you are converting a single profile, it is easiest to use the functionality in this class through the convert_V_profile function. If you are converting multiple profiles, then it is easiest to access the functionality in this class through the MultiConversion class.

Parameters:

  • data (dict) – Seismic profile data containing depth, velocity, and density. This is the output from the read_data method of the Convert class.

  • parameters (numpy.ndarray) – Instance of the Constants, VpConstants or VsConstants class. If T_dependence is True, then this must be an instance of Constants class and must contain material_constants attribute.

  • profile (str) – Profile type, “Vp” or “Vs.”

  • constant_depth (float) – Depth for constant density for uppermost x km.

  • constant_density (float) – Constant density value for uppermost x km.

  • T_dependence (bool) – True for temperature dependence inclusion (default True).

  • geotherm (list) – Instance of the Geotherm class, must be set if T_dependence is True (default None).

data

Seismic profile data containing depth, velocity, and density. This is the output from the read_data method of the Convert class.

Type:

dict

parameters

Instance of the Constants, VpConstants or VsConstants class. If T_dependence is True, then this must be an instance of Constants class and must contain material_constants attribute.

Type:

numpy.ndarray

profile

Profile type, “Vp” or “Vs.”

Type:

str

constant_depth

Depth for constant density for uppermost x km.

Type:

float

constant_density

Constant density value for uppermost x km.

Type:

float

T_dependence

True for temperature dependence inclusion (default True).

Type:

bool

geotherm

Instance of the Geotherm class, must be set if T_dependence is True (default None).

Type:

list

calculate_density_profile(dz=0.1)[source]

Calculate the density profile.

This method calculates the density profile based on the provided seismic data and parameters, taking into account temperature dependence if specified.

Parameters:

dz (float, optional) – Depth interval for calculations (default is 0.1 km).

Returns:

data – Seismic profile data including density. Also returns as an attribute of the class.

Return type:

dict

density_functions.V2rho_stephenson(data, parameters)[source]

Convert seismic wave velocity (Vp or Vs) into density at standard temperature and pressure (s.t.p.) given pressure. This function is based on mineral physics and accounts for pressure effects on velocity and density. WARNING: this function returns density at standard temperature and pressure. Please correct for thermal expansion and compression by lithostatic pressure to obtain rho(z).

NB. This function can be used for both Vp and Vs conversions, but the parameter values will differ. Ensure you use the correct parameters for the velocity type.

Parameters:

  • data (array-like) – An array of pressure and velocity data where data[0] is pressure (in GPa) and data[1] is velocity (in km/s).

  • parameters (class instance) –

    Instance of the VpConstants or VsConstants classes containing the attributes listed below:

    • v0 (float): Intercept velocity.

    • b (float): Velocity gradient with respect to density.

    • d0 (float): Partial derivative of velocity gradient (dvdr) with

      respect to pressure (dp).

    • dp (float): Velocity gradient with respect to pressure (dvdpr).

    • c (float): Amplitude of the velocity drop-off at low pressure.

    • k (float): Lengthscale of the velocity drop-off.

Returns:

Density at standard temperature and pressure (s.t.p.) in g/cm³.

Return type:

float or numpy.ndarray

density_functions.Vp2rho_brocher(Vp)[source]

Convert compressional wave velocity (Vp) to density (g/cm³) using the Nafe-Drake relationship.

Parameters:

Vp (float or numpy.ndarray) – Compressional wave velocity in km/s.

Returns:

Density in g/cm³.

Return type:

float or numpy.ndarray

density_functions.Vs2Vp_brocher(Vs)[source]

Convert shear wave velocity (Vs) to compressional wave velocity (Vp) using Brocher’s (2005) regression fit.

Parameters:

Vs (float or numpy.ndarray) – Shear wave velocity in km/s.

Returns:

Compressional wave velocity (Vp) in km/s.

Return type:

float or numpy.ndarray

density_functions.av_profile(profiles, base_depth, average_type='median')[source]

Calculate the average profile and bulk property profile given a family of profiles.

This function takes a list of profiles, where each profile is represented as a NumPy array with columns [depth, property]. It calculates the average profile and the bulk property profile as a function of depth.

Parameters:

  • profiles (list of np.ndarray) – List of profiles (e.g., Vp, Vs, P, T, rho, etc.). Each profile must contain at least two columns: depth and property.

  • base_depth (float) – Depth to which the bulk property should be averaged.

  • average_type (str, optional) – Method of averaging (‘median’ or ‘mean’). Default is ‘median’.

Returns:

  • av_z_f (scipy.interpolate.CubicSpline) – Average profile as a function of depth.

  • bulk_base_depth_f (scipy.interpolate.CubicSpline) – Bulk property profile as a function of base depth.

Examples

>>> profiles = [
...     np.array([
...         [0.0, 1.5],
...         [5.0, 2.0],
...         [10.0, 2.5]
...     ]),
...     np.array([
...         [0.0, 1.6],
...         [5.0, 2.1],
...         [10.0, 2.6]
...     ])
... ]
>>> av_z_f, bulk_base_depth_f = av_profile(
...     profiles,
...     base_depth=15.0,
...     average_type='mean'
... )
density_functions.check_arguments(T_dependence, constant_depth, constant_density, approach, parameters, geotherm=None)[source]

Check if required arguments are provided and raise errors if not.

This function checks the provided arguments for the conversion process and raises errors if any required information is missing based on the selected conversion approach. Prevents conflicting options.

Parameters:

  • T_dependence (bool) – Flag indicating whether to correct for temperature and pressure.

  • constant_depth (float) – Depth range for constant density from the surface.

  • constant_density (float) – Constant density value for the uppermost x kilometers.

  • approach (str) – Density conversion scheme, “brocher” or “stephenson”.

  • parameters (str) – Constants class instance if using the Stephenson scheme.

  • geotherm (Geotherm, optional) – Geotherm class instance object for temperature dependence. Defaults to None.

Raises:

ValueError – If required arguments are missing or incompatible with the selected approach.

density_functions.convert_V_profile(file, profile_type, write_data=False, path=None, approach='stephenson', location=None, parameters=None, constant_depth=None, constant_density=None, T_dependence=False, geotherm=None, print_working_file=False)[source]

Convert a single Vp or Vs velocity profile to a density profile.

Parameters:

  • file (str) – Path to the input velocity profile file.

  • profile_type (str) – Type of velocity profile, ‘Vp’ or ‘Vs’.

  • write_data (bool, optional) – If True, writes the converted data to files. Default is False.

  • path (str, optional) – Path to write the results. Default is None.

  • approach (str, optional) – Density conversion scheme to use, “brocher” or “stephenson”. Default is “stephenson”.

  • location (str, optional) – Regional location of the profile. Default is None.

  • parameters (np.ndarray, optional) – Instance of the Constants, VpConstants or VsConstants classes. Must be instance of the Constants class with a material_constants object if T_dependence is set to True.

  • constant_depth (float, optional) – Depth range for constant density from the surface. Default is None.

  • constant_density (float, optional) – Constant density value for the uppermost x kilometers. Default is None.

  • T_dependence (bool, optional) – If True, corrects velocity and density for temperature and pressure. Default is False.

  • geotherm (object, optional) – Instance of the Geotherm class for temperature dependence. Must be provided if T_dependence is set to True. Default is None.

  • print_working_file (bool, optional) – If True, prints the file that is being converted. Default is False.

Returns:

DataFrame containing the converted velocity and density data.

Return type:

pd.DataFrame

Examples

>>> convert_V_profile(
...     'velocity_profile.dat',
...     'Vp',
...     write_data=True,
...     path='output/',
...     approach='stephenson',
...     location='Region1',
...     parameters=constants_object,
...     constant_depth=10.0,
...     constant_density=2.7,
...     T_dependence=True,
...     geotherm=geotherm_object,
...     print_working_file=True
... )
density_functions.convert_to_same_depth_intervals(profile1, profile2)[source]

Convert profiles to have the same depth intervals. Takes depth as first column and any arbitrary y value as second. profiles 1 and 2 need not have the same y field.

This function bins velocity profiles so that they have the same depth ranges. It is suitable for profiles with distinct layers rather than those with gradients, although it is designed to produce a reasonable solution in both cases.

Parameters:

  • profile1 (np.ndarray) – The first profile as a NumPy array with columns [depth, velocity].

  • profile2 (np.ndarray) – The second profile as a NumPy array with columns [depth, velocity].

Returns:

A binned velocity profile aligned with the depth intervals of the lowest-resolution profile.

Return type:

np.ndarray

Examples

>>> binned_profile = convert_to_same_depth_intervals(profile1, profile2)
>>> print(binned_profile)
array([[ 0. , 1600. ],
       [ 2.5, 1600. ],
       [ 2.5 , 2050. ],
       [10. , 2050. ]])

Notes

  • This function finds the lowest-resolution profile and interpolates the

    other profile accordingly.

  • It returns a binned profile with the same depth intervals as the

    lowest-resolution profile.

  • It works best when lowest-resolution profile has discrete layers

    rather than gradients, although it will return reasonable values in both cases.

density_functions.ensure_dir(filename)[source]

Ensure that the directory for the given filename exists. If the directory does not exist, it will be created.

Parameters:

filename (str) – The path to the file.

Return type:

None

density_functions.p(rho, h, g=9.81)[source]

Calculate lithostatic pressure given average density of overburden and thickness of the overburden.

Parameters:

  • rho (float) – Average density of the overburden in Mg/m³.

  • h (float) – Thickness of the overburden in kilometers.

  • g (float, optional) – Acceleration due to gravity in m/s². Default is 9.81 m/s².

Returns:

Lithostatic pressure in Pascals (Pa).

Return type:

float

Notes

The lithostatic pressure (P) is calculated using the formula: P = rho * g * h * 1e6

density_functions.progress_bar(iterable, length=50)[source]

Create a text-based progress bar to track the progress of an iterable.

Parameters:

  • iterable (iterable) – The iterable object (e.g., list, range) to iterate through.

  • length (int, optional) – The length of the progress bar. Defaults to 50.

Yields:

item – The current item from the iterable.

Examples

>>> for item in progress_bar(range(100)):
...     # Simulate some processing
...     time.sleep(0.1)
...
[=====================>                   ] 50%
density_functions.read_file(filename, delimiter=None)[source]

Read the contents of a file and return them as a list of lines.

Parameters:

  • filename (str) – The path to the file.

  • delimiter (str, optional) – The delimiter used to split lines into elements (default is None, which means lines are not split).

Returns:

A list of lines from the file.

Return type:

list

density_functions.read_no_profile_data(data_file, region_name=None)[source]

Read crustal data from a database of crustal thickness estimates and vp/vs ratios.

Parameters:

  • data_file (str) – Path to the data file to be read.

  • region_name (str, optional) – Optional regional location of the profile.

Returns:

A Pandas DataFrame containing the following columns:

  • ’region’ (str): Regional location of the profile (if provided,

    otherwise None).

  • ’lon’ (float): Longitude of the station.

  • ’lat’ (float): Latitude of the station.

  • ’moho’ (float): Crustal thickness estimate.

  • ’vp_vs’ (float): Ratio of seismic velocities (vp/vs).

Return type:

pd.DataFrame

Examples

>>> df = read_no_profile_data(
...         'crustal_data.csv',
...         region_name='Example_Region')
>>> print(df.head())
    region     lon     lat   moho   vp_vs
0  Example Region  -75.56   40.21  36.58  1.73
1  Example Region  -72.89   42.17  38.12  1.68
2  Example Region  -80.12   35.69  33.45  1.80
3  Example Region  -76.55   39.29  40.27  1.64
4  Example Region  -74.21   41.54  35.78  1.75
density_functions.save_bulk_profiles(outpath, filename, data_function, start, stop, step)[source]

Save bulk profiles data to a file.

This function generates x-values using numpy.arange() based on the specified start, stop, and step values. It computes the data values using the provided data_function and saves the data to a file in the specified outpath.

Parameters:

  • outpath (str) – The output directory where the data file will be saved.

  • filename (str) – The name of the output data file.

  • data_function (function) – A function that computes the average property at depth z.

  • start (float) – The start value for the x-values.

  • stop (float) – The stop value for the x-values.

  • step (float) – The step size for generating x-values.

Return type:

None

Examples

>>> save_bulk_profiles(
...     '/output/directory/',
...     'bulk_data.txt',
...     function,
...     0,
...     50.5,
...     0.5
... )
density_functions.write_profile_to_file(data, filename)[source]

Write data to a file.

Parameters:

  • data (list) – A list of data to be written to the file.

  • filename (str) – The path to the output file.

Return type:

None