HMTK Geological Tools Demonstration

This notepad demonstrates the use of the HMTK geological tools for preparing fault source models for input into OpenQuake

Construction of the Geological Input File

An active fault model input file contains two sections:

1) A tectonic regionalisation - this can provide a container for a set of properties that may be assigned to multiple faults by virtue of a common tectonic region

2) A set of active faults

Tectonic Regionalisation Representation in the Fault Source File

#***************************************************************************** #FAULT FILE IN YAML (Yet Another Markup Language) FORMAT #***************************************************************************** #\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\ tectonic_regionalisation: - Name: Active Shallow Crust Code: 001 # Magnitude scaling relation (see http://docs.openquake.org/oq-hazardlib) #for currently available choices! Magnitude_Scaling_Relation: { Value: [WC1994], Weight: [1.0]} # Shear Modulus (in gigapascals, GPa) Shear_Modulus: { Value: [30.0], Weight: [1.0]} # Fault displacement to length ratio Displacement_Length_Ratio: { Value: [1.25E-5], Weight: [1.0]}

In the tectonic regionalisation information each of the three properties can be represented according to a set of weighted values. For example, in the case below faults in an arbitrarily named tectonic region (called here "GEM Region 1") will share the same set of magnitude scaling relations and shear moduli, unless over-written by the specific fault. Those faults assigned to "GEM Region 2" will have the magnitude scaling relation fixed as WC1994 and the shear modulus of 30 GPa

tectonic_regionalisation: - Name: GEM Region 1 Code: 001 # Magnitude scaling relation (see http://docs.openquake.org/oq-hazardlib) #for currently available choices! Magnitude_Scaling_Relation: { Value: [WC1994, PeerMSR], Weight: [0.7, 0.3]} # Shear Modulus (in gigapascals, GPa) Shear_Modulus: { Value: [25.0, 30.0, 35.0], Weight: [0.2, 0.6, 0.2]} # Fault displacement to length ratio Displacement_Length_Ratio: { Value: [1.25E-5], Weight: [1.0]} - Name: GEM Region 2 Code: 002 # Magnitude scaling relation (see http://docs.openquake.org/oq-hazardlib) #for currently available choices! Magnitude_Scaling_Relation: { Value: [WC1994], Weight: [1.0]} # Shear Modulus (in gigapascals, GPa) Shear_Modulus: { Value: [30.0], Weight: [1.0]} # Fault displacement to length ratio Displacement_Length_Ratio: { Value: [1.25E-5], Weight: [1.0]}

Active Fault Model

A set of active faults will be defined with a common ID and name.

Fault_Model_ID: 001 Fault_Model_Name: Template Simple Fault

An active fault set containing a single fault is shown below:

#\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\ Fault_Model_ID: 001 Fault_Model_Name: Template Simple Fault Fault_Model: - ID: 001 Tectonic_Region: Active Shallow Crust Fault_Name: A Simple Fault Fault_Geometry: { Fault_Typology: Simple, # For simple typology, defines the trace in terms of Long., Lat. Fault_Trace: [30.0, 30.0, 30.0, 31.5], # Upper Seismogenic Depth (km) Upper_Depth: 0.0, # Lower Seismogenic Depth (km) Lower_Depth: 20.0, Strike: , # Dip (degrees) Dip: 60.} Rake: -90.0 Slip_Type: Thrust Slip_Completeness_Factor: 1 #slip [value_1, value_2, ... value_n] # [weight_1, weight_2, ... weight_n] Slip: { Value: [18., 20.0, 23.], Weight: [0.3, 0.5, 0.2]} #Aseismic Slip Factor Aseismic: 0.0 MFD_Model: # Example of constructor for characteristic earthquake - Model_Name: Characteristic # Spacing (magnitude units) of the magnitude frequency distribution MFD_spacing: 0.1 # Weight of the model Model_Weight: 0.7 # Magnitude of the Characteristic Earthquake Maximum_Magnitude: # Uncertainty on Characteristic Magnitude (in magnitude units) Sigma: 0.12 # Lower bound truncation (in number of standard deviations) Lower_Bound: -3.0 # Upper bound truncation (in number of standard deviations) Upper_Bound: 3.0 #################################################### - Model_Name: AndersonLucoArbitrary # Example constructor of the Anderson & Luco (1983) - Arbitrary Exponential # Type - chooses between type 1 ('First'), type 2 ('Second') or type 3 ('Third') Model_Type: First MFD_spacing: 0.1 Model_Weight: 0.3 # Maximum Magnitude of the exponential distribution Maximum_Magnitude: Minimum_Magnitude: 4.5 # b-value of the exponential distribution as [expected, uncertainty] b_value: [0.8, 0.05] Megazone: Shear_Modulus: { Value: [30., 35.0], Weight: [0.8, 0.2]} Magnitude_Scaling_Relation: { Value: [WC1994], Weight: [1.0]} Scaling_Relation_Sigma: { Value: [-1.5, 0.0, 1.5], Weight: [0.15, 0.7, 0.15]} Aspect_Ratio: 1.5 Displacement_Length_Ratio: { Value: [1.25E-5, 1.5E-5], Weight:[0.5, 0.5]}

Fault Geometry Representations - Example 1: Simple Fault

Fault_Geometry: { Fault_Typology: Simple, # For simple typology, defines the trace in terms of Long., Lat. Fault_Trace: [30.0, 30.0, 30.0, 31.5], # Upper Seismogenic Depth (km) Upper_Depth: 0.0, # Lower Seismogenic Depth (km) Lower_Depth: 20.0, Strike: , # Dip (degrees) Dip: 60.}

Fault Geometry Representations - Example 2: Complex Fault

Fault_Geometry: { Fault_Typology: Complex, # For complex fault geometry define a list of edgees from top to # bottom, where each edges is [Long_1, Lat_1, Depth_1, Long_2, Lat_2, # Depth_2, ... Long_N, Lat_N, Depth_N Fault_Trace: [ # Top Edge [30.0, 30.0, 0.0, 30.0, 31.5, 0.0], # Intermediate Edge [30.2, 29.8, 30.0, 30.2, 31.7, 31.0], # Bottom Edge [30.3, 29.6, 58.0, 30.3, 31.9, 60.0] ]}

Rupture Properties

The rupture requires characterisation of the rake (using the Aki & Richards 2002 convention), the slip-type, the slip completeness factor (an integer constraining the quality of the slip information with 1 being the hights quality), the range of slip values and their corresponding weights, and the aseismic slip coefficient (the proportion of slip released aseismically, 1.0 - coupling coefficient)

Rake: -90.0 Slip_Type: Thrust Slip_Completeness_Factor: 1 #slip [value_1, value_2, ... value_n] # [weight_1, weight_2, ... weight_n] Slip: { Value: [18., 20.0, 23.], Weight: [0.3, 0.5, 0.2]} #Aseismic Slip Factor Aseismic: 0.0

The Magnitude Frequency Distributions



In [ ]:

    
%load_ext autoreload
%autoreload 2
import warnings; warnings.filterwarnings('ignore')



In [ ]:

    
#Import tools
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from openquake.hmtk.plotting.faults.geology_mfd_plot import plot_recurrence_models
from openquake.hazardlib.scalerel.wc1994 import WC1994  # In all the following examples the Wells & Coppersmith (1994) Scaling Relation is Used

The following examples refer to a fault with the following properties:

Length (Along-strike) = 100 km, Width (Down-Dip) = 20 km, Slip = 10.0 mm/yr, Rake = 0. (Strike Slip), Magnitude Scaling Relation = Wells & Coppersmith (1994), Shear Modulus = 30.0 GPa



In [ ]:

    
# Set up fault parameters
slip = 10.0  # Slip rate in mm/yr

# Area = along-strike length (km) * down-dip with (km)
area = 100.0 * 20.0

# Rake = 0.
rake = 0.

# Magnitude Scaling Relation
msr = WC1994()

Anderson & Luco (Arbitrary)

This describes a set of distributons where the maximum magnitude is assumed to rupture the whole fault surface



In [ ]:

    
#Magnitude Frequency Distribution Example

anderson_luco_config1 = {'Model_Name': 'AndersonLucoArbitrary',
                         'Model_Type': 'First',
                         'Model_Weight': 1.0,  # Weight is a required key - normally weights should sum to 1.0 - current example is simply illustrative! 
                         'MFD_spacing': 0.1,
                         'Maximum_Magnitude': None,
                         'Minimum_Magnitude': 4.5,
                         'b_value': [0.8, 0.05]}
anderson_luco_config2 = {'Model_Name': 'AndersonLucoArbitrary',
                         'Model_Type': 'Second',
                         'Model_Weight': 1.0,
                         'MFD_spacing': 0.1,
                         'Maximum_Magnitude': None,
                         'Minimum_Magnitude': 4.5,
                         'b_value': [0.8, 0.05]}
anderson_luco_config3 = {'Model_Name': 'AndersonLucoArbitrary',
                         'Model_Type': 'Third',
                         'Model_Weight': 1.0,   
                         'MFD_spacing': 0.1,
                         'Maximum_Magnitude': None,
                         'Minimum_Magnitude': 4.5,
                         'b_value': [0.8, 0.05]}
# Create a list of the configurations
anderson_luco_arb = [anderson_luco_config1, anderson_luco_config2, anderson_luco_config3]

# View the corresponding magnitude recurrence model
plot_recurrence_models(anderson_luco_arb, area, slip, msr, rake, msr_sigma=0.0)

Anderson & Luco (Area - MMax)

This describes a set of distributons where the maximum rupture extent is limited to only part of the fault surface



In [ ]:

    
anderson_luco_config1 = {'Model_Name': 'AndersonLucoAreaMmax',
                         'Model_Type': 'First',
                         'Model_Weight': 1.0,  # Weight is a required key - normally weights should sum to 1.0 - current example is simply illustrative! 
                         'MFD_spacing': 0.1,
                         'Maximum_Magnitude': None,
                         'Minimum_Magnitude': 4.5,
                         'b_value': [0.8, 0.05]}
anderson_luco_config2 = {'Model_Name': 'AndersonLucoAreaMmax',
                         'Model_Type': 'Second',
                         'Model_Weight': 1.0,
                         'MFD_spacing': 0.1,
                         'Maximum_Magnitude': None,
                         'Minimum_Magnitude': 4.5,
                         'b_value': [0.8, 0.05]}
anderson_luco_config3 = {'Model_Name': 'AndersonLucoAreaMmax',
                         'Model_Type': 'Third',
                         'Model_Weight': 1.0,   
                         'MFD_spacing': 0.1,
                         'Maximum_Magnitude': None,
                         'Minimum_Magnitude': 4.5,
                         'b_value': [0.8, 0.05]}

# For these models a displacement to length ratio is needed
disp_length_ratio = 1.25E-5

# Create a list of the configurations
anderson_luco_area_mmax = [anderson_luco_config1, anderson_luco_config2, anderson_luco_config3]

# View the corresponding magnitude recurrence model
plot_recurrence_models(anderson_luco_area_mmax, area, slip, msr, rake, msr_sigma=0.0)

Characteristic Earthquake

The following example illustrates a "Characteristic" Model, represented by a Truncated Gaussian Distribution



In [ ]:

    
characteristic = [{'Model_Name': 'Characteristic',
                   'MFD_spacing': 0.05,
                   'Model_Weight': 1.0,
                   'Maximum_Magnitude': None,
                   'Sigma': 0.15,  # Standard Deviation of Distribution (in Magnitude Units) - omit for fixed value
                   'Lower_Bound': -3.0,   # Bounds of the distribution correspond to the number of sigma for truncation
                   'Upper_Bound': 3.0}]

# View the corresponding magnitude recurrence model
plot_recurrence_models(characteristic, area, slip, msr, rake, msr_sigma=0.0)

Youngs & Coppersmith (1985) Models

The following describes the recurrence from two distributions presented by Youngs & Coppersmith (1985): 1) Exponential Distribution, 2) Hybrid Exponential-Characteristic Distribution



In [ ]:

    
exponential = {'Model_Name': 'YoungsCoppersmithExponential',
               'MFD_spacing': 0.1,
               'Maximum_Magnitude': None,
               'Maximum_Magnitude_Uncertainty': None,
               'Minimum_Magnitude': 5.0,
               'Model_Weight': 1.0,
               'b_value': [0.8, 0.1]}

hybrid = {'Model_Name': 'YoungsCoppersmithCharacteristic',
          'MFD_spacing': 0.1,
          'Maximum_Magnitude': None,
          'Maximum_Magnitude_Uncertainty': None,
          'Minimum_Magnitude': 5.0,
          'Model_Weight': 1.0,
          'b_value': [0.8, 0.1],
          'delta_m': None}

youngs_coppersmith = [exponential, hybrid]

# View the corresponding magnitude recurrence model
plot_recurrence_models(youngs_coppersmith, area, slip, msr, rake, msr_sigma=0.0)

Epistemic Uncertainty Examples

At present the toolkit can support, directly, epistemic uncertainty in six of the primary attributes within a single workflow: 1) Slip Rate 2) Choice of MAgnitude Scaling Relation 3) Shear Modulus 4) Displacement to Length Ratio (Applies to Anderson & Luco Area-MMax relation only) 5) Standard Deviation of scaling relation 6) Magnitude Frequency Distribution Configuration

This example considers the fault defined at the top of the page. This fault defines two values of slip rate and two different magnitude frequency distributions



In [ ]:

    
def show_file_contents(filename):
    """
    Shows the file contents
    """
    fid = open(filename, 'r')
    for row in fid.readlines():
        print row
    fid.close()

input_file = 'input_data/simple_fault_example_4branch.yml'
show_file_contents(input_file)

Example 1 - Full Enumeration

In this example each individual MFD for each branch is determined. In the resulting file the fault is duplicated n_branches number of times, with the corresponding MFD multiplied by the end-branch weight



In [ ]:

    
# Import the Parser
from openquake.hmtk.parsers.faults.fault_yaml_parser import FaultYmltoSource

# Fault mesh discretization step
mesh_spacing = 1.0 # (km)

# Read in the fault model
reader = FaultYmltoSource(input_file)
fault_model, tectonic_region = reader.read_file(mesh_spacing)

# Construct the fault source model (this is really running the MFD calculation code)
fault_model.build_fault_model()

# Write to an output NRML file
output_file_1 = 'output_data/fault_example_enumerated.xml'
fault_model.source_model.serialise_to_nrml(output_file_1)

show_file_contents(output_file_1)

Example 2: Collapsed Branches

In the following example we implement the same model, this time collapsing the branched. This means that the MFD is discretised and the incremental rate in each magnitude bin is the weighted sum of the rates in that bin from all the end branches of the logic tree.

When collapsing the branches, however, it is necessary to define a single Magnitude Scaling Relation that will need to be assigned to the fault for use in OpenQuake.



In [ ]:

    
# Read in the fault model
reader = FaultYmltoSource(input_file)
fault_model, tectonic_region = reader.read_file(mesh_spacing)

# Scaling relation for export
output_msr = WC1994()

# Construct the fault source model - collapsing the branches
fault_model.build_fault_model(collapse=True, rendered_msr=output_msr)


# Write to an output NRML file
output_file_2 = 'output_data/fault_example_collapsed.xml'
fault_model.source_model.serialise_to_nrml(output_file_2)

show_file_contents(output_file_2)

Example 3: Gulf of Corinth Faults

In the last example we will take some real fault data from the Gulf of Corinth



In [ ]:

    
# For plotting
from openquake.hmtk.plotting.mapping import HMTKBaseMap
corinth_fault_file = "input_data/GulfOfCorinth_Faults.yml"

reader = FaultYmltoSource(corinth_fault_file)
corinth_faults, tectonic_region = reader.read_file(mesh_spacing)

print "           ID                         Name   Area (km^2)"
for i, fault in enumerate(corinth_faults.faults):
    print "%2g    %s    %25s    %10.3f" % (i, fault.id, fault.name, fault.area)

For convenience we will collapse the branches again



In [ ]:

    
# Scaling relation for export
output_msr = WC1994()

# Construct the fault source model - collapsing the branches
corinth_faults.build_fault_model(collapse=True, rendered_msr=output_msr)

Let's take a look at the fault sources



In [ ]:

    
map_config = {'min_lon': 20.9, 'max_lon': 24.0,
              'min_lat': 37.5, 'max_lat': 39.0, 'resolution':'h'}

# Map the Source
src_basemap = HMTKBaseMap(map_config, "Gulf of Corinth Faults")
src_basemap.add_source_model(corinth_faults.source_model)

Find and Show Events within 30 km of the Xylocastro Fault



In [ ]:

    
from openquake.hmtk.parsers.catalogue.csv_catalogue_parser import CsvCatalogueParser
from openquake.hmtk.seismicity.selector import CatalogueSelector

# Load in the catalogue
catalogue_file = "input_data/Aegean_ExtendedCat1.csv"
parser = CsvCatalogueParser(catalogue_file)
catalogue = parser.read_file()
print catalogue.get_number_events()

#Create a selector object
selector = CatalogueSelector(catalogue)
xylocastro_fault = corinth_faults.faults[15]

# Select earthquakes within 30 km (rupture distance) of the fault
xylocastro_fault.select_catalogue(selector, 30.0, distance_metric="rupture")
print "%g events with 30 km of the Xylocastro Fault" % xylocastro_fault.catalogue.get_number_events()

Show the events



In [ ]:

    
map_config = {'min_lon': 20.9, 'max_lon': 24.0,
              'min_lat': 37.5, 'max_lat': 39.0, 'resolution':'h'}

# Map the Source
src_basemap = HMTKBaseMap(map_config, "Gulf of Corinth Faults")
src_basemap.add_source_model(corinth_faults.source_model, overlay=True)
src_basemap.add_catalogue(xylocastro_fault.catalogue)

Export output to XML File



In [ ]:

    
corinth_output_file = "output_data/New_Gulf_Of_Corinth_Faults.xml"
corinth_faults.source_model.serialise_to_nrml(corinth_output_file)
show_file_contents(corinth_output_file)