.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "python_module/12_external_data/pet_SST_collocation.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. or to run this example in your browser via Binder .. rst-class:: sphx-glr-example-title .. _sphx_glr_python_module_12_external_data_pet_SST_collocation.py: Collocating external data ========================= Script will use py-eddy-tracker methods to upload external data (sea surface temperature, SST) in a common structure with altimetry. Figures higlights the different steps. .. GENERATED FROM PYTHON SOURCE LINES 10-30 .. code-block:: Python from datetime import datetime from matplotlib import pyplot as plt from py_eddy_tracker import data from py_eddy_tracker.dataset.grid import RegularGridDataset date = datetime(2016, 7, 7) filename_alt = data.get_demo_path( f"dt_blacksea_allsat_phy_l4_{date:%Y%m%d}_20200801.nc" ) filename_sst = data.get_demo_path( f"{date:%Y%m%d}000000-GOS-L4_GHRSST-SSTfnd-OISST_HR_REP-BLK-v02.0-fv01.0.nc" ) var_name_sst = "analysed_sst" extent = [27, 42, 40.5, 47] .. GENERATED FROM PYTHON SOURCE LINES 31-33 Loading data ------------ .. GENERATED FROM PYTHON SOURCE LINES 33-42 .. code-block:: Python sst = RegularGridDataset(filename=filename_sst, x_name="lon", y_name="lat") alti = RegularGridDataset( data.get_demo_path(filename_alt), x_name="longitude", y_name="latitude" ) # We can use `Grid` tools to interpolate ADT on the sst grid sst.regrid(alti, "sla") sst.add_uv("sla") .. GENERATED FROM PYTHON SOURCE LINES 43-44 Functions to initiate figure axes .. GENERATED FROM PYTHON SOURCE LINES 44-61 .. code-block:: Python def start_axes(title, extent=extent): fig = plt.figure(figsize=(13, 6), dpi=120) ax = fig.add_axes([0.03, 0.05, 0.89, 0.91]) ax.set_xlim(extent[0], extent[1]) ax.set_ylim(extent[2], extent[3]) ax.set_title(title) ax.set_aspect("equal") return ax def update_axes(ax, mappable=None, unit=""): ax.grid() if mappable: cax = ax.figure.add_axes([0.93, 0.05, 0.01, 0.9], title=unit) plt.colorbar(mappable, cax=cax) .. GENERATED FROM PYTHON SOURCE LINES 62-64 ADT first display ----------------- .. GENERATED FROM PYTHON SOURCE LINES 64-68 .. code-block:: Python ax = start_axes("SLA", extent=extent) m = sst.display(ax, "sla", vmin=0.05, vmax=0.35) update_axes(ax, m, unit="[m]") .. image-sg:: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_001.png :alt: SLA, [m] :srcset: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 69-71 SST first display ----------------- .. GENERATED FROM PYTHON SOURCE LINES 73-74 We can now plot SST from `sst` .. GENERATED FROM PYTHON SOURCE LINES 74-78 .. code-block:: Python ax = start_axes("SST") m = sst.display(ax, "analysed_sst", vmin=295, vmax=300) update_axes(ax, m, unit="[°K]") .. image-sg:: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_002.png :alt: SST, [°K] :srcset: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 79-85 .. code-block:: Python ax = start_axes("SST") m = sst.display(ax, "analysed_sst", vmin=295, vmax=300) u, v = sst.grid("u").T, sst.grid("v").T ax.quiver(sst.x_c[::3], sst.y_c[::3], u[::3, ::3], v[::3, ::3], scale=10) update_axes(ax, m, unit="[°K]") .. image-sg:: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_003.png :alt: SST, [°K] :srcset: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 86-87 Now, with eddy contours, and displaying SST anomaly .. GENERATED FROM PYTHON SOURCE LINES 87-89 .. code-block:: Python sst.bessel_high_filter("analysed_sst", 400) .. GENERATED FROM PYTHON SOURCE LINES 90-91 Eddy detection .. GENERATED FROM PYTHON SOURCE LINES 91-98 .. code-block:: Python sst.bessel_high_filter("sla", 400) # ADT filtered ax = start_axes("SLA", extent=extent) m = sst.display(ax, "sla", vmin=-0.1, vmax=0.1) update_axes(ax, m, unit="[m]") a, c = sst.eddy_identification("sla", "u", "v", date, 0.002) .. image-sg:: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_004.png :alt: SLA, [m] :srcset: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_004.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none /home/docs/checkouts/readthedocs.org/user_builds/py-eddy-tracker/conda/latest/lib/python3.12/site-packages/py_eddy_tracker/dataset/grid.py:1915: RuntimeWarning: invalid value encountered in sqrt self._speed_ev = sqrt(u * u + v * v) /home/docs/checkouts/readthedocs.org/user_builds/py-eddy-tracker/conda/latest/lib/python3.12/site-packages/numpy/lib/function_base.py:4824: UserWarning: Warning: 'partition' will ignore the 'mask' of the MaskedArray. arr.partition( .. GENERATED FROM PYTHON SOURCE LINES 99-107 .. code-block:: Python kwargs_a = dict(lw=2, label="Anticyclonic", ref=-10, color="b") kwargs_c = dict(lw=2, label="Cyclonic", ref=-10, color="r") ax = start_axes("SST anomaly") m = sst.display(ax, "analysed_sst", vmin=-1, vmax=1) a.display(ax, **kwargs_a), c.display(ax, **kwargs_c) ax.legend() update_axes(ax, m, unit="[°K]") .. image-sg:: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_005.png :alt: SST anomaly, [°K] :srcset: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_005.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 108-111 Example of post-processing -------------------------- Get mean of sst anomaly_high in each internal contour .. GENERATED FROM PYTHON SOURCE LINES 111-114 .. code-block:: Python anom_a = a.interp_grid(sst, "analysed_sst", method="mean", intern=True) anom_c = c.interp_grid(sst, "analysed_sst", method="mean", intern=True) .. GENERATED FROM PYTHON SOURCE LINES 115-116 Are cyclonic (resp. anticyclonic) eddies generally associated with positive (resp. negative) SST anomaly ? .. GENERATED FROM PYTHON SOURCE LINES 116-127 .. code-block:: Python fig = plt.figure(figsize=(7, 5)) ax = fig.add_axes([0.05, 0.05, 0.90, 0.90]) ax.set_xlabel("SST anomaly") ax.set_xlim([-1, 1]) ax.set_title("Histograms of SST anomalies") ax.hist( anom_a, 5, alpha=0.5, color="b", label="Anticyclonic (mean:%s)" % (anom_a.mean()) ) ax.hist(anom_c, 5, alpha=0.5, color="r", label="Cyclonic (mean:%s)" % (anom_c.mean())) ax.legend() .. image-sg:: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_006.png :alt: Histograms of SST anomalies :srcset: /python_module/12_external_data/images/sphx_glr_pet_SST_collocation_006.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 128-129 Not clearly so in that case .. .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 4.530 seconds) .. _sphx_glr_download_python_module_12_external_data_pet_SST_collocation.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: binder-badge .. image:: images/binder_badge_logo.svg :target: https://mybinder.org/v2/gh/AntSimi/py-eddy-tracker/master?urlpath=lab/tree/notebooks/python_module/12_external_data/pet_SST_collocation.ipynb :alt: Launch binder :width: 150 px .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: pet_SST_collocation.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: pet_SST_collocation.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: pet_SST_collocation.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_