SMV2rho: Tutorial 1

This notebook will show you how to load a velocity profile into SMV2rho.

Before starting this tutorial, please refer to the README for installation instructions.

The installation of SMV2rho should have checked that all dependencies are installed. If errors arise relating to missing modules such as numpy, then please refere to README.md for further information.

Importing the relevant modules

For this tutorial we need to import the plotting and density_functions modules from the SMV2rho package.

from SMV2rho import plotting as smplt
from SMV2rho import density_functions as smd

Intel MKL WARNING: Support of Intel(R) Streaming SIMD Extensions 4.2 (Intel(R) SSE4.2) enabled only processors has been deprecated. Intel oneAPI Math Kernel Library 2025.0 will require Intel(R) Advanced Vector Extensions (Intel(R) AVX) instructions.
Intel MKL WARNING: Support of Intel(R) Streaming SIMD Extensions 4.2 (Intel(R) SSE4.2) enabled only processors has been deprecated. Intel oneAPI Math Kernel Library 2025.0 will require Intel(R) Advanced Vector Extensions (Intel(R) AVX) instructions.

File layout

First let’s take a look at the necessary file structure required by this family of programs.

The required file format is strict because the scripts will read the file in a particular order. If the file is not presented in this format then errors will occur that could lead to erroneous results. Please note that this requirement means that continuous 3-D velocity models must first be converted into 1-D profiles before using the software. Run the following cell to see the required file layout. Note only ascii text files with the .dat suffix are acceptable.

Depth must be in the second column and negative. At this stage it does not matter, but when coverting profiles later on, depth must be in km/s and depth in negative km.

First we will use a local file stored in the TEST_DATA directory that is provided with this distribution.

# plot an example data file
print("\nExample file:\n \
       ID \n \
       lon lat \n \
       crustal_thickness \n \
       V  -depth \n \
       .   . \n \
       .   . \n \
       .   . \n")

Example file:
        ID
        lon lat
        crustal_thickness
        V  -depth
        .   .
        .   .
        .   .

# path to test velocity file
#  - this file comes with the distribution so there is no need to change this path
vp_file = "../TEST_DATA/EUROPE/Vp/RECEIVER_FUNCTION/DATA/M19_AQU_Vp.dat"

# print the file format
print("Real file...\n")
smplt.print_text_file(vp_file)

Real file...

M19_AQU_Vp
13.48   42.34
37.2
4.84865  0
4.84865  -2.5
7.23144  -2.5
7.23144  -16.2
6.42768  -16.2
6.42768  -37.2
6.55704  -37.2
6.55704  -65.9

Data Structure

While arbitrary file structures are supported, it is highly advisable that the file structure is organised as follows. A typical velocity profile file path will look like this:

Vp/METHOD/DATA/my_profile_name.dat

Using this file structure will allow us to automatically detect metadata like profile location etc. from the file path string. We can draw a file tree using the draw_fie_tree function in the plotting module to look at the recommended file structure

# draw a file tree
smplt.draw_file_tree('../TEST_DATA', include_files=False,
               suppress_pycache=True, suppress_hidden=True)

TEST_DATA
|- TEST_DATA
|  |- EUROPE
|  |  `- Vp
|  |     `- RECEIVER_FUNCTION
|  |        `- DATA
|  `- HUDSON_BAY
|     `- Vs
|        `- RECEIVER_FUNCTION
|           `- DATA

Loading files

We will now load this velocity profile into the program using the Convert class. Note that if the profile_type variable is not set, you will receive a ValueError asking for profile_type to be set. This error ensures that we keep a log of whether the profile is a \(V_P\) or \(V_S\) profile throughout the workflow.

Here we will also set the region_name and seismic_method_name variables. Note that if the recommended file structure is used then these arguments will be read from the file path string automatically and do not need to be provided. Since the profiles in TEST_DATA follow this recommended layout, these arguments are not sctrictly necessary here. You can play with seting these values to None and seeing that the relevant information is read from the file path string. If these parameters are not known, then a placeholder such as "unknown" or "na" can be used.

# load a profile into the Convert class
profile = smd.Convert(vp_file, profile_type = "Vp",
                      region_name = "EUROPE",
                      seismic_method_name = "RF")

Documentation for all objects (e.g. modules, classes, functions) is availiable by typing the object of interest followed by a question mark. Let’s take a look at the docstring for the Convert calss…

# manual for Convert class
smd.Convert?

Init signature:
smd.Convert(
    profile,
    profile_type=None,
    region_name=None,
    seismic_method_name=None,
    geotherm=None,
)
Docstring:
Convert seismic velocity profiles to various parameters using different
approaches.

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

Args:
    profile (str): The file path to the seismic profile data.
    profile_type (str): The type of the seismic profile, either "Vp" or
        "Vs."
    region_name (str, optional): The geographic location of the profile,
        e.g., "MADAGASCAR." (default is None)
    seismic_method_name (str): The method used to acquire the vlocity
        profile.  e.g. "RECEIVER_FUNCTION".  If set to None, the read_data
        metod 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.

Attributes:
    data (dict): A dictionary containing parsed seismic profile data.
    moho (float): The Moho depth parsed from the profile data.

Methods:
    read_data: Read in data file and parse it into a data dictionary.
    convert_profile_brocher: Convert Vs profile to Vp profile using
        Brocher's (2005) approach.
    Vp_to_density_brocher: Convert Vp profile to density using Brocher
        (2005) method.
    V_to_density_stephenson: Convert Vp profile to density using the
        Stephenson method described in the study.
    write_data: Write the converted data to appropriate file locations
        based on the specified conversion approach and temperature
        dependence settings.
File:           ~/Work/SMV2rho/src/SMV2rho/density_functions.py
Type:           type
Subclasses:

Reading the file

We will now read the profile and generate a "Vp" attribute within the profile object. The read_data method will unpack the information contained in the file and calculate the average velocity for the profile.

The data, profile.data, will be stored as a disctionary.

# read in the file
profile.read_data()

# print the data dictionary
profile.data

{'station': 'M19_AQU_Vp',
 'Vp_file': '../TEST_DATA/EUROPE/Vp/RECEIVER_FUNCTION/DATA/M19_AQU_Vp.dat',
 'region': 'EUROPE',
 'moho': 37.2,
 'location': array([13.48, 42.34]),
 'av_Vp': 6.617358544354839,
 'Vp': array([[  0.     ,   4.84865],
        [ -2.5    ,   4.84865],
        [ -2.51   ,   7.23144],
        [-16.2    ,   7.23144],
        [-16.21   ,   6.42768],
        [-37.2    ,   6.42768]]),
 'type': 'Vp',
 'method': 'RF',
 'geotherm': None}

We can inspect the different data types that are stored using…

# show velocity profile dictionary keys
profile.data.keys()

dict_keys(['station', 'Vp_file', 'region', 'moho', 'location', 'av_Vp', 'Vp', 'type', 'method'])

We can access different variables, such as the average, bulk velocity by querying this dictionary.

# print average velocity and moho depth
print(profile.data['type'])
print(profile.data['method'])
print(profile.data['av_Vp'])
print(profile.data['moho'])

Vp
RF
6.617358544354839
37.2

Summary

We have now imported modules, loaded a velocity profile and interrogated some of their attributes. Let’s plot the velocity profile to check what it looks like.

# plot velocity profile
data1 = [{'x': profile.data["Vp"][:,1], 'y': profile.data["Vp"][:,0]}]

# Call the plot_panels function
smplt.plot_panels([data1], plot_type='line',
            cmap=None, titles=None,
            xlabels=[r'${V_P}$ (km/s)'], ylabels=['Depth (km)'],
            z_values=None, figure_scale=0.7,
            save_path=None)