Tutorial
========

Implicit modelling
^^^^^^^^^^^^^^^^^^

In this tutorial, we use `Geoscience ANALYST <https://www.mirageoscience.com/mining-industry-software/geoscience-analyst/pro-geophysics/>`_
and ``gempy-drivers`` to create a physical property model within an Octree mesh.  This
model can then be used in geophysical forward or inverse modelling. The tutorial will
use data from the Flin Flon VMS deposit and is available to `download here <https://github.com/MiraGeoscience/gempy-drivers/raw/refs/heads/release/0.1.0/gempy_drivers-assets/demo.geoh5?download=>`_.

.. contents:: Table of Contents
    :depth: 3

Getting started
---------------

We will be starting with a pre-existing project in `ANALYST <https://www.mirageoscience.com/mining-industry-software/geoscience-analyst/pro-geophysics/>`_.
If you are starting from scratch, the first step will be to import some data from
which to make the observations that will constrain the implicit modelling.  This may
include drillhole data, surface mapping or point data.  Whatever the data, it must
contain at least some orientation measurements and some contact points to satisfy the
requirements of Gempy (see :ref:`gempy rules <gempy rules>`).

Before starting to create any of the required data input for gempy-drivers package, it
is a good idea to identify what are the main geological units and events to model and
what are their age relationships.  For this tutorial, we will focus our attention to
the Blue Lagoon, Millrock, Hidden, and Mississippi formations, as well as a thrust fault
that strikes Northeast across the area, and offsets both the Millrock and Hidden
formations.

.. _history_table:

.. table:: Geological units and events

    +-----------+-------------+
    | Age       | Unit/Event  |
    | Index     |             |
    +===========+=============+
    | 1         | Fault       |
    +-----------+-------------+
    | 2         | Mississippi |
    +-----------+-------------+
    | 3         | Hidden      |
    +-----------+-------------+
    | 4         | Millrock    |
    +-----------+-------------+
    | 5         | Blue Lagoon |
    +-----------+-------------+


Adding contacts
---------------

In our Flin Flon project, we have a geological map from which we can pick the base of
the main geological units.  We begin by creating some points (``P+left-clicks``)
along the boundary of the Blue Lagoon and Millrock units.  The picks we make here will
be tagged as Millrock since the younging direction is to the East and Gempy treats
contact points as the base of a unit.

.. figure:: /images/tutorial/1-create_points.png
    :align: center
    :width: 80%

    *Picking contact points at the base of the Millrock unit using P+left-click.*

This will result in a new ``Points`` object in our project tree.  It is good practise
to rename this object for easy referencing later.

.. figure:: /images/tutorial/2-rename_observations.png
    :align: center
    :width: 80%

    *Renaming the 'Points' object to 'Observations' for easy referencing.*

Notice that we have started with the *second* oldest unit that we identified in the
:ref:`units and events table <history_table>`.  This is because ``Gempy`` will
automatically add a basement layer that fills the remaining space below the oldest
unit.  So, by starting at the *Millrock*, ``Gempy`` will add our *Blue Lagoon* unit
for us and we simply need to rename it after ``gempy-drivers`` has run.

In order to tag the points as Millrock, we need to create a new data channel in our
*observations* object.  This is done by right-clicking on the *observations* object
and selecting the Script menu.  We will assign the integer value corresponding to the
*Millrock* unit as planned in the :ref:`units and events table <history_table>`.

.. figure:: /images/tutorial/4-history_script.png
    :align: center
    :width: 60%

    *Creating a new integer history channel in the 'Observations' object.*

We can now pick and mark all of the other geological contacts with their respective unit ids
in the :ref:`units and events table <history_table>`.

.. figure:: /images/tutorial/13-update_events.png
    :align: center
    :width: 40%

    *Assigning unit ids to the picked contact points.*

Once these are finished we can transfer the *History* data to a ``ReferenceData`` object by
right-clicking on the object and selecting the transfer option.

.. figure:: /images/tutorial/5-transfer_to_reference.png
    :align: center
    :width: 80%

    *Transferring the history data to a ReferenceData object.*

We can assign the units names to the integer values by loading a lookup table from text file.
The file must contain integer values and their corresponding names, but may also optionally
contain RGB values for visualization.

.. figure:: /images/tutorial/6-reference_table.png
    :align: center
    :width: 80%

    *Editing a text file to map units indices to unit names and optionally assign rgb ratios to control color assignment.*

The unit names can now be updated by loading from file and selecting the import by index option.

.. figure:: /images/tutorial/10-load_edited_history.png
    :align: center
    :width: 70%

    *Loading the edited history data to update the unit names.*

This will assign the unit names and any provided colors to the indexed contact points.

.. figure:: /images/tutorial/11-load_history.png
    :align: center
    :width: 100%

    *Updated contact classification.*


Adding orientations
-------------------

``Gempy`` also requires orientation data.  We can provide either ``Strike & dip`` or ``Dip
direction & dip`` data depending on the convention of the orientation markers provided.
First we need to create two new data channels for the observations strike
(or dip direction) and dip.

.. figure:: /images/tutorial/16-create_strike_channel.png
    :align: center
    :width: 80%

    *Creating new data channels for the strike and dip. Using the scripting function ndv() here
    gives a good default since most points won't have any orientation.*

Then we can add them to a group of the appropriate type (``Strike & dip`` in our case).

.. figure:: /images/tutorial/grouping_orientation_data.png
    :align: center
    :width: 100%

    *Grouping the orientation data into a Strike & dip' group.*


The orientation data can now be added to the existing picks by editing the observations and
entering values in the *Strike* and *Dip* columns.  For structural markers lying within the
middle of a geological unit, we must indicate to gempy that it doesn't belong to a contact
surface.  This is done by setting a :ref:`contact <contacts_and_faults>` channel,

.. figure:: /images/tutorial/22a-create-contact_channel.png
    :align: center
    :width: 80%

    *Creating a new boolean data channel for the contact points.*

and marking the observation with as a non-contact.

.. figure:: /images/tutorial/add_orientations.png
    :align: center
    :width: 100%

    *Adding a point to represent an orientation field measurement.*


Adding faults
-------------

If you have a fault in your area of interest, you can add it to the model by creating a new data
channel to mark fault observations.

.. figure:: /images/tutorial/24-create_fault_channel.png
    :align: center
    :width: 80%

    *Creating a new data channel for the fault observations.*

Since faults are not considered conformable, they will always end up in their own group,
so they must contain at least one orientation.  They must also all be marked as a
:ref:`contact <contacts_and_faults>` point.

.. figure:: /images/tutorial/25-update_faults.png
    :align: center
    :width: 80%

    *Add fault observations with the fault channel and contact data columns set to true.*


Running gempy-drivers
---------------------

Now that all the data objects have been created, we can edit the options and run
gempy-drivers from the contextual menu (right+click) of the application group.

.. figure:: /images/tutorial/26-edit_options.png
    :align: center
    :width: 80%

    *Opening the gempy-drivers dialogue.*

The inputs to the gempy-drivers program are described in the :ref:`gempy rules section <gempy-drivers input>`.
We can now select the *observations* object, the required and optional data channels,
a mesh to interpolate the geological model to, and a optionally a topography surface
to cut the model.

.. figure:: /images/tutorial/ui_json_menu.png
    :align: center
    :width: 60%

    *Choosing the input objects and data to run gempy-drivers.*

To run gempy-drivers with the selected data, right-click the gempy-drivers group and
choose the ``Run locally`` option.

.. figure:: /images/tutorial/27-run_local.png
    :align: center
    :width: 80%

    *Running gempy-drivers locally.*

This will spawn a new terminal and run gempy-drivers with some feedback on the progress.

.. figure:: /images/tutorial/28-gempy_drivers_logging.png
    :align: center
    :width: 80%

    *Feedback on the structural stack generated by Gempy from our input data and the condition number of the interpolation matrix.*

This includes a rundown of the structural stack generated by gempy from the input data
and the condition number in the matrix responsible for the interpolation.  This
information may be useful for quality control and debugging models that may have gone
wrong due to bad input data.

When the results are ready, they will be stored in a new group containing all the
generated surfaces as well as a model containing the geological unit reference values.

.. figure:: /images/tutorial/29-load_results.png
    :align: center
    :width: 70%

    *Loading the results of the gempy-drivers run.*

.. figure:: /images/tutorial/30-results_group.png
    :align: center
    :width: 70%

    *The results group contains the generated surfaces and geological model.*

Forward modelling
^^^^^^^^^^^^^^^^^

Physical property model
-----------------------

In order to generate synthetic data over this model we must first populate the model
with physical property values.  This can be done from the *History* reference data
map.  Select the *History* in the project tree and click on the *add a column* button
in the *Data Colors* tab to create a physical property column and populate it with
physical property values.  This will create a new data channel in the provided mesh
that contains the physical property values assigned to each unit.


.. figure:: /images/tutorial/31-add_column_to_references.png
    :align: center
    :width: 80%

    *Adding a new column to the reference data map to enter the physical properties.*

.. figure:: /images/tutorial/32-add_physical_properties.png
    :align: center
    :width: 50%

    *Assigning physical property values to the geological units.*

In order to use the physical property data you will also need to transfer the geometric
data into a ``FloatData`` object.

.. figure:: /images/tutorial/33-transfer_fx_to_float.png
    :align: center
    :width: 70%

    *Transferring the geometric data object to float data.*

Simulation
----------

At this point we have everything we need to run a forward modelling of our geological
(physical property) model.  This can be done easily with the `simpeg-driver <https://github.com/MiraGeoscience/simpeg-drivers>`_ package.
To learn more about this, you may refer to the `documentation <https://mirageoscience-simpeg-drivers.readthedocs-hosted.com/en/stable/intro.html>`_
for the simpeg-drivers package.  Pro Geophysics users of `Geoscience ANALYST <https://www.mirageoscience.com/mining-industry-software/geoscience-analyst/pro-geophysics/>`_
can access this functionality from the *Geophysics -> Simpeg Python Interface* sub-menu.

.. figure:: /images/tutorial/35-simpeg_sub_menu.png
    :align: center
    :width: 80%

    *Accessing the Gravity forward modelling from within ANALYST.*

The end result of the forward modelling will be accessible from a SimPEG group and may
be visualized over the created gempy model.

.. figure:: /images/tutorial/34-simulated_data.png
    :align: center
    :width: 100%

    *Simulated gravity over the geology model populated with density values.*