.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "auto_examples/howto_meshlib/plot_sample_mesh_line.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

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

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

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

.. _sphx_glr_auto_examples_howto_meshlib_plot_sample_mesh_line.py:


Extract a 1D profile from 2D and plot it
========================================

.. sectionauthor:: Feliks Kiszkurno (Helmholtz Centre for Environmental Research GmbH - UFZ)

For this example we load a 2D meshseries from within the ``meshplotlib``
examples. In the ``meshplotlib.setup`` we can provide a dictionary to map names
to material ids. First, let's plot the material ids (cell_data). Per default in
the setup, this will automatically show the element edges.

.. GENERATED FROM PYTHON SOURCE LINES 14-16

0. Setup and imports
--------------------

.. GENERATED FROM PYTHON SOURCE LINES 18-32

.. code-block:: Python

    import matplotlib.pyplot as plt
    import numpy as np

    import ogstools.meshplotlib as mpl
    from ogstools import examples
    from ogstools.meshlib import sample_polyline
    from ogstools.meshplotlib import lineplot, plot_profile
    from ogstools.propertylib import Scalar, properties

    mpl.setup.reset()
    mpl.setup.combined_colorbar = False
    plt.rcdefaults()
    # mpl.setup.length.output_unit = "km"








.. GENERATED FROM PYTHON SOURCE LINES 33-37

1. Single fracture
------------------
Define profile line by providing a list of points in an X, Y, Z coordinates
and loading an example data set:

.. GENERATED FROM PYTHON SOURCE LINES 39-43

.. code-block:: Python

    ms_HT = examples.load_meshseries_HT_2D_XDMF()

    profile_HT = np.array([[4, 2, 0], [4, 18, 0]])








.. GENERATED FROM PYTHON SOURCE LINES 44-50

.. code-block:: Python

    ms_HT_sp, ms_HT_kp = sample_polyline(
        ms_HT.read(-1),
        ["pressure", "temperature"],
        profile_HT,
    )








.. GENERATED FROM PYTHON SOURCE LINES 51-54

It has returned a DataFrame containing all information about profile and
numpy array with position of the "knot-points".
Let's investigate the DataFrame first:

.. GENERATED FROM PYTHON SOURCE LINES 56-58

.. code-block:: Python

    ms_HT_sp.head(10)






.. raw:: html

    <div class="output_subarea output_html rendered_html output_result">
    <div>
    <style scoped>
        .dataframe tbody tr th:only-of-type {
            vertical-align: middle;
        }

        .dataframe tbody tr th {
            vertical-align: top;
        }

        .dataframe thead th {
            text-align: right;
        }
    </style>
    <table border="1" class="dataframe">
      <thead>
        <tr style="text-align: right;">
          <th></th>
          <th>x</th>
          <th>y</th>
          <th>z</th>
          <th>pressure</th>
          <th>temperature</th>
          <th>dist</th>
          <th>dist_in_segment</th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>0</th>
          <td>4.0</td>
          <td>2.00</td>
          <td>0.0</td>
          <td>0.498440</td>
          <td>79.850025</td>
          <td>0.00</td>
          <td>0.00</td>
        </tr>
        <tr>
          <th>1</th>
          <td>4.0</td>
          <td>2.16</td>
          <td>0.0</td>
          <td>0.498725</td>
          <td>79.850031</td>
          <td>0.16</td>
          <td>0.16</td>
        </tr>
        <tr>
          <th>2</th>
          <td>4.0</td>
          <td>2.32</td>
          <td>0.0</td>
          <td>0.499150</td>
          <td>79.850021</td>
          <td>0.32</td>
          <td>0.32</td>
        </tr>
        <tr>
          <th>3</th>
          <td>4.0</td>
          <td>2.48</td>
          <td>0.0</td>
          <td>0.499663</td>
          <td>79.849999</td>
          <td>0.48</td>
          <td>0.48</td>
        </tr>
        <tr>
          <th>...</th>
          <td>...</td>
          <td>...</td>
          <td>...</td>
          <td>...</td>
          <td>...</td>
          <td>...</td>
          <td>...</td>
        </tr>
        <tr>
          <th>6</th>
          <td>4.0</td>
          <td>2.96</td>
          <td>0.0</td>
          <td>0.501203</td>
          <td>79.849934</td>
          <td>0.96</td>
          <td>0.96</td>
        </tr>
        <tr>
          <th>7</th>
          <td>4.0</td>
          <td>3.12</td>
          <td>0.0</td>
          <td>0.501717</td>
          <td>79.849913</td>
          <td>1.12</td>
          <td>1.12</td>
        </tr>
        <tr>
          <th>8</th>
          <td>4.0</td>
          <td>3.28</td>
          <td>0.0</td>
          <td>0.502230</td>
          <td>79.849891</td>
          <td>1.28</td>
          <td>1.28</td>
        </tr>
        <tr>
          <th>9</th>
          <td>4.0</td>
          <td>3.44</td>
          <td>0.0</td>
          <td>0.502936</td>
          <td>79.849931</td>
          <td>1.44</td>
          <td>1.44</td>
        </tr>
      </tbody>
    </table>
    <p>10 rows × 7 columns</p>
    </div>
    </div>
    <br />
    <br />

.. GENERATED FROM PYTHON SOURCE LINES 59-68

We can see the spatial coordinates of points on the profile ("x", "y", "z"
- columns), distances from the beginning of the profile ("dist") and within
current segment ("dist_in_segment"). Note, that since we defined our profile
on only two points, there is only one segment, hence in this special case
columns dist and dist_in_segment are identical. At the end of the DataFrame
we can can find two columns with the properties that we are interested in:
"temperature" and "pressure". Each occupies one column, as those are scalar
values. Using columns "dist", "pressure" and "temperature" we can easily
plot the data:

.. GENERATED FROM PYTHON SOURCE LINES 70-83

.. code-block:: Python

    fig, ax = plt.subplots(1, 1, figsize=(7, 5))
    ax = lineplot(
        x="dist",
        y=["pressure", "temperature"],
        mesh=ms_HT.read(-1),
        profile_points=profile_HT,
        ax=ax,
        twinx=True,
        fontsize=15,
    )
    fig.tight_layout()





.. image-sg:: /auto_examples/howto_meshlib/images/sphx_glr_plot_sample_mesh_line_001.png
   :alt: plot sample mesh line
   :srcset: /auto_examples/howto_meshlib/images/sphx_glr_plot_sample_mesh_line_001.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 84-86

What happens when we are interested in vector property? We can see it in the
following example using Darcy velocity:

.. GENERATED FROM PYTHON SOURCE LINES 88-94

.. code-block:: Python

    ms_HT_sp, ms_HT_kp = sample_polyline(
        ms_HT.read(-1),
        "darcy_velocity",
        profile_HT,
    )








.. GENERATED FROM PYTHON SOURCE LINES 95-97

.. code-block:: Python

    ms_HT_sp.head(5)






.. raw:: html

    <div class="output_subarea output_html rendered_html output_result">
    <div>
    <style scoped>
        .dataframe tbody tr th:only-of-type {
            vertical-align: middle;
        }

        .dataframe tbody tr th {
            vertical-align: top;
        }

        .dataframe thead th {
            text-align: right;
        }
    </style>
    <table border="1" class="dataframe">
      <thead>
        <tr style="text-align: right;">
          <th></th>
          <th>x</th>
          <th>y</th>
          <th>z</th>
          <th>darcy_velocity_0</th>
          <th>darcy_velocity_1</th>
          <th>dist</th>
          <th>dist_in_segment</th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>0</th>
          <td>4.0</td>
          <td>2.00</td>
          <td>0.0</td>
          <td>1.169737e-93</td>
          <td>-2.231183e-94</td>
          <td>0.00</td>
          <td>0.00</td>
        </tr>
        <tr>
          <th>1</th>
          <td>4.0</td>
          <td>2.16</td>
          <td>0.0</td>
          <td>1.196810e-93</td>
          <td>-2.421097e-94</td>
          <td>0.16</td>
          <td>0.16</td>
        </tr>
        <tr>
          <th>2</th>
          <td>4.0</td>
          <td>2.32</td>
          <td>0.0</td>
          <td>1.237164e-93</td>
          <td>-2.637041e-94</td>
          <td>0.32</td>
          <td>0.32</td>
        </tr>
        <tr>
          <th>3</th>
          <td>4.0</td>
          <td>2.48</td>
          <td>0.0</td>
          <td>1.285968e-93</td>
          <td>-2.869550e-94</td>
          <td>0.48</td>
          <td>0.48</td>
        </tr>
        <tr>
          <th>4</th>
          <td>4.0</td>
          <td>2.64</td>
          <td>0.0</td>
          <td>1.334773e-93</td>
          <td>-3.102060e-94</td>
          <td>0.64</td>
          <td>0.64</td>
        </tr>
      </tbody>
    </table>
    </div>
    </div>
    <br />
    <br />

.. GENERATED FROM PYTHON SOURCE LINES 98-106

Now we have two columns for the property. Darcy velocity is a vector,
therefore "sample_over_polyline" has split it into two columns and appended
the property name with increasing integer. Note, that this suffix has no
physical meaning and only indicates order. It is up to user to interpret it
in a meaningful way. By the
`OpenGeoSys conventions <https://www.opengeosys.org/docs/userguide/basics/conventions/#a-namesymmetric-tensorsa--symmetric-tensors-and-kelvin-mapping>`_,
"darcy_velocity_0" will be in the x-direction and "darcy_velocity_1" in
y-direction.

.. GENERATED FROM PYTHON SOURCE LINES 109-114

2. Elder benchmark
------------------
In this example we will provide a Scalar object from the propertylib to
define the property of interest. "sample_over_polyline" will automatically
convert the data from "data_unit" to "output_unit":

.. GENERATED FROM PYTHON SOURCE LINES 116-129

.. code-block:: Python

    profile_CT = np.array(
        [
            [47.0, 1.17, 72.0],  # Point A
            [-4.5, 1.17, -59.0],  # Point B
        ]
    )

    ms_CT = examples.load_meshseries_CT_2D_XDMF()

    si = Scalar(
        data_name="Si", data_unit="", output_unit="%", output_name="Saturation"
    )








.. GENERATED FROM PYTHON SOURCE LINES 130-136

.. code-block:: Python

    ms_CT_sp, ms_CT_kp = sample_polyline(
        ms_CT.read(11),
        si,
        profile_CT,
    )








.. GENERATED FROM PYTHON SOURCE LINES 137-139

As before we can see the profile parameters and properties values in a
DataFrame:

.. GENERATED FROM PYTHON SOURCE LINES 141-144

.. code-block:: Python

    ms_CT_sp.head(5)







.. raw:: html

    <div class="output_subarea output_html rendered_html output_result">
    <div>
    <style scoped>
        .dataframe tbody tr th:only-of-type {
            vertical-align: middle;
        }

        .dataframe tbody tr th {
            vertical-align: top;
        }

        .dataframe thead th {
            text-align: right;
        }
    </style>
    <table border="1" class="dataframe">
      <thead>
        <tr style="text-align: right;">
          <th></th>
          <th>x</th>
          <th>y</th>
          <th>z</th>
          <th>Si</th>
          <th>dist</th>
          <th>dist_in_segment</th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>0</th>
          <td>47.000</td>
          <td>1.17</td>
          <td>72.00</td>
          <td>79.601252</td>
          <td>0.000000</td>
          <td>0.000000</td>
        </tr>
        <tr>
          <th>1</th>
          <td>46.485</td>
          <td>1.17</td>
          <td>70.69</td>
          <td>73.478046</td>
          <td>1.407595</td>
          <td>1.407595</td>
        </tr>
        <tr>
          <th>2</th>
          <td>45.970</td>
          <td>1.17</td>
          <td>69.38</td>
          <td>70.034378</td>
          <td>2.815191</td>
          <td>2.815191</td>
        </tr>
        <tr>
          <th>3</th>
          <td>45.455</td>
          <td>1.17</td>
          <td>68.07</td>
          <td>67.833706</td>
          <td>4.222786</td>
          <td>4.222786</td>
        </tr>
        <tr>
          <th>4</th>
          <td>44.940</td>
          <td>1.17</td>
          <td>66.76</td>
          <td>67.245709</td>
          <td>5.630382</td>
          <td>5.630382</td>
        </tr>
      </tbody>
    </table>
    </div>
    </div>
    <br />
    <br />

.. GENERATED FROM PYTHON SOURCE LINES 145-147

This time we will prepare more complicated plot showing both mesh data and
the profile.

.. GENERATED FROM PYTHON SOURCE LINES 149-160

.. code-block:: Python

    fig, ax = plot_profile(
        ms_CT.read(11),
        si,
        profile_CT,
        resolution=100,
        profile_plane=[0, 2],  # This profile is in XZ plane, not XY!
    )
    fig, ax = mpl.update_font_sizes(label_axes="none", fig=fig)
    fig.tight_layout()





.. image-sg:: /auto_examples/howto_meshlib/images/sphx_glr_plot_sample_mesh_line_002.png
   :alt: plot sample mesh line
   :srcset: /auto_examples/howto_meshlib/images/sphx_glr_plot_sample_mesh_line_002.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 161-166

3. THM
------
It is also possible to obtain more than one property at the same time using
more complex profiles. They can be constructed by providing more than
2 points. With those points:

.. GENERATED FROM PYTHON SOURCE LINES 168-185

.. code-block:: Python

    profile_THM = np.array(
        [
            [-1000.0, -175.0, 6700.0],  # Point A
            [-600.0, -600.0, 6700.0],  # Point B
            [100.0, -300.0, 6700.0],  # Point C
            [910.0, -590.0, 6700.0],  # Point D
        ]
    )

    profile_THM = np.array(
        [
            [-1000.0, -175.0, 6700.0],  # Point A
            [-600.0, -600.0, 6700.0],  # Point B
            [100.0, -300.0, 6700.0],  # Point C
            [3500, -900.0, 6700.0],  # Point D
        ]
    )







.. GENERATED FROM PYTHON SOURCE LINES 186-196

the profile will run as follows:

.. math::

  \text{AB} \rightarrow \text{BC} \rightarrow \text{CD}

Point B will at the same time be the last point in the first segment AB
and first one in second segment BC, however in the returned array,
it will occur only once.
For this example we will use a different dataset:

.. GENERATED FROM PYTHON SOURCE LINES 198-200

.. code-block:: Python

    meshseries_THM = examples.load_meshseries_THM_2D_PVD()








.. GENERATED FROM PYTHON SOURCE LINES 201-208

.. code-block:: Python

    ms_THM_sp, dist_at_knot = sample_polyline(
        meshseries_THM.read(-1),
        [properties.pressure, properties.temperature],
        profile_THM,
        resolution=100,
    )








.. GENERATED FROM PYTHON SOURCE LINES 209-211

Again, we can investigate the returned DataFrame, but this time we will
have a look at its beginning:

.. GENERATED FROM PYTHON SOURCE LINES 213-215

.. code-block:: Python

    ms_THM_sp.head(5)






.. raw:: html

    <div class="output_subarea output_html rendered_html output_result">
    <div>
    <style scoped>
        .dataframe tbody tr th:only-of-type {
            vertical-align: middle;
        }

        .dataframe tbody tr th {
            vertical-align: top;
        }

        .dataframe thead th {
            text-align: right;
        }
    </style>
    <table border="1" class="dataframe">
      <thead>
        <tr style="text-align: right;">
          <th></th>
          <th>x</th>
          <th>y</th>
          <th>z</th>
          <th>pressure</th>
          <th>temperature</th>
          <th>dist</th>
          <th>dist_in_segment</th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>0</th>
          <td>-1000.000000</td>
          <td>-175.000000</td>
          <td>6700.0</td>
          <td>2.230181</td>
          <td>8.832194</td>
          <td>0.000000</td>
          <td>0.000000</td>
        </tr>
        <tr>
          <th>1</th>
          <td>-969.230769</td>
          <td>-207.692308</td>
          <td>6700.0</td>
          <td>2.551761</td>
          <td>8.828883</td>
          <td>44.894683</td>
          <td>44.894683</td>
        </tr>
        <tr>
          <th>2</th>
          <td>-938.461538</td>
          <td>-240.384615</td>
          <td>6700.0</td>
          <td>2.873396</td>
          <td>8.825926</td>
          <td>89.789366</td>
          <td>89.789366</td>
        </tr>
        <tr>
          <th>3</th>
          <td>-907.692308</td>
          <td>-273.076923</td>
          <td>6700.0</td>
          <td>3.195135</td>
          <td>8.825316</td>
          <td>134.684048</td>
          <td>134.684048</td>
        </tr>
        <tr>
          <th>4</th>
          <td>-876.923077</td>
          <td>-305.769231</td>
          <td>6700.0</td>
          <td>3.516919</td>
          <td>8.821331</td>
          <td>179.578731</td>
          <td>179.578731</td>
        </tr>
      </tbody>
    </table>
    </div>
    </div>
    <br />
    <br />

.. GENERATED FROM PYTHON SOURCE LINES 216-217

and end:

.. GENERATED FROM PYTHON SOURCE LINES 219-221

.. code-block:: Python

    ms_THM_sp.tail(10)






.. raw:: html

    <div class="output_subarea output_html rendered_html output_result">
    <div>
    <style scoped>
        .dataframe tbody tr th:only-of-type {
            vertical-align: middle;
        }

        .dataframe tbody tr th {
            vertical-align: top;
        }

        .dataframe thead th {
            text-align: right;
        }
    </style>
    <table border="1" class="dataframe">
      <thead>
        <tr style="text-align: right;">
          <th></th>
          <th>x</th>
          <th>y</th>
          <th>z</th>
          <th>pressure</th>
          <th>temperature</th>
          <th>dist</th>
          <th>dist_in_segment</th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>92</th>
          <td>3075.000000</td>
          <td>-825.000000</td>
          <td>6700.0</td>
          <td>0.0</td>
          <td>32.327025</td>
          <td>4366.176575</td>
          <td>3020.968388</td>
        </tr>
        <tr>
          <th>93</th>
          <td>3122.222222</td>
          <td>-833.333333</td>
          <td>6700.0</td>
          <td>0.0</td>
          <td>35.663869</td>
          <td>4414.128454</td>
          <td>3068.920267</td>
        </tr>
        <tr>
          <th>94</th>
          <td>3169.444444</td>
          <td>-841.666667</td>
          <td>6700.0</td>
          <td>0.0</td>
          <td>40.358542</td>
          <td>4462.080333</td>
          <td>3116.872146</td>
        </tr>
        <tr>
          <th>95</th>
          <td>3216.666667</td>
          <td>-850.000000</td>
          <td>6700.0</td>
          <td>0.0</td>
          <td>44.984874</td>
          <td>4510.032212</td>
          <td>3164.824025</td>
        </tr>
        <tr>
          <th>...</th>
          <td>...</td>
          <td>...</td>
          <td>...</td>
          <td>...</td>
          <td>...</td>
          <td>...</td>
          <td>...</td>
        </tr>
        <tr>
          <th>98</th>
          <td>3358.333333</td>
          <td>-875.000000</td>
          <td>6700.0</td>
          <td>0.0</td>
          <td>34.045579</td>
          <td>4653.887850</td>
          <td>3308.679663</td>
        </tr>
        <tr>
          <th>99</th>
          <td>3405.555556</td>
          <td>-883.333333</td>
          <td>6700.0</td>
          <td>0.0</td>
          <td>31.436554</td>
          <td>4701.839729</td>
          <td>3356.631542</td>
        </tr>
        <tr>
          <th>100</th>
          <td>3452.777778</td>
          <td>-891.666667</td>
          <td>6700.0</td>
          <td>0.0</td>
          <td>28.867132</td>
          <td>4749.791608</td>
          <td>3404.583421</td>
        </tr>
        <tr>
          <th>101</th>
          <td>3500.000000</td>
          <td>-900.000000</td>
          <td>6700.0</td>
          <td>0.0</td>
          <td>26.706150</td>
          <td>4797.743487</td>
          <td>3452.535300</td>
        </tr>
      </tbody>
    </table>
    <p>10 rows × 7 columns</p>
    </div>
    </div>
    <br />
    <br />

.. GENERATED FROM PYTHON SOURCE LINES 222-225

Note, that unlike in the first example, here the columns "dist" and
"dist_in_segment" are not identical, as this time profile consists of
multiple segments. Following figure illustrates the difference:

.. GENERATED FROM PYTHON SOURCE LINES 225-234

.. code-block:: Python

    plt.rcdefaults()
    fig, ax = plt.subplots(1, 1, figsize=(7, 3))
    ax.plot(ms_THM_sp["dist"], label="dist")
    ax.plot(ms_THM_sp["dist_in_segment"], label="dist_in_segment")
    ax.set_xlabel("Point ID / -")
    ax.set_ylabel("Distance / m")
    ax.legend()
    fig.tight_layout()




.. image-sg:: /auto_examples/howto_meshlib/images/sphx_glr_plot_sample_mesh_line_003.png
   :alt: plot sample mesh line
   :srcset: /auto_examples/howto_meshlib/images/sphx_glr_plot_sample_mesh_line_003.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 235-240

Orange line returns to 0 twice. It is because of how the overlap of nodal
points between segments is handled. Nodal point always belongs to the
segment it starts: point B is included in segment BC but not AB and point
C in CD but not in in BC.
Following figure shows the profile on the mesh:

.. GENERATED FROM PYTHON SOURCE LINES 242-243

plt.rcdefaults()

.. GENERATED FROM PYTHON SOURCE LINES 243-251

.. code-block:: Python

    fig, ax = plot_profile(
        meshseries_THM.read(-1),
        [properties.pressure, properties.temperature],
        profile_THM,
        resolution=100,
    )
    fig, ax = mpl.update_font_sizes(label_axes="none", fig=fig)
    fig.tight_layout()



.. image-sg:: /auto_examples/howto_meshlib/images/sphx_glr_plot_sample_mesh_line_004.png
   :alt: plot sample mesh line
   :srcset: /auto_examples/howto_meshlib/images/sphx_glr_plot_sample_mesh_line_004.png
   :class: sphx-glr-single-img






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

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


.. _sphx_glr_download_auto_examples_howto_meshlib_plot_sample_mesh_line.py:

.. only:: html

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

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

      :download:`Download Jupyter notebook: plot_sample_mesh_line.ipynb <plot_sample_mesh_line.ipynb>`

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

      :download:`Download Python source code: plot_sample_mesh_line.py <plot_sample_mesh_line.py>`