Two Layers¶

This example is devoted to demonstration of interface geometry optimization (interface inversion) applied to a simple synthetic model consisting of two layers.

Description¶

Original model¶

The original (or true) model contains two bodies (layers) of identical sizes: an upper body with zero density (an air layer) and a lower body (a subsurface layer) with density of 200 kg/m $^3$ , separated by a horizontal plane at a depth of 5 km (see Figure below).

Both bodies are 20 $\times$ 20 km horizontally and 5 km thick.

Input data¶

The original model was used to create a gravity dataset calculated for a set of 6561 stations regularly placed on the plane at a zero depth level (see Figure below).

Bouguer anomaly modelled above a simple two-layer model measured with a regularly spaced set of stations placed at a zero depth level.

Stations cover the area 40 $\times$ 40 km, exceeding the model edges by 10 km in each direction.
For simplicity, there is no extension of the bodies beyond the model boundaries, and the edge effect is clearly visible in the dataset.

The described dataset was used as an input measured gravity data for inversion.

Starting model¶

In order to have a starting model which would differ from the original model, additional vertices have been added on the boundary between the two bodies and then manually moved up or down (see Figure below).

From a 3D view of the starting model one can see that the boundary is modified at each of the 5 working sections (see Figure below).

The 3D view of the starting model showing modified boundary at 5 working sections. The measured gravity is shown on the top.

Download¶

The input data required for this example as well as the two models are available for download here.

The share contains 3 files:

Two_Layers_Measured_Gravity.csv: the input measured gravity data in CSV format
Two_Layers_Original.zip: the original IGMAS+ model of two layers separated with a horizontal plane used to calculate the input gravity data
Two_Layers_Inversion.zip: the IGMAS+ model used for inversion consisting of two timeline steps:
- starting model: a model where the boundary between the two layers was distorted
- final model - a model after application of the interface inversion to the starting model

The input gravity data file Two_Layers_Measured_Gravity.csv is in CSV format and has 4 data columns: "x" "y" "z" "measured z component" and 6561 data rows corresponding to stations.
The values in columns are delimited with space.

The two models are zip archives with IGMAS+ projects. Simply unpack and load projects in IGMAS+.

Modelling¶

Goal

The goal of this modelling example is to apply interface geometry optimization to the starting model in order to demonstrate how it can be iteratively inverted to a final model close to the true one based on the input gravity data, assuming that densities are known.

The interface geometry optimization functionality is provided by the Interface Inversion plugin. With this plugin it is possible to iteratively optimize the geometry of an interface using the misfit between the measured and calculated (at each iteration) field.

Open the starting model¶

Start IGMAS+

Select File → Open Project
Select the folder Two_Layers_Inversion
In the Timeline window select the earliest timeline - this is the starting model (see Figure on the right)

The starting model will open in the 3D view
However, there are no stations - for that you need to load the input gravity data.

Load input gravity data¶

Select File → Import → Stations
Navigate to the folder with Two_Layers_Measured_Gravity.csv file
Change the file type to [csv|xyz] - Comma Separated Values to see the file

Now you can see the files with .csv and .xyz extensions
Click on the Two_Layers_Measured_Gravity.csv file
Information on the units that will be used for input is visible on the right
Click Open

The preview window will show the data from the file distributed in columns
Click Finish

The loaded data will be displayed on top of the model together with stations
The data will also appear as measured data ("meas Gz") in the Object Tree under "Fields" → "Gravity: z-component"
Now it is necessary to calculate the gravity field for the loaded.

Calculate gravity field¶

First, triangulate the model using Edit → Model Triangulation or by clicking
Then calculate the anomaly of the model using Tools → Calculate Anomalies or by clicking
- Select "Calc Gz"
- Click Finish

Open 2D View using Add View → 2D View or by clicking
The curves above the vertical cross section shows the misfit between the measured and calculated gravity anomalies for the starting model
In the Object Tree under "Fields" select "Gravity: z-component"
Click on the Property Editor Tab and uncheck "Auto Shift" if you want to run inversion without the autoshift.
Now, before starting the inversion, it is necessary to add the inversion category.

Add inversion category¶

In the Object Tree, click right mouse button on "Interfaces" and select Add Category
In the dropdown list choose "Inversion" and click OK
A new "Inversion" entry will appear in "Interfaces"

Add an interface to the "Inversion" category
- In the [Object Tree] under "Interfaces" find the interface between the upper and lower bodies called "Unten <> Oben"
- Pick this interface with left mouse button and drag it to the "Inversion"
The interface "Unten <> Oben" will be copied under the "Inversion" (see Figure on the right)
Now you are ready to start the interface inversion.

Tip

The model at this stage can be opened from the middle timeline of the Two_Layers_Inversion model.

Start interface inversion¶

Open the interface inversion wizard using the start icon

Interface inversion settings window will open (see Figure on the right)
Select the "Optimizer" to be "Optimizer no-check"
Adjust the "Standard Deviation" to be in the range from 0.15 to 0.05 - this is the initial standard deviation for variation of depth coordinates of the interface vertices. The more the value, the larger is the initial variation. We recommend to keep 0.1 here.
Adjust the "Stop-Quality SD" to be in the range from 0.1 to 0.002. The less the value is, the longer the inversion will last and the better the fit will be in the end. A recommended value for the "Stop-Quality SD" to reach an optimal accuracy in a reasonable time is 0.05.
Make sure the "Use Triangle Effect" is checked: it involves calculation of gravity effect for the triangulated bodies
Make sure the "Use Voxel Effect" is unchecked: there are no voxel cubes and we don't need to involve it
Once ready, click Next

Area settings window will open:
The default area is the area covered by stations
Click Finish

Interface inversion process will start
A Population Quality window will open automatically and will dynamically show the statistics on the optimization process
Usually it takes about 20 to 30 minutes to reach the quality of 0.002 for this model
To stop the inversion process before it reaches the stop quality threshold, use the stop icon
Once inversion is done, the stop icon will change back to the start icon

The animation of the inversion process is shown on Figure below (video version is here):

Final model¶

The final result, as explained earlier depends on the stop quality and the initial standard deviation, and, besides, is not possible to reproduce because of the random nature of the evolution strategy.

The result obtained for the stop quality of 0.002 and the initial standard deviation of 0.1 shows that the original flat interface is almost perfectly reconstructed with minimal residual gravity field:

Result of interface inversion: vertical cross-section through the final model and the misfit between the measured and calculated gravity fields.

Result of interface inversion: map view on the measured and calculated gravity fields and the residual.

Result of interface inversion: 3D view of the final model.

Tip

The final model can be opened from the latest timeline of the Two_Layers_Inversion model.