27#include "matpack_complex.h"
640 ArrayOfMatrix& iy_aux,
641 ArrayOfTensor3& diy_dx,
651 Tensor4& ppvar_trans_cumulat,
652 Tensor4& ppvar_trans_partial,
653 const Index& stokes_dim,
654 const Vector& f_grid,
655 const Index& atmosphere_dim,
656 const Vector& p_grid,
657 const Tensor3& t_field,
659 const Tensor4& vmr_field,
661 const Tensor3& wind_u_field,
662 const Tensor3& wind_v_field,
663 const Tensor3& wind_w_field,
664 const Tensor3& mag_u_field,
665 const Tensor3& mag_v_field,
666 const Tensor3& mag_w_field,
667 const Index& cloudbox_on,
669 const Index& gas_scattering_do,
670 const Tensor4& pnd_field,
671 const ArrayOfTensor4& dpnd_field_dx,
675 const Index& jacobian_do,
678 const Matrix& iy_transmitter,
679 const Agenda& propmat_clearsky_agenda,
680 const Agenda& water_p_eq_agenda,
681 const Agenda& gas_scattering_agenda,
682 const Index& iy_agenda_call1,
683 const Tensor3& iy_transmittance,
684 const Numeric& rte_alonglos_v,
687 const Index j_analytical_do = jacobian_do ? do_analytical_jacobian<1>(jacobian_quantities) : 0;
690 const Index nf = f_grid.nelem();
691 const Index ns = stokes_dim;
692 const Index np = ppath.
np;
693 const Index nq = j_analytical_do ? jacobian_quantities.
nelem() : 0;
700 if (!iy_agenda_call1)
702 "Recursive usage not possible (iy_agenda_call1 must be 1)");
703 if (!iy_transmittance.empty())
704 throw runtime_error(
"*iy_transmittance* must be empty");
705 if (rbi < 1 || rbi > 9)
707 "ppath.background is invalid. Check your "
708 "calculation of *ppath*?");
710 if (dpnd_field_dx.nelem() != jacobian_quantities.
nelem())
712 "*dpnd_field_dx* not properly initialized:\n"
713 "Number of elements in dpnd_field_dx must be equal number of jacobian"
714 " quantities.\n(Note: jacobians have to be defined BEFORE *pnd_field*"
715 " is calculated/set.");
719Jacobian calculation are not supported when gas scattering or suns are included.
720This feature will be added in a future version.
727 if (iy_transmitter.ncols() != ns || iy_transmitter.nrows() != nf) {
729 os <<
"The size of *iy_transmitter* returned from *iy_transmitter_agenda* is\n"
731 <<
" expected size = [" << nf <<
"," << stokes_dim <<
"]\n"
732 <<
" size of iy_transmitter = [" << iy_transmitter.nrows() <<
"," << iy_transmitter.ncols() <<
"]\n";
733 throw runtime_error(os.str());
737 ArrayOfTensor3 diy_dpath = j_analytical_do ?
get_standard_diy_dpath(jacobian_quantities, np, nf, ns,
false) : ArrayOfTensor3(0);
749 const Index naux = iy_aux_vars.
nelem();
752 Index auxOptDepth = -1;
754 for (Index i = 0; i < naux; i++) {
755 iy_aux[i].resize(nf, ns);
758 if (iy_aux_vars[i] ==
"Radiative background")
759 iy_aux[i](joker, 0) = (Numeric)
min((Index)2, rbi - 1);
760 else if (iy_aux_vars[i] ==
"Optical depth")
764 os <<
"The only allowed strings in *iy_aux_vars* are:\n"
765 <<
" \"Radiative background\"\n"
766 <<
" \"Optical depth\"\n"
767 <<
"but you have selected: \"" << iy_aux_vars[i] <<
"\"";
768 throw runtime_error(os.str());
774 ppvar_trans_cumulat.resize(np, nf, ns, ns);
775 ppvar_trans_partial.resize(np, nf, ns, ns);
789 if (np == 1 && rbi == 1) {
792 ppvar_vmr.resize(0, 0);
793 ppvar_wind.resize(0, 0);
794 ppvar_mag.resize(0, 0);
795 ppvar_pnd.resize(0, 0);
796 ppvar_f.resize(0, 0);
797 ppvar_iy.resize(0, 0, 0);
798 ppvar_trans_cumulat = 0;
799 ppvar_trans_partial = 0;
800 for (Index iv = 0; iv < nf; iv++) {
801 for (Index is = 0; is < ns; is++) {
802 ppvar_trans_cumulat(0,iv,is,is) = 1;
803 ppvar_trans_partial(0,iv,is,is) = 1;
809 ppvar_iy.resize(nf, ns, np);
810 ppvar_iy(joker, joker, np - 1) = iy_transmitter;
833 ppvar_f, ppath, f_grid, atmosphere_dim, rte_alonglos_v, ppvar_wind);
836 ArrayOfMatrix ppvar_dpnd_dx(0);
848 clear2cloudy.resize(np);
849 for (Index ip = 0; ip < np; ip++) clear2cloudy[ip] = -1;
856 PropagationMatrix K_this(nf, ns), K_past(nf, ns), Kp(nf, ns);
857 StokesVector
a(nf, ns), S(nf, ns), Sp(nf, ns);
862 Vector sca_fct_dummy;
863 PropagationMatrix K_sca;
864 if (gas_scattering_do) {
865 K_sca = PropagationMatrix(nf, ns);
868 Vector gas_scattering_los_in, gas_scattering_los_out;
872 ArrayOfPropagationMatrix dK_this_dx(nq), dK_past_dx(nq), dKp_dx(nq);
873 ArrayOfStokesVector da_dx(nq), dS_dx(nq), dSp_dx(nq);
875 Index temperature_derivative_position = -1;
878 if (j_analytical_do) {
881 dK_past_dx[iq] = PropagationMatrix(nf, ns);
882 dKp_dx[iq] = PropagationMatrix(nf, ns);
883 da_dx[iq] = StokesVector(nf, ns);
884 dS_dx[iq] = StokesVector(nf, ns);
885 dSp_dx[iq] = StokesVector(nf, ns);
886 if (jacobian_quantities[iq] == Jacobian::Atm::Temperature) {
887 temperature_derivative_position = iq;
888 do_hse = jacobian_quantities[iq].Subtag() ==
894 for (Index ip = 0; ip < np; ip++) {
901 propmat_clearsky_agenda,
903 Vector{ppvar_f(joker, ip)},
904 Vector{ppvar_mag(joker, ip)},
905 Vector{ppath.
los(ip, joker)},
907 Vector{ppvar_vmr(joker, ip)},
913 if (gas_scattering_do) {
921 Vector{ppvar_vmr(joker, ip)},
922 gas_scattering_los_in,
923 gas_scattering_los_out,
925 gas_scattering_agenda);
930 if (j_analytical_do) {
935 ppath.
los(ip, joker),
941 if (clear2cloudy[ip] + 1) {
947 ppvar_pnd(joker, Range(ip, 1)),
951 ppath.
los(ip, joker),
952 ppvar_t[Range(ip, 1)],
957 if (j_analytical_do) {
964 const Numeric dr_dT_past =
965 do_hse ? ppath.
lstep[ip - 1] / (2.0 * ppvar_t[ip - 1]) : 0;
966 const Numeric dr_dT_this =
967 do_hse ? ppath.
lstep[ip - 1] / (2.0 * ppvar_t[ip]) : 0;
978 temperature_derivative_position);
981 swap(K_past, K_this);
982 swap(dK_past_dx, dK_this_dx);
990 if (auxOptDepth >= 0) {
991 for (Index iv = 0; iv < nf; iv++)
992 iy_aux[auxOptDepth](iv, 0) = -std::log(tot_tra[np - 1](iv, 0, 0));
995 lvl_rad[np - 1] = iy_transmitter;
998 for (Index ip = np - 2; ip >= 0; ip--) {
999 lvl_rad[ip] = lvl_rad[ip + 1];
1009 dlyr_tra_above[ip + 1],
1010 dlyr_tra_below[ip + 1],
1011 PropagationMatrix(),
1012 PropagationMatrix(),
1013 ArrayOfPropagationMatrix(),
1014 ArrayOfPropagationMatrix(),
1027 for (Index ip = 0; ip < lvl_rad.
nelem(); ip++) {
1028 ppvar_trans_cumulat(ip, joker, joker, joker) = tot_tra[ip];
1029 ppvar_trans_partial(ip, joker, joker, joker) = lyr_tra[ip];
1030 ppvar_iy(joker, joker, ip) = lvl_rad[ip];
1031 if (j_analytical_do)
1038 if (j_analytical_do) {
1053 jacobian_quantities,
1060 const Index& stokes_dim,
1061 const Vector& f_grid,
1064 const Index nf = f_grid.nelem();
1066 if (instrument_pol.
nelem() != nf)
1067 throw runtime_error(
1068 "The length of *f_grid* and the number of elements "
1069 "in *instrument_pol* must be equal.");
1071 iy_transmitter.resize(nf, stokes_dim);
1073 for (Index i = 0; i < nf; i++) {
1074 stokes2pol(iy_transmitter(i, joker), stokes_dim, instrument_pol[i], 1);
1080 const Index& stokes_dim,
1081 const Vector& f_grid,
1084 const Index nf = f_grid.nelem();
1086 if (instrument_pol.
nelem() != 1)
1087 throw runtime_error(
1088 "The number of elements in *instrument_pol* must be 1.");
1090 iy_transmitter.resize(nf, stokes_dim);
1092 stokes2pol(iy_transmitter(0, joker), stokes_dim, instrument_pol[0], 1);
1094 for (Index i = 1; i < nf; i++) {
1095 iy_transmitter(i, joker) = iy_transmitter(0, joker);
Array< Index > ArrayOfIndex
An array of Index.
base min(const Array< base > &x)
Min function.
The global header file for ARTS.
Constants of physical expressions as constexpr.
void gas_scattering_agendaExecute(Workspace &ws, PropagationMatrix &gas_scattering_coef, TransmissionMatrix &gas_scattering_mat, Vector &gas_scattering_fct_legendre, const Vector &f_grid, const Numeric rtp_pressure, const Numeric rtp_temperature, const Vector &rtp_vmr, const Vector &gas_scattering_los_in, const Vector &gas_scattering_los_out, const Index gas_scattering_output_type, const Agenda &input_agenda)
Index nelem() const ARTS_NOEXCEPT
#define ARTS_USER_ERROR_IF(condition,...)
ArrayOfIndex get_pointers_for_analytical_species(const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfArrayOfSpeciesTag &abs_species)
Help function for analytical jacobian calculations.
ArrayOfTensor3 get_standard_starting_diy_dx(const ArrayOfRetrievalQuantity &jacobian_quantities, Index np, Index nf, Index ns, bool active)
Help function for analytical jacobian calculations.
ArrayOfTensor3 get_standard_diy_dpath(const ArrayOfRetrievalQuantity &jacobian_quantities, Index np, Index nf, Index ns, bool active)
Help function for analytical jacobian calculations.
ArrayOfIndex get_pointers_for_scat_species(const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfString &scat_species, const bool cloudbox_on)
Help function for analytical jacobian calculations.
Routines for setting up the jacobian.
#define FOR_ANALYTICAL_JACOBIANS_DO(what_to_do)
void iy_transmitterMultiplePol(Matrix &iy_transmitter, const Index &stokes_dim, const Vector &f_grid, const ArrayOfIndex &instrument_pol, const Verbosity &)
WORKSPACE METHOD: iy_transmitterMultiplePol.
constexpr Numeric SPEED_OF_LIGHT
constexpr Numeric DEG2RAD
constexpr Numeric RAD2DEG
void iy_transmitterSinglePol(Matrix &iy_transmitter, const Index &stokes_dim, const Vector &f_grid, const ArrayOfIndex &instrument_pol, const Verbosity &)
WORKSPACE METHOD: iy_transmitterSinglePol.
void iyTransmissionStandard(Workspace &ws, Matrix &iy, ArrayOfMatrix &iy_aux, ArrayOfTensor3 &diy_dx, Vector &ppvar_p, Vector &ppvar_t, EnergyLevelMap &ppvar_nlte, Matrix &ppvar_vmr, Matrix &ppvar_wind, Matrix &ppvar_mag, Matrix &ppvar_pnd, Matrix &ppvar_f, Tensor3 &ppvar_iy, Tensor4 &ppvar_trans_cumulat, Tensor4 &ppvar_trans_partial, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Tensor3 &t_field, const EnergyLevelMap &nlte_field, const Tensor4 &vmr_field, const ArrayOfArrayOfSpeciesTag &abs_species, const Tensor3 &wind_u_field, const Tensor3 &wind_v_field, const Tensor3 &wind_w_field, const Tensor3 &mag_u_field, const Tensor3 &mag_v_field, const Tensor3 &mag_w_field, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Index &gas_scattering_do, const Tensor4 &pnd_field, const ArrayOfTensor4 &dpnd_field_dx, const ArrayOfString &scat_species, const ArrayOfArrayOfSingleScatteringData &scat_data, const ArrayOfString &iy_aux_vars, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Ppath &ppath, const Matrix &iy_transmitter, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &gas_scattering_agenda, const Index &iy_agenda_call1, const Tensor3 &iy_transmittance, const Numeric &rte_alonglos_v, const Verbosity &)
WORKSPACE METHOD: iyTransmissionStandard.
Declarations having to do with the four output streams.
constexpr Numeric pi
The following mathematical constants are generated in python Decimal package by the code:
constexpr Numeric speed_of_light
Speed of light [m/s] From: https://en.wikipedia.org/wiki/2019_redefinition_of_SI_base_units Date: 201...
constexpr auto deg2rad(auto x) noexcept
Converts degrees to radians.
constexpr auto rad2deg(auto x) noexcept
Converts radians to degrees.
Index ppath_what_background(const Ppath &ppath)
Returns the case number for the radiative background.
void get_stepwise_scattersky_propmat(StokesVector &ap, PropagationMatrix &Kp, ArrayOfStokesVector &dap_dx, ArrayOfPropagationMatrix &dKp_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, const ConstMatrixView &ppath_1p_pnd, const ArrayOfMatrix &ppath_dpnd_dx, const Index ppath_1p_id, const ArrayOfArrayOfSingleScatteringData &scat_data, const ConstVectorView &ppath_line_of_sight, const ConstVectorView &ppath_temperature, const Index &atmosphere_dim, const bool &jacobian_do)
Computes the contribution by scattering at propagation path point.
void get_ppath_atmvars(Vector &ppath_p, Vector &ppath_t, EnergyLevelMap &ppath_nlte, Matrix &ppath_vmr, Matrix &ppath_wind, Matrix &ppath_mag, const Ppath &ppath, const Index &atmosphere_dim, const ConstVectorView &p_grid, const ConstTensor3View &t_field, const EnergyLevelMap &nlte_field, const ConstTensor4View &vmr_field, const ConstTensor3View &wind_u_field, const ConstTensor3View &wind_v_field, const ConstTensor3View &wind_w_field, const ConstTensor3View &mag_u_field, const ConstTensor3View &mag_v_field, const ConstTensor3View &mag_w_field)
Determines pressure, temperature, VMR, winds and magnetic field for each propgataion path point.
void get_ppath_f(Matrix &ppath_f, const Ppath &ppath, const ConstVectorView &f_grid, const Index &atmosphere_dim, const Numeric &rte_alonglos_v, const ConstMatrixView &ppath_wind)
Determines the Doppler shifted frequencies along the propagation path.
void rtmethods_jacobian_finalisation(Workspace &ws, ArrayOfTensor3 &diy_dx, ArrayOfTensor3 &diy_dpath, const Index &ns, const Index &nf, const Index &np, const Index &atmosphere_dim, const Ppath &ppath, const Vector &ppvar_p, const Vector &ppvar_t, const Matrix &ppvar_vmr, const Index &iy_agenda_call1, const Tensor3 &iy_transmittance, const Agenda &water_p_eq_agenda, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfIndex &jac_species_i)
This function fixes the last steps to made on the Jacobian in some radiative transfer WSMs.
void get_stepwise_clearsky_propmat(Workspace &ws, PropagationMatrix &K, StokesVector &S, Index <e, ArrayOfPropagationMatrix &dK_dx, ArrayOfStokesVector &dS_dx, const Agenda &propmat_clearsky_agenda, const ArrayOfRetrievalQuantity &jacobian_quantities, const Vector &ppath_f_grid, const Vector &ppath_magnetic_field, const Vector &ppath_line_of_sight, const EnergyLevelMap &ppath_nlte, const Vector &ppath_vmrs, const Numeric &ppath_temperature, const Numeric &ppath_pressure, const bool &jacobian_do)
Gets the clearsky propgation matrix and NLTE contributions.
void adapt_stepwise_partial_derivatives(ArrayOfPropagationMatrix &dK_dx, ArrayOfStokesVector &dS_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, const ConstVectorView &ppath_f_grid, const ConstVectorView &ppath_line_of_sight, const Index <e, const Index &atmosphere_dim, const bool &jacobian_do)
Adapts clearsky partial derivatives.
void get_ppath_cloudvars(ArrayOfIndex &clear2cloudy, Matrix &ppath_pnd, ArrayOfMatrix &ppath_dpnd_dx, const Ppath &ppath, const Index &atmosphere_dim, const ArrayOfIndex &cloudbox_limits, const Tensor4 &pnd_field, const ArrayOfTensor4 &dpnd_field_dx)
Determines the particle fields along a propagation path.
Declaration of functions in rte.cc.
void stokes2pol(VectorView w, const Index &stokes_dim, const Index &ipol_1based, const Numeric nv)
stokes2pol
Sensor modelling functions.
The structure to describe a propagation path and releated quantities.
Matrix los
Line-of-sight at each ppath point.
Index np
Number of points describing the ppath.
Vector lstep
The length between ppath points.
Radiation Vector for Stokes dimension 1-4.
Class to keep track of Transmission Matrices for Stokes Dim 1-4.
ArrayOfTransmissionMatrix cumulative_transmission(const ArrayOfTransmissionMatrix &T, const CumulativeTransmission type)
Accumulate the transmission matrix over all layers.
void update_radiation_vector(RadiationVector &I, ArrayOfRadiationVector &dI1, ArrayOfRadiationVector &dI2, const RadiationVector &J1, const RadiationVector &J2, const ArrayOfRadiationVector &dJ1, const ArrayOfRadiationVector &dJ2, const TransmissionMatrix &T, const TransmissionMatrix &PiT, const ArrayOfTransmissionMatrix &dT1, const ArrayOfTransmissionMatrix &dT2, const PropagationMatrix &K1, const PropagationMatrix &K2, const ArrayOfPropagationMatrix &dK1, const ArrayOfPropagationMatrix &dK2, const Numeric r, const Vector &dr1, const Vector &dr2, const Index ia, const Index iz, const RadiativeTransferSolver solver)
Update the Radiation Vector.
void stepwise_transmission(TransmissionMatrix &T, ArrayOfTransmissionMatrix &dT1, ArrayOfTransmissionMatrix &dT2, const PropagationMatrix &K1, const PropagationMatrix &K2, const ArrayOfPropagationMatrix &dK1, const ArrayOfPropagationMatrix &dK2, const Numeric &r, const Numeric &dr_dtemp1, const Numeric &dr_dtemp2, const Index temp_deriv_pos)
Set the stepwise transmission matrix.
Array< RadiationVector > ArrayOfRadiationVector
Array< TransmissionMatrix > ArrayOfTransmissionMatrix