.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "examples/11_plurigaussian/08_multi_field_pgs.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

    .. note::
        :class: sphx-glr-download-link-note

        :ref:`Go to the end <sphx_glr_download_examples_11_plurigaussian_08_multi_field_pgs.py>`
        to download the full example code.

.. rst-class:: sphx-glr-example-title

.. _sphx_glr_examples_11_plurigaussian_08_multi_field_pgs.py:


From bigaussian to plurigaussian simulation
--------------------------------------------------

In many PGS implementations, the dimensions of the simulation domain often
matches the number of fields that are supplied. However, this is not a
requirement of the PGS algorithm. In fact, it is possible to use multiple
spatial random fields in PGS, which can be useful for more complex lithotype
definitions. In this example, we will demonstrate how to use multiple SRFs in
PGS. In GSTools, this is done by utilising the tree based architecture.

Typically, PGS in two dimensions is carried out as a bigaussian simulation,
where two random fields are used. Here, we will employ four. We begin by
setting up the simulation domain and generating the necessary random fields,
where the length scales of two of the fields are much larger than the other two.

.. GENERATED FROM PYTHON SOURCE LINES 17-39

.. code-block:: Python


    import matplotlib.pyplot as plt
    import numpy as np

    import gstools as gs

    dim = 2

    N = [200, 200]

    x = np.arange(N[0])
    y = np.arange(N[1])

    model = gs.Gaussian(dim=dim, var=1, len_scale=15)
    srf = gs.SRF(model)
    field1 = srf.structured([x, y], seed=201519)
    field2 = srf.structured([x, y], seed=199221)
    model = gs.Gaussian(dim=dim, var=1, len_scale=3)
    srf = gs.SRF(model)
    field3 = srf.structured([x, y], seed=1345134)
    field4 = srf.structured([x, y], seed=1351455)








.. GENERATED FROM PYTHON SOURCE LINES 40-41

As in the previous example, an ellipse is used as the decision boundary.

.. GENERATED FROM PYTHON SOURCE LINES 41-54

.. code-block:: Python



    def ellipse(data, key1, key2, c1, c2, s1, s2, angle=0):
        x, y = data[key1] - c1, data[key2] - c2

        if angle:
            theta = np.deg2rad(angle)
            c, s = np.cos(theta), np.sin(theta)
            x, y = x * c + y * s, -x * s + y * c

        return (x / s1) ** 2 + (y / s2) ** 2 <= 1









.. GENERATED FROM PYTHON SOURCE LINES 55-60

The configuration dictionary for the decision tree is defined as before, but
this time we pass the additional keys `Z3` and `Z4`, which refer to the
additional fields `field3` and `field4`. The decision tree is structured in a
way that the first decision node is based on the first two fields, and the
second decision node is based on the last two fields.

.. GENERATED FROM PYTHON SOURCE LINES 60-96

.. code-block:: Python


    config = {
        "root": {
            "type": "decision",
            "func": ellipse,
            "args": {
                "key1": "Z1",
                "key2": "Z2",
                "c1": 0.7,
                "c2": 0.7,
                "s1": 1.3,
                "s2": 1.3,
            },
            "yes_branch": "phase1",
            "no_branch": "node1",
        },
        "node1": {
            "type": "decision",
            "func": ellipse,
            "args": {
                "key1": "Z3",
                "key2": "Z4",
                "c1": -0.7,
                "c2": -0.7,
                "s1": 1.3,
                "s2": 1.3,
                "angle": 30,
            },
            "yes_branch": "phase2",
            "no_branch": "phase0",
        },
        "phase2": {"type": "leaf", "action": 2},
        "phase1": {"type": "leaf", "action": 1},
        "phase0": {"type": "leaf", "action": 0},
    }








.. GENERATED FROM PYTHON SOURCE LINES 97-99

Next, we create the PGS object, passing in all four fields, and generate the
plurigaussian field `P`

.. GENERATED FROM PYTHON SOURCE LINES 99-104

.. code-block:: Python


    pgs = gs.PGS(dim, [field1, field2, field3, field4])

    P = pgs(tree=config)








.. GENERATED FROM PYTHON SOURCE LINES 105-109

Finally, we plot `P`

NB: In the current implementation, the calculation of the equivalent spatial
lithotype `L` is not supported for multiple fields

.. GENERATED FROM PYTHON SOURCE LINES 109-113

.. code-block:: Python


    plt.figure(figsize=(8, 6))
    plt.imshow(P, cmap="copper", origin="lower")




.. image-sg:: /examples/11_plurigaussian/images/sphx_glr_08_multi_field_pgs_001.png
   :alt: 08 multi field pgs
   :srcset: /examples/11_plurigaussian/images/sphx_glr_08_multi_field_pgs_001.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 114-117

.. image:: ../../pics/2d_multi_tree_pgs.png
   :width: 400px
   :align: center


.. rst-class:: sphx-glr-timing

   **Total running time of the script:** (0 minutes 4.791 seconds)


.. _sphx_glr_download_examples_11_plurigaussian_08_multi_field_pgs.py:

.. only:: html

  .. container:: sphx-glr-footer sphx-glr-footer-example

    .. container:: sphx-glr-download sphx-glr-download-jupyter

      :download:`Download Jupyter notebook: 08_multi_field_pgs.ipynb <08_multi_field_pgs.ipynb>`

    .. container:: sphx-glr-download sphx-glr-download-python

      :download:`Download Python source code: 08_multi_field_pgs.py <08_multi_field_pgs.py>`

    .. container:: sphx-glr-download sphx-glr-download-zip

      :download:`Download zipped: 08_multi_field_pgs.zip <08_multi_field_pgs.zip>`


.. only:: html

 .. rst-class:: sphx-glr-signature

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_