.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples/plot_mercator_argo.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_examples_plot_mercator_argo.py: Compare Mercator to ARGO ======================== In this notebook, we show: * how to easily rename variables and coordinates in an ARGO and a Mercator datasets, * how to plot a T-S diagram, * how to stack geographical coordinates, * how to compute a geographical distance between two datasets with coordinates, * how to interpolate a vertical profile with different methods. .. GENERATED FROM PYTHON SOURCE LINES 18-20 Initialisations --------------- .. GENERATED FROM PYTHON SOURCE LINES 20-38 .. code-block:: Python import xarray as xr import matplotlib as mpl import matplotlib.pyplot as plt import cartopy.crs as ccrs import gsw import cmocean import xoa import xoa.grid as xgrid import xoa.regrid as xregrid import xoa.meta as xmeta import xoa.coords as xcoords import xoa.geo as xgeo import xoa.plot as xplot mpl.rc("axes", grid=True) xr.set_options(display_style="text") .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 39-40 Register the :ref:`xoa ` accessors : .. GENERATED FROM PYTHON SOURCE LINES 40-42 .. code-block:: Python xoa.register_accessors() .. GENERATED FROM PYTHON SOURCE LINES 43-44 Register the Mercator and ARGO naming specifications .. GENERATED FROM PYTHON SOURCE LINES 44-49 .. code-block:: Python mercator_cfg = xoa.get_meta_config_file("mercator") argo_cfg = xoa.get_meta_config_file("argo") xmeta.register_meta_specs(mercator_cfg) xmeta.register_meta_specs(argo_cfg) .. GENERATED FROM PYTHON SOURCE LINES 50-51 Here is what these CF specifications contain .. GENERATED FROM PYTHON SOURCE LINES 51-56 .. code-block:: Python with open(mercator_cfg) as f: print(f.read()) with open(argo_cfg) as f: print(f.read()) .. rst-class:: sphx-glr-script-out .. code-block:: none [register] name=mercator [[attrs]] source=IBI-MFC* [data_vars] [[temp]] name=thetao [[sal]] name=so [[u]] name=uo [[v]] name=vo [[ssh]] name=zos [register] name=argo [[attrs]] DATA_ID=ARGO source=drifting subsurface profiling float [sglocator] name_format=False [data_vars] [[temp]] name=TEMP search_order=n [[temp_adjusted]] name=TEMP_ADJUSTED search_order=n [[sal]] name=PSAL search_order=n [[sal_adjusted]] name=PSAL_ADJUSTED search_order=n [[pres]] name=PRES search_order=n [[pres_adjusted]] search_order=n name=PRES_ADJUSTED [coords] [[lon]] name=LONGITUDE [[lat]] name=LATITUDE [[level]] name=N_LEVELS [[time]] name=TIME .. GENERATED FROM PYTHON SOURCE LINES 57-59 Read and decode the datasets ---------------------------- .. GENERATED FROM PYTHON SOURCE LINES 61-73 ARGO profiles ~~~~~~~~~~~~~ The ARGO data were downloaded with the `argopy `_ package with a piece of code close to the following:: from argopy import DataFetcher adf = DataFetcher().float(7900573) ds = adf.to_xarray().argo.point2profile().isel(N_PROF=slice(-10, None)) ds.to_netcdf("argo-7900573.nc") Here is the web page of this float: https://fleetmonitoring.euro-argo.eu/float/7900573 .. GENERATED FROM PYTHON SOURCE LINES 73-75 .. code-block:: Python ds_argo = xoa.open_data_sample("OBS/ARGO/argo-7900573.nc") .. GENERATED FROM PYTHON SOURCE LINES 76-78 We decode the names to make them more generic thanks to the CF specifications and select only the variables of interest. .. GENERATED FROM PYTHON SOURCE LINES 78-80 .. code-block:: Python ds_argo = ds_argo.xoa.decode()[["temp", "sal", "pres"]] .. GENERATED FROM PYTHON SOURCE LINES 81-83 We compute depths with the :mod:`gsw` package and assign them as coordinates of the ARGO dataset. .. GENERATED FROM PYTHON SOURCE LINES 83-91 .. code-block:: Python alat = ds_argo.lat.broadcast_like(ds_argo.pres) adepth = gsw.z_from_p(ds_argo.pres, alat) adepth = adepth.where(adepth.notnull(), adepth.min()) adepth.attrs.update({"long_name": "Depth", "units": "m"}) ds_argo = ds_argo.assign_coords(depth=adepth) ds_argo = ds_argo.isel(level=slice(None, None, -1)) ds_argo = ds_argo.set_index(N_PROF="time").rename(N_PROF="time") .. GENERATED FROM PYTHON SOURCE LINES 92-93 Quick plot of the salinity that highligths the mediterranean water. .. GENERATED FROM PYTHON SOURCE LINES 93-106 .. code-block:: Python fig, axs = plt.subplots(ncols=2, figsize=(14, 6)) ds_argo.sal.plot.contourf("time", "depth", cmap="cmo.haline", levels=20, ax=axs[0]) ds_argo.sal.plot.contour("time", "depth", levels=[35.65], linewidths=1, colors="k", ax=axs[0]) xplot.plot_ts( ds_argo.temp, ds_argo.sal, potential=False, scatter_c=ds_argo.depth, scatter_cmap="cmo.deep", axes=axs[1], ) axs[1].axvline(x=35.65, ls="--", color="k") .. image-sg:: /examples/images/sphx_glr_plot_mercator_argo_001.png :alt: plot mercator argo :srcset: /examples/images/sphx_glr_plot_mercator_argo_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none /home/docs/checkouts/readthedocs.org/user_builds/xoa/conda/develop/lib/python3.14/site-packages/xarray/plot/accessor.py:663: FutureWarning: Using positional arguments is deprecated for plot methods, use keyword arguments instead. return dataarray_plot.contourf(self._da, *args, **kwargs) /home/docs/checkouts/readthedocs.org/user_builds/xoa/conda/develop/lib/python3.14/site-packages/xarray/plot/accessor.py:542: FutureWarning: Using positional arguments is deprecated for plot methods, use keyword arguments instead. return dataarray_plot.contour(self._da, *args, **kwargs) .. GENERATED FROM PYTHON SOURCE LINES 107-109 Here we select the last ARGO profile as a reference profile since it better samples de mediterranean water vein. .. GENERATED FROM PYTHON SOURCE LINES 109-111 .. code-block:: Python ds_argo_prof = ds_argo.isel(time=-1).squeeze(drop=True) .. GENERATED FROM PYTHON SOURCE LINES 112-120 Mercator ~~~~~~~~ The Mercator data come from the IBI configuration and were downloaded from the CMEMS site: https://resources.marine.copernicus.eu/product-detail/IBI_ANALYSISFORECAST_PHY_005_001/INFORMATION The dataset has been undersampled by a factor 3 in longitude and latitude to limit disk usage. .. GENERATED FROM PYTHON SOURCE LINES 120-125 .. code-block:: Python ds_merc = xoa.open_data_sample("MODELS/CMEMS-IBI/ibi-argo-7900573.nc").xoa.decode() ds_merc = ds_merc.isel(depth=slice(None, None, -1)) ds_merc = ds_merc.assign_coords(depth=-ds_merc.depth) print(ds_merc) .. rst-class:: sphx-glr-script-out .. code-block:: none Size: 134kB Dimensions: (time: 2, depth: 40, lat: 8, lon: 13) Coordinates: * time (time) datetime64[ns] 16B 2022-02-06T12:00:00 2022-02-07T12:00:00 * depth (depth) float32 160B -1.942e+03 -1.684e+03 ... -1.541 -0.494 * lat (lat) float32 32B 43.61 43.69 43.78 43.86 43.94 44.03 44.11 44.19 * lon (lon) float32 52B -10.31 -10.22 -10.14 ... -9.472 -9.389 -9.306 Data variables: sal (time, depth, lat, lon) float32 33kB ... temp (time, depth, lat, lon) float32 33kB ... u (time, depth, lat, lon) float32 33kB ... v (time, depth, lat, lon) float32 33kB ... ssh (time, lat, lon) float32 832B ... Attributes: (12/25) Conventions: CF-1.0 source: IBI-MFC (PdE Production Center) institution: Puertos del Estado (PdE) title: Ocean 3D daily mean fields for the Ib... easting: longitude northing: latitude ... ... z_max: 5727.917f contact: mailto: servicedesk.cmems@mercator-oc... _CoordSysBuilder: ucar.nc2.dataset.conv.CF1Convention comment: history: Thu Feb 17 14:44:06 2022: ncks -d lon... NCO: netCDF Operators version 4.7.9 (Homep... .. GENERATED FROM PYTHON SOURCE LINES 126-130 By analogy with the ensemble data assimilation, this block of data can be viewed as a ensemble of profiles that emulates the uncertainties of the model, so we stack geographical coordinates in a `member` dimension thanks to :func:`~xoa.coords.geo_stack`. .. GENERATED FROM PYTHON SOURCE LINES 130-133 .. code-block:: Python ds_merc_ens_rect = xcoords.geo_stack(ds_merc, "member") print(ds_merc_ens_rect) .. rst-class:: sphx-glr-script-out .. code-block:: none Size: 135kB Dimensions: (time: 2, depth: 40, member: 104) Coordinates: * time (time) datetime64[ns] 16B 2022-02-06T12:00:00 2022-02-07T12:00:00 * depth (depth) float32 160B -1.942e+03 -1.684e+03 ... -1.541 -0.494 lat (member) float32 416B 43.61 43.61 43.61 43.61 ... 44.19 44.19 44.19 lon (member) float32 416B -10.31 -10.22 -10.14 ... -9.472 -9.389 -9.306 Dimensions without coordinates: member Data variables: sal (time, depth, member) float32 33kB 35.03 35.04 ... 35.56 35.58 temp (time, depth, member) float32 33kB 3.735 3.772 ... 13.12 13.27 u (time, depth, member) float32 33kB -0.072 -0.051 ... -0.08 -0.049 v (time, depth, member) float32 33kB 0.026 0.066 ... -0.021 0.013 ssh (time, member) float32 832B -0.391 -0.384 -0.379 ... -0.374 -0.375 Attributes: (12/25) Conventions: CF-1.0 source: IBI-MFC (PdE Production Center) institution: Puertos del Estado (PdE) title: Ocean 3D daily mean fields for the Ib... easting: longitude northing: latitude ... ... z_max: 5727.917f contact: mailto: servicedesk.cmems@mercator-oc... _CoordSysBuilder: ucar.nc2.dataset.conv.CF1Convention comment: history: Thu Feb 17 14:44:06 2022: ncks -d lon... NCO: netCDF Operators version 4.7.9 (Homep... .. GENERATED FROM PYTHON SOURCE LINES 134-139 This is a rectangular selection and we find it more appropriate to retain only profiles that are in a given radius of proximity, that is chosen close to the surface radius of deformation. We compute the distance from each model profile to the selected ARGO profile with :func:`~xoa.geo.get_distances`. .. GENERATED FROM PYTHON SOURCE LINES 139-149 .. code-block:: Python radius = 35e3 # 35 km dist_merc2argo = ( xgeo.get_distances(ds_merc_ens_rect, ds_argo_prof) .squeeze() .rename(npts0="member") .assign_coords(member=ds_merc_ens_rect.member) ) ds_merc_ens = ds_merc_ens_rect.where(dist_merc2argo < radius, drop=True) .. GENERATED FROM PYTHON SOURCE LINES 150-153 Compare Mercator and ARGO ------------------------- .. GENERATED FROM PYTHON SOURCE LINES 155-156 Let's plot the situation. .. GENERATED FROM PYTHON SOURCE LINES 156-168 .. code-block:: Python pmerc = ccrs.Mercator() pcar = ccrs.PlateCarree() fig, ax = plt.subplots(subplot_kw={"projection": pmerc}, figsize=(8, 8)) ax.set_extent(xgeo.get_extent(ds_merc, margin=1, square=True)) ax.gridlines(draw_labels=True) ax.add_wms("https://ows.emodnet-bathymetry.eu/wms", "emodnet:mean_atlas_land", alpha=0.5) kws = dict(c="C3", s=15, marker="s", transform=pcar) ax.scatter(ds_merc_ens_rect.lon, ds_merc_ens_rect.lat, label="Rejected", alpha=0.15, **kws) ax.scatter(ds_merc_ens.lon, ds_merc_ens.lat, label='Mercator', **kws) ax.scatter(ds_argo_prof.lon, ds_argo_prof.lat, s=100, c="C2", transform=pcar, label='ARGO') plt.legend() .. image-sg:: /examples/images/sphx_glr_plot_mercator_argo_002.png :alt: plot mercator argo :srcset: /examples/images/sphx_glr_plot_mercator_argo_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 169-171 We used the :func:`xoa.geo.get_extent` to compute a square geographical extent based on the coordinates of a dataset with an added margin. .. GENERATED FROM PYTHON SOURCE LINES 173-177 Now interpolate Mercator onto ARGO time and position for both the ensemble and the reference profile. All interpolations can be performed with the :meth:`~xarray.Dataset.interp` method since we are dealing with axis (1D) coordinates. .. GENERATED FROM PYTHON SOURCE LINES 177-181 .. code-block:: Python ds_merc_ens = ds_merc_ens.interp(time=ds_argo_prof.time.values) ds_merc_ens_std = ds_merc_ens.std("member") ds_merc_prof = ds_merc.interp(time=ds_argo_prof.time, lon=ds_argo_prof.lon, lat=ds_argo_prof.lat) .. GENERATED FROM PYTHON SOURCE LINES 182-183 We plot now the model, its simulated uncertainty and the observations. .. GENERATED FROM PYTHON SOURCE LINES 183-202 .. code-block:: Python plt.figure(figsize=(4, 7)) plt.fill_betweenx( ds_merc_prof.sal.depth, ds_merc_prof.sal - 1.96 * ds_merc_ens_std.sal, ds_merc_prof.sal + 1.96 * ds_merc_ens_std.sal, fc="C1", alpha=0.2, ) ds_merc_prof.sal.plot(y="depth", color="C1", label="Mercator") plt.plot(ds_merc_prof.sal.values, ds_merc_prof.depth, "o-", color="C1") ds_argo_prof.sal.plot(y="depth", label="ARGO", color="C0") plt.title("Uncertain profile") plt.legend() xplot.plot_double_minimap( ds_argo_prof, regional_ax="left", global_ax=(0.88, 0.88, 0.11), ) .. image-sg:: /examples/images/sphx_glr_plot_mercator_argo_003.png :alt: Uncertain profile :srcset: /examples/images/sphx_glr_plot_mercator_argo_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none (, ) .. GENERATED FROM PYTHON SOURCE LINES 203-206 Therefore, if we accept an uncertainty in the positioning of the ocean structure in the model, it differences with (assimilated) observations become acceptables. .. GENERATED FROM PYTHON SOURCE LINES 208-220 As you can see on the previous figures, the ARGO profiles come with a very high vertical resolution at all depths, whereas the model comes with a crude vertical resolution far from the surface. Each model thermodynamical model value is representative of the cube that defines its XYZ cell. This means that despite the model does simulate the small scale variability, the cell average may be accurate. So we propose here to perform comparisons in the model and the observational spaces. We use the function :func:`xoa.regrid.regrid1d` that supports the linear (:attr:`~xoa.regrid.regrid1d_methods.linear` method) interpolation to higher resolutions and the cell averaging (:attr:`~xoa.regrid.regrid1d_methods.cellave` method) to lower resolutions. .. GENERATED FROM PYTHON SOURCE LINES 220-225 .. code-block:: Python sal_merc_ens_o = xregrid.regrid1d(ds_merc_ens.sal, ds_argo_prof.depth, method="linear") sal_merc_prof_o = xregrid.regrid1d(ds_merc_prof.sal, ds_argo_prof.depth, method="linear") sal_argo_prof_m = xregrid.regrid1d(ds_argo_prof.sal, ds_merc_prof.depth, method="cellave") merc_prof_edges = xgrid.get_edges(ds_merc_prof.depth, "depth") .. GENERATED FROM PYTHON SOURCE LINES 226-227 Let's plot them. .. GENERATED FROM PYTHON SOURCE LINES 227-239 .. code-block:: Python fig, axes = plt.subplots(ncols=2, sharex=True, sharey=True, figsize=(10, 5)) sal_merc_prof_o.plot(y="depth", color="C1", label="Mercator", ax=axes[0]) ds_argo_prof.sal.plot(y="depth", label="ARGO", color="C0", ax=axes[0]) axes[0].set_title("On ARGO depths") axis = axes[0].axis() ds_merc_prof.sal.plot(y="depth", color="C1", label="Mercator", ax=axes[1]) plt.stairs(ds_merc_prof.sal, merc_prof_edges, orientation="horizontal", color="C1", lw=0.5) sal_argo_prof_m.plot(y="depth", label="ARGO", color="C0", ax=axes[1]) plt.stairs(sal_argo_prof_m, merc_prof_edges, orientation="horizontal", color="C0", lw=0.5) axes[1].set_title("On Mercator depths") axes[1].axis(axis) .. image-sg:: /examples/images/sphx_glr_plot_mercator_argo_004.png :alt: On ARGO depths, On Mercator depths :srcset: /examples/images/sphx_glr_plot_mercator_argo_004.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none (np.float64(34.96304988861084), np.float64(36.083950996398926), np.float64(-2072.132029431184), np.float64(92.43771385860644)) .. GENERATED FROM PYTHON SOURCE LINES 240-245 These plots show that the model is too haline, except at the mediterranean water depth where the model is very close to the observations once an appropriate vertical regridding is performed. However, the vertical resolution is not sufficient enough to simulate the maximum of salinity. .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 11.948 seconds) .. _sphx_glr_download_examples_plot_mercator_argo.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/shom-fr/xoa/master?urlpath=lab/tree/notebooks/examples/plot_mercator_argo.ipynb :alt: Launch binder :width: 150 px .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_mercator_argo.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_mercator_argo.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_mercator_argo.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_