radis.lbl.factory module

Contains the SpectrumFactory class, which is the core of the RADIS Line-by-Line module.

Examples

Calculate a CO Spectrum, fetching the lines from HITRAN

# This is how you get a spectrum (see spectrum.py for front-end functions
# that do just that)
sf = SpectrumFactory(2125, 2249.9,
                     molecule='CO',
                     isotope=1,
                     cutoff=1e-30,   # for faster calculations. See
                                     # `plot_linestrength_hist` for more details
                     **kwargs)
sf.fetch_databank()        # autodownload from HITRAN
s = sf.eq_spectrum(Tgas=300)
s.plot('abscoeff')          # opacity

# Here we get some extra informations:
s.plot('radiance', wunit='nm',
                   Iunit='µW/cm2/sr/nm',  # Iunit is arbitrary. Use whatever makes sense
                   show_points=True)  # show_points to have an idea of the resolution

See radis.test.lbl.test_factory for more examples

Routine Listing

(the following sections are marked on the slider bar in the Spyder IDE with the # XXX symbol)

PUBLIC METHODS

Most methods are written in inherited class with the following inheritance scheme:

DatabankLoader > BaseFactory > BroadenFactory > BandFactory > SpectrumFactory

Inheritance diagram of radis.lbl.factory.SpectrumFactory

Notes

Debugging:

Remember one can access all locals variable in a function either by running %debug after a crash in IPython. Or by adding the line globals().update(locals()) at any point

Performance:

Fast version: iterate over chunks of dataframe Note that we can’t use the full dataframe because it’s too big and takes too much memory See Performance for more details.

for Developers:

  • To implement new database formats, see the databases parsers in cdsd.py / hitran.py, and the partition function interpolators / calculators methods of SpectrumFactory: _build_partition_function_calculator() and _build_partition_function_interpolator()


class SpectrumFactory(wmin=None, wmax=None, wunit=<radis.misc.utils.Default object>, wavenum_min=None, wavenum_max=None, wavelength_min=None, wavelength_max=None, Tref=296, pressure=1.01325, mole_fraction=1, path_length=1, wstep=0.01, molecule=None, isotope='all', medium='air', broadening_max_width=10, pseudo_continuum_threshold=0, self_absorption=True, chunksize=None, optimization='simple', folding_thresh=1e-06, zero_padding=-1, broadening_method=<radis.misc.utils.Default object>, cutoff=0, verbose=True, warnings=True, save_memory=False, export_populations=None, export_lines=False, **kwargs)[source]

Bases: radis.lbl.bands.BandFactory

A class to put together all functions related to loading CDSD / HITRAN databases, calculating the broadenings, and summing over all the lines.

Parameters
  • wmin, wmax (float or Quantity) – a hybrid parameter which can stand for minimum (maximum) wavenumber or minimum (maximum) wavelength depending upon the unit accompanying it. If dimensionless, wunit is considered as the accompanying unit.

  • wunit (string) – the unit accompanying wmin and wmax. Can only be passed with wmin and wmax. Default is cm-1.

  • wavenum_min, wavenum_max (float(cm^-1) or Quantity) – minimum (maximum) wavenumber to be processed in \(cm^{-1}\). use astropy.units to specify arbitrary inverse-length units.

  • wavelength_min, wavelength_max (float(nm) or Quantity) – minimum (maximum) wavelength to be processed in \(nm\). This wavelength can be in 'air' or 'vacuum' depending on the value of the parameter medium=. use astropy.units to specify arbitrary length units.

  • pressure (float(bar) or Quantity) – partial pressure of gas in bar. Default 1.01325 (1 atm). use astropy.units to specify arbitrary pressure units. For example, 1013.25 * u.mbar.

  • mole_fraction (N/D) – species mole fraction. Default 1. Note that the rest of the gas is considered to be air for collisional broadening.

  • path_length (float(cm) or Quantity) – path length in cm. Default 1. use astropy.units to specify arbitrary length units.

  • molecule (int, str, or None) – molecule id (HITRAN format) or name. If None, the molecule can be infered from the database files being loaded. See the list of supported molecules in MOLECULES_LIST_EQUILIBRIUM and MOLECULES_LIST_NONEQUILIBRIUM. Default None.

  • isotope (int, list, str of the form ‘1,2’, or ‘all’) – isotope id (sorted by relative density: (eg: 1: CO2-626, 2: CO2-636 for CO2). See HITRAN documentation for isotope list for all species. If ‘all’, all isotopes in database are used (this may result in larger computation times!). Default ‘all’

  • medium ('air', 'vacuum') – propagating medium when giving inputs with 'wavenum_min', 'wavenum_max'. Does not change anything when giving inputs in wavenumber. Default 'air'

Other Parameters
  • *Computation parameters (see :py:attr:`~radis.lbl.loader.DatabankLoader.params`)*

  • Tref (K) – Reference temperature for calculations (linestrength temperature correction). HITRAN database uses 296 Kelvin. Default 296 K

  • self_absorption (boolean) – Compute self absorption. If False, spectra are optically thin. Default True.

  • broadening_max_width (float (cm-1)) – Full width over which to compute the broadening. Large values will create a huge performance drop (scales as ~broadening_width^2 without DLM) The calculated spectral range is increased (by broadening_max_width/2 on each side) to take into account overlaps from out-of-range lines. Default 10 cm-1.

  • wstep (float (cm-1) or 'auto') – Resolution of wavenumber grid. Default 0.01 cm-1. If 'auto', it is ensured that there are slightly more or less than GRIDPOINTS_PER_LINEWIDTH_WARN_THRESHOLD points for each linewidth.

    Note

    wstep = ‘auto’ is optimized for performances while ensuring accuracy, but is still experimental in 0.9.30. Feedback welcome!

  • cutoff (float (~ unit of Linestrength: cm-1/(#.cm-2))) – discard linestrengths that are lower that this, to reduce calculation times. 1e-27 is what is generally used to generate databases such as CDSD. If 0, no cutoff. Default 1e-27.

  • pseudo_continuum_threshold (float) – if not 0, first calculate a rough approximation of the spectrum, then moves all lines whose linestrength intensity is less than this threshold of the maximum in a semi-continuum. Values above 0.01 can yield significant errors, mostly in highly populated areas. 80% of the lines can typically be moved in a continuum, resulting in 5 times faster spectra. If 0, no semi-continuum is used. Default 0.

  • save_memory (boolean) – if True, removes databases calculated by intermediate functions (for instance, delete the full database once the linestrength cutoff criteria was applied). This saves some memory but requires to reload the database & recalculate the linestrength for each new parameter. Default False.

  • export_populations ('vib', 'rovib', None) – if not None, store populations in Spectrum. Either store vibrational populations (‘vib’) or rovibrational populations (‘rovib’). Default None

  • export_lines (boolean) – if True, saves details of all calculated lines in Spectrum. This is necessary to later use line_survey(), but can take some space. Default False.

  • chunksize (int, or None) – Splits the lines database in several chuncks during calculation, else the multiplication of lines over all spectral range takes too much memory and slows the system down. Chunksize let you change the default chunck size. If None, all lines are processed directly. Usually faster but can create memory problems. Default None

  • optimization ("simple", "min-RMS", None) – If either "simple" or "min-RMS" DLM optimization for lineshape calculation is used: - "min-RMS" : weights optimized by analytical minimization of the RMS-error (See: [Spectral-Synthesis-Algorithm]) - "simple" : weights equal to their relative position in the grid

    If using the DLM optimization, broadening method is automatically set to 'fft'. If None, no lineshape interpolation is performed and the lineshape of all lines is calculated.

    Refer to [Spectral-Synthesis-Algorithm] for more explanation on the DLM method for lineshape interpolation.

    Default "min-RMS"

  • folding_thresh (float) – Folding is a correction procedure thet is applied when the lineshape is calculated with the fft broadening method and the linewidth is comparable to wstep, that prevents sinc(v) modulation of the lineshape. Folding continues until the lineshape intensity is below folding_threshold. Setting to 1 or higher effectively disables folding correction.

    Range: 0.0 < folding_thresh <= 1.0 Default: 1e-6

  • zero_padding (int) – Zero padding is used in conjunction with the fft broadening method to prevent circular convolution at the cost of performance. When set to -1, padding is set equal to the spectrum length, which guarantees a linear convolution.

    Range: 0 <= zero_padding <= len(w), or zero_padding = -1 Default: -1

  • broadening_method ("voigt", "convolve", "fft") – Calculates broadening with a direct voigt approximation (‘voigt’) or by convoluting independantly calculated Doppler and collisional broadening (‘convolve’). First is much faster, 2nd can be used to compare results. This SpectrumFactory parameter can be manually adjusted a posteriori with:

    sf = SpectrumFactory(...)
    sf.params.broadening_method = 'voigt'
    

    Fast fourier transform 'fft' is only available if using the DLM lineshape calculation optimization. Because the DLM convolves all lines at the same time, and thus operates on large arrays, 'fft' becomes more appropriate than convolutions in real space ('voit', 'convolve' )

    By default, use "fft" for any optimization, and "voigt" if optimization is None .

  • warnings (bool, or one of ['warn', 'error', 'ignore'], dict) – If one of ['warn', 'error', 'ignore'], set the default behaviour for all warnings. Can also be a dictionary to set specific warnings only. Example:

    warnings = {'MissingSelfBroadeningWarning':'ignore',
                'NegativeEnergiesWarning':'ignore',
                'HighTemperatureWarning':'ignore'}
    

    See default_warning_status for more information.

  • verbose (boolean, or int) – If False, stays quiet. If True, tells what is going on. If >=2, gives more detailed messages (for instance, details of calculation times). Default True.

Examples

An example using SpectrumFactory, load_databank(), the Spectrum methods, and units

from radis import SpectrumFactory
from astropy import units as u
sf = SpectrumFactory(wavelength_min=4165 * u.nm,
                     wavelength_max=4200 * u.nm,
                     isotope='1,2',
                     broadening_max_width=10,  # cm-1
                     medium='vacuum',
                     verbose=1,    # more for more details
                     )
sf.load_databank('HITRAN-CO2-TEST')        # predefined in ~/radis.json
s = sf.eq_spectrum(Tgas=300 * u.K, path_length=1 * u.cm)
s.rescale_path_length(0.01)    # cm
s.plot('radiance_noslit', Iunit='µW/cm2/sr/nm')

Refer to the online Examples for more cases.

See also

Alternative, calc_spectrum(), Main, load_databank(), eq_spectrum(), non_eq_spectrum(), For, fetch_databank(), init_databank(), init_database(), eq_bands(), non_eq_bands(), Inputs

input

physical input

params

computational parameters

misc

miscallenous parameters (don’t change output)

eq_spectrum(Tgas, mole_fraction=None, path_length=None, pressure=None, name=None)[source]

Generate a spectrum at equilibrium.

Parameters
  • Tgas (float or Quantity) – Gas temperature (K)

  • mole_fraction (float) – database species mole fraction. If None, Factory mole fraction is used.

  • path_length (float or Quantity) – slab size (cm). If None, the default Factory path_length is used.

  • pressure (float or Quantity) – pressure (bar). If None, the default Factory pressure is used.

  • name (str) – output Spectrum name (useful in batch)

Returns

s

Returns a Spectrum object

Use the get() method to get something among ['radiance', 'radiance_noslit', 'absorbance', etc...]

Or directly the plot() method to plot it. See [1]_ to get an overview of all Spectrum methods

Return type

Spectrum

References

1

RADIS doc: Spectrum how to?

eq_spectrum_gpu(Tgas, mole_fraction=None, path_length=None, pressure=None, name=None)[source]

Generate a spectrum at equilibrium with calculation of lineshapes and broadening done on the GPU.

Note

This method requires CUDA compatible hardware to execute. For more information on how to setup your system to run GPU-accelerated methods using CUDA and Cython, check GPU Spectrum Calculation on RADIS

Parameters
  • Tgas (float or Quantity) – Gas temperature (K)

  • mole_fraction (float) – database species mole fraction. If None, Factory mole fraction is used.

  • path_length (float or Quantity) – slab size (cm). If None, the default Factory path_length is used.

  • pressure (float or Quantity) – pressure (bar). If None, the default Factory pressure is used.

  • name (str) – output Spectrum name (useful in batch)

Returns

s

Returns a Spectrum object

Use the get() method to get something among ['radiance', 'radiance_noslit', 'absorbance', etc...]

Or directly the plot() method to plot it. See [1]_ to get an overview of all Spectrum methods

Return type

Spectrum

See also

eq_spectrum()

non_eq_spectrum(Tvib, Trot, Ttrans=None, mole_fraction=None, path_length=None, pressure=None, vib_distribution='boltzmann', rot_distribution='boltzmann', overpopulation=None, name=None)[source]

Calculate emission spectrum in non-equilibrium case. Calculates absorption with broadened linestrength and emission with broadened Einstein coefficient.

Parameters
  • Tvib (float) – vibrational temperature [K] can be a tuple of float for the special case of more-than-diatomic molecules (e.g: CO2)

  • Trot (float) – rotational temperature [K]

  • Ttrans (float) – translational temperature [K]. If None, translational temperature is taken as rotational temperature (valid at 1 atm for times above ~ 2ns which is the RT characteristic time)

  • mole_fraction (float) – database species mole fraction. If None, Factory mole fraction is used.

  • path_length (float or Quantity) – slab size (cm). If None, the default Factory path_length is used.

  • pressure (float or Quantity) – pressure (bar). If None, the default Factory pressure is used.

Other Parameters
  • vib_distribution ('boltzmann', 'treanor') – vibrational distribution

  • rot_distribution ('boltzmann') – rotational distribution

  • overpopulation (dict, or None) –

    add overpopulation factors for given levels:

    {level:overpopulation_factor}
    
  • name (str) – output Spectrum name (useful in batch)

Returns

s

Returns a Spectrum object

Use the get() method to get something among ['radiance', 'radiance_noslit', 'absorbance', etc...]

Or directly the plot() method to plot it. See [1]_ to get an overview of all Spectrum methods

Return type

Spectrum

References

1

RADIS doc: Spectrum how to?

optically_thin_power(Tgas=None, Tvib=None, Trot=None, Ttrans=None, vib_distribution='boltzmann', rot_distribution='boltzmann', mole_fraction=None, path_length=None, unit='mW/cm2/sr')[source]

Calculate total power emitted in equilibrium or non-equilibrium case in the optically thin approximation: it sums all emission integral over the total spectral range.

Warning

this is a fast implementation that doesnt take into account the contribution of lines outside the given spectral range. It is valid for spectral ranges surrounded by no lines, and spectral ranges much broaded than the typical line broadening (~ 1-10 cm-1 in the infrared)

If what you’re looking for is an accurate simulation on a narrow spectral range you better calculate the spectrum (that does take all of that into account) and integrate it with get_power()

Parameters
  • Tgas (float) – equilibrium temperature [K] If doing a non equilibrium case it should be None. Use Ttrans for translational temperature

  • Tvib (float) – vibrational temperature [K]

  • Trot (float) – rotational temperature [K]

  • Ttrans (float) – translational temperature [K]. If None, translational temperature is taken as rotational temperature (valid at 1 atm for times above ~ 2ns which is the RT characteristic time)

  • mole_fraction (float) – database species mole fraction. If None, Factory mole fraction is used.

  • path_length (float) – slab size (cm). If None, Factory mole fraction is used.

  • unit (str) – output unit. Default 'mW/cm2/sr'

Returns

float – see unit=.

Return type

Returns total power density in mW/cm2/sr (unless different unit is chosen),