.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/02_resampling_using_1d_interpolation/02_01_interpolate_1d_samples.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 JupyterLite or Binder .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_02_resampling_using_1d_interpolation_02_01_interpolate_1d_samples.py: Interpolate 1D samples ====================== Interpolate 1D samples with standard interpolation. 1. Assume that a user-provided 1D list of samples :math:`f[k]` has been obtained by sampling a spline on a unit grid. 2. From the samples, recover the continuously defined spline :math:`f(x)`. .. GENERATED FROM PYTHON SOURCE LINES 17-19 Imports ------- .. GENERATED FROM PYTHON SOURCE LINES 19-32 .. code-block:: Python import numpy as np import matplotlib.pyplot as plt from splineops.interpolate.tensorspline import TensorSpline plt.rcParams.update({ "font.size": 14, # Base font size "axes.titlesize": 18, # Title font size "axes.labelsize": 16, # Label font size "xtick.labelsize": 14, "ytick.labelsize": 14 }) .. GENERATED FROM PYTHON SOURCE LINES 33-41 Initial 1D Samples ------------------ We generate 1D samples and treat them as discrete signal points. Let :math:`\mathbf{f} = (f[0], f[1], f[2], \dots, f[K-1])` be a 1D array of data. These are the input samples that we are going to interpolate. .. GENERATED FROM PYTHON SOURCE LINES 41-72 .. code-block:: Python number_of_samples = 27 f_support = np.arange(number_of_samples) f_support_length = len(f_support) # It's equal to number_of_samples f_samples = np.array([ -0.657391, -0.641319, -0.613081, -0.518523, -0.453829, -0.385138, -0.270688, -0.179849, -0.11805, -0.0243016, 0.0130667, 0.0355389, 0.0901577, 0.219599, 0.374669, 0.384896, 0.301386, 0.128646, -0.00811776, 0.0153119, 0.106126, 0.21688, 0.347629, 0.419532, 0.50695, 0.544767, 0.555373 ]) plt.figure(figsize=(10, 4)) plt.title("f[k] samples") plt.stem(f_support, f_samples, basefmt=" ") # Add a black horizontal line at y=0: plt.axhline( y=0, color="black", linewidth=1, zorder=0 ) plt.xlabel("k") plt.ylabel("f[k]") plt.grid(True) plt.tight_layout() plt.show() .. image-sg:: /auto_examples/02_resampling_using_1d_interpolation/images/sphx_glr_02_01_interpolate_1d_samples_001.png :alt: f[k] samples :srcset: /auto_examples/02_resampling_using_1d_interpolation/images/sphx_glr_02_01_interpolate_1d_samples_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 73-89 Interpolate the Samples with a Spline ------------------------------------- We interpolate the 1D samples with a spline to obtain the continuously defined function .. math:: f(x) = \sum_{k\in{\mathbb{Z}}}\,c[k]\,\beta^{n}(x-k), where - the B-spline of degree :math:`n` is :math:`\beta^n`; - the spline coefficients :math:`c[k]` are determined from the input samples, such that :math:`f(k) = f[k]`. Let us now plot :math:`f`. .. GENERATED FROM PYTHON SOURCE LINES 89-97 .. code-block:: Python # Plot points plot_points_per_unit = 12 # Interpolated signal base = "bspline3" mode = "mirror" .. GENERATED FROM PYTHON SOURCE LINES 98-102 TensorSpline ~~~~~~~~~~~~ Here is one way to perform the standard interpolation. .. GENERATED FROM PYTHON SOURCE LINES 102-111 .. code-block:: Python f = TensorSpline(data=f_samples, coordinates=f_support, bases=base, modes=mode) f_coords = np.array([q / plot_points_per_unit for q in range(plot_points_per_unit * f_support_length)]) # Syntax hint: pass (plot_coords,) not plot_coords f_data = f(coordinates=(f_coords,), grid=False) .. GENERATED FROM PYTHON SOURCE LINES 112-116 Resize Method ~~~~~~~~~~~~~ The resize method with standard interpolation yields the same result. .. GENERATED FROM PYTHON SOURCE LINES 116-140 .. code-block:: Python from splineops.resize.resize import resize # We'll produce the same number of output samples as in f_coords desired_length = plot_points_per_unit * f_support_length # IMPORTANT: We explicitly define a coordinate array from 0..(f_support_length - 1) # with `desired_length` points. This matches the domain and size that the `resize` # function will produce below, ensuring the two outputs are sampled at the exact # same x-positions, and thus comparable point-by-point. f_coords_resize = np.linspace(0, f_support_length - 1, desired_length) f_data_resize = resize( data=f_samples, # 1D input output_size=(desired_length,), method="cubic" # ensures TensorSpline standard interpolation, not least-squares or oblique ) # Ensure both arrays have identical shapes f_data_spline = f(coordinates=(f_coords_resize,), grid=False) assert f_data_spline.shape == f_data_resize.shape, "Arrays must match in shape." mse_diff = np.mean((f_data_spline - f_data_resize)**2) print(f"MSE between TensorSpline result and resize result = {mse_diff:.6e}") .. rst-class:: sphx-glr-script-out .. code-block:: none MSE between TensorSpline result and resize result = 0.000000e+00 .. GENERATED FROM PYTHON SOURCE LINES 141-143 Plot of the Spline f ~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 143-162 .. code-block:: Python plt.figure(figsize=(10, 4)) plt.title("f[k] samples with interpolated f spline") plt.stem(f_support, f_samples, basefmt=" ", label="f[k] samples") # Add a black horizontal line at y=0: plt.axhline( y=0, color="black", linewidth=1, zorder=0 # draw behind other plot elements ) plt.plot(f_coords_resize, f_data_resize, color="green", linewidth=2, label="f spline") plt.xlabel("k") plt.ylabel("f") plt.legend() plt.grid(True) plt.tight_layout() plt.show() .. image-sg:: /auto_examples/02_resampling_using_1d_interpolation/images/sphx_glr_02_01_interpolate_1d_samples_002.png :alt: f[k] samples with interpolated f spline :srcset: /auto_examples/02_resampling_using_1d_interpolation/images/sphx_glr_02_01_interpolate_1d_samples_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.214 seconds) .. _sphx_glr_download_auto_examples_02_resampling_using_1d_interpolation_02_01_interpolate_1d_samples.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/splineops/splineops.github.io/main?urlpath=lab/tree/notebooks_binder/auto_examples/02_resampling_using_1d_interpolation/02_01_interpolate_1d_samples.ipynb :alt: Launch binder :width: 150 px .. container:: lite-badge .. image:: images/jupyterlite_badge_logo.svg :target: ../../lite/lab/index.html?path=auto_examples/02_resampling_using_1d_interpolation/02_01_interpolate_1d_samples.ipynb :alt: Launch JupyterLite :width: 150 px .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: 02_01_interpolate_1d_samples.ipynb <02_01_interpolate_1d_samples.ipynb>` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: 02_01_interpolate_1d_samples.py <02_01_interpolate_1d_samples.py>` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: 02_01_interpolate_1d_samples.zip <02_01_interpolate_1d_samples.zip>` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_