{ "cells": [ { "cell_type": "markdown", "id": "024a62fd-9975-4107-9f3c-1b97601168dd", "metadata": {}, "source": [ "# Make Structured 3D Mesh\n", "\n", "A common input into 3D inversions is a structured mesh with equal cells within the station area and some padding cells. There are some tools in MTpy that can help design a mesh based on station locations, add topography or bathymetry, write to file specific for the inversion code. " ] }, { "cell_type": "code", "execution_count": 1, "id": "95c96fc1-fb21-455f-89e9-0cf9f6ec3e28", "metadata": {}, "outputs": [], "source": [ "%matplotlib widget" ] }, { "cell_type": "code", "execution_count": 2, "id": "3bc3b093-3499-4953-8ba9-a3867eb461eb", "metadata": {}, "outputs": [], "source": [ "from pathlib import Path\n", "from mtpy import MTData\n", "from mtpy.modeling import StructuredGrid3D" ] }, { "cell_type": "markdown", "id": "920c3158-dd9f-4735-bfed-c05d44bc134e", "metadata": {}, "source": [ "## 1. Read Data File\n", "\n", "In the previous notebook we created a data file, lets read that in and use the station locations to create a structured mesh for finite element codes like ModEM, WS3DINV, JIF3D." ] }, { "cell_type": "code", "execution_count": 3, "id": "7e03ab2f-c781-4d2b-8395-ba01659df049", "metadata": {}, "outputs": [], "source": [ "mtd = MTData()\n", "mtd.from_modem(r\"example_modem_data_file.dat\")" ] }, { "cell_type": "markdown", "id": "b7b5a323-6ced-495b-a1d9-deb2b2b859b9", "metadata": {}, "source": [ "## 2. Create a Model\n", "\n", "Now that we have station locations we can create a 3D mesh." ] }, { "cell_type": "code", "execution_count": 4, "id": "a4d2571d-949c-46ba-9c37-fd36002ac186", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "MT Location: \n", "--------------------\n", " Latitude (deg): 38.873786\n", " Longitude (deg): -118.196245\n", " Elevation (m): 0.0000\n", " Datum crs: EPSG:4326\n", "\n", " Easting (m): 396229.653\n", " Northing (m): 4303450.561\n", " UTM crs: EPSG:32611\n", "\n", " Model Easting (m): 396229.653\n", " Model Northing (m): 4303450.561\n", " Model Elevation (m): 0.000\n", " Profile Offset (m): 0.000\n" ] } ], "source": [ "mesh = StructuredGrid3D()\n", "mesh.station_locations = mtd.station_locations\n", "mesh.center_point = mtd.center_point\n", "print(mesh.center_point)" ] }, { "cell_type": "code", "execution_count": null, "id": "9f345aa2", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "482b6383-0790-42ac-b4a4-aaf41fc67b81", "metadata": {}, "source": [ "### Horizontal Set Cell Sizes\n", "\n", "We can set the lateral cell sizes in the east and north direction. There is a trade-off between cell size, computation speed, and resolution. Ideally you would like a few cells between each station to allow the inversion to find a better solution, but you don't want too many that it will take ages to calculate. \n", "\n", "The example data was collected at about 8 km station spacing, so lets try 2km cell sizes within the station area." ] }, { "cell_type": "code", "execution_count": 5, "id": "1a9af124-7730-4193-838e-39d256185fe4", "metadata": {}, "outputs": [], "source": [ "mesh.cell_size_east = 2000\n", "mesh.cell_size_north = 2000" ] }, { "cell_type": "markdown", "id": "f7957e38-ebbe-476f-af52-d0219bf6893f", "metadata": {}, "source": [ "#### Horizontal Padding cells\n", "\n", "There are a few different padding cells. \n", "\n", "- `pad_num` is the number of padding cells of the same size as the station cells that extent outside the station area (default is 3)\n", "- `pad_east` and `pad_north` describe the number of padding cells between the edge of the station area and the desired extent, these are calculated on an exponential distance. " ] }, { "cell_type": "code", "execution_count": 6, "id": "9267ffc3-7949-4252-b240-8efeb2b7e700", "metadata": {}, "outputs": [], "source": [ "mesh.pad_east = 10\n", "mesh.pad_north = 10" ] }, { "cell_type": "code", "execution_count": 7, "id": "4f60633c-3563-4b36-bb79-e6ae9b1fa573", "metadata": {}, "outputs": [], "source": [ "mesh.ew_ext = 250000\n", "mesh.ns_ext = 250000" ] }, { "cell_type": "markdown", "id": "942dd378-ee68-4d85-9ce1-a7eeccaff84b", "metadata": {}, "source": [ "#### Vertical Cells\n", "\n", "Because MT is a diffusive method we want to have a vertical grid that increases in size as a function of depth. There are a few ways to do this. Here we provide the first layer thickness `z1_layer`, the number of layers `n_layers` and the target depth `z_target_depth`. Then the cells will increase geometrically from the `z1_layer` to `z_target_depth`. These data are meant to image regional structures using broadband data, so the first layer should be thin and target depth around 50 km." ] }, { "cell_type": "code", "execution_count": 8, "id": "eaf99162-4e90-4515-9f35-a2849d9723d2", "metadata": {}, "outputs": [], "source": [ "mesh.z1_layer = 30\n", "mesh.n_layers = 50\n", "mesh.z_target_depth = 50000" ] }, { "cell_type": "markdown", "id": "b511dc81-b818-4a6f-adf1-4867ff5468a7", "metadata": {}, "source": [ "#### Vertical Padding cells\n", "\n", "We also need some padding cells in the vertical direction. We can set the geometric factor if needed." ] }, { "cell_type": "code", "execution_count": 9, "id": "7c83995d-c6e2-45e9-a35d-cff1f410a128", "metadata": {}, "outputs": [], "source": [ "mesh.pad_z = 5\n", "mesh.pad_stretch_v = 1.5" ] }, { "cell_type": "markdown", "id": "04a3256e-8fe0-484f-ba76-776a0643ea98", "metadata": {}, "source": [ "### Make Mesh \n", "\n", "Now we can create the mesh and see what it looks like." ] }, { "cell_type": "code", "execution_count": 10, "id": "69ca3a04-65b4-47ff-b215-ff05c31b8131", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Structured3DMesh Model Object:\n", "--------------------\n", "\tNumber of stations = 59\n", "\tMesh Parameter: \n", "\t\tcell_size_east: 2000\n", "\t\tcell_size_north: 2000\n", "\t\tpad_east: 10\n", "\t\tpad_north: 10\n", "\t\tpad_num: 3\n", "\t\tz1_layer: 30\n", "\t\tz_target_depth: 50000\n", "\t\tn_layers: 50\n", "\t\tn_air_layers: 0\n", "\t\tres_initial_value: 100.0\n", "\tDimensions: \n", "\t\te-w: 61\n", "\t\tn-s: 62\n", "\t\tz: 51 (without 7 air layers)\n", "\tExtensions: \n", "\t\te-w: 252000.0 (m)\n", "\t\tn-s: 250000.0 (m)\n", "\t\t0-z: 161492.0 (m)\n", "--------------------\n" ] } ], "source": [ "mesh.make_mesh()" ] }, { "cell_type": "markdown", "id": "398d486c-eff7-4990-afaa-3eb4834d62a0", "metadata": {}, "source": [ "
| \n", " | survey | \n", "station | \n", "latitude | \n", "longitude | \n", "elevation | \n", "datum_epsg | \n", "east | \n", "north | \n", "utm_epsg | \n", "model_east | \n", "model_north | \n", "model_elevation | \n", "profile_offset | \n", "
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | \n", "data | \n", "gv100 | \n", "38.611 | \n", "-118.535 | \n", "-1437.40 | \n", "4326 | \n", "0.0 | \n", "0.0 | \n", "None | \n", "-29000.0 | \n", "-28000.0 | \n", "-1545.999 | \n", "0.0 | \n", "
| 1 | \n", "data | \n", "gv101 | \n", "38.595 | \n", "-118.351 | \n", "-1540.55 | \n", "4326 | \n", "0.0 | \n", "0.0 | \n", "None | \n", "-13020.0 | \n", "-30020.0 | \n", "-1725.999 | \n", "0.0 | \n", "
| 2 | \n", "data | \n", "gv102 | \n", "38.594 | \n", "-118.277 | \n", "-1554.80 | \n", "4326 | \n", "0.0 | \n", "0.0 | \n", "None | \n", "-7000.0 | \n", "-30020.0 | \n", "-1635.999 | \n", "0.0 | \n", "
| 3 | \n", "data | \n", "gv103 | \n", "38.585 | \n", "-118.202 | \n", "-1543.90 | \n", "4326 | \n", "0.0 | \n", "0.0 | \n", "None | \n", "-1000.0 | \n", "-32020.0 | \n", "-1575.999 | \n", "0.0 | \n", "
| 4 | \n", "data | \n", "gv104 | \n", "38.596 | \n", "-118.137 | \n", "-1801.80 | \n", "4326 | \n", "0.0 | \n", "0.0 | \n", "None | \n", "5000.0 | \n", "-30020.0 | \n", "-1935.999 | \n", "0.0 | \n", "