.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/2_Experimental_spectra/plot_chain_spectrum_edition.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note You can download :ref:`below ` the full example code and run it with 🔬 `Radis-Lab `__, .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_2_Experimental_spectra_plot_chain_spectrum_edition.py: ============================================== Chain Editing and Lineshape Fitting a Spectrum ============================================== RADIS includes a powerful chaining syntax to edit (crop/offset/multiply/etc) calculated and experimental spectra. Some examples are given below. We also do some lineshape fitting using the :py:mod:`specutils` :py:func:`specutils.fitting.fit_lines` routine and some :py:mod:`astropy` models, among :py:class:`astropy.modeling.functional_models.Gaussian1D`, :py:class:`astropy.modeling.functional_models.Lorentz1D` or :py:class:`astropy.modeling.functional_models.Voigt1D` See more loading and post-processing functions on the :ref:`Spectrum page `. .. GENERATED FROM PYTHON SOURCE LINES 18-21 .. code-block:: Python from radis import Spectrum from radis.test.utils import getTestFile .. GENERATED FROM PYTHON SOURCE LINES 22-44 .. code-block:: Python s = Spectrum.from_mat( getTestFile("trimmed_1857_VoigtCO_Minesi.mat"), "absorbance", wunit="cm-1", unit="", index=10, ) # Plot default Spectrum: s.plot( nfig=1, # nfig=1 and show=False to plot on the same figure after show=False, # needed when using `inline` ploting (e.g. default in Spyder) ) # Now crop to the range we want to study s.crop(2010.66, 2010.80).plot( nfig=1, show=True, # show = True is default color="r", lw=2, # plot accepts matlplotlib args ) .. image-sg:: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_001.png :alt: plot chain spectrum edition :srcset: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none /home/docs/checkouts/readthedocs.org/user_builds/radis/checkouts/latest/radis/spectrum/spectrum.py:4921: UserWarning: Wavespace is not evenly spaced (0.000%) for absorbance. This may create problems if later convolving with slit function (`s.apply_slit()`). You can use `s.resample_even()` .. GENERATED FROM PYTHON SOURCE LINES 45-47 It seems the baseline is slightly offset (negative). Fix it with algebraic operations and show the difference .. GENERATED FROM PYTHON SOURCE LINES 47-51 .. code-block:: Python s.plot(nfig=2, show=False) s += 0.003 # could have used an Astropy unit here too ! s.plot(nfig=2, color="r") .. image-sg:: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_002.png :alt: plot chain spectrum edition :srcset: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 52-59 Some operations, such as crop, happen "in-place" by default, i.e. they modify the Spectrum. If you don't want to modify the Spectrum make sure you specify ``inplace=False`` For instance we compare the current Spectrum to a mock-up new experimental spectrum obtained by adding some random noise (15%) on top of the the previous one, with some offset and a coarser grid. Note the use of ``len(s)`` to get the number of spectral points .. GENERATED FROM PYTHON SOURCE LINES 59-78 .. code-block:: Python import matplotlib.pyplot as plt import numpy as np from radis.phys.units import Unit s.normalize(inplace=False).plot(lw=1.5, show=False) noise_array = 0.15 * ( np.random.rand(len(s)) * Unit("") ) # must be dimensioned to be multipled to spectra noise_offset = np.random.rand(1) * 0.01 s2 = (-1 * (s.normalize(inplace=False) + noise_array)).offset( noise_offset, "cm-1", inplace=False ) # Resample the spectrum on a coarser (1 every 3 points) grid s2.resample(s2.get_wavenumber()[::3], inplace=True, energy_threshold=None) plt.axhline(0) # draw horizontal line s2.plot(nfig="same", lw=1.5) .. image-sg:: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_003.png :alt: plot chain spectrum edition :srcset: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none /home/docs/checkouts/readthedocs.org/user_builds/radis/checkouts/latest/radis/misc/warning.py:427: UnevenWaverangeWarning: When resampling the spectrum, the new waverange had unequal spacing. .. GENERATED FROM PYTHON SOURCE LINES 79-81 Back to our original spectrum, we get the line positions using the :py:meth:`~radis.spectrum.spectrum.Spectrum.max` or :py:meth:`~radis.spectrum.spectrum.Spectrum.argmax` function .. GENERATED FROM PYTHON SOURCE LINES 81-84 .. code-block:: Python print("Maximum identified at ", s.argmax()) .. rst-class:: sphx-glr-script-out .. code-block:: none Maximum identified at 2010.7269449336343 1 / cm .. GENERATED FROM PYTHON SOURCE LINES 85-88 For a more accurate measurement of the line position, we fit a Lorentzian lineshape using Astropy & Specutil models (with straightforward Radis ==> Specutil conversion) and compare the fitted line center. .. GENERATED FROM PYTHON SOURCE LINES 88-111 .. code-block:: Python from astropy import units as u from astropy.modeling import models from specutils.fitting import fit_lines # Fit the spectrum and calculate the fitted flux values (``y_fit``) g_init = models.Lorentz1D( amplitude=s.max(), x_0=s.argmax(), fwhm=0.2 ) # see also models.Voigt1D g_fit = fit_lines(s.to_specutils(), g_init) w_fit = s.get_wavenumber() / u.cm y_fit = g_fit(w_fit) print("Fitted Line Center : ", g_fit.x_0.value) # Plot the original spectrum and the fitted. import matplotlib.pyplot as plt s.plot(lw=2, show=False) plt.plot(w_fit, y_fit, label="Fit result") plt.grid(True) plt.legend() .. image-sg:: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_004.png :alt: plot chain spectrum edition :srcset: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_004.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Fitted Line Center : 2010.7261919753003 .. GENERATED FROM PYTHON SOURCE LINES 112-113 We can also mesure the area under the line : .. GENERATED FROM PYTHON SOURCE LINES 113-122 .. code-block:: Python import numpy as np print("Absorbance of original line : ", s.get_integral("absorbance")) print( "Absorbance of fitted line :", np.trapz(y_fit, w_fit) ) # negative because w_fit in descending order # .. rst-class:: sphx-glr-script-out .. code-block:: none Absorbance of original line : 0.001629699288083762 Absorbance of fitted line : -0.0016606547425965664 1 / cm .. GENERATED FROM PYTHON SOURCE LINES 123-125 Finally note that the fitting routine can be achieved directly using the :py:meth:`~radis.spectrum.spectrum.Spectrum.fit_model` function: .. GENERATED FROM PYTHON SOURCE LINES 125-128 .. code-block:: Python from astropy.modeling import models s.fit_model(models.Lorentz1D(), plot=True) .. image-sg:: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_005.png :alt: plot chain spectrum edition :srcset: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_005.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none ([], [, , 0.00037279456823682285]) .. GENERATED FROM PYTHON SOURCE LINES 129-130 We see that fitting a Voigt profile yields substantially better results .. GENERATED FROM PYTHON SOURCE LINES 130-133 .. code-block:: Python from astropy.modeling import models s.fit_model(models.Voigt1D(), plot=True) .. image-sg:: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_006.png :alt: plot chain spectrum edition :srcset: /auto_examples/2_Experimental_spectra/images/sphx_glr_plot_chain_spectrum_edition_006.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none ([], [, , 0.000891958921859789, 0.0009052129810098868]) .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 5.905 seconds) .. _sphx_glr_download_auto_examples_2_Experimental_spectra_plot_chain_spectrum_edition.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_chain_spectrum_edition.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_chain_spectrum_edition.py `