21 PropagationMatrix& propmat_clearsky,
22 StokesVector& nlte_source,
23 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
24 ArrayOfStokesVector& dnlte_source_dx,
31 const Vector& rtp_vmr,
33 const Vector& rtp_mag,
34 const Vector& rtp_los,
35 const Numeric& rtp_pressure,
36 const Numeric& rtp_temperature,
38 const Index& manual_tag,
40 const Numeric& theta0,
41 const Numeric& eta0) {
44 const Index nf = f_grid.nelem();
45 const Index nq = jacobian_quantities.
nelem();
46 const Index ns = abs_species.
nelem();
51 "Only for 3D *rtp_mag* or a manual magnetic field")
53 "*rtp_vmr* must match *abs_species*")
55 "*f_grid* must match *propmat_clearsky*")
57 "*propmat_clearsky* must have *stokes_dim* 4")
59 "*f_grid* must match *nlte_source*")
61 "*nlte_source* must have *stokes_dim* 4")
63 "*dpropmat_clearsky_dx* must match derived form of *jacobian_quantities*")
65 "*dpropmat_clearsky_dx* must have Stokes dim 4 and frequency dim same as *f_grid*")
67 "*dnlte_source_dx* must match derived form of *jacobian_quantities* when non-LTE is on")
69 "*dnlte_source_dx* must have Stokes dim 4 and frequency dim same as *f_grid* when non-LTE is on")
90 const auto polarization_scale_dtheta_data =
92 const auto polarization_scale_deta_data =
96 const Vector f_grid_sparse(0);
97 const Numeric sparse_limit = 0;
113 for (Index ispecies = 0; ispecies < ns; ispecies++) {
115 if (not abs_species[ispecies].nelem() or
116 not abs_species[ispecies].
Zeeman() or
117 not abs_lines_per_species[ispecies].nelem())
119 if (select_abs_species.
nelem() and
120 select_abs_species not_eq abs_species[ispecies])
123 for (
auto& band : abs_lines_per_species[ispecies]) {
125 band, jacobian_quantities, rtp_nlte,
126 band.BroadeningSpeciesVMR(rtp_vmr, abs_species), abs_species[ispecies], rtp_vmr[ispecies],
127 isotopologue_ratios[band.Isotopologue()], rtp_pressure, rtp_temperature, X.H, sparse_limit,
128 polar, Options::LblSpeedup::None,
false);
137 for (Index j=0; j<nq; j++) {
138 auto& deriv = jacobian_quantities[j];
140 if (not deriv.propmattype())
continue;
142 if (deriv == Jacobian::Atm::MagneticU) {
144 pol, dpol_dtheta, dpol_deta,
145 X.dH_du, X.dtheta_du, X.deta_du);
146 }
else if (deriv == Jacobian::Atm::MagneticV) {
148 pol, dpol_dtheta, dpol_deta,
149 X.dH_dv, X.dtheta_dv, X.deta_dv);
150 }
else if (deriv == Jacobian::Atm::MagneticW) {
152 pol, dpol_dtheta, dpol_deta,
153 X.dH_dw, X.dtheta_dw, X.deta_dw);
164 for (Index j=0; j<nq; j++) {
165 auto& deriv = jacobian_quantities[j];
167 if (not deriv.propmattype())
continue;
169 if (deriv == Jacobian::Atm::MagneticU) {
171 pol, dpol_dtheta, dpol_deta,
172 X.dH_du, X.dtheta_du, X.deta_du,
false);
173 }
else if (deriv == Jacobian::Atm::MagneticV) {
175 pol, dpol_dtheta, dpol_deta,
176 X.dH_dv, X.dtheta_dv, X.deta_dv,
false);
177 }
else if (deriv == Jacobian::Atm::MagneticW) {
179 pol, dpol_dtheta, dpol_deta,
180 X.dH_dw, X.dtheta_dw, X.deta_dw,
false);
182 Zeeman::sum(dnlte_source_dx[j], com.
dN(joker, j), pol,
false);
Index nelem() const ARTS_NOEXCEPT
#define ARTS_USER_ERROR_IF(condition,...)
Constains various line scaling functions.
constexpr bool any_negative(const MatpackType &var) noexcept
Checks for negative values.
constexpr auto deg2rad(auto x) noexcept
Converts degrees to radians.
void compute(ComputeData &com, ComputeData &sparse_com, const AbsorptionLines &band, const ArrayOfRetrievalQuantity &jacobian_quantities, const EnergyLevelMap &nlte, const Vector &vmrs, const ArrayOfSpeciesTag &self_tag, const Numeric &self_vmr, const Numeric &isot_ratio, const Numeric &P, const Numeric &T, const Numeric &H, const Numeric &sparse_lim, const Zeeman::Polarization zeeman_polarization, const Options::LblSpeedup speedup_type, const bool robust) ARTS_NOEXCEPT
Compute the absorption of an absorption band.
Implements Zeeman modeling.
const PolarizationVector & SelectPolarization(const AllPolarizationVectors &data, Polarization type) noexcept
Selects the polarization vector depending on polarization type.
constexpr Derived FromPreDerived(Numeric H, Numeric theta, Numeric eta) noexcept
Sets Derived from predefined Derived parameters.
Derived FromGrids(Numeric u, Numeric v, Numeric w, Numeric z, Numeric a) noexcept
Computes the derived plane from ARTS grids.
AllPolarizationVectors AllPolarization_deta(Numeric theta, Numeric eta) noexcept
The derivative of AllPolarization wrt eta.
AllPolarizationVectors AllPolarization_dtheta(Numeric theta, const Numeric eta) noexcept
The derivative of AllPolarization wrt theta.
void dsum(PropagationMatrix &pm, const ComplexVectorView &abs, const ComplexVectorView &dabs, const PolarizationVector &polvec, const PolarizationVector &dpolvec_dtheta, const PolarizationVector &dpolvec_deta, const Numeric dH, const Numeric dt, const Numeric de, const bool do_phase) ARTS_NOEXCEPT
Sums the Zeeman components derivatives into a propagation matrix.
void sum(PropagationMatrix &pm, const ComplexVectorView &abs, const PolarizationVector &polvec, const bool do_phase) ARTS_NOEXCEPT
Sums the Zeeman components into a propagation matrix.
AllPolarizationVectors AllPolarization(Numeric theta, Numeric eta) noexcept
Computes the polarization of each polarization type.
Some molecular constants.
Main computational data for the line shape and strength calculations.
void zeeman_on_the_fly(PropagationMatrix &propmat_clearsky, StokesVector &nlte_source, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, ArrayOfStokesVector &dnlte_source_dx, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfSpeciesTag &select_abs_species, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const SpeciesIsotopologueRatios &isotopologue_ratios, const Vector &f_grid, const Vector &rtp_vmr, const EnergyLevelMap &rtp_nlte, const Vector &rtp_mag, const Vector &rtp_los, const Numeric &rtp_pressure, const Numeric &rtp_temperature, const Index &nlte_do, const Index &manual_tag, const Numeric &H0, const Numeric &theta0, const Numeric &eta0)
Main and only way to compute Zeeman effect.