ARTS
2.4.0(git:4fb77825)
|
Line functions related to line shapes and line strength. More...
Classes | |
class | InternalData |
Functions | |
constexpr Index | ExpectedDataSize () |
Size required for data buffer. More... | |
void | set_lineshape (Eigen::Ref< Eigen::VectorXcd > F, const Eigen::Ref< const Eigen::VectorXd > f_grid, const Absorption::SingleLine &line, const Numeric &temperature, const Numeric &zeeman_df, const Numeric &magnetic_magnitude, const Numeric &doppler_constant, const LineShape::Output &lso, const LineShape::Type lineshape_type, const Absorption::MirroringType mirroring_type, const Absorption::NormalizationType norm_type) |
Sets the lineshape normalized to unity. More... | |
void | set_lorentz (Eigen::Ref< Eigen::VectorXcd > F, Eigen::Ref< Eigen::MatrixXcd > dF, Eigen::Ref< Eigen::Matrix< Complex, Eigen::Dynamic, ExpectedDataSize()>> data, const Eigen::Ref< const Eigen::VectorXd > f_grid, const Numeric &zeeman_df, const Numeric &magnetic_magnitude, const Numeric &F0_noshift, const LineShape::Output &lso, const AbsorptionLines &band=AbsorptionLines(), const Index &line_ind=0, const ArrayOfRetrievalQuantity &derivatives_data=ArrayOfRetrievalQuantity(), const ArrayOfIndex &derivatives_data_position=ArrayOfIndex(), const LineShape::Output &dT={0, 0, 0, 0, 0, 0, 0, 0, 0}, const LineShape::Output &dVMR={0, 0, 0, 0, 0, 0, 0, 0, 0}) |
Sets the Lorentz line shape. More... | |
void | set_htp (Eigen::Ref< Eigen::VectorXcd > F, Eigen::Ref< Eigen::MatrixXcd > dF, const Eigen::Ref< const Eigen::VectorXd > f_grid, const Numeric &zeeman_df, const Numeric &magnetic_magnitude, const Numeric &F0_noshift, const Numeric &GD_div_F0, const LineShape::Output &lso, const AbsorptionLines &band=AbsorptionLines(), const Index &line_ind=0, const ArrayOfRetrievalQuantity &derivatives_data=ArrayOfRetrievalQuantity(), const ArrayOfIndex &derivatives_data_position=ArrayOfIndex(), const Numeric &dGD_div_F0_dT=0.0, const LineShape::Output &dT={0, 0, 0, 0, 0, 0, 0, 0, 0}, const LineShape::Output &dVMR={0, 0, 0, 0, 0, 0, 0, 0, 0}) |
Sets the HTP line shape. More... | |
void | set_voigt (Eigen::Ref< Eigen::VectorXcd > F, Eigen::Ref< Eigen::MatrixXcd > dF, Eigen::Ref< Eigen::Matrix< Complex, Eigen::Dynamic, ExpectedDataSize()>> data, const Eigen::Ref< const Eigen::VectorXd > f_grid, const Numeric &zeeman_df, const Numeric &magnetic_magnitude, const Numeric &F0_noshift, const Numeric &GD_div_F0, const LineShape::Output &lso, const AbsorptionLines &band=AbsorptionLines(), const Index &line_ind=0, const ArrayOfRetrievalQuantity &derivatives_data=ArrayOfRetrievalQuantity(), const ArrayOfIndex &derivatives_data_position=ArrayOfIndex(), const Numeric &dGD_div_F0_dT=0.0, const LineShape::Output &dT={0, 0, 0, 0, 0, 0, 0, 0, 0}, const LineShape::Output &dVMR={0, 0, 0, 0, 0, 0, 0, 0, 0}) |
Sets the Voigt line shape. More... | |
void | set_doppler (Eigen::Ref< Eigen::VectorXcd > F, Eigen::Ref< Eigen::MatrixXcd > dF, Eigen::Ref< Eigen::Matrix< Complex, Eigen::Dynamic, ExpectedDataSize()>> data, const Eigen::Ref< const Eigen::VectorXd > f_grid, const Numeric &zeeman_df, const Numeric &magnetic_magnitude, const Numeric &F0_noshift, const Numeric &GD_div_F0, const AbsorptionLines &band=AbsorptionLines(), const Index &line_ind=0, const ArrayOfRetrievalQuantity &derivatives_data=ArrayOfRetrievalQuantity(), const ArrayOfIndex &derivatives_data_position=ArrayOfIndex(), const Numeric &dGD_div_F0_dT=0.0) |
Sets the Doppler line shape. More... | |
void | apply_linemixing_scaling_and_mirroring (Eigen::Ref< Eigen::VectorXcd > F, Eigen::Ref< Eigen::MatrixXcd > dF, const Eigen::Ref< Eigen::VectorXcd > Fm, const Eigen::Ref< Eigen::MatrixXcd > dFm, const LineShape::Output &lso, const bool with_mirroring, const AbsorptionLines &band=AbsorptionLines(), const Index &line_ind=0, const ArrayOfRetrievalQuantity &derivatives_data=ArrayOfRetrievalQuantity(), const ArrayOfIndex &derivatives_data_position=ArrayOfIndex(), const LineShape::Output &dT={0, 0, 0, 0, 0, 0, 0, 0, 0}, const LineShape::Output &dVMR={0, 0, 0, 0, 0, 0, 0, 0, 0}) |
Applies line mixing scaling to already set lineshape and line mirror. More... | |
void | apply_rosenkranz_quadratic_scaling (Eigen::Ref< Eigen::VectorXcd > F, Eigen::Ref< Eigen::MatrixXcd > dF, const Eigen::Ref< const Eigen::VectorXd > f_grid, const Numeric &F0, const Numeric &T, const AbsorptionLines &band=AbsorptionLines(), const Index &line_ind=0, const ArrayOfRetrievalQuantity &derivatives_data=ArrayOfRetrievalQuantity(), const ArrayOfIndex &derivatives_data_position=ArrayOfIndex()) |
Applies Rosenkranz quadratic normalization to already set line shape. More... | |
void | apply_VVH_scaling (Eigen::Ref< Eigen::VectorXcd > F, Eigen::Ref< Eigen::MatrixXcd > dF, Eigen::Ref< Eigen::Matrix< Complex, Eigen::Dynamic, ExpectedDataSize()>> data, const Eigen::Ref< const Eigen::VectorXd > f_grid, const Numeric &F0, const Numeric &T, const AbsorptionLines &band=AbsorptionLines(), const Index &line_ind=0, const ArrayOfRetrievalQuantity &derivatives_data=ArrayOfRetrievalQuantity(), const ArrayOfIndex &derivatives_data_position=ArrayOfIndex()) |
Applies Van Vleck and Huber normalization to already set line shape. More... | |
void | apply_VVW_scaling (Eigen::Ref< Eigen::VectorXcd > F, Eigen::Ref< Eigen::MatrixXcd > dF, const Eigen::Ref< const Eigen::VectorXd > f_grid, const Numeric &F0, const AbsorptionLines &band=AbsorptionLines(), const Index &line_ind=0, const ArrayOfRetrievalQuantity &derivatives_data=ArrayOfRetrievalQuantity(), const ArrayOfIndex &derivatives_data_position=ArrayOfIndex()) |
Applies Van Vleck and Weiskopf normalization to already set line shape. More... | |
Numeric | lte_linestrength (Numeric S0, Numeric E0, Numeric F0, Numeric QT0, Numeric T0, Numeric QT, Numeric T) |
Gets the local thermodynamic equilibrium line strength. More... | |
void | apply_linestrength_scaling_by_lte (Eigen::Ref< Eigen::VectorXcd > F, Eigen::Ref< Eigen::MatrixXcd > dF, Eigen::Ref< Eigen::VectorXcd > N, Eigen::Ref< Eigen::MatrixXcd > dN, const Absorption::SingleLine &line, const Numeric &T, const Numeric &T0, const Numeric &isotopic_ratio, const Numeric &QT, const Numeric &QT0, const AbsorptionLines &band=AbsorptionLines(), const Index &line_ind=0, const ArrayOfRetrievalQuantity &derivatives_data=ArrayOfRetrievalQuantity(), const ArrayOfIndex &derivatives_data_position=ArrayOfIndex(), const Numeric &dQT_dT=0.0) |
Applies linestrength to already set line shape by LTE population type. More... | |
void | apply_linestrength_scaling_by_vibrational_nlte (Eigen::Ref< Eigen::VectorXcd > F, Eigen::Ref< Eigen::MatrixXcd > dF, Eigen::Ref< Eigen::VectorXcd > N, Eigen::Ref< Eigen::MatrixXcd > dN, const Absorption::SingleLine &line, const Numeric &T, const Numeric &T0, const Numeric &Tu, const Numeric &Tl, const Numeric &Evu, const Numeric &Evl, const Numeric &isotopic_ratio, const Numeric &QT, const Numeric &QT0, const AbsorptionLines &band=AbsorptionLines(), const Index &line_ind=0, const ArrayOfRetrievalQuantity &derivatives_data=ArrayOfRetrievalQuantity(), const ArrayOfIndex &derivatives_data_position=ArrayOfIndex(), const Numeric &dQT_dT=0.0) |
Applies linestrength to already set line shape by vibrational level temperatures. More... | |
void | apply_lineshapemodel_jacobian_scaling (Eigen::Ref< Eigen::MatrixXcd > dF, const AbsorptionLines &band, const Index &line_ind, const ArrayOfRetrievalQuantity &derivatives_data, const ArrayOfIndex &derivatives_data_position, const Numeric &T, const Numeric &P, const Vector &vmrs) |
Applies the line-by-line pressure broadening jacobian for the matching lines. More... | |
Numeric | DopplerConstant (Numeric T, Numeric mass) |
Returns the frequency-independent part of the Doppler broadening. More... | |
Numeric | dDopplerConstant_dT (const Numeric &T, const Numeric &dc) |
Returns the temperature derivative of the frequency-independent part of the Doppler broadening. More... | |
void | find_cutoff_ranges (Index &start_cutoff, Index &nelem_cutoff, const Eigen::Ref< const Eigen::VectorXd > f_grid, const Numeric &fmin, const Numeric &fmax) |
Sets cutoff frequency indices. More... | |
void | apply_linestrength_from_nlte_level_distributions (Eigen::Ref< Eigen::VectorXcd > F, Eigen::Ref< Eigen::MatrixXcd > dF, Eigen::Ref< Eigen::VectorXcd > N, Eigen::Ref< Eigen::MatrixXcd > dN, const Numeric &r1, const Numeric &r2, const Numeric &g1, const Numeric &g2, const Numeric &A21, const Numeric &F0, const Numeric &T, const AbsorptionLines &band=AbsorptionLines(), const Index &line_ind=0, const ArrayOfRetrievalQuantity &derivatives_data=ArrayOfRetrievalQuantity(), const ArrayOfIndex &derivatives_data_position=ArrayOfIndex()) |
Applies non-lte linestrength to already set line shape. More... | |
void | set_cross_section_of_band (InternalData &scratch, InternalData &sum, const ConstVectorView f_grid, const AbsorptionLines &band, const ArrayOfRetrievalQuantity &derivatives_data, const ArrayOfIndex &derivatives_data_active, const Vector &vmrs, const EnergyLevelMap &nlte, const Numeric &P, const Numeric &T, const Numeric &isot_ratio, const Numeric &H, const Numeric &DC, const Numeric &dDCdT, const Numeric &QT, const Numeric &dQTdT, const Numeric &QT0, const bool no_negatives=false, const bool zeeman=false, const Zeeman::Polarization zeeman_polarization=Zeeman::Polarization::Pi) |
Computes the cross-section of an absorption band. More... | |
Line functions related to line shapes and line strength.
void Linefunctions::apply_linemixing_scaling_and_mirroring | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
Eigen::Ref< Eigen::MatrixXcd > | dF, | ||
const Eigen::Ref< Eigen::VectorXcd > | Fm, | ||
const Eigen::Ref< Eigen::MatrixXcd > | dFm, | ||
const LineShape::Output & | lso, | ||
const bool | with_mirroring, | ||
const AbsorptionLines & | band = AbsorptionLines() , |
||
const Index & | line_ind = 0 , |
||
const ArrayOfRetrievalQuantity & | derivatives_data = ArrayOfRetrievalQuantity() , |
||
const ArrayOfIndex & | derivatives_data_position = ArrayOfIndex() , |
||
const LineShape::Output & | dT = {0, 0, 0, 0, 0, 0, 0, 0, 0} , |
||
const LineShape::Output & | dVMR = {0, 0, 0, 0, 0, 0, 0, 0, 0} |
||
) |
Applies line mixing scaling to already set lineshape and line mirror.
Equation: with_mirroring: F := (1+G-iY) * F + (1+G+iY) * Fm else: F := (1+G-iY) * F
and appropriate derivatives
[in,out] | F | Lineshape. Must be right size |
[in,out] | dF | Lineshape derivative. Must be right size |
[in] | Fm | Mirrored lineshape. Must be right size |
[in] | dFm | Mirrored lineshape derivative. Must be right size |
[in] | lso | Line shape parameters |
[in] | with_mirroring | Mirror lineshape check |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
[in] | dT | Temperature derivatives of line shape parameters |
[in] | dVMR | VMR derivatives of line shape parameters |
Definition at line 420 of file linefunctions.cc.
References conj(), LineShape::Output::G, Array< base >::nelem(), and LineShape::Output::Y.
Referenced by set_cross_section_of_band().
void Linefunctions::apply_lineshapemodel_jacobian_scaling | ( | Eigen::Ref< Eigen::MatrixXcd > | dF, |
const AbsorptionLines & | band, | ||
const Index & | line_ind, | ||
const ArrayOfRetrievalQuantity & | derivatives_data, | ||
const ArrayOfIndex & | derivatives_data_position, | ||
const Numeric & | T, | ||
const Numeric & | P, | ||
const Vector & | vmrs | ||
) |
Applies the line-by-line pressure broadening jacobian for the matching lines.
[in,out] | dF | Lineshape derivative. Must be right size |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
[in] | T | Atmospheric temperature |
[in] | P | Atmospheric pressure |
[in] | vmrs | VMRs for line shape broadeners |
Definition at line 787 of file linefunctions.cc.
References Absorption::id_in_line(), is_lineshape_parameter(), Array< base >::nelem(), RetrievalQuantity::QuantumIdentity(), Absorption::Lines::ShapeParameter_dInternal(), and LineShape::vmrs().
Referenced by set_cross_section_of_band().
void Linefunctions::apply_linestrength_from_nlte_level_distributions | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
Eigen::Ref< Eigen::MatrixXcd > | dF, | ||
Eigen::Ref< Eigen::VectorXcd > | N, | ||
Eigen::Ref< Eigen::MatrixXcd > | dN, | ||
const Numeric & | r1, | ||
const Numeric & | r2, | ||
const Numeric & | g1, | ||
const Numeric & | g2, | ||
const Numeric & | A21, | ||
const Numeric & | F0, | ||
const Numeric & | T, | ||
const AbsorptionLines & | band = AbsorptionLines() , |
||
const Index & | line_ind = 0 , |
||
const ArrayOfRetrievalQuantity & | derivatives_data = ArrayOfRetrievalQuantity() , |
||
const ArrayOfIndex & | derivatives_data_position = ArrayOfIndex() |
||
) |
Applies non-lte linestrength to already set line shape.
Works on ratio-inputs, meaning that the total distribution does not have to be known
Cannot support partial derivatives at this point due to ARTS not possessing its own NLTE ratio calculation agenda
[in,out] | F | Lineshape. Must be right size |
[in,out] | dF | Lineshape derivative. Must be right size |
[in,out] | N | Source lineshape |
[in,out] | dN | Source lineshape derivative |
[in] | r1 | Ratio of molecules at energy level 1 |
[in] | r2 | Ratio of molecules at energy level 2 |
[in] | g1 | Statistical weight of energy level 1 |
[in] | g2 | Statistical weight of energy level 2 |
[in] | A21 | Einstein coefficient for the transition from energy level 2 to energy level 1 |
[in] | F0 | Central frequency without any shifts |
[in] | T | Atmospheric temperature |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
Definition at line 844 of file linefunctions.cc.
References F0, Array< base >::nelem(), Constant::pow2(), and ARTS::Var::x().
Referenced by set_cross_section_of_band().
void Linefunctions::apply_linestrength_scaling_by_lte | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
Eigen::Ref< Eigen::MatrixXcd > | dF, | ||
Eigen::Ref< Eigen::VectorXcd > | N, | ||
Eigen::Ref< Eigen::MatrixXcd > | dN, | ||
const Absorption::SingleLine & | line, | ||
const Numeric & | T, | ||
const Numeric & | T0, | ||
const Numeric & | isotopic_ratio, | ||
const Numeric & | QT, | ||
const Numeric & | QT0, | ||
const AbsorptionLines & | band = AbsorptionLines() , |
||
const Index & | line_ind = 0 , |
||
const ArrayOfRetrievalQuantity & | derivatives_data = ArrayOfRetrievalQuantity() , |
||
const ArrayOfIndex & | derivatives_data_position = ArrayOfIndex() , |
||
const Numeric & | dQT_dT = 0.0 |
||
) |
Applies linestrength to already set line shape by LTE population type.
[in,out] | F | Lineshape. Must be right size |
[in,out] | dF | Lineshape derivative. Must be right size |
[in,out] | N | Source lineshape |
[in,out] | dN | Source lineshape derivative |
[in] | line | The absorption line |
[in] | T | The atmospheric temperature |
[in] | T0 | The reference temperature |
[in] | isotopic_ratio | The ratio of the isotopologue in the atmosphere |
[in] | QT | Partition function at atmospheric temperature of level |
[in] | QT0 | Partition function at reference temperature |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
[in] | dQT_dT | Temperature derivative of QT |
Definition at line 632 of file linefunctions.cc.
Referenced by set_cross_section_of_band().
void Linefunctions::apply_linestrength_scaling_by_vibrational_nlte | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
Eigen::Ref< Eigen::MatrixXcd > | dF, | ||
Eigen::Ref< Eigen::VectorXcd > | N, | ||
Eigen::Ref< Eigen::MatrixXcd > | dN, | ||
const Absorption::SingleLine & | line, | ||
const Numeric & | T, | ||
const Numeric & | T0, | ||
const Numeric & | Tu, | ||
const Numeric & | Tl, | ||
const Numeric & | Evu, | ||
const Numeric & | Evl, | ||
const Numeric & | isotopic_ratio, | ||
const Numeric & | QT, | ||
const Numeric & | QT0, | ||
const AbsorptionLines & | band = AbsorptionLines() , |
||
const Index & | line_ind = 0 , |
||
const ArrayOfRetrievalQuantity & | derivatives_data = ArrayOfRetrievalQuantity() , |
||
const ArrayOfIndex & | derivatives_data_position = ArrayOfIndex() , |
||
const Numeric & | dQT_dT = 0.0 |
||
) |
Applies linestrength to already set line shape by vibrational level temperatures.
[in,out] | F | Lineshape. Must be right size |
[in,out] | dF | Lineshape derivative. Must be right size |
[in,out] | N | Source lineshape |
[in,out] | dN | Source lineshape derivative |
[in] | line | The absorption line |
[in] | T | The atmospheric temperature |
[in] | T0 | The reference temperature |
[in] | Tu | The upper state vibrational temperature; must be T if level is LTE |
[in] | Tl | The lower state vibrational temperature; must be T if level is LTE |
[in] | Evu | The upper state vibrational energy; if set funny, yields funny results |
[in] | Evl | The lower state vibrational energy; if set funny, yields funny results |
[in] | isotopic_ratio | The ratio of the isotopologue in the atmosphere |
[in] | QT | Partition function at atmospheric temperature of level |
[in] | QT0 | Partition function at reference temperature |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
[in] | dQT_dT | Temperature derivative of QT |
Definition at line 684 of file linefunctions.cc.
Referenced by set_cross_section_of_band().
void Linefunctions::apply_rosenkranz_quadratic_scaling | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
Eigen::Ref< Eigen::MatrixXcd > | dF, | ||
const Eigen::Ref< const Eigen::VectorXd > | f_grid, | ||
const Numeric & | F0, | ||
const Numeric & | T, | ||
const AbsorptionLines & | band = AbsorptionLines() , |
||
const Index & | line_ind = 0 , |
||
const ArrayOfRetrievalQuantity & | derivatives_data = ArrayOfRetrievalQuantity() , |
||
const ArrayOfIndex & | derivatives_data_position = ArrayOfIndex() |
||
) |
Applies Rosenkranz quadratic normalization to already set line shape.
[in,out] | F | Lineshape. Must be right size |
[in,out] | dF | Lineshape derivative. Must be right size |
[in] | f_grid | Frequency grid of computations |
[in] | F0 | Central frequency without any shifts |
[in] | T | Atmospheric temperature at level |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
Definition at line 481 of file linefunctions.cc.
References do_temperature_jacobian(), F0, ARTS::Var::f_grid(), Absorption::id_in_line(), is_frequency_parameter(), Array< base >::nelem(), and Absorption::Lines::QuantumIdentity().
Referenced by set_cross_section_of_band(), and set_lineshape().
void Linefunctions::apply_VVH_scaling | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
Eigen::Ref< Eigen::MatrixXcd > | dF, | ||
Eigen::Ref< Eigen::Matrix< Complex, Eigen::Dynamic, ExpectedDataSize()>> | data, | ||
const Eigen::Ref< const Eigen::VectorXd > | f_grid, | ||
const Numeric & | F0, | ||
const Numeric & | T, | ||
const AbsorptionLines & | band = AbsorptionLines() , |
||
const Index & | line_ind = 0 , |
||
const ArrayOfRetrievalQuantity & | derivatives_data = ArrayOfRetrievalQuantity() , |
||
const ArrayOfIndex & | derivatives_data_position = ArrayOfIndex() |
||
) |
Applies Van Vleck and Huber normalization to already set line shape.
[in,out] | F | Lineshape. Must be right size |
[in,out] | dF | Lineshape derivative. Must be right size |
[in,out] | data | Block of allocated memory. Output nonsensical |
[in] | f_grid | Frequency grid of computations |
[in] | F0 | Central frequency without any shifts |
[in] | T | Atmospheric temperature at level |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
Definition at line 531 of file linefunctions.cc.
References data, F0, ARTS::Var::f_grid(), Absorption::id_in_line(), is_frequency_parameter(), Array< base >::nelem(), and Absorption::Lines::QuantumIdentity().
Referenced by set_cross_section_of_band(), and set_lineshape().
void Linefunctions::apply_VVW_scaling | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
Eigen::Ref< Eigen::MatrixXcd > | dF, | ||
const Eigen::Ref< const Eigen::VectorXd > | f_grid, | ||
const Numeric & | F0, | ||
const AbsorptionLines & | band = AbsorptionLines() , |
||
const Index & | line_ind = 0 , |
||
const ArrayOfRetrievalQuantity & | derivatives_data = ArrayOfRetrievalQuantity() , |
||
const ArrayOfIndex & | derivatives_data_position = ArrayOfIndex() |
||
) |
Applies Van Vleck and Weiskopf normalization to already set line shape.
[in,out] | F | Lineshape. Must be right size |
[in,out] | dF | Lineshape derivative. Must be right size |
[in] | f_grid | Frequency grid of computations |
[in] | F0 | Central frequency without any shifts |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
Definition at line 580 of file linefunctions.cc.
References F0, ARTS::Var::f_grid(), fac(), Absorption::id_in_line(), is_frequency_parameter(), Array< base >::nelem(), and Absorption::Lines::QuantumIdentity().
Referenced by set_cross_section_of_band(), and set_lineshape().
Returns the temperature derivative of the frequency-independent part of the Doppler broadening.
[in] | T | Atmospheric temperature at level |
[in] | dc | Output of Linefunctions::DopplerConstant(T, mass) |
Definition at line 811 of file linefunctions.cc.
Referenced by xsec_species(), and zeeman_on_the_fly().
Returns the frequency-independent part of the Doppler broadening.
[in] | T | Atmospheric temperature at level |
[in] | mass | Mass of molecule under consideration |
Definition at line 807 of file linefunctions.cc.
References sqrt().
Referenced by Absorption::PredefinedModel::makarov2020_o2_lines_mpm(), xsec_species(), and zeeman_on_the_fly().
|
constexpr |
Size required for data buffer.
Definition at line 42 of file linefunctions.h.
Referenced by Linefunctions::InternalData::InternalData(), and set_lineshape().
void Linefunctions::find_cutoff_ranges | ( | Index & | start_cutoff, |
Index & | nelem_cutoff, | ||
const Eigen::Ref< const Eigen::VectorXd > | f_grid, | ||
const Numeric & | fmin, | ||
const Numeric & | fmax | ||
) |
Sets cutoff frequency indices.
[out] | start_cutoff | Start pos of cutoff frequency |
[out] | end_cutoff | End pos of cutoff frequency |
[in] | f_grid | Frequency grid of computations |
[in] | fmin | Minimum frequency |
[in] | fmax | Maximum frequency |
Definition at line 816 of file linefunctions.cc.
References ARTS::Var::f_grid().
Referenced by set_cross_section_of_band().
Numeric Linefunctions::lte_linestrength | ( | Numeric | S0, |
Numeric | E0, | ||
Numeric | F0, | ||
Numeric | QT0, | ||
Numeric | T0, | ||
Numeric | QT, | ||
Numeric | T | ||
) |
Gets the local thermodynamic equilibrium line strength.
[in] | S0 | Reference linestrength |
[in] | E0 | Reference line lower energy |
[in] | F0 | Central frequency without any shifts |
[in] | QT0 | Partition function at reference temperature |
[in] | T0 | Reference temperature |
[in] | QT | Partition function at atmospheric temperature |
[in] | T | Atmospheric temperature at level |
Definition at line 618 of file linefunctions.cc.
void Linefunctions::set_cross_section_of_band | ( | InternalData & | scratch, |
InternalData & | sum, | ||
const ConstVectorView | f_grid, | ||
const AbsorptionLines & | band, | ||
const ArrayOfRetrievalQuantity & | derivatives_data, | ||
const ArrayOfIndex & | derivatives_data_active, | ||
const Vector & | vmrs, | ||
const EnergyLevelMap & | nlte, | ||
const Numeric & | P, | ||
const Numeric & | T, | ||
const Numeric & | isot_ratio, | ||
const Numeric & | H, | ||
const Numeric & | DC, | ||
const Numeric & | dDCdT, | ||
const Numeric & | QT, | ||
const Numeric & | dQTdT, | ||
const Numeric & | QT0, | ||
const bool | no_negatives = false , |
||
const bool | zeeman = false , |
||
const Zeeman::Polarization | zeeman_polarization = Zeeman::Polarization::Pi |
||
) |
Computes the cross-section of an absorption band.
[in,out] | scratch | Data that is overwritten by every line |
[in,out] | sun | Data that is set to zero then added onto by every line |
[in] | f_grid | As WSV |
[in] | band | The absorption band |
[in] | derivatives_data | Derivatives |
[in] | derivatives_data_active | Derivatives that are active |
[in] | vmrs | The VMRs of this band's broadening species |
[in] | nlte | A map of NLTE energy levels |
[in] | P | The pressure |
[in] | T | The temperature |
[in] | isot_ratio | The band isotopic ratio |
[in] | H | The strength of the magnetic field |
[in] | DC | As per DopplerConstant |
[in] | dDCdT | Temperature derivative of DC |
[in] | QT | The partition function at the temperature |
[in] | dQTdT | Temperature derivative of QT |
[in] | QT0 | The partition function at the band reference temperature |
[in] | no_negatives | Check sum.F before output of any real negative values, and removes them if present |
[in] | zeeman | Attempts adding up the fine Zeeman lines |
[in] | zeeman_polarization | The polarization of Zeeman model (to know how many Zeeman lines there will be) |
Definition at line 1291 of file linefunctions.cc.
References Absorption::Lines::A(), apply_linemixing_scaling_and_mirroring(), apply_lineshapemodel_jacobian_scaling(), apply_linestrength_from_nlte_level_distributions(), apply_linestrength_scaling_by_lte(), apply_linestrength_scaling_by_vibrational_nlte(), apply_rosenkranz_quadratic_scaling(), apply_VVH_scaling(), apply_VVW_scaling(), Absorption::ByLTE, Absorption::Lines::Cutoff(), Absorption::Lines::CutoffFreq(), Absorption::Lines::CutoffFreqMinus(), data, Linefunctions::InternalData::data, Linefunctions::InternalData::datac, Linefunctions::InternalData::dF, Linefunctions::InternalData::dFc, Linefunctions::InternalData::dN, Linefunctions::InternalData::dNc, do_temperature_jacobian(), do_vmr_jacobian(), LineShape::DP, Linefunctions::InternalData::F, Absorption::Lines::F0(), ARTS::Var::f_grid(), Absorption::Lines::F_mean(), Linefunctions::InternalData::Fc, find_cutoff_ranges(), Absorption::Lines::g_low(), Absorption::Lines::g_upp(), EnergyLevelMap::get_ratio_params(), EnergyLevelMap::get_vibtemp_params(), Absorption::Lines::Line(), Absorption::Lines::LinemixingLimit(), Absorption::Lines::LineShapeType(), MapToEigen(), LineShape::mirroredOutput(), Absorption::Lines::Mirroring(), Linefunctions::InternalData::N, N, Linefunctions::InternalData::Nc, Array< base >::nelem(), Absorption::nelem(), Absorption::None, Absorption::Lines::Normalization(), Absorption::Lines::NumLines(), Absorption::Lines::Population(), Absorption::Lines::QuantumIdentity(), Absorption::relaxationtype_relmat(), set_doppler(), set_htp(), set_lorentz(), set_voigt(), Linefunctions::InternalData::SetZero(), Absorption::Lines::ShapeParameters(), Absorption::Lines::ShapeParameters_dT(), Absorption::Lines::ShapeParameters_dVMR(), Zeeman::start(), Absorption::Lines::T0(), LineShape::vmrs(), Absorption::Lines::ZeemanCount(), Absorption::Lines::ZeemanSplitting(), and Absorption::Lines::ZeemanStrength().
Referenced by xsec_species(), and zeeman_on_the_fly().
void Linefunctions::set_doppler | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
Eigen::Ref< Eigen::MatrixXcd > | dF, | ||
Eigen::Ref< Eigen::Matrix< Complex, Eigen::Dynamic, ExpectedDataSize()>> | data, | ||
const Eigen::Ref< const Eigen::VectorXd > | f_grid, | ||
const Numeric & | zeeman_df, | ||
const Numeric & | magnetic_magnitude, | ||
const Numeric & | F0_noshift, | ||
const Numeric & | GD_div_F0, | ||
const AbsorptionLines & | band = AbsorptionLines() , |
||
const Index & | line_ind = 0 , |
||
const ArrayOfRetrievalQuantity & | derivatives_data = ArrayOfRetrievalQuantity() , |
||
const ArrayOfIndex & | derivatives_data_position = ArrayOfIndex() , |
||
const Numeric & | dGD_div_F0_dT = 0.0 |
||
) |
Sets the Doppler line shape.
Normalization is unity.
[in,out] | F | Lineshape. Must be right size |
[in,out] | dF | Lineshape derivative. Must be right size |
[in,out] | data | Block of allocated memory. Output nonsensical |
[in] | f_grid | Frequency grid of computations |
[in] | zeeman_df | Zeeman shift parameter for the line |
[in] | magnetic_magnitude | Absolute strength of the magnetic field |
[in] | F0_noshift | Central frequency without any shifts |
[in] | GD_div_F0 | Frequency-independent part of the Doppler broadening |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
[in] | dGD_div_F0_dT | Temperature derivative of GD_div_F0 |
Definition at line 375 of file linefunctions.cc.
References data, F0, ARTS::Var::f_grid(), Absorption::id_in_line(), is_frequency_parameter(), Array< base >::nelem(), Absorption::Lines::QuantumIdentity(), and ARTS::Var::x().
Referenced by set_cross_section_of_band(), and set_lineshape().
void Linefunctions::set_htp | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
Eigen::Ref< Eigen::MatrixXcd > | dF, | ||
const Eigen::Ref< const Eigen::VectorXd > | f_grid, | ||
const Numeric & | zeeman_df, | ||
const Numeric & | magnetic_magnitude, | ||
const Numeric & | F0_noshift, | ||
const Numeric & | GD_div_F0, | ||
const LineShape::Output & | lso, | ||
const AbsorptionLines & | band = AbsorptionLines() , |
||
const Index & | line_ind = 0 , |
||
const ArrayOfRetrievalQuantity & | derivatives_data = ArrayOfRetrievalQuantity() , |
||
const ArrayOfIndex & | derivatives_data_position = ArrayOfIndex() , |
||
const Numeric & | dGD_div_F0_dT = 0.0 , |
||
const LineShape::Output & | dT = {0, 0, 0, 0, 0, 0, 0, 0, 0} , |
||
const LineShape::Output & | dVMR = {0, 0, 0, 0, 0, 0, 0, 0, 0} |
||
) |
Sets the HTP line shape.
Normalization is unity.
Note: We are not experienced with this line shape and cannot tell what parameters depends on what input. It is therefore likely that this function will have to be adapted in the future
[in,out] | F | Lineshape. Must be right size |
[in,out] | dF | Lineshape derivative. Must be right size |
[in] | f_grid | Frequency grid of computations |
[in] | zeeman_df | Zeeman shift parameter for the line |
[in] | magnetic_magnitude | Absolute strength of the magnetic field |
[in] | F0_noshift | Central frequency without any shifts |
[in] | GD_div_F0 | Frequency-independent part of the Doppler broadening |
[in] | lso | Line shape parameters |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
[in] | dGD_div_F0_dT | Temperature derivative of GD_div_F0 |
[in] | dT | Temperature derivatives of line shape parameters |
[in] | dVMR | VMR derivatives of line shape parameters |
Definition at line 931 of file linefunctions.cc.
References abs, Conversion::freq2kaycm(), Constant::pow2(), Constant::pow3(), Constant::pow4(), and sqrt().
Referenced by set_cross_section_of_band(), and set_lineshape().
void Linefunctions::set_lineshape | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
const Eigen::Ref< const Eigen::VectorXd > | f_grid, | ||
const Absorption::SingleLine & | line, | ||
const Numeric & | temperature, | ||
const Numeric & | zeeman_df, | ||
const Numeric & | magnetic_magnitude, | ||
const Numeric & | doppler_constant, | ||
const LineShape::Output & | lso, | ||
const LineShape::Type | lineshape_type, | ||
const Absorption::MirroringType | mirroring_type, | ||
const Absorption::NormalizationType | norm_type | ||
) |
Sets the lineshape normalized to unity.
No line mixing or linestrength is computed.
[in,out] | F | Lineshape. Must be right size |
[in] | f_grid | Frequency grid of computations |
[in] | line | The absortion line |
[in] | temperature | Atmospheric temperature |
[in] | zeeman_df | Zeeman splitting coefficient |
[in] | magnetic_magnitude | Strength of local magnetic field |
[in] | lso | Line shape parameters |
[in] | lineshape_type | Line shape scheme |
[in] | mirroring_type | Mirroring scheme |
[in] | norm_type | Normalization scheme |
Definition at line 87 of file linefunctions.cc.
References apply_rosenkranz_quadratic_scaling(), apply_VVH_scaling(), apply_VVW_scaling(), data, LineShape::DP, ExpectedDataSize(), Absorption::SingleLine::F0(), ARTS::Var::f_grid(), LineShape::mirroredOutput(), Absorption::None, set_doppler(), set_htp(), set_lorentz(), and set_voigt().
void Linefunctions::set_lorentz | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
Eigen::Ref< Eigen::MatrixXcd > | dF, | ||
Eigen::Ref< Eigen::Matrix< Complex, Eigen::Dynamic, ExpectedDataSize()>> | data, | ||
const Eigen::Ref< const Eigen::VectorXd > | f_grid, | ||
const Numeric & | zeeman_df, | ||
const Numeric & | magnetic_magnitude, | ||
const Numeric & | F0_noshift, | ||
const LineShape::Output & | lso, | ||
const AbsorptionLines & | band = AbsorptionLines() , |
||
const Index & | line_ind = 0 , |
||
const ArrayOfRetrievalQuantity & | derivatives_data = ArrayOfRetrievalQuantity() , |
||
const ArrayOfIndex & | derivatives_data_position = ArrayOfIndex() , |
||
const LineShape::Output & | dT = {0, 0, 0, 0, 0, 0, 0, 0, 0} , |
||
const LineShape::Output & | dVMR = {0, 0, 0, 0, 0, 0, 0, 0, 0} |
||
) |
Sets the Lorentz line shape.
Normalization is unity.
[in,out] | F | Lineshape. Must be right size |
[in,out] | dF | Lineshape derivative. Must be right size |
[in,out] | data | Block of allocated memory. Output nonsensical |
[in] | f_grid | Frequency grid of computations |
[in] | zeeman_df | Zeeman shift parameter for the line |
[in] | magnetic_magnitude | Absolute strength of the magnetic field |
[in] | F0_noshift | Central frequency without any shifts |
[in] | lso | Line shape parameters |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
[in] | dT | Temperature derivatives of line shape parameters |
[in] | dVMR | VMR derivatives of line shape parameters |
Definition at line 234 of file linefunctions.cc.
References LineShape::Output::D0, data, LineShape::Output::DV, dw(), F0, ARTS::Var::f_grid(), LineShape::Output::G0, Absorption::id_in_line(), is_frequency_parameter(), is_magnetic_parameter(), is_pressure_broadening_D0(), is_pressure_broadening_DV(), is_pressure_broadening_G0(), Array< base >::nelem(), and Absorption::Lines::QuantumIdentity().
Referenced by set_cross_section_of_band(), and set_lineshape().
void Linefunctions::set_voigt | ( | Eigen::Ref< Eigen::VectorXcd > | F, |
Eigen::Ref< Eigen::MatrixXcd > | dF, | ||
Eigen::Ref< Eigen::Matrix< Complex, Eigen::Dynamic, ExpectedDataSize()>> | data, | ||
const Eigen::Ref< const Eigen::VectorXd > | f_grid, | ||
const Numeric & | zeeman_df, | ||
const Numeric & | magnetic_magnitude, | ||
const Numeric & | F0_noshift, | ||
const Numeric & | GD_div_F0, | ||
const LineShape::Output & | lso, | ||
const AbsorptionLines & | band = AbsorptionLines() , |
||
const Index & | line_ind = 0 , |
||
const ArrayOfRetrievalQuantity & | derivatives_data = ArrayOfRetrievalQuantity() , |
||
const ArrayOfIndex & | derivatives_data_position = ArrayOfIndex() , |
||
const Numeric & | dGD_div_F0_dT = 0.0 , |
||
const LineShape::Output & | dT = {0, 0, 0, 0, 0, 0, 0, 0, 0} , |
||
const LineShape::Output & | dVMR = {0, 0, 0, 0, 0, 0, 0, 0, 0} |
||
) |
Sets the Voigt line shape.
Normalization is unity.
[in,out] | F | Lineshape. Must be right size |
[in,out] | dF | Lineshape derivative. Must be right size |
[in,out] | data | Block of allocated memory. Output nonsensical |
[in] | f_grid | Frequency grid of computations |
[in] | zeeman_df | Zeeman shift parameter for the line |
[in] | magnetic_magnitude | Absolute strength of the magnetic field |
[in] | F0_noshift | Central frequency without any shifts |
[in] | GD_div_F0 | Frequency-independent part of the Doppler broadening |
[in] | lso | Line shape parameters |
[in] | band | The absorption lines |
[in] | line_ind | The current line's ID |
[in] | derivatives_data | The derivatives in dF |
[in] | derivatives_data_position | The derivatives positions in dF |
[in] | dGD_div_F0_dT | Temperature derivative of GD_div_F0 |
[in] | dT | Temperature derivatives of line shape parameters |
[in] | dVMR | VMR derivatives of line shape parameters |
Definition at line 298 of file linefunctions.cc.
References LineShape::Output::D0, data, LineShape::Output::DV, dw(), F0, ARTS::Var::f_grid(), fac(), LineShape::Output::G0, Absorption::id_in_line(), is_frequency_parameter(), is_magnetic_parameter(), is_pressure_broadening_D0(), is_pressure_broadening_DV(), is_pressure_broadening_G0(), Array< base >::nelem(), Absorption::Lines::QuantumIdentity(), and w().
Referenced by set_cross_section_of_band(), and set_lineshape().