.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/grid/ex_2d.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_2d.py: .. _example_2d_interpolation: 2D Interpolation ================ This example illustrates how to perform 2D interpolation of a variable on a regular grid. The pyinterp library provides several interpolation methods, and this guide will walk you through bivariate and bicubic interpolation. .. 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-33 Bivariate Interpolation ----------------------- Bivariate interpolation is a common method for estimating the value of a variable at a new point based on its surrounding grid points. In this section, we will perform bivariate interpolation using pyinterp. First, we load the data and create the interpolator object. The constructor will automatically detect the longitude and latitude axes. If it fails, you can specify them by setting the ``units`` attribute to ``degrees_east`` and ``degrees_north``. If your grid is not geodetic, set the ``geodetic`` parameter to ``False``. .. GENERATED FROM PYTHON SOURCE LINES 33-36 .. code-block:: Python ds = pyinterp.tests.load_grid2d() interpolator = pyinterp.backends.xarray.Grid2D(ds.mss) .. GENERATED FROM PYTHON SOURCE LINES 37-40 Next, we define the coordinates where we want to interpolate the grid. To avoid interpolating at the grid points themselves, we shift the coordinates slightly. .. GENERATED FROM PYTHON SOURCE LINES 40-44 .. code-block:: Python mx, my = numpy.meshgrid(numpy.arange(-180, 180, 1) + 1 / 3.0, numpy.arange(-89, 89, 1) + 1 / 3.0, indexing='ij') .. GENERATED FROM PYTHON SOURCE LINES 45-47 Now, we interpolate the grid to the new coordinates using the :py:meth:`bivariate ` method. .. GENERATED FROM PYTHON SOURCE LINES 47-52 .. code-block:: Python mss = interpolator.bivariate(coords={ 'lon': mx.ravel(), 'lat': my.ravel() }).reshape(mx.shape) .. GENERATED FROM PYTHON SOURCE LINES 53-55 To visualize the results, we can plot the original grid and the interpolated grid. .. GENERATED FROM PYTHON SOURCE LINES 55-85 .. 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)) lons, lats = numpy.meshgrid(ds.lon, ds.lat, indexing='ij') pcm = ax1.pcolormesh(lons, lats, ds.mss.T, cmap='jet', shading='auto', transform=cartopy.crs.PlateCarree(), vmin=-0.1, vmax=0.1) ax1.set_extent([-180, 180, -90, 90], crs=cartopy.crs.PlateCarree()) ax1.coastlines() ax1.set_title('Original MSS') ax2 = fig.add_subplot(212, projection=cartopy.crs.PlateCarree()) pcm = ax2.pcolormesh(mx, my, mss, cmap='jet', shading='auto', transform=cartopy.crs.PlateCarree(), vmin=-0.1, vmax=0.1) ax2.set_extent([-180, 180, -90, 90], crs=cartopy.crs.PlateCarree()) ax2.coastlines() ax2.set_title('Bilinear Interpolated MSS') fig.colorbar(pcm, ax=[ax1, ax2], shrink=0.8) .. image-sg:: /auto_examples/grid/images/sphx_glr_ex_2d_001.png :alt: Original MSS, Bilinear Interpolated MSS :srcset: /auto_examples/grid/images/sphx_glr_ex_2d_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 86-106 The bivariate method supports several interpolation techniques, including bilinear, nearest neighbor, and inverse distance weighting. Distance calculations are performed using the `Haversine formula `_. Bicubic Interpolation --------------------- Bicubic interpolation provides a smoother result compared to bivariate interpolation by considering a 4x4 neighborhood of grid points. .. warning:: When using this interpolator, pay attention to NaN values. If the calculation window contains even a single NaN, the result of the interpolation will also be NaN, due to NaN propagation in arithmetic operations. This means the masked region effectively grows during interpolation. To avoid this behavior, you should :doc:`pre-process ` the grid to replace or remove NaN values. .. GENERATED FROM PYTHON SOURCE LINES 108-109 The following code performs bicubic interpolation on the same grid. .. GENERATED FROM PYTHON SOURCE LINES 109-114 .. code-block:: Python mss_bicubic = interpolator.bicubic(coords={ 'lon': mx.ravel(), 'lat': my.ravel() }).reshape(mx.shape) .. GENERATED FROM PYTHON SOURCE LINES 115-117 Let's visualize the result of the bicubic interpolation and compare it with the original data. .. GENERATED FROM PYTHON SOURCE LINES 117-146 .. 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)) pcm = ax1.pcolormesh(lons, lats, ds.mss.T, cmap='jet', shading='auto', transform=cartopy.crs.PlateCarree(), vmin=-0.1, vmax=0.1) ax1.set_extent([-180, 180, -90, 90], crs=cartopy.crs.PlateCarree()) ax1.coastlines() ax1.set_title('Original MSS') ax2 = fig.add_subplot(212, projection=cartopy.crs.PlateCarree()) pcm = ax2.pcolormesh(mx, my, mss_bicubic, cmap='jet', shading='auto', transform=cartopy.crs.PlateCarree(), vmin=-0.1, vmax=0.1) ax2.set_extent([-180, 180, -90, 90], crs=cartopy.crs.PlateCarree()) ax2.coastlines() ax2.set_title('Bicubic Interpolated MSS') fig.colorbar(pcm, ax=[ax1, ax2], shrink=0.8) .. image-sg:: /auto_examples/grid/images/sphx_glr_ex_2d_002.png :alt: Original MSS, Bicubic Interpolated MSS :srcset: /auto_examples/grid/images/sphx_glr_ex_2d_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 147-150 The interpolation :py:meth:`bicubic ` function has more parameters to define the data frame used by the spline functions and how to process the edges of the regional grids: .. GENERATED FROM PYTHON SOURCE LINES 150-159 .. code-block:: Python mss = interpolator.bicubic( coords={ 'lon': mx.ravel(), 'lat': my.ravel() }, nx=3, ny=3, ).reshape(mx.shape) .. GENERATED FROM PYTHON SOURCE LINES 160-171 .. warning:: The grid provided must have strictly increasing axes to meet the specifications of the interpolation. When building the grid, specify the ``increasing_axes`` option to flip the decreasing axes and the grid automatically. For example: .. code:: python interpolator = pyinterp.backends.xarray.Grid2D( ds.mss, increasing_axes=True) .. GENERATED FROM PYTHON SOURCE LINES 171-198 .. 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)) pcm = ax1.pcolormesh(lons, lats, ds.mss.T, cmap='jet', shading='auto', transform=cartopy.crs.PlateCarree(), vmin=-0.1, vmax=0.1) ax1.set_extent([-180, 180, -90, 90], crs=cartopy.crs.PlateCarree()) ax1.coastlines() ax1.set_title('Original MSS') ax2 = fig.add_subplot(212, projection=cartopy.crs.PlateCarree()) pcm = ax2.pcolormesh(mx, my, mss, cmap='jet', shading='auto', transform=cartopy.crs.PlateCarree(), vmin=-0.1, vmax=0.1) ax2.set_extent([-180, 180, -90, 90], crs=cartopy.crs.PlateCarree()) ax2.coastlines() ax2.set_title('Bicubic Interpolated MSS') .. image-sg:: /auto_examples/grid/images/sphx_glr_ex_2d_003.png :alt: Original MSS, Bicubic Interpolated MSS :srcset: /auto_examples/grid/images/sphx_glr_ex_2d_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Text(0.5, 1.0, 'Bicubic Interpolated MSS') .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 7.155 seconds) .. _sphx_glr_download_auto_examples_grid_ex_2d.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_2d.ipynb :alt: Launch binder :width: 150 px .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: ex_2d.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: ex_2d.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: ex_2d.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_