radis.spectrum.spectrum module

Summary

Spectrum class holder

Contains both the emission and absorption features calculated by a SpectrumFactory, or generated by another spectral code or an experiment and imported in RADIS. Allows use post-processing functions such as applying an instrumental slit, saving to disk, plot all spectral quantities with arbitrary units.

See the The Spectrum object for the list of all post-processing functions, how to save a Spectrum, or how to generate a Spectrum from text.

Examples

Typical use:

from radis calculated_spectrum
s = calculated_spectrum(w, I, conditions={'case':'previously calculated by ##'})
s.plot('radiance_noslit')
s.apply_slit(0.5, shape='triangular')
s.plot('radiance')

Spectrum objects can be modified, stored, resampled, rescaled, or retrieved after they have been created, with store(), rescale_path_length(), rescale_mole_fraction(), resample(), store(), load_spec()

from radis import load_spec
s = load_spec('co_calculation.spec')
s.rescale_path_length(0.5)                  # calculate for new path_length
s.rescale_mole_fraction(0.02)   # calculate for new mole fraction
s.resample(w_new)               # resample on new wavespace
s.store('co_calculation2.spec')

More in The Spectrum object.


class Spectrum(quantities, units=None, conditions=None, cond_units=None, populations=None, lines=None, waveunit=None, name=None, warnings=True)[source]

Bases: object

This class holds results calculated with the SpectrumFactory calculation, with other radiative codes, or experimental data. It can be used to plot different quantities a posteriori, or manipulate output units (for instance convert a spectral radiance per wavelength units to a spectral radiance per wavenumber).

See more information on how to generate, edit or combine Spectrum objects on the Spectrum object guide.

Parameters
  • quantities (dict of tuples {‘quantity’:(wavenum, quantity)}) – where quantities are spectral quantities (absorbance, radiance, etc.) and wavenum is in \(cm^{-1}\) example:

    {'radiance_noslit':(wavenum, radiance_noslit),
    

    ‘absorbance’:(wavenum, absorbance)}

  • units (dict) – units for quantities

Other Parameters
  • conditions (dict) – physical conditions and calculation parameters

  • cond_units (dict) – units for conditions

  • populations (dict) – a dictionary of all species, and levels. Should be compatible with other radiative codes such as Specair output. Suggested format: {molecules: {isotopes: {elec state: rovib levels}}} e.g:

    {'CO2':{1: 'X': df}}   # with df a Pandas Dataframe
    
  • lines (pandas Dataframe) – all lines in databank (necessary for using line_survey()). Warning if you want to play with the lines content: The signification of columns in lines may be specific to a database format. Plus, some additional columns may have been added by the calculation (e.g: Ei and S for emission integral and linestrength in SpectrumFactory). Refer to the code to know what they mean (and their units)

  • wavespace ('nm', 'cm-1', 'nm_vac' or None) – wavelength in air ('nm'), wavenumber ('cm-1'), or wavelength in vacuum ('nm_vac'). Quantities should be evenly distributed along this space for fast convolution with the slit function If None, 'wavespace' must be defined in conditions. (non-uniform slit function is not implemented anyway… ) Defaults None (but raises an error if wavespace is not defined in conditions neither)

Other Parameters
  • name (str, or None) – Give a name to this Spectrum object (helps debugging in multislab configurations). Default None

  • warnings (boolean) – if True, test if inputs are valid, e.g, spectra are evenly distributed in wavelength, and raise a warning if not. Note that this take ~ 3.5 ms for a 20k points spectrum, when the rest of the creation process is only ~ 1.8ms (makes it 3 times longer, and can be a problem if hundreds of spectra are created in a row). Default True

Examples

Manipulate a Spectrum calculated by RADIS:

s = calc_spectrum(2125, 2300, Tgas=2000, databank='CDSD')
s.print_conditions()
s.plot('absorbance')
s.line_survey(overlay='absorbance')
s.plot('radiance_noslit', wunits='cm-1', Iunits='W/m2/sr/cm-1')
s.apply_slit(5)
s.plot('radiance')
w, t = s.get('transmittance_noslit')  # for use in multi-slabs configs

Any tuple of numpy arrays (w, I) can also be converted into a Spectrum object from the Spectrum class directly, or using the calculated_spectrum() function. All the following methods are equivalent:

from radis import Spectrum, calculated_spectrum
s1 = calculated_spectrum(w, I, wunit='nm', Iunit='mW/cm2/sr/nm')
s2 = Spectrum.from_array(w, I, 'radiance_noslit',
                       waveunit='nm', unit='mW/cm2/sr/nm')
s3 = Spectrum({'radiance_noslit': (w, I)},
              units={'radiance_noslit':'mW/cm2/sr/nm'},
              waveunit='nm')

See more examples in the [Spectrum] page.

Spectrum objects can be stored, retrieved, rescaled, resampled:

from radis import load_spec
s = load_spec('co_calculation.spec')
s.rescale_path_length(0.5)                  # calculate for new path_length
s.rescale_mole_fraction(0.02)   # calculate for new mole fraction
s.resample(w_new)               # resample on new wavespace
s.store('co_calculation2.spec')

Notes

Implementation:

quantities are stored in self._q and self._q_conv dictionaries. They are better accessed with the get() method that deals with units and wavespace

Wavebase:

quantites are stored either in wavenum or wavelength base, but this doesnt matter as they are retrieved / plotted with the get() and plot() methods which have units as input arguments

conditions[source]

Stores computation / measurement conditions

Type

dict

populations[source]

Stores molecules, isotopes, electronic states and vibrational or rovibrational populations

Type

dict

References

Spectrum

See the Spectrum object page

apply_slit()[source]

Apply an instrumental slit function to all quantities in Spectrum. Slit function can be generated with usual shapes (see shape=) or imported from an experimental slit function (path to a text file or numpy array of shape n*2). Convoluted spectra are cut on the edge compared to non-convoluted spectra, to remove side effects. See mode= to change this behaviour.

Warning with units: read about 'unit' and 'return_unit' parameters.

Parameters
  • slit_function (float or str or array) –

    If float:

    generate slit function with FWHM of slit function (in nm or cm-1 depending on unit=). A (top, base) tuple of (float, float) is required when asking for a trapezoidal slit function.

    If .txt:

    import experimental slit function from .txt file: format must be 2-columns with wavelengths and intensity (doesn’t have to be normalized)

    If array:

    format must be 2-columns with wavelengths and intensity (doesn’t have to be normalized)

  • unit ('nm' or 'cm-1') – unit of slit_function (FWHM, or imported file)

  • shape ('triangular', 'trapezoidal', 'gaussian', or any of SLIT_SHAPES) –

    which shape to use when generating a slit. Will call,

    respectively, triangular_slit(), trapezoidal_slit(), gaussian_slit(). Default ‘triangular’

  • center_wavespace (float, or None) – center of slit when generated (in unit). Not used if slit is imported.

  • norm_by ('area', 'max') – normalisation type: - 'area' normalizes the slit function to an area

    of 1. It conserves energy, and keeps the same units.

    • 'max' normalizes the slit function to a maximum of 1. The convoluted spectrum units change (they are multiplied by the spectrum waveunit, e.g: a radiance non convoluted in mW/cm2/sr/nm on a wavelength (nm). range will yield a convoluted radiance in mW/cm2/sr. Note that the slit is set to 1 in the Spectrum wavespace (i.e: a Spectrum calculated in cm-1 will have a slit set to 1 in cm-1).

    Default 'area'

  • mode ('valid', 'same') –

    'same' returns output of same length as initial spectra,

    but boundary effects are still visible. 'valid' returns output of length len(spectra) - len(slit) + 1, for which lines outside of the calculated range have no impact. Default 'valid'.

Other Parameters
  • auto_recenter_crop (bool) – if True, recenter slit and crop zeros on the side when importing an experimental slit. Default True. See recenter_slit(), crop_slit()

  • plot_slit (boolean) – if True, plot slit

  • store (boolean) – if True, store slit in the Spectrum object so it can be retrieved with get_slit() and plot with plot_slit(). Default True

  • slit_dispersion (func of (lambda, in 'nm'), or None) – spectrometer reciprocal function : dλ/dx(λ) (in nm) If not None, then the slit_dispersion function is used to correct the slit function for the whole range. Can be important if slit function was measured far from the measured spectrum (e.g: a slit function measured at 632.8 nm will look broader at 350 nm because the spectrometer dispersion is higher at 350 nm. Therefore it should be corrected) Default None

    Warning

    slit dispersion function is assumed to be given in nm if your spectrum is stored in cm-1 the wavenumbers are converted to wavelengths before being forwarded to the dispersion function

    See test_auto_correct_dispersion() for an example of the slit dispersion effect.

    A Python implementation of the slit dispersion:

    >>> def f(lbd):
    >>>    return  w/(2*f)*(tan(Φ)+sqrt((2*d/m/(w*1e-9)*cos(Φ))^2-1))
    

    Theoretical / References:

    >>> /dx ~ d/mf    # at first order
    >>> /dx = w/(2*f)*(tan(Φ)+sqrt((2*d/m/(w)*cos(Φ))^2-1))  # cf
    

    with:

    • Φ: spectrometer angle (°)

    • f: focal length (mm)

    • m: order of dispersion

    • d: grooves spacing (mm) = 1/gr with gr in (gr/mm)

    See Laux 1999 “Experimental study and modeling of infrared air plasma radiation” for more information

  • slit_dispersion_warning_threshold (float) – if not None, check that slit dispersion is about constant (< threshold change) on the calculated range. Default 0.01 (1%). See offset_dilate_slit_function()

  • *args, **kwargs – are forwarded to slit generation or import function

  • verbose (boolean) – print stuff

  • energy_threshold (float) – tolerance fraction when resampling. Default 1e-3 (0.1%) If areas before and after resampling differ by more than that an error is raised.

Notes

Units:

the slit function is first converted to the wavespace (wavelength/wavenumber) that the Spectrum is stored in, and applied to the spectral quantities in their native wavespace.

Implementation:

convolve_with_slit() is applied to all quantities in get_vars() that ends with _noslit. Generate a triangular instrumental slit function (or any other shape depending of shape=) with base slit_function_base (Uses the central wavelength of the spectrum for the slit function generation)

We deal with several special cases (which makes the code a little heavy, but the method very versatile):

  • when slit unit and spectrum unit arent the same

  • when spectrum is not evenly spaced

Examples

s.apply_slit(1.2, 'nm')

To manually apply the slit to a particular quantity use:

wavenum, quantity = s['quantity']
s['convolved_quantity'] = convolve_slit(wavenum, quantity,
    slit_function_base)

See convolve_with_slit() for more details on Units and Normalization

The slit is made considering the “center wavelength” which is the mean wavelength of the full spectrum you are applying it to.

Examples using apply_slit :

compare_with(other, spectra_only=False, plot=True, wunit='default', verbose=True, rtol=1e-05, ignore_nan=False, ignore_outliers=False, normalize=False, **kwargs)[source]

Compare Spectrum with another Spectrum object.

Parameters
  • other (type Spectrum) – another Spectrum to compare with

  • spectra_only (boolean, or str) – if True, only compares spectral quantities (in the same waveunit) and not lines or conditions. If str, compare a particular quantity name. If False, compare everything (including lines and conditions and populations). Default False

  • plot (boolean) – if True, use plot_diff to plot all quantities for the 2 spectra and the difference between them. Default True.

  • wunit ("nm", "cm-1", "default") – in which wavespace to compare (and plot). If "default", natural wavespace of first Spectrum is taken.

  • rtol (float) – relative difference to use for spectral quantities comparison

  • ignore_nan (boolean) – if True, nans are ignored when comparing spectral quantities

  • ignore_outliers (boolean, or float) –

    if not False, outliers are discarded. i.e, output is determined by:

    out = (~np.isclose(I, Ie, rtol=rtol, atol=0)).sum()/len(I) < ignore_outliers
    
  • normalize (bool) – Normalize the spectra to be ploted

Other Parameters

kwargs (dict) – arguments are forwarded to plot_diff()

Returns

equals – return True if spectra are equal (respective to tolerance defined by rtol and other input conditions)

Return type

boolean

Examples

Compare two Spectrum objects, or specifically the transmittance:

s1.compare_with(s2)
s1.compare_with(s2, 'transmittance')

Note that you can also simply use s1 == s2, that uses compare_with() internally:

s1 == s2       # will return True or False
cond_units[source]
conditions[source]
copy(copy_lines=True, quantity='all')[source]

Returns a copy of this Spectrum object (performs a smart deepcopy)

Parameters
  • copy_lines (bool) – default True

  • quantity (‘all’, or one of ‘radiance_noslit’, ‘absorbance’, etc.) – if not ‘all’, copy only one quantity. Default 'all'

crop(wmin=None, wmax=None, wunit='default', inplace=True)[source]

Crop spectrum to wmin-wmax range in wunit (inplace)

Parameters
  • wmin, wmax (float, or None) – boundaries of spectral range (in wunit)

  • wunit ('nm', 'cm-1', 'nm_vac') – which waveunit to use for wmin, wmax. If default: use the default Spectrum wavespace defined with get_waveunit().

Other Parameters

inplace (bool) – if True, modifies the Spectrum object directly. Else, returns a copy. Default True.

Returns

s – Cropped Spectrum. If inplace=True, Spectrum has been updated directly anyway.

Return type

Spectrum

Examples

Crop to experimental Spectrum, and compare:

from radis import calc_spectrum, load_spec, plot_diff
s = calc_spectrum(...)
s_exp = load_spec('typical_result.spec')
s.crop(s_exp.get_wavelength.min(), s_exp.get_wavelength.max(), 'nm')
plot_diff(s_exp, s)

See also

radis.spectrum.operations.crop()

MergeSlabs()

if used with ``resample=’full’,

out, this, to

file[source]
classmethod from_array(w, I, quantity, waveunit, unit, *args, **kwargs)[source]

Construct Spectrum from 2 arrays.

Parameters
  • w, I (array) – waverange and vector

  • quantity (str) – spectral quantity name

  • waveunit ('nm', 'cm-1', 'nm_vac') – unit of waverange: wavelength in air ('nm'), wavenumber ('cm-1'), or wavelength in vacuum ('nm_vac').

  • unit (str) – spectral quantity unit (arbitrary). Ex: ‘mW/cm2/sr/nm’ for radiance_noslit

  • *args, **kwargs – see Spectrum doc

Other Parameters
  • conditions (dict) – physical conditions and calculation parameters

  • cond_units (dict) – units for conditions

  • populations (dict) – a dictionary of all species, and levels. Should be compatible with other radiative codes such as Specair output. Suggested format: {molecules: {isotopes: {elec state: rovib levels}}} e.g:

    {'CO2':{1: 'X': df}}   # with df a Pandas Dataframe
    
  • lines (pandas Dataframe) – all lines in databank (necessary for using line_survey()). Warning if you want to play with the lines content: The signification of columns in lines may be specific to a database format. Plus, some additional columns may have been added by the calculation (e.g: Ei and S for emission integral and linestrength in SpectrumFactory). Refer to the code to know what they mean (and their units)

Returns

s – creates a Spectrum object

Return type

Spectrum

Examples

Create a spectrum:

from radis import Spectrum
s = Spectrum.from_array(w, I, 'radiance_noslit',
                       waveunit='nm', unit='mW/cm2/sr/nm')

To create a spectrum with absorption and emission components (e.g: radiance_noslit and transmittance_noslit, or emisscoeff and abscoeff) call the Spectrum class directly. Ex:

from radis import Spectrum
s = Spectrum({'abscoeff': (w, A), 'emisscoeff': (w, E)},
             units={'abscoeff': 'cm-1', 'emisscoeff':'W/cm2/sr/nm'},
             waveunit='nm')
classmethod from_txt(file, quantity, waveunit, unit, *args, **kwargs)[source]

Construct Spectrum from txt file.

Parameters
  • file (str) – file name

  • quantity (str) – spectral quantity name

  • waveunit ('nm', 'cm-1', 'nm_vac') – unit of waverange: wavelength in air ('nm'), wavenumber ('cm-1'), or wavelength in vacuum ('nm_vac').

  • unit (str) – spectral quantity unit

  • *args, **kwargs – the following inputs are forwarded to loadtxt: 'delimiter', 'skiprows' The rest if forwarded to Spectrum. see Spectrum doc

Other Parameters
  • delimiter (',', etc.) – see numpy.loadtxt()

  • skiprows (int) – see numpy.loadtxt()

  • argsort (bool) – sorts the arrays in file by wavespace. Convenient way to load a file where points have been manually added at the end. Default False.

  • *Optional Spectrum parameters*

  • conditions (dict) – physical conditions and calculation parameters

  • cond_units (dict) – units for conditions

  • populations (dict) – a dictionary of all species, and levels. Should be compatible with other radiative codes such as Specair output. Suggested format: {molecules: {isotopes: {elec state: rovib levels}}} e.g:

    {'CO2':{1: 'X': df}}   # with df a Pandas Dataframe
    
  • lines (pandas Dataframe) – all lines in databank (necessary for using line_survey()). Warning if you want to play with the lines content: The signification of columns in lines may be specific to a database format. Plus, some additional columns may have been added by the calculation (e.g: Ei and S for emission integral and linestrength in SpectrumFactory). Refer to the code to know what they mean (and their units)

Returns

s – creates a Spectrum object

Return type

Spectrum

Examples

Generate an experimental spectrum from txt. In that example the delimiter key is forwarded to loadtxt():

from radis import Spectrum
s = Spectrum.from_txt('spectrum.csv', 'radiance', waveunit='nm',
                          unit='W/cm2/sr/nm', delimiter=',')

To create a spectrum with absorption and emission components (e.g: radiance_noslit and transmittance_noslit, or emisscoeff and abscoeff) call the Spectrum class directly. Ex:

from radis import Spectrum
s = Spectrum({'abscoeff': (w, A), 'emisscoeff': (w, E)},
             units={'abscoeff': 'cm-1', 'emisscoeff':'W/cm2/sr/nm'},
             waveunit='nm')

Notes

Internally, the numpy loadtxt() function is used and transposed:

w, I = np.loadtxt(file).T

You can use 'delimiter' and ‘skiprows' as arguments.

get()[source]

Retrieve a spectral quantity from a Spectrum object. You can select wavespace unit, intensity unit, or propagation medium.

Parameters
  • var (variable (‘absorbance’, ‘transmittance’, etc.)) – Should be a defined quantity among CONVOLUTED_QUANTITIES or NON_CONVOLUTED_QUANTITIES. To get the full list of quantities defined in this Spectrum object use the get_vars() method.

  • wunit ('nm', 'cm', 'nm_vac'.) – wavespace unit: wavelength in air ('nm'), wavenumber ('cm-1'), or wavelength in vacuum ('nm_vac'). if "default", default unit for waveunit is used. See get_waveunit().

  • Iunit (unit for variable var) – if "default", default unit for quantity var is used. See the units attribute. For var="radiance", one can use per wavelength (~ ‘W/m2/sr/nm’) or per wavenumber (~ ‘W/m2/sr/cm-1’) units

Other Parameters

copy (boolean) – if True, returns a copy of the stored quantity (modifying it wont change the Spectrum object). Default True.

Returns

w, I – wavespace, quantity (ex: wavelength, radiance_noslit). For numpy users, note that these are copies (values) of the Spectrum quantity and not a view (reference): if you modify them the Spectrum is not changed

Return type

array-like

Examples

Get transmittance in cm-1:

w, I = s.get('transmittance_noslit', wunit='cm-1')

Get radiance (in wavelength in air):

_, R = s.get('radiance_noslit', wunit='nm', Iunit='W/cm2/sr/nm')
get_conditions()[source]

Get all physical / computational parameters.

See also

print_conditions(), the Spectrum page

get_integral(var, wunit='default', Iunit='default', **kwargs)[source]

Returns integral of variable ‘var’ over waverange.

Parameters
  • var (str) – spectral quantity to integate

  • wunit (str) – over which waverange to integrated. If default, use the default Spectrum wavespace defined with get_waveunit().

  • Iunit (str) – default 'default'

    Warning

    this is the unit of the quantity, not the unit of the integral. Don’t forget to multiply by wunit

Other Parameters

kwargs (**dict) – forwarded to get()

Returns

integral – integral in [Iunit]*[wunit]

Return type

float

get_name()[source]

Return Spectrum name.

If not defined, returns either the file name if Spectrum was loaded from a file, or the 'spectrum{id}' with the Python id object

get_populations(molecule=None, isotope=None, electronic_state=None)[source]

Return populations that are featured in the spectrum, either as upper or lower levels.

Parameters
  • molecule (str, or None) – if None, only one molecule must be defined. Else, an error is raised

  • isotope (int, or None) – isotope number. if None, only one isotope must be defined. Else, an error is raised

  • electronic_state (str) – if None, only one electronic state must be defined. Else, an error is raised

Returns

  • pandas dataframe of levels, where levels are the index,

  • and ‘Evib’ and ‘nvib’ are featured

Notes

Structure:

{molecule: {isotope: {electronic_state: {‘vib’: pandas Dataframe, # (copy of) vib levels

‘rovib’: pandas Dataframe, # (copy of) rovib levels ‘Ia’: float # isotopic abundance }}}}

(If Spectrum generated with RADIS, structure should match that of SpectrumFactory.get_populations())

get_power(unit='mW/cm2/sr')[source]

Returns integrated radiance (no slit) power density.

Parameters

Iunit (str) – power unit.

Returns

P – radiated power in unit

Return type

float

get_quantities(which='any')[source]

Returns all spectral quantities stored in this object (convoluted or non convoluted). Wrapper to get_vars()

Parameters

which (‘any’, ‘convoluted’, ‘non convoluted’)

get_radiance(Iunit='mW/cm2/sr/nm', copy=True)[source]

Return radiance in whatever unit, and can even convert from ~1/nm to ~1/cm-1 (and the other way round)

Other Parameters

copy (boolean) – if True, returns a copy of the stored waverange (modifying it wont change the Spectrum object). Default True.

get_radiance_noslit(Iunit='mW/cm2/sr/nm', copy=True)[source]

Return radiance (non convoluted) in whatever unit, and can even convert from ~1/nm to ~1/cm-1 (and the other way round)

Other Parameters

copy (boolean) – if True, returns a copy of the stored waverange (modifying it wont change the Spectrum object). Default True.

get_rovib_levels(molecule=None, isotope=None, electronic_state=None, first=None)[source]

Return rovibrational levels calculated in the spectrum (energies, populations)

Parameters
  • molecule (str, or None) – if None, only one molecule must be defined. Else, an error is raised

  • isotope (int, or None) – isotope number. if None, only one isotope must be defined. Else, an error is raised

  • electronic_state (str) – if None, only one electronic state must be defined. Else, an error is raised

  • first (int, or ‘all’ or None) – only show the first N levels. If None or ‘all’, all levels are shown

Returns

out – pandas dataframe of levels, where levels are the index, and ‘Evib’ and ‘nvib’ are featured

Return type

pandas DataFrame

get_slit(unit='same')[source]

Get slit function that was applied to the Spectrum.

Returns

wslit, Islit – slit function with wslit in Spectrum waveunit. See get_waveunit()

Return type

array

get_vars(which='any')[source]

Returns all spectral quantities stored in this object (convoluted or non convoluted)

Parameters

which (‘any’, ‘convoluted’, ‘non convoluted’)

get_vib_levels(molecule=None, isotope=None, electronic_state=None, first=None)[source]

Return vibrational levels in the spectrum (energies, populations)

Parameters
  • molecule (str, or None) – if None, only one molecule must be defined. Else, an error is raised

  • isotope (int, or None) – isotope number. if None, only one isotope must be defined. Else, an error is raised

  • electronic_state (str) – if None, only one electronic state must be defined. Else, an error is raised

  • first (int, or ‘all’ or None) – only show the first N levels. If None or ‘all’, all levels are shown

Returns

out – pandas dataframe of levels, where levels are the index, and ‘Evib’ and ‘nvib’ are featured

Return type

pandas DataFrame

get_wavelength(medium='air', which='any', copy=True)[source]

Return wavelength in defined medium.

Parameters
  • which (‘convoluted’, ‘non_convoluted’, ‘any’) – return wavelength for convoluted quantities, non convoluted quantities, or any. If any and both are defined, they have to be the same else an error is raised. Default any.

  • medium ('air', 'vacuum') – returns wavelength as seen in air, or vacuum. Default 'air'. See vacuum2air(), air2vacuum()

Other Parameters

copy (boolean) – if True, returns a copy of the stored waverange (modifying it wont change the Spectrum object). Default True.

Returns

w – (a copy of) spectrum wavelength for convoluted or non convoluted quantities

Return type

array_like

get_wavenumber(which='any', copy=True)[source]

Return wavenumber (if the same for all quantities)

Parameters

which (‘convoluted’, ‘non_convoluted’, ‘any’) – return wavenumber for convoluted quantities, non convoluted quantities, or any. If any and both are defined, they have to be the same else an error is raised. Default any.

Other Parameters

copy (boolean) – if True, returns a copy of the stored waverange (modifying it wont change the Spectrum object). Default True.

Returns

w – (a copy of) spectrum wavenumber for convoluted or non convoluted quantities

Return type

array_like

get_waveunit()[source]

Returns whether this spectrum is defined in wavelength (nm) or wavenumber (cm-1)

is_at_equilibrium(check='warn', verbose=False)[source]

Returns whether this spectrum is at (thermal) equilibrium. Reads the thermal_equilibrium key in Spectrum conditions. It does not imply chemical equilibrium (mole fractions are still arbitrary)

If they are defined, also check that the following assertions are True:

Tvib = Trot = Tgas self_absorption = True overpopulation = None

If they are not, still trust the value in Spectrum conditions, but raise a warning.

Other Parameters
  • check ('warn', 'error', 'ignore') – what to do if Spectrum conditions dont match the given equilibrium state: raise a warning, raise an error, or just ignore and dont even check. Default 'warn'.

  • verbose (bool) – if True, print why is the spectrum is not at equilibrium, if applicable.

is_optically_thin()[source]

Returns whether the spectrum is optically thin, based on the value on the self_absorption key in conditions.

If not given, raises an error

line_survey(overlay=None, wunit='default', writefile=None, cutoff=None, *args, **kwargs)[source]

Plot Line Survey (all linestrengths used for calculation) Output in Plotly (html)

Parameters
  • spec (Spectrum) – result from SpectrumFactory calculation (see spectrum.py)

  • overlay (‘absorbance’, ‘transmittance’, ‘radiance’, etc… or list of the above, or None) – overlay Linestrength with specified variable calculated in spec. Get the full list with the get_vars() method. Default None.

  • wunit ('default', 'nm', 'cm-1', 'nm_vac',) – wavelength air, wavenumber, or wavelength vacuum. If 'default', Spectrum get_waveunit() is used.

  • medium ({‘air’, ‘vacuum’, ‘default’}) – Choose whether wavelength are shown in air or vacuum. If 'default' lines are shown as stored in the spectrum.

Other Parameters
  • writefile (str) – if not None, a valid filename to save the plot under .html format. If None, use the fig object returned to show the plot.

  • kwargs:: dict – Other inputs are passed to LineSurvey(). Example below (see LineSurvey() documentation for more details):

  • Iunit (hitran, splot) – Linestrength output units:

    • hitran: (cm-1/(molecule/cm-2))

    • splot: (cm-1/atm) (Spectraplot units 2)

    Note: if not None, cutoff criteria is applied in this unit. Not used if plot is not ‘S’

  • barwidth (float) – With of bars in LineSurvey. Default 0.07

Returns

  • fig (a Plotly figure.) – If using a Jupyter Notebook, the plot will appear. Else, use writefile to export to an html file.

  • Plot in Plotly. See Output in [1]_

Examples

An example using the SpectrumFactory to generate a spectrum:

from radis import SpectrumFactory
sf = SpectrumFactory(
                     wavenum_min=2380,
                     wavenum_max=2400,
                     mole_fraction=400e-6,
                     path_length=100,  # cm
                     isotope=[1],
                     export_lines=True,    # required for LineSurvey!
                     db_use_cached=True)
sf.load_databank('HITRAN-CO2-TEST')
s = sf.eq_spectrum(Tgas=1500)
s.apply_slit(0.5)
s.line_survey(overlay='radiance_noslit', barwidth=0.01)

See the output in Examples

References

1

RADIS Online Documentation (LineSurvey)

2

SpectraPlot

lines[source]
max()[source]

Maximum of the Spectrum, if only one spectral quantity is available:

s.max()

Else, use Radiance(), Radiance_noslit(), Transmittance() or Transmittance_noslit()

Radiance(s).max()
min()[source]

Minimum of the Spectrum, if only one spectral quantity is available

s.min()

Else, use Radiance(), Radiance_noslit(), Transmittance() or Transmittance_noslit()

Radiance(s).min()
name[source]
normalize(normalize_how='max', wrange=(), wunit=None, inplace=False, force=False, verbose=True)[source]

Normalise the Spectrum, if only one spectral quantity is available.

Parameters
  • normalize_how ('max', 'area', 'mean') – how to normalize. 'max' is the default but may not be suited for very noisy experimental spectra. 'area' will normalize the integral to 1. 'mean' will normalize by the mean amplitude value

  • wrange (tuple) – if not empty, normalize on this range

  • wunit ("nm", "cm-1", "nm_vac") – unit of the normalisation range above. If None, use the spectrum default waveunit.

  • inplace (bool) – if True, changes the Spectrum.

Other Parameters

force (boolean) – By default, normalizing some parametres such as transmittance is forbidden because considered non-physical. Use force=True if you really want to.

Examples

s.normalize(“max”, (4200, 4800), inplace=True)

offset(offset, unit, inplace=True)[source]

Offset the spectrum by a wavelength or wavenumber (inplace)

Parameters
  • offset (float) – Constant to add to all quantities in the Spectrum.

  • unit (‘nm’ or ‘cm-1’) – unit for offset

Other Parameters

inplace (bool) – if True, modifies the Spectrum object directly. Else, returns a copy. Default True.

Returns

s – Offset Spectrum. If inplace=True, Spectrum has been updated directly anyway.

Return type

Spectrum

plot()[source]

Plot a Spectrum object.

Parameters
  • var (variable (absorbance, transmittance, transmittance_noslit, etc.)) – For full list see get_vars(). If None, plot the first thing in the Spectrum. Default None.

  • wunit ('default', 'nm', 'cm-1', 'nm_vac',) – wavelength air, wavenumber, or wavelength vacuum. If 'default', Spectrum get_waveunit() is used.

  • Iunit (unit for variable) – if default, default unit for quantity var is used. for radiance, one can use per wavelength (~ W/m2/sr/nm) or per wavenumber (~ W/m2/sr/cm-1) units

Other Parameters
  • show_points (boolean) – show calculated points. Default True.

  • nfig (int, None, or ‘same’) – plot on a particular figure. ‘same’ plots on current figure. For instance:

    s1.plot()
    s2.plot(nfig='same')
    
  • show_medium (bool, 'vacuum_only') – if True and wunit are wavelengths, plot the propagation medium in the xaxis label ([air] or [vacuum]). If 'vacuum_only', plot only if wunit=='nm_vac'. Default 'vacuum_only' (prevents from inadvertently plotting spectra with different propagation medium on the same graph).

  • yscale (‘linear’, ‘log’) – plot yscale

  • normalize (boolean, or tuple.) – option to normalize quantity to 1 (ex: for radiance). Default False

  • plot_by_parts (bool) – if True, look for discontinuities in the wavespace and plot the different parts without connecting lines. Useful for experimental spectra produced by overlapping step-and-glue. Additional parameters read from kwargs : split_threshold and cutwings. See more in split_and_plot_by_parts().

  • force (bool) – plotting on an existing figure is forbidden if labels are not the same. Use force=True to ignore that.

  • show_ruler (bool) – if True, add a ruler tool to the Matplotlib toolbar.

    Warning

    still experimental in 0.9.30 ! Try it, feedback welcome !

  • **kwargs (**dict) – kwargs forwarded as argument to plot (e.g: lineshape attributes: lw=3, color='r')

Returns

line plot

Return type

line

Examples

Plot an experimental_spectrum() in arbitrary units:

s = experimental_spectrum(..., Iunit='mW/cm2/sr/nm')
s.plot(Iunit='W/cm2/sr/cm-1')

See more examples in the plot Spectral quantities page.

plot_populations(what=None, nunit='', correct_for_abundance=False, **kwargs)[source]

Plots vib populations if given and format is valid.

Parameters
  • what (‘vib’, ‘rovib’, None) – if None plot everything

  • nunit (‘’, ‘cm-3’) – plot either in a fraction of vibrational levels, or a molecule number in in cm-3

  • correct_for_abundance (boolean) – if True, multiplies each population by the isotopic abundance (as it is done during the calculation of emission integral)

  • kwargs (**dict) – are forwarded to the plot

plot_slit(wunit=None)[source]

Plot slit function that was applied to the Spectrum.

If dispersion was used (see apply_slit()) the different slits are built again and plotted too (dotted).

Parameters

wunit ('nm', 'cm-1', or None) – plot slit in wavelength or wavenumber. If None, use the unit the slit in which the slit function was given. Default None

Returns

fix, ax – figure and ax

Return type

matplotlib objects

populations[source]
print_conditions(**kwargs)[source]

Prints all physical / computational parameters. You can also simply print the Spectrum object directly:

print(s)
Parameters

kwargs (dict) – refer to print_conditions()

resample(w_new, unit='same', out_of_bounds='nan', if_conflict_drop='error', energy_threshold=0.001, print_conservation=False, inplace=True, **kwargs)[source]

Resample spectrum over a new wavelength. Fills with transparent medium when out of bound (transmittance 1, radiance 0)

Warning

This may result in information loss. Resampling is done with oversampling and spline interpolation. These parameters can be adjusted, and energy conservation ensured with the appropriate parameters.

Uses the radis.misc.signal.resample() function.

Parameters
  • w_new (array, or Spectrum) – new wavespace to resample the spectrum on. Must be inclosed in the current wavespace (we won’t extrapolate) One can also give a Spectrum directly:

    s1.resample(s2.get_wavenumber())
    s1.resample(s2)            # also valid
    
  • unit ('same', 'nm', 'cm-1', 'nm_vac') – unit of new wavespace. It 'same' it is assumed to be the current waveunit. Default 'same'. The spectrum waveunit is changed to this unit after resampling (i.e: a spectrum calculated and stored in cm-1 but resampled in nm will be stored in nm from now on). If 'nm', wavelength in air. If 'nm_vac', wavelength in vacuum.

  • out_of_bounds (‘transparent’, ‘nan’, ‘error’) – what to do if resampling is out of bounds. ‘transparent’: fills with transparent medium. ‘nan’: fill with nan. ‘error’: raises an error. Default ‘nan’

  • if_conflict_drop (‘error’, ‘convoluted’, ‘non_convoluted’) – There is a problem if both convoluted and non convoluted (*no_slit) quantities coexists, as they aren’t scaled on the same wavespace grid. If ‘error’ an error is raised. If ‘convoluted’, convoluted quantities will be dropped. If ‘non_convoluted’ non convoluted quantities are dropped. Default ‘error’

  • medium (‘air’, ‘vacuum’, or ‘default’) – in which medium is the new waverange is calculated if it is given in ‘nm’. Ignored if unit=’cm-1’

Other Parameters
  • *Inputs forwarded to :func:`radis.misc.signal.resample`*

  • energy_threshold (float) – if energy conservation (integrals) is above this threshold, raise an error.

  • print_conservation (boolean) – if True, prints energy conservation. Default False.

  • inplace (boolean) – if True, modifies the Spectrum object directly. Else, returns a copy. Default True.

  • **kwargs (**dict) – all other arguments are sent to resample()

Returns

s – resampled Spectrum object. If using inplace=True, the Spectrum object has been modified anyway.

Return type

Spectrum

See also

resample()

rescale_mole_fraction(new_mole_fraction, old_mole_fraction=None, inplace=True, ignore_warnings=False, force=False, verbose=True)[source]

Update spectrum with new molar fraction Convoluted values (with slit) are dropped in the process.

Parameters
  • new_mole_fraction (float) – new mole fraction

  • old_mole_fraction (float, or None) – if None, current mole fraction (conditions[‘mole_fraction’]) is used

Other Parameters
  • inplace (boolean) – if True, modifies the Spectrum object directly. Else, returns a copy. Default True.

  • force (boolean) – if False, won’t allow rescaling to 0 (not to loose information). Default False

Returns

s – Cropped Spectrum. If inplace=True, Spectrum has been updated directly anyway.

Return type

Spectrum

Notes

Implementation:

similar to rescale_path_length() but we have to scale abscoeff & emisscoeff Note that this is valid only for small changes in mole fractions. Then, the change in line broadening becomes significant

rescale_path_length(new_path_length, old_path_length=None, inplace=True, force=False)[source]

Rescale spectrum to new path length. Starts from absorption coefficient and emission coefficient, and solves the RTE again for the new path length Convoluted values (with slit) are dropped in the process.

Parameters
  • new_path_length (float) – new path length

  • old_path_length (float, or None) – if None, current path length (conditions[‘path_length’]) is used

Other Parameters
  • inplace (boolean) – if True, modifies the Spectrum object directly. Else, returns a copy. Default True.

  • force (boolean) – if False, won’t allow rescaling to 0 (not to loose information). Default False

Returns

s – Cropped Spectrum. If inplace=True, Spectrum has been updated directly anyway.

Return type

Spectrum

Notes

Implementation:

To deal with all the input cases, we first make a list of what has to be recomputed, and what has to be recalculated

save(*args, **kwargs)[source]

Alias to Spectrum.store.

See Spectrum.store for documentation

savetxt(filename, var, wunit='default', Iunit='default')[source]

Export spectral quantity var to filename.

(note that by doing this you will loose additional information, such

as the calculation conditions or the units. You better save a Spectrum object under a .spec file with store() and load it afterwards with load_spec())

Parameters
  • filename (str) – file name

  • var (str) – which spectral variable ot export

Other Parameters

wunit, Iunit, medium (str) – see get() for more information

Notes

Export variable as:

np.savetxt(filename, np.vstack(self.get(var, wunit=wunit, Iunit=Iunit,
                                    medium=medium)).T, header=header)
store(path, discard=['lines', 'populations'], compress=True, add_info=None, add_date=None, if_exists_then='error', verbose=True)[source]

Save a Spectrum object in JSON format. Object can be recovered with load_spec(). If many Spectrum are saved in a same folder you can view their properties with the SpecDatabase structure.

Parameters
  • path (path to folder (database) or file) – if a folder, file is saved to database and name is generated automatically. if a file name, then Spectrum is saved to this file and the later formatting options dont apply

  • file (str) – explicitely give a filename to save

  • compress (boolean) – if False, save under text format, readable with any editor. if True, saves under binary format. Faster and takes less space. If 2, removes all quantities that can be regenerated with s.update(), e.g, transmittance if abscoeff and path length are given, radiance if emisscoeff and abscoeff are given in non-optically thin case, etc. Default True.

  • add_info (list) – append these parameters and their values if they are in conditions example:

    add_info = ['Tvib', 'Trot']
    
  • discard (list of str) – parameters to exclude, for instance to save some memory for instance Default [lines, populations]: retrieved Spectrum looses the line_survey() ability, and plot_populations() (but it saves tons of memory!)

  • if_exists_then ('increment', 'replace', 'error', 'ignore') – what to do if file already exists. If increment an incremental digit is added. If replace file is replaced (yeah). If 'ignore' no file is created. If error (or anything else) an error is raised. Default error

Returns

Return type

Returns filename used

Notes

If many spectra are stored in a folder, it may be time to set up a SpecDatabase structure to easily see all Spectrum conditions and get Spectrum that suits specific parameters.

Implementation:

Shouldnt rely on a Database. One may just want to store/load a Spectrum once.

Examples

Store a spectrum in compressed mode, regenerate quantities after loading:

from radis import load_spec
s.store('test.spec', compress=True)   # s is a Spectrum
s2 = load_spec('test.spec')
s2.update()                           # regenerate missing quantities
take(var)[source]
Parameters

var (str) – spectral quantity

Returns

s – same Spectrum with only the var spectral quantity

Return type

Spectrum

Examples

Use it to chain other commands

s.take('radiance').normalize().plot()
units[source]

units for spectral quantities.

Type

dict

update(quantity='all', optically_thin='default', verbose=True)[source]

Calculate missing quantities: ex: if path_length and emisscoeff are given, recalculate radiance_noslit.

Parameters
  • spec (Spectrum)

  • quantity (str) – name of the spectral quantity to recompute. If ‘same’, only the quantities in the Spectrum are recomputed. If ‘all’, then all quantities that can be derived are recomputed. Default ‘all’.

  • optically_thin (True, False, or ‘default’) – determines whether to calculate radiance with or without self absorption. If ‘default’, the value is determined from the self_absorption key in Spectrum.conditions. If not given, False is taken. Default ‘default’ Also updates the self_absorption value in conditions (creates it if doesnt exist)

is_spectrum(a)[source]

Returns whether a is a Spectrum object.

Parameters

a (anything) – a Python object

Returns

bool

Return type

True if a is a Spectrum object

Notes

is_spectrum compares the object class name (str): in some cases the Spectrum class gets imported twice (when databases are involved, mostly), and a purely isinstance() comparison fails