.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/grid/ex_4d.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_auto_examples_grid_ex_4d.py: .. _example_4d_interpolation: 4D Interpolation ================ This example demonstrates how to perform 4D interpolation on a regular grid. The pyinterp library supports both quadrivariate and bicubic interpolation for 4D data. .. GENERATED FROM PYTHON SOURCE LINES 11-20 .. code-block:: Python import cartopy.crs import matplotlib import matplotlib.pyplot import numpy import pyinterp import pyinterp.backends.xarray import pyinterp.tests .. GENERATED FROM PYTHON SOURCE LINES 21-30 Quadrivariate Interpolation --------------------------- Quadrivariate interpolation extends trivariate interpolation to four dimensions. It performs bilinear interpolation on the 2D spatial plane (longitude and latitude) and then linear interpolation on the third and fourth dimensions. First, we load the 4D dataset and create the interpolator object. .. GENERATED FROM PYTHON SOURCE LINES 30-33 .. code-block:: Python ds = pyinterp.tests.load_grid4d() interpolator = pyinterp.backends.xarray.Grid4D(ds.pressure) .. GENERATED FROM PYTHON SOURCE LINES 34-42 Next, we define the coordinates for interpolation. .. warning:: When using a time axis, ensure that the date units are consistent between the grid and the interpolation coordinates. The :py:meth:`pyinterp.TemporalAxis.safe_cast` method can help manage date conversions and prevent inconsistencies. .. GENERATED FROM PYTHON SOURCE LINES 42-48 .. code-block:: Python mx, my, mz, mu = numpy.meshgrid(numpy.arange(-125, -70, 0.5), numpy.arange(25, 50, 0.5), numpy.datetime64('2000-01-01T12:00'), 0.5, indexing='ij') .. GENERATED FROM PYTHON SOURCE LINES 49-50 Now, we perform the quadrivariate interpolation. .. GENERATED FROM PYTHON SOURCE LINES 50-57 .. code-block:: Python quadrivariate = interpolator.quadrivariate({ 'longitude': mx.ravel(), 'latitude': my.ravel(), 'time': mz.ravel(), 'level': mu.ravel() }).reshape(mx.shape) .. GENERATED FROM PYTHON SOURCE LINES 58-69 Bicubic Interpolation on a 4D Grid ---------------------------------- For smoother results, you can use bicubic interpolation for the spatial dimensions, followed by linear interpolation on the other dimensions. .. note:: Bicubic interpolation requires that the grid axes are strictly increasing. If your latitudes are in descending order, you can set the `increasing_axes` parameter to ``True`` to automatically flip them. .. GENERATED FROM PYTHON SOURCE LINES 69-72 .. code-block:: Python interpolator = pyinterp.backends.xarray.Grid4D(ds.pressure, increasing_axes=True) .. GENERATED FROM PYTHON SOURCE LINES 73-74 We then perform the bicubic interpolation. .. GENERATED FROM PYTHON SOURCE LINES 74-84 .. code-block:: Python bicubic = interpolator.bicubic( { 'longitude': mx.ravel(), 'latitude': my.ravel(), 'time': mz.ravel(), 'level': mu.ravel() }, nx=2, ny=2).reshape(mx.shape) .. GENERATED FROM PYTHON SOURCE LINES 85-87 To visualize the results, we reshape the output arrays and extract the longitude and latitude coordinates. .. GENERATED FROM PYTHON SOURCE LINES 87-92 .. code-block:: Python quadrivariate = quadrivariate.squeeze(axis=(2, 3)) bicubic = bicubic.squeeze(axis=(2, 3)) lons = mx[:, 0].squeeze() lats = my[0, :].squeeze() .. GENERATED FROM PYTHON SOURCE LINES 93-101 Finally, let's plot the results of both quadrivariate and bicubic interpolation. .. note:: The resolution of the example grid is low (one pixel per degree), so the bicubic interpolation may not find enough pixels at the edges, resulting in undefined values. .. GENERATED FROM PYTHON SOURCE LINES 101-130 .. code-block:: Python fig = matplotlib.pyplot.figure(figsize=(10, 8)) fig.subplots_adjust(left=0.05, right=0.95, top=0.95, bottom=0.05, hspace=0.25) ax1 = fig.add_subplot( 211, projection=cartopy.crs.PlateCarree(central_longitude=180)) ax1.set_extent([lons.min(), lons.max(), lats.min(), lats.max()], crs=cartopy.crs.PlateCarree()) pcm = ax1.pcolormesh(lons, lats, quadrivariate.T, cmap='jet', shading='auto', transform=cartopy.crs.PlateCarree()) ax1.coastlines() ax1.set_title('Quadrivariate Interpolation') ax2 = fig.add_subplot( 212, projection=cartopy.crs.PlateCarree(central_longitude=180)) ax2.set_extent([lons.min(), lons.max(), lats.min(), lats.max()], crs=cartopy.crs.PlateCarree()) pcm = ax2.pcolormesh(lons, lats, bicubic.T, cmap='jet', shading='auto', transform=cartopy.crs.PlateCarree()) ax2.coastlines() ax2.set_title('Bicubic Interpolation') fig.colorbar(pcm, ax=[ax1, ax2], shrink=0.8) .. image-sg:: /auto_examples/grid/images/sphx_glr_ex_4d_001.png :alt: Quadrivariate Interpolation, Bicubic Interpolation :srcset: /auto_examples/grid/images/sphx_glr_ex_4d_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.953 seconds) .. _sphx_glr_download_auto_examples_grid_ex_4d.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/CNES/pangeo-pyinterp/master?urlpath=lab/tree/notebooks/auto_examples/grid/ex_4d.ipynb :alt: Launch binder :width: 150 px .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: ex_4d.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: ex_4d.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: ex_4d.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_