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

.. only:: html

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

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

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

.. _sphx_glr_examples_Creating_GS_Files_plot_tifs.py:


GeoTIFFs (2D & 3D)
------------------

In this example, we demonstrate the workflow for creating a GS file from the GeoTIFF (.tif/.tiff) file format. This includes adding individual TIF files as single 2-D variables, as well as how to create a 3-D variable by stacking multiple TIF files along a specified dimension.

The GS standard requires a single set of ``x``, ``y``, ``z``, ``t`` coordinate variables per data leaf group. Therefore, this example also shows how tifs with differing x-y grids need to be added to separate groups, and all variables in a group should have matching coordinates:

1. Raster Dataset #1
    - 2-D magnetic grid, original x-y discretization (600 m cell size)
2. Raster Dataset #2
    - 2-D magnetic grid, aligned to match the x-y dimensions of the resistivity layers (1000 m cell size)
    - 3-D resistivity grid

Lastly, GSPy provides a "to_tif()" method to export raster data as GeoTIFF. This example demonstrates how to use this method for both 2D and 3D variables.

Source References: 

Minsley, B.J., James, S.R., Bedrosian, P.A., Pace, M.D., Hoogenboom, B.E., and Burton, B.L., 2021, Airborne electromagnetic, magnetic, and radiometric survey of the Mississippi Alluvial Plain, November 2019 - March 2020: U.S. Geological Survey data release, https://doi.org/10.5066/P9E44CTQ.

James, S.R., and Minsley, B.J., 2021, Combined results and derivative products of hydrogeologic structure and properties from airborne electromagnetic surveys in the Mississippi Alluvial Plain: U.S. Geological Survey data release, https://doi.org/10.5066/P9382RCI.

.. GENERATED FROM PYTHON SOURCE LINES 25-30

.. code-block:: Python

    import matplotlib.pyplot as plt
    from os.path import join
    import gspy
    from gspy import Survey








.. GENERATED FROM PYTHON SOURCE LINES 31-33

Initialize the Survey
^^^^^^^^^^^^^^^^^^^^^

.. GENERATED FROM PYTHON SOURCE LINES 33-43

.. code-block:: Python


    # Path to example files
    data_path = "..//data_files//tempest_aseg"

    # Survey metadata file
    metadata = join(data_path, "data//Tempest_survey_md.yml")

    # Establish the Survey
    survey = Survey.from_dict(metadata)








.. GENERATED FROM PYTHON SOURCE LINES 44-49

.. literalinclude:: /../../examples/data_files/tempest_aseg/data/Tempest_survey_md.yml
   :language: yaml
   :linenos:
   :caption: Survey YAML file


.. GENERATED FROM PYTHON SOURCE LINES 52-54

Create a branch for all maps
^^^^^^^^^^^^^^^^^^^^^^^^^^^^

.. GENERATED FROM PYTHON SOURCE LINES 54-59

.. code-block:: Python


    container = survey.gs.add_container('derived_products',
                                        **dict(content = "container of gridded maps of magnetic and electrical resistivity values",
                                               comment = "Magnetic map is contractor-derived, resistivity maps are USGS-derived"))








.. GENERATED FROM PYTHON SOURCE LINES 60-62

Attach the 2D Magnetic Raster Dataset (600 m cell size)
"""""""""""""""""""""""""""""""""""""""""""""""""""""""

.. GENERATED FROM PYTHON SOURCE LINES 64-66

Define input metadata file, which contains the TIF filename linked 
with desired variable name and info. 

.. GENERATED FROM PYTHON SOURCE LINES 66-72

.. code-block:: Python


    d_supp1 = join(data_path, 'data//Tempest_raster_md.yml')

    # Attach the magnetic map to the container
    mm = container.gs.add(key="mag_map", metadata_file=d_supp1)








.. GENERATED FROM PYTHON SOURCE LINES 73-78

.. literalinclude:: /../../examples/data_files/tempest_aseg/data/Tempest_raster_md.yml
   :language: yaml
   :linenos:
   :caption: Magnetic raster (600 m cell size) YAML file


.. GENERATED FROM PYTHON SOURCE LINES 81-83

Attach 3D Resistivity Grids + Aligned 2D Magnetic Raster (1000 m cell size)
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""

.. GENERATED FROM PYTHON SOURCE LINES 85-86

Import both 3-D resistivity and 2-D magnetic data, aligned onto a common 1000 m x 1000 m grid

.. GENERATED FROM PYTHON SOURCE LINES 86-91

.. code-block:: Python

    d_supp2 = join(data_path, 'data//Tempest_rasters_md.yml')

    # Attach rasters to the container
    rm = container.gs.add(key="all_maps", metadata_file=d_supp2)








.. GENERATED FROM PYTHON SOURCE LINES 92-97

.. literalinclude:: /../../examples/data_files/tempest_aseg/data/Tempest_rasters_md.yml
   :language: yaml
   :linenos:
   :caption: Resistivity & Magnetic grids (1000 m cell size) YAML file


.. GENERATED FROM PYTHON SOURCE LINES 100-105

Note: the stack dimension is defined as "depth" but for the resistivity variable 
the dimensions are listed as [x, y, z] ... this is because under coordinates
we are linking the `z` coordinate dimension to the incoming depth dimension. 
Therefore the name "depth" ends up not getting used because that dimension becomes
"z" instead.

.. GENERATED FROM PYTHON SOURCE LINES 107-109

Inspect the data tree
^^^^^^^^^^^^^^^^^^^^^

.. GENERATED FROM PYTHON SOURCE LINES 109-111

.. code-block:: Python

    print(survey)





.. rst-class:: sphx-glr-script-out

 .. code-block:: none

    <xarray.DataTree 'survey'>
    Group: /survey
    │   Dimensions:                 ()
    │   Coordinates:
    │       spatial_ref             float64 8B 0.0
    │   Data variables:
    │       survey_information      float64 8B nan
    │       survey_units            float64 8B nan
    │       flightline_information  float64 8B nan
    │       survey_equipment        float64 8B nan
    │   Attributes:
    │       type:         survey
    │       title:        Example Tempest Airborne Electromagnetic (AEM) Dataset
    │       institution:  USGS Geology, Geophysics, & Geochemistry Science Center
    │       source:       Contractor provided ASEG-GDF2 formatted data
    │       history:      This example dataset includes the raw AEM data and gridded ...
    │       references:   Minsley, B.J., James, S.R., Bedrosian, P.A., Pace, M.D., Ho...
    │       comment:      This dataset is incomplete and has been downsampled for the...
    │       content:      Tempest AEM Survey from the Mississippi Alluvial Plain
    └── Group: /survey/derived_products
        │   Dimensions:      ()
        │   Data variables:
        │       spatial_ref  float64 8B 0.0
        │   Attributes:
        │       content:  container of gridded maps of magnetic and electrical resistivit...
        │       comment:  Magnetic map is contractor-derived, resistivity maps are USGS-d...
        │       type:     container
        ├── Group: /survey/derived_products/mag_map
        │   │   Dimensions:       (x: 599, nv: 2, y: 1212)
        │   │   Coordinates:
        │   │     * x             (x) float64 5kB 2.928e+05 2.934e+05 ... 6.51e+05 6.516e+05
        │   │     * nv            (nv) int64 16B 0 1
        │   │     * y             (y) float64 10kB 1.607e+06 1.606e+06 ... 8.808e+05 8.802e+05
        │   │       spatial_ref   float64 8B 0.0
        │   │   Data variables:
        │   │       x_bnds        (x, nv) float64 10kB 2.925e+05 2.931e+05 ... 6.519e+05
        │   │       y_bnds        (y, nv) float64 19kB 1.607e+06 1.606e+06 ... 8.799e+05
        │   │       magnetic_tmi  (y, x) float64 6MB nan nan nan nan nan ... nan nan nan nan nan
        │   │   Attributes:
        │   │       comment:     contractor-derived product
        │   │       content:     gridded map of total magnetic intensity
        │   │       type:        data
        │   │       mode:        airborne
        │   │       method:      magnetic
        │   │       instrument:  Scintrex CS-3 cesium vapor magnetometer
        │   │       structure:   raster
        │   └── Group: /survey/derived_products/mag_map/magnetic_system
        │           Dimensions:                             (base_mag_locations: 6, nv: 2,
        │                                                    n_transmitter: 1, n_receiver: 1,
        │                                                    n_couplet: 1, n_base_magnetometer: 1)
        │           Coordinates:
        │             * base_mag_locations                  (base_mag_locations) int64 48B 1 2 ... 6
        │             * n_transmitter                       (n_transmitter) int64 8B 0
        │             * n_receiver                          (n_receiver) int64 8B 0
        │             * n_couplet                           (n_couplet) int64 8B 0
        │             * n_base_magnetometer                 (n_base_magnetometer) int64 8B 0
        │           Data variables: (12/35)
        │               base_mag_locations_bnds             (base_mag_locations, nv) float64 96B ...
        │               transmitter_label                   (n_transmitter) <U7 28B 'passive'
        │               transmitter_description             (n_transmitter) <U107 428B "No artifi...
        │               receiver_label                      (n_receiver) <U19 76B 'scalar_magneto...
        │               receiver_sensor_type                (n_receiver) <U12 48B 'cesium_vapor'
        │               receiver_sensor_model               (n_receiver) <U4 16B 'CS-3'
        │               ...                                  ...
        │               diurnal_correction                  <U142 568B 'Diurnal base values appli...
        │               igrf_model_date                     <U10 40B '2019-12-01'
        │               igrf_model_height                   <U7 28B '167.9 m'
        │               igrf_removed_model_epoch            <U6 24B '2015.0'
        │               tieline_levelling                   <U137 548B 'RMI levelled to prior RES...
        │               deliverables                        <U63 252B 'Total Magnetic Intensity (...
        │           Attributes:
        │               type:        system
        │               mode:        airborne
        │               method:      magnetic
        │               instrument:  Scintrex CS-3 c cesium-vapor magnetometer
        │               name:        magnetic_system
        └── Group: /survey/derived_products/all_maps
                Dimensions:       (z: 5, nv: 2, x: 363, y: 770)
                Coordinates:
                  * z             (z) float64 40B 0.0 5.0 10.0 15.0 20.0
                  * nv            (nv) int64 16B 0 1
                  * x             (x) float64 3kB 2.915e+05 2.925e+05 ... 6.525e+05 6.535e+05
                  * y             (y) float64 6kB 1.648e+06 1.647e+06 ... 8.798e+05 8.788e+05
                    spatial_ref   float64 8B 0.0
                Data variables:
                    z_bnds        (z, nv) float64 80B -2.5 2.5 2.5 7.5 ... 12.5 17.5 17.5 22.5
                    x_bnds        (x, nv) float64 6kB 2.91e+05 2.92e+05 ... 6.53e+05 6.54e+05
                    y_bnds        (y, nv) float64 12kB 1.648e+06 1.647e+06 ... 8.783e+05
                    resistivity   (z, y, x) float64 11MB nan nan nan nan nan ... nan nan nan nan
                    magnetic_tmi  (y, x) float32 1MB nan nan nan nan nan ... nan nan nan nan nan
                Attributes:
                    comment:     resistivity models are USGS-derived, the total magnetic inte...
                    content:     interpolated resistivity models (3-D depth grid) and total m...
                    type:        data, models
                    method:      ['electromagnetic, time domain', 'magnetic']
                    instrument:  30Hz Tempest
                    mode:        airborne
                    property:    electrical resistivity, total magnetic intensity
                    structure:   raster




.. GENERATED FROM PYTHON SOURCE LINES 112-114

Save to NetCDF file
^^^^^^^^^^^^^^^^^^^

.. GENERATED FROM PYTHON SOURCE LINES 114-117

.. code-block:: Python

    d_out = join(data_path, 'tifs.nc')
    survey.gs.to_netcdf(d_out)








.. GENERATED FROM PYTHON SOURCE LINES 118-121

Option 1: 

Pass a variable name to export just that variable

.. GENERATED FROM PYTHON SOURCE LINES 121-123

.. code-block:: Python

    survey['derived_products']["all_maps"].gs.to_tif('magnetic_tmi')





.. rst-class:: sphx-glr-script-out

 .. code-block:: none

    magnetic_tmi.tif




.. GENERATED FROM PYTHON SOURCE LINES 124-129

Option 2: 

Export all the variables by NOT passing any variable names,
but need to specify a slice dimension for the 3D resistivity variable.
Can optionally pass a directory path to export tiffs to.

.. GENERATED FROM PYTHON SOURCE LINES 129-131

.. code-block:: Python

    survey['derived_products']["all_maps"].gs.to_tif(slice_dim='z', out_dir=data_path)





.. rst-class:: sphx-glr-script-out

 .. code-block:: none

    ../data_files/tempest_aseg/resistivity_0.0.tif
    ../data_files/tempest_aseg/resistivity_5.0.tif
    ../data_files/tempest_aseg/resistivity_10.0.tif
    ../data_files/tempest_aseg/resistivity_15.0.tif
    ../data_files/tempest_aseg/resistivity_20.0.tif
    ../data_files/tempest_aseg/magnetic_tmi.tif




.. GENERATED FROM PYTHON SOURCE LINES 132-133

Reading back in the GS NetCDF file

.. GENERATED FROM PYTHON SOURCE LINES 133-135

.. code-block:: Python

    new_survey = gspy.open_datatree(d_out)['survey']








.. GENERATED FROM PYTHON SOURCE LINES 136-137

Plotting

.. GENERATED FROM PYTHON SOURCE LINES 139-141

Make a map-view plot of a specific data variable, using Xarray's plotter
In this case, we slice the 3-D resistivity variable along the depth dimension

.. GENERATED FROM PYTHON SOURCE LINES 141-144

.. code-block:: Python


    r_plot = new_survey['derived_products']["all_maps"]['resistivity'].plot(col='z', vmax=3, cmap='jet', robust=True)




.. image-sg:: /examples/Creating_GS_Files/images/sphx_glr_plot_tifs_001.png
   :alt: z = 0.0, z = 5.0, z = 10.0, z = 15.0, z = 20.0
   :srcset: /examples/Creating_GS_Files/images/sphx_glr_plot_tifs_001.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 145-146

Make a map-view plot comparing the different x-y discretization of the two magnetic variables, using Xarray's plotter

.. GENERATED FROM PYTHON SOURCE LINES 146-153

.. code-block:: Python

    plt.figure()
    ax=plt.gca()
    plot1 = new_survey['derived_products']["all_maps"]['magnetic_tmi'].plot(ax=ax, cmap='jet', vmin=0, vmax=1000, robust=True)
    plot2 = new_survey['derived_products']["mag_map"]['magnetic_tmi'].plot(ax=ax, cmap='Greys', cbar_kwargs={'label': ''}, robust=True)
    plt.ylim([1.20556e6, 1.21476e6])
    plt.xlim([3.5201e5, 3.6396e5])
    plt.show()



.. image-sg:: /examples/Creating_GS_Files/images/sphx_glr_plot_tifs_002.png
   :alt: spatial_ref = 0.0
   :srcset: /examples/Creating_GS_Files/images/sphx_glr_plot_tifs_002.png
   :class: sphx-glr-single-img






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

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


.. _sphx_glr_download_examples_Creating_GS_Files_plot_tifs.py:

.. only:: html

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

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

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

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

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

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

      :download:`Download zipped: plot_tifs.zip <plot_tifs.zip>`


.. only:: html

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

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