Defining connected fault systems
This notebook demonstrates how to identify which fault segments in a shapefile should be connected up into large fault systems.
We’ll begin by importing relevant modules
# Import modules
from fault_mesh.faults.leapfrog import LeapfrogMultiFault
import os
import numpy as np
import geopandas as gpd
Setting useful parameters
There are a few parameters that are worth setting here:
Coordinate system via EPSG code (2193 for New Zealand) — this is optional but useful if you wish to visualize contours in GIS software.
Trimming gradient (\(\alpha_{trim}\)) for clipping depth contours — see Section 3.3 of Howell et al. (in review)
Dip multiplier and strike multiplier for calculating \(\Theta_{change}\) — see Section 3.3 of Howell et al. (in review)
# Set coordinate system (optional) EPSG code
# If not desired, set to None
epsg = 2193
# Trimming gradient (alpha) and strike and dip multipliers for trimming depth contours of multi-segment faults
trimming_gradient = 1.
dip_multiplier = 1.
strike_multiplier = 0.5
Reading in faults
First, you need to read in your GIS representation of faults. In this example, we use a subset of faults from the New Zealand Community Fault Model (Seebeck et al., 2022).
# Read in fault data from shapefile
fault_data = LeapfrogMultiFault.from_shp("tutorial_gis/central_nz_minimal_data.shp", remove_colons=True, epsg=epsg, trimming_gradient=trimming_gradient,
dip_multiplier=dip_multiplier, strike_multiplier=strike_multiplier)
---------------------------------------------------------------------------
AssertionError Traceback (most recent call last)
Cell In[3], line 2
1 # Read in fault data from shapefile
----> 2 fault_data = LeapfrogMultiFault.from_shp("tutorial_gis/central_nz_minimal_data.shp", remove_colons=True, epsg=epsg, trimming_gradient=trimming_gradient,
3 dip_multiplier=dip_multiplier, strike_multiplier=strike_multiplier)
File ~/miniconda3/envs/leapfrog-fault-models/lib/python3.10/site-packages/fault_mesh/faults/leapfrog.py:290, in LeapfrogMultiFault.from_shp(cls, filename, remove_colons, sort_sr, epsg, dip_choice, trimming_gradient, dip_multiplier, strike_multiplier, check_optional_fields)
286 @classmethod
287 def from_shp(cls, filename: str, remove_colons: bool = False, sort_sr: bool = False, epsg: int = None,
288 dip_choice: str = "pref", trimming_gradient: float = 1.,
289 dip_multiplier: float = 1., strike_multiplier: float = 0.5, check_optional_fields: bool = False):
--> 290 assert os.path.exists(filename)
291 trimmed_fault_gdf = cls.gdf_from_nz_cfm_shp(filename=filename, exclude_region_polygons=None, depth_type=None,
292 exclude_region_min_sr=1.8, include_names=None,
293 exclude_aus=False,
294 exclude_zero=False, )
295 multi_fault = cls(trimmed_fault_gdf, sort_sr=sort_sr, remove_colons=remove_colons, epsg=epsg,
296 dip_choice=dip_choice, trimming_gradient=trimming_gradient,
297 dip_multiplier=dip_multiplier, strike_multiplier=strike_multiplier,
298 check_optional_fields=check_optional_fields)
AssertionError:
Finding connections between segments
Now your data are read in, you need to set the distance tolerance. This tolerance is the minimum horizontal distance between two fault traces that is allowed to count as a connection.
dist_tolerance = 200.
The next cell uses python module networkx to find segment traces that are within the specified distance tolerance of each other.
fault_data.find_connections(verbose=False)
Found 156 connections
Found 142 connections between segment ends
It is now necessary to write out these connections for manual editing and review. The file will be written out into the same directory as this Jupyter notebook. It will have a prefix supplied by you and the suffix “_suggested.csv”.
fault_data.suggest_fault_systems("central_gt1_5_connected")
This will create a CSV file that looks like this:
The name of the combined fault system is in the first column, and names of the faults that make up the connected system are in subsequent columns.
Making and incorporating manual edits
The automatically-generated fault system suggestions will (by design) include hyper-connected fault systems that need to be broken up. At this stage, the best way to break up these networks into smaller fault systems is to do it manually by editing the CSV file. An example below shows a new line added to the CSV representing the Hope Fault system – The Hope Fault is grouped with the Alpine and Kekerengu-Needles fault systems in the automatically-generated connections CSV. Make sure you save this new file with a different name to avoid overwriting it!
Once you have made the necessary edits, read your new CSV to overwite the automatically-generated connected fault systems:
fault_data.read_fault_systems("./define_connections_data/central_gt1_5_connected_edited.csv")
fault_data.generate_curated_faults()
Defining a cutting hierarchy
Now you have defined your fault system (which you will later use to create meshes), it is necessary to specify which faults terminate against other faults. For example, it seems highly likely that the western end of the Hope Fault is truncated at depth by the Alpine fault. This complex mesh cutting is best achieved using dedicated software such as leapfrog (see below), but for cutting to be done automatically, it is best to first specify a cutting hierarchy.
Suggesting a hierarchy based on slip rate
A first pass at a hierarchy can be generated based purely on fault slip rate. Assuming that the fault data you have already read in have slip rates associated with them, this first pass is easy to make:
fault_data.suggest_cutting_hierarchy("central_gt1_5_hierarchy")
This operation simply orders the faults (or fault systems) in your model in descending order of slip rate. For connected fault systems that have segments with different slip rates, the maximum slip rate of any segment in that fault system is used to place the fault system in the cutting hierarchy.
Editing the cutting hierarchy
You can then edit this hierarchy by switching the order of lines in the file. For any pair of faults/systems that intersect, the fault closer to the bottom of the file will terminate against the fault closer to the top of the list.
An example of a situation where editing is desirable is illustrated below. The maximum slip rate of the Jordan-Kekerengu-Needles Fault System (23 mm/yr) is faster than the corresponding maximum for the Hope Fault (15.8 mm/yr), but we wish to create a fault model where the Jordan Fault terminates against the Hope Fault. We effect this termination by moving Hope combined
above Jordan - Kek - Needles combined
in the CSV file. For similar reasons, we move the Hanmer
Fault above Hope Hanmer NW
.
We read in this new hierarchy as follows:
fault_data.read_cutting_hierarchy("./define_connections_data/central_gt1_5_hierarchy_edited.csv")
Create shapefiles for mesh creation
The final pre-meshing step is the creation of files that can be combined with meshing software to create triangular mesh representations of faults. Although it is possible to construct these triangular surfaces in multiple software packages (for example, MOVE 3D), the following discussion is geared towards use with Leapfrog Geo software.
Create directories to hold shapefiles
For organisational reasons, it helps to have the different shapefiles in different directories
for dir_name in ["depth_contours", "traces", "footprints"]:
if not os.path.exists(dir_name):
os.mkdir(dir_name)
Write out shapefiles
for fault in fault_data.curated_faults:
# Generate depth contours
fault.generate_depth_contours(np.arange(2000, 32000., 2000.), smoothing=False)
# Write contours to files
fault.contours.to_file(f"depth_contours/{fault.name}_contours.shp")
# Write traces
fault.nztm_trace_geoseries.to_file(f"traces/{fault.name}_trace.shp")
# Write fault footprints
for fault in reversed(fault_data.curated_faults):
fault.adjust_footprint()
fault.footprint_geoseries.to_file(f"footprints/{fault.name}_footprint.shp")