39 ArrayOfPropagationMatrix& dK_dx,
40 ArrayOfStokesVector& dS_dx,
42 const ConstVectorView& ppath_f_grid,
43 const ConstVectorView& ppath_line_of_sight,
45 const Index& atmosphere_dim,
46 const bool& jacobian_do) {
47 if (not jacobian_do)
return;
50 const Index nq = jacobian_quantities.
nelem();
55 for (Index i = 0; i < nq; i++) {
56 if (jacobian_quantities[i] == Jacobian::Type::Sensor or jacobian_quantities[i] == Jacobian::Special::SurfaceString)
continue;
58 if (jacobian_quantities[i].is_wind()) {
59 const auto scale =
get_stepwise_f_partials(ppath_line_of_sight, ppath_f_grid, jacobian_quantities[i].Target().atm, atmosphere_dim);
61 if (not lte) dS_dx[i] *= scale;
66void adjust_los(VectorView los,
const Index& atmosphere_dim) {
67 if (atmosphere_dim == 1) {
70 }
else if (los[0] > 180) {
71 los[0] = 360 - los[0];
73 }
else if (atmosphere_dim == 2) {
75 los[0] = los[0] + 360;
76 }
else if (los[0] > 180) {
77 los[0] = los[0] - 360;
81 if (abs(los[0] - 90) > 90 || abs(los[1]) > 180) {
91 const ConstVectorView& f_grid,
97 const Index nf = iy.nrows();
98 const Index ns = iy.ncols();
103 if (iy_unit ==
"1") {
109 else if (iy_unit ==
"RJBT") {
110 for (Index iv = 0; iv < nf; iv++) {
111 const Numeric scfac =
invrayjean(1, f_grid[iv]);
112 for (Index is = 0; is < ns; is++) {
118 iy(iv, is) *= 2 * scfac;
124 else if (iy_unit ==
"PlanckBT") {
125 for (Index iv = 0; iv < nf; iv++) {
126 for (Index is = ns - 1; is >= 0; is--)
128 if (i_pol[is] == 1) {
129 iy(iv, is) =
invplanck(iy(iv, is), f_grid[iv]);
130 }
else if (i_pol[is] < 5) {
132 iy(iv, is) =
invplanck(0.5 * (iy(iv, 0) + iy(iv, is)), f_grid[iv]) -
133 invplanck(0.5 * (iy(iv, 0) - iy(iv, is)), f_grid[iv]);
135 iy(iv, is) =
invplanck(2 * iy(iv, is), f_grid[iv]);
141 else if (iy_unit ==
"W/(m^2 m sr)") {
142 for (Index iv = 0; iv < nf; iv++) {
143 const Numeric scfac = n * n * f_grid[iv] * (f_grid[iv] /
SPEED_OF_LIGHT);
144 for (Index is = 0; is < ns; is++) {
150 else if (iy_unit ==
"W/(m^2 m-1 sr)") {
156 "Unknown option: iy_unit = \"", iy_unit,
"\"\n"
157 "Recognised choices are: \"1\", \"RJBT\", \"PlanckBT\""
158 "\"W/(m^2 m sr)\" and \"W/(m^2 m-1 sr)\"")
163 const ConstMatrixView& iy,
165 const ConstVectorView& f_grid,
171 const Index nf = iy.nrows();
172 const Index ns = iy.ncols();
173 const Index np = J.npages();
180 if (iy_unit ==
"1") {
186 else if (iy_unit ==
"RJBT") {
187 for (Index iv = 0; iv < nf; iv++) {
188 const Numeric scfac =
invrayjean(1, f_grid[iv]);
189 for (Index is = 0; is < ns; is++) {
192 for (Index ip = 0; ip < np; ip++) {
193 J(ip, iv, is) *= scfac;
197 for (Index ip = 0; ip < np; ip++) {
198 J(ip, iv, is) *= 2 * scfac;
205 else if (iy_unit ==
"PlanckBT") {
206 for (Index iv = 0; iv < f_grid.nelem(); iv++) {
207 for (Index is = ns - 1; is >= 0; is--) {
209 if (i_pol[is] == 1) {
211 }
else if (i_pol[is] < 5) {
213 scfac =
dinvplanckdI(0.5 * (iy(iv, 0) + iy(iv, is)), f_grid[iv]) +
214 dinvplanckdI(0.5 * (iy(iv, 0) - iy(iv, is)), f_grid[iv]);
219 for (Index ip = 0; ip < np; ip++) {
220 J(ip, iv, is) *= scfac;
226 else if (iy_unit ==
"W/(m^2 m sr)") {
227 for (Index iv = 0; iv < nf; iv++) {
228 const Numeric scfac = n * n * f_grid[iv] * (f_grid[iv] /
SPEED_OF_LIGHT);
229 for (Index ip = 0; ip < np; ip++) {
230 for (Index is = 0; is < ns; is++) {
231 J(ip, iv, is) *= scfac;
237 else if (iy_unit ==
"W/(m^2 m-1 sr)") {
243 "Unknown option: iy_unit = \"", iy_unit,
"\"\n"
244 "Recognised choices are: \"1\", \"RJBT\", \"PlanckBT\""
245 "\"W/(m^2 m sr)\" and \"W/(m^2 m-1 sr)\"")
298 const Vector& rte_pos,
300 const Agenda& ppath_step_agenda,
301 const Numeric& ppath_lmax,
302 const Numeric& ppath_lraytrace,
303 const Index& atmosphere_dim,
304 const Vector& p_grid,
305 const Vector& lat_grid,
306 const Vector& lon_grid,
307 const Tensor3& z_field,
308 const Vector& f_grid,
309 const Vector& refellipsoid,
310 const Matrix& z_surface,
314 bool invert_lat =
false;
315 if (atmosphere_dim == 1 && (rte_los[0] < 0 || rte_los[0] > 180)) {
349 Index ilast = ppx.
np - 1;
351 for (Index i = 1; i <= ilast; i++) {
353 lox[i - 1] + ppx.
lstep[i - 1] * (ppx.
nreal[i - 1] + ppx.
nreal[i]) / 2.0;
356 pos.resize(
max(Index(2), atmosphere_dim));
359 if (lox[ilast] < lo0) {
360 const Numeric dl = lo0 - lox[ilast];
361 if (atmosphere_dim < 3) {
362 Numeric x, z, dx, dz;
364 x, z, dx, dz, ppx.
r[ilast], ppx.
pos(ilast, 1), ppx.
los(ilast, 0));
372 Numeric x, y, z, dx, dy, dz;
404 pos[0] =
interp(itw, ppx.
r, gp);
405 pos[1] =
interp(itw, ppx.
pos(joker, 1), gp);
406 if (atmosphere_dim == 3) {
407 pos[2] =
interp(itw, ppx.
pos(joker, 2), gp);
418 const Agenda& ppath_step_agenda,
419 const Index& atmosphere_dim,
420 const Vector& p_grid,
421 const Vector& lat_grid,
422 const Vector& lon_grid,
423 const Tensor3& z_field,
424 const Vector& f_grid,
425 const Vector& refellipsoid,
426 const Matrix& z_surface,
428 const Numeric& ppath_lmax,
429 const Numeric& ppath_lraytrace,
435 for (Index i = 0; i <= ppath.
np - 2; i++) {
436 lp += ppath.
lstep[i];
440 const Vector rte_pos{ppath.
start_pos[Range(0, atmosphere_dim)]};
442 Vector rte_los0(
max(Index(1), atmosphere_dim - 1)), rte_los;
444 rte_los0 = rte_los[Range(0,
max(Index(1), atmosphere_dim - 1))];
501 if (backg1 == backg2) {
503 if (atmosphere_dim < 3) {
504 distance2D(l12, pos1[0], pos1[1], pos2[0], pos2[1]);
506 distance3D(l12, pos1[0], pos1[1], pos1[2], pos2[0], pos2[1], pos2[2]);
514 if (atmosphere_dim == 1) {
515 const Numeric r = refellipsoid[0];
517 }
else if (atmosphere_dim == 2) {
537 const Agenda& ppath_step_agenda,
538 const Index& atmosphere_dim,
539 const Vector& p_grid,
540 const Vector& lat_grid,
541 const Vector& lon_grid,
542 const Tensor3& z_field,
543 const Vector& f_grid,
544 const Vector& refellipsoid,
545 const Matrix& z_surface,
547 const Numeric& ppath_lmax,
548 const Numeric& ppath_lraytrace,
552 "The function *defocusing_sat2sat* can only be used "
553 "for limb sounding geometry.");
562 for (Index i = it; i <= ppath.
np - 2; i++) {
563 lt += ppath.
lstep[i];
565 for (Index i = 0; i < it; i++) {
566 lr += ppath.
lstep[i];
576 const Numeric lf = lr * lt / (lr + lt);
577 const Numeric alt = 1 / (1 - alpha0 * lf / refellipsoid[0]);
580 Numeric alpha1, a1, alpha2, a2, dada;
582 Vector rte_pos{ppath.
end_pos[Range(0, atmosphere_dim)]};
610 rte_los[0] += 2 * dza;
639 dada = (alpha2 - alpha1) / (a2 - a1);
641 dada = (alpha2 - alpha0) / (a2 - a0);
645 const Numeric zlt = 1 / (1 - dada * lf);
655 const Index& atmosphere_dim) {
657 const Numeric f = sqrt(
u *
u +
v *
v +
w *
w);
660 const Numeric za_f = acos(
w / f);
661 const Numeric aa_f = atan2(
u,
v);
669 return f * (cos(za_f) * cos(za_p) + sin(za_f) * sin(za_p) * cos(aa_f - aa_p));
674 const Index& cloudbox_on,
675 const Vector& f_grid,
677 const Vector& rte_pos,
678 const Vector& rte_los,
679 const Vector& rte_pos2,
681 const Agenda& iy_main_agenda) {
682 ArrayOfTensor3 diy_dx;
683 ArrayOfMatrix iy_aux;
686 Tensor3 iy_transmittance(0, 0, 0);
687 const Index iy_agenda_call1 = 1;
689 const Index iy_id = 0;
690 const Index jacobian_do = 0;
715 ArrayOfTensor3& diy_dx,
716 const Tensor3& iy_transmittance,
718 const Index& jacobian_do,
721 const Vector& rte_pos2,
722 const Index& atmosphere_dim,
724 const Index& cloudbox_on,
725 const Index& stokes_dim,
726 const Vector& f_grid,
728 const Tensor3& surface_props_data,
729 const Agenda& iy_main_agenda,
730 const Agenda& iy_space_agenda,
731 const Agenda& iy_surface_agenda,
732 const Agenda& iy_cloudbox_agenda,
733 const Index& iy_agenda_call1,
738 const Index nf = f_grid.nelem();
739 const Index np = ppath.
np;
747 Vector rtp_pos, rtp_los;
748 rtp_pos.resize(atmosphere_dim);
749 rtp_pos = ppath.
pos(np - 1, Range(0, atmosphere_dim));
750 rtp_los.resize(ppath.
los.ncols());
751 rtp_los = ppath.
los(np - 1, joker);
753 out3 <<
"Radiative background: " << ppath.
background <<
"\n";
762 agenda_name =
"iy_space_agenda";
769 agenda_name =
"iy_surface_agenda";
772 const Index los_id = iy_id % (Index)1000;
773 Index iy_id_new = iy_id + (Index)9 * los_id;
777 if (jacobian_do && iy_agenda_call1) {
778 for (Index i = 0; i < jacobian_quantities.
nelem(); i++) {
779 if (jacobian_quantities[i] == Jacobian::Special::SurfaceString) {
780 dsurface_names.push_back(jacobian_quantities[i].Subtag());
784 ArrayOfTensor4 dsurface_rmatrix_dx(dsurface_names.nelem());
785 ArrayOfMatrix dsurface_emission_dx(dsurface_names.nelem());
791 dsurface_emission_dx,
810 agenda_name =
"iy_cloudbox_agenda";
813 ws, iy, f_grid, rtp_pos, rtp_los, iy_cloudbox_agenda);
822 "The size of *iy* returned from *", agenda_name,
"* is\n"
824 " expected size = [", nf,
",", stokes_dim,
"]\n"
825 " size of iy = [", iy.nrows(),
",", iy.ncols(),
"]\n")
835 const Index& atmosphere_dim,
836 const ConstVectorView& p_grid,
837 const ConstTensor3View& t_field,
839 const ConstTensor4View& vmr_field,
840 const ConstTensor3View& wind_u_field,
841 const ConstTensor3View& wind_v_field,
842 const ConstTensor3View& wind_w_field,
843 const ConstTensor3View& mag_u_field,
844 const ConstTensor3View& mag_v_field,
845 const ConstTensor3View& mag_w_field) {
846 const Index np = ppath.
np;
851 itw2p(ppath_p, p_grid, ppath.
gp_p, itw_p);
867 const Index ns = vmr_field.nbooks();
868 ppath_vmr.resize(ns, np);
869 for (Index is = 0; is < ns; is++) {
872 vmr_field(is, joker, joker, joker),
883 ppath_wind.resize(3, np);
886 if (wind_u_field.npages() > 0) {
895 if (wind_v_field.npages() > 0) {
904 if (wind_w_field.npages() > 0) {
915 ppath_mag.resize(3, np);
918 if (mag_u_field.npages() > 0) {
927 if (mag_v_field.npages() > 0) {
936 if (mag_w_field.npages() > 0) {
949 ArrayOfMatrix& ppath_dpnd_dx,
951 const Index& atmosphere_dim,
953 const Tensor4& pnd_field,
954 const ArrayOfTensor4& dpnd_field_dx) {
955 const Index np = ppath.
np;
958 ppath_pnd.resize(pnd_field.nbooks(), np);
960 ppath_dpnd_dx.resize(dpnd_field_dx.nelem());
962 bool any_dpnd =
false;
963 for (Index iq = 0; iq < dpnd_field_dx.nelem(); iq++) {
964 if (dpnd_field_dx[iq].empty()) {
965 ppath_dpnd_dx[iq].resize(0, 0);
968 ppath_dpnd_dx[iq].resize(pnd_field.nbooks(), np);
974 clear2cloudy.resize(np);
978 for (Index ip = 0; ip < np; ip++)
980 Matrix itw(1, Index(pow(2.0, Numeric(atmosphere_dim))));
984 if (atmosphere_dim >= 2) {
987 if (atmosphere_dim == 3) {
1006 for (Index i = 0; i < pnd_field.nbooks(); i++) {
1009 pnd_field(i, joker, joker, joker),
1015 bool any_ppath_dpnd =
false;
1017 for (Index iq = 0; iq < dpnd_field_dx.nelem();
1020 if (!dpnd_field_dx[iq].empty()) {
1021 for (Index i = 0; i < pnd_field.nbooks();
1026 dpnd_field_dx[iq](i, joker, joker, joker),
1032 if (
max(ppath_dpnd_dx[iq](joker, ip)) > 0. ||
1033 min(ppath_dpnd_dx[iq](joker, ip)) < 0.)
1034 any_ppath_dpnd =
true;
1038 if (
max(ppath_pnd(joker, ip)) > 0. ||
min(ppath_pnd(joker, ip)) < 0. ||
1040 clear2cloudy[ip] = nin;
1043 clear2cloudy[ip] = -1;
1046 clear2cloudy[ip] = -1;
1053 const ConstVectorView& f_grid,
1054 const Index& atmosphere_dim,
1055 const Numeric& rte_alonglos_v,
1056 const ConstMatrixView& ppath_wind) {
1058 const Index nf = f_grid.nelem();
1059 const Index np = ppath.
np;
1061 ppath_f.resize(nf, np);
1065 for (Index ip = 0; ip < np; ip++) {
1067 Numeric v_doppler = rte_alonglos_v;
1070 if (ppath_wind(1, ip) != 0 || ppath_wind(0, ip) != 0 ||
1071 ppath_wind(2, ip) != 0) {
1082 if (v_doppler == 0) {
1083 ppath_f(joker, ip) = f_grid;
1088 for (Index iv = 0; iv < nf; iv++) {
1089 ppath_f(iv, ip) =
a * f_grid[iv];
1096 const Index& mblock_index) {
1097 const Index n1y = sensor_response.nrows();
1098 return Range(n1y * mblock_index, n1y);
1103 const ConstVectorView& ppath_f_grid,
1104 const Numeric& ppath_temperature,
1105 const bool& do_temperature_derivative) {
1106 const Index nf = ppath_f_grid.nelem();
1108 for (Index i = 0; i < nf; i++)
1109 B[i] =
planck(ppath_f_grid[i], ppath_temperature);
1111 if (do_temperature_derivative)
1112 for (Index i = 0; i < nf; i++)
1113 dB_dT[i] =
dplanck_dt(ppath_f_grid[i], ppath_temperature);
1118 PropagationMatrix& K,
1121 ArrayOfPropagationMatrix& dK_dx,
1122 ArrayOfStokesVector& dS_dx,
1123 const Agenda& propmat_clearsky_agenda,
1125 const Vector& ppath_f_grid,
1126 const Vector& ppath_magnetic_field,
1127 const Vector& ppath_line_of_sight,
1129 const Vector& ppath_vmrs,
1130 const Numeric& ppath_temperature,
1131 const Numeric& ppath_pressure,
1132 const bool& jacobian_do) {
1142 ppath_magnetic_field,
1143 ppath_line_of_sight,
1148 propmat_clearsky_agenda);
1151 lte = S.allZeroes();
1155 const ConstVectorView& f_grid,
1156 const Jacobian::Atm wind_type,
1157 const Index& atmosphere_dim) {
1159 Numeric dv_doppler_dx = 0.0;
1161 Vector deriv(f_grid);
1163 switch (wind_type) {
1164 case Jacobian::Atm::WindMagnitude:
1165 dv_doppler_dx = 1.0;
1167 case Jacobian::Atm::WindU:
1171 case Jacobian::Atm::WindV:
1175 case Jacobian::Atm::WindW:
1180 ARTS_ASSERT(
false,
"Not allowed to call this function without a wind parameter as wind_type");
1190 PropagationMatrix& Kp,
1191 ArrayOfStokesVector& dap_dx,
1192 ArrayOfPropagationMatrix& dKp_dx,
1194 const ConstMatrixView& ppath_1p_pnd,
1195 const ArrayOfMatrix&
1197 const Index ppath_1p_id,
1199 const ConstVectorView& ppath_line_of_sight,
1200 const ConstVectorView& ppath_temperature,
1201 const Index& atmosphere_dim,
1202 const bool& jacobian_do) {
1203 const Index nf = Kp.NumberOfFrequencies(), stokes_dim = Kp.StokesDimensions();
1213 mirror_los(dir, ppath_line_of_sight, atmosphere_dim);
1214 Matrix dir_array(1, 2, 0.);
1215 dir_array(0, joker) = dir;
1217 ArrayOfArrayOfTensor5 ext_mat_Nse;
1218 ArrayOfArrayOfTensor4 abs_vec_Nse;
1221 ArrayOfTensor5 ext_mat_ssbulk;
1222 ArrayOfTensor4 abs_vec_ssbulk;
1224 Tensor5 ext_mat_bulk;
1225 Tensor4 abs_vec_bulk;
1237 Vector{ppath_temperature},
1256 const Index nf_ssd = abs_vec_bulk.nbooks();
1260 for (Index iv = 0; iv < nf; iv++) {
1262 ap.SetAtPosition(abs_vec_bulk(iv, 0, 0, joker), iv);
1263 Kp.SetAtPosition(ext_mat_bulk(iv, 0, 0, joker, joker), iv);
1265 ap.SetAtPosition(abs_vec_bulk(0, 0, 0, joker), iv);
1266 Kp.SetAtPosition(ext_mat_bulk(0, 0, 0, joker, joker), iv);
1272 if (ppath_dpnd_dx[iq].empty()) {
1273 dap_dx[iq].SetZero();
1274 dKp_dx[iq].SetZero();
1284 ppath_dpnd_dx[iq](joker, Range(ppath_1p_id, 1)),
1292 for (Index iv = 0; iv < nf; iv++) {
1294 dap_dx[iq].SetAtPosition(abs_vec_bulk(iv, 0, 0, joker), iv);
1295 dKp_dx[iq].SetAtPosition(ext_mat_bulk(iv, 0, 0, joker, joker),
1298 dap_dx[iq].SetAtPosition(abs_vec_bulk(0, 0, 0, joker), iv);
1299 dKp_dx[iq].SetAtPosition(ext_mat_bulk(0, 0, 0, joker, joker), iv);
1307 ArrayOfStokesVector& dSp_dx,
1309 const ConstVectorView& ppath_1p_pnd,
1310 const ArrayOfMatrix&
1312 const Index ppath_1p_id,
1314 const ConstTensor7View& cloudbox_field,
1315 const ConstVectorView& za_grid,
1316 const ConstVectorView& aa_grid,
1317 const ConstMatrixView& ppath_line_of_sight,
1318 const GridPos& ppath_pressure,
1319 const Vector& temperature,
1320 const Index& atmosphere_dim,
1321 const bool& jacobian_do,
1322 const Index& t_interp_order) {
1324 "This function handles so far only 1D atmospheres.");
1326 const Index nf = Sp.NumberOfFrequencies();
1327 const Index stokes_dim = Sp.StokesDimensions();
1328 const Index ne = ppath_1p_pnd.nelem();
1330 const Index nza = za_grid.nelem();
1331 const Index naa = aa_grid.nelem();
1332 const Index nq = jacobian_do ? jacobian_quantities.
nelem() : 0;
1340 Tensor3 inc_field(nf, nza, stokes_dim, 0.);
1341 for (Index iv = 0; iv < nf; iv++) {
1342 for (Index iza = 0; iza < nza; iza++) {
1343 for (Index i = 0; i < stokes_dim; i++) {
1344 inc_field(iv, iza, i) =
1345 interp(itw_p, cloudbox_field(iv, joker, 0, 0, iza, 0, i), gp_p);
1352 Matrix idir(nza * naa, 2);
1354 for (Index iza = 0; iza < nza; iza++) {
1355 for (Index iaa = 0; iaa < naa; iaa++) {
1356 idir(ia, 0) = za_grid[iza];
1357 idir(ia, 1) = aa_grid[iaa];
1368 pdir(0, 0) = ppath_line_of_sight(0, 0);
1373 Index nf_ssd = scat_data[0][0].pha_mat_data.nlibraries();
1374 Index duplicate_freqs = ((nf == nf_ssd) ? 0 : 1);
1375 Tensor6 pha_mat_1se(nf_ssd, 1, 1, nza * naa, stokes_dim, stokes_dim);
1378 Tensor3 scat_source_1se(ne, nf, stokes_dim, 0.);
1381 for (Index i_ss = 0; i_ss < scat_data.
nelem(); i_ss++) {
1382 for (Index i_se = 0; i_se < scat_data[i_ss].
nelem(); i_se++) {
1386 if (ppath_1p_pnd[ise_flat] != 0) {
1388 }
else if (jacobian_do) {
1389 for (Index iq = 0; (!val_pnd) && (iq < nq); iq++) {
1390 if ((not(jacobian_quantities[iq] == Jacobian::Type::Sensor) and not(jacobian_quantities[iq] == Jacobian::Special::SurfaceString)) &&
1391 !ppath_dpnd_dx[iq].empty() &&
1392 ppath_dpnd_dx[iq](ise_flat, ppath_1p_id) != 0) {
1402 scat_data[i_ss][i_se],
1409 "Interpolation error for (flat-array) scattering "
1410 "element #", ise_flat,
"\n"
1411 "at location/temperature point #", ppath_1p_id,
"\n")
1414 for (Index iv = 0; iv < nf; iv++) {
1415 if (!duplicate_freqs) {
1418 Tensor3 product_fields(nza, naa, stokes_dim, 0.);
1421 for (Index iza = 0; iza < nza; iza++) {
1422 for (Index iaa = 0; iaa < naa; iaa++) {
1423 for (Index i = 0; i < stokes_dim; i++) {
1424 for (Index j = 0; j < stokes_dim; j++) {
1425 product_fields(iza, iaa, i) +=
1426 pha_mat_1se(this_iv, 0, 0, ia, i, j) *
1427 inc_field(iv, iza, j);
1434 for (Index i = 0; i < stokes_dim; i++) {
1436 product_fields(joker, joker, i), za_grid, aa_grid);
1446 for (Index iv = 0; iv < nf; iv++) {
1447 Vector scat_source(stokes_dim, 0.);
1448 for (ise_flat = 0; ise_flat < ne; ise_flat++) {
1449 for (Index i = 0; i < stokes_dim; i++) {
1451 scat_source_1se(ise_flat, iv, i) * ppath_1p_pnd[ise_flat];
1455 Sp.SetAtPosition(scat_source, iv);
1459 if (ppath_dpnd_dx[iq].empty()) { dSp_dx[iq].SetZero(); }
else {
1461 for (ise_flat = 0; ise_flat < ne; ise_flat++) {
1462 for (Index i = 0; i < stokes_dim; i++) {
1463 scat_source[i] += scat_source_1se(ise_flat, iv, i) *
1464 ppath_dpnd_dx[iq](ise_flat, ppath_1p_id);
1465 dSp_dx[iq].SetAtPosition(scat_source, iv);
1475 ArrayOfArrayOfMatrix& iy_aux_array,
1479 ArrayOfMatrix& diyb_dx,
1481 const Index& mblock_index,
1482 const Index& atmosphere_dim,
1484 const Index& cloudbox_on,
1485 const Index& stokes_dim,
1486 const Matrix& sensor_pos,
1487 const Matrix& sensor_los,
1488 const Matrix& transmitter_pos,
1489 const Matrix& mblock_dlos,
1491 const Agenda& iy_main_agenda,
1492 const Index& j_analytical_do,
1495 const Vector& f_grid,
1504 Vector los(sensor_los.ncols());
1506 los = sensor_los(mblock_index, joker);
1507 if (mblock_dlos.ncols() == 1) {
1508 los[0] += mblock_dlos(ilos, 0);
1515 mblock_dlos(ilos, 0),
1516 mblock_dlos(ilos, 1));
1521 Vector rtp_pos, rtp_pos2(0);
1523 rtp_pos = sensor_pos(mblock_index, joker);
1524 if (!transmitter_pos.empty()) {
1525 rtp_pos2 = transmitter_pos(mblock_index, joker);
1531 ArrayOfTensor3 diy_dx;
1532 Tensor3 iy_transmittance(0, 0, 0);
1533 const Index iy_agenda_call1 = 1;
1535 (Index)1e6 * (mblock_index + 1) + (Index)1e3 * (ilos + 1);
1559 const Index row0 = ilos * nf * stokes_dim;
1563 if (j_analytical_do) {
1566 ip < jacobian_indices[iq][1] - jacobian_indices[iq][0] + 1;
1568 for (Index is = 0; is < stokes_dim; is++) {
1569 diyb_dx[iq](Range(row0 + is, nf, stokes_dim), ip) =
1570 diy_dx[iq](ip, joker, is);
1576 for (Index is = 0; is < stokes_dim; is++) {
1577 iyb[Range(row0 + is, nf, stokes_dim)] = iy(joker, is);
1582 catch (
const std::exception& e) {
1583#pragma omp critical(iyb_calc_fail)
1585 fail_msg = e.what();
1593 ArrayOfVector& iyb_aux,
1594 ArrayOfMatrix& diyb_dx,
1595 Matrix& geo_pos_matrix,
1596 const Index& mblock_index,
1597 const Index& atmosphere_dim,
1599 const Index& cloudbox_on,
1600 const Index& stokes_dim,
1601 const Vector& f_grid,
1602 const Matrix& sensor_pos,
1603 const Matrix& sensor_los,
1604 const Matrix& transmitter_pos,
1605 const Matrix& mblock_dlos,
1607 const Agenda& iy_main_agenda,
1608 const Index& j_analytical_do,
1616 const Index nf = f_grid.nelem();
1617 const Index nlos = mblock_dlos.nrows();
1618 const Index niyb = nf * nlos * stokes_dim;
1623 if (j_analytical_do) {
1624 diyb_dx.resize(jacobian_indices.
nelem());
1626 niyb, jacobian_indices[iq][1] - jacobian_indices[iq][0] + 1);)
1631 geo_pos_matrix.resize(nlos, 5);
1632 geo_pos_matrix = NAN;
1636 ArrayOfArrayOfMatrix iy_aux_array(nlos);
1639 bool failed =
false;
1641 out3 <<
" Parallelizing los loop (" << nlos <<
" iterations, " << nf
1642 <<
" frequencies)\n";
1647#pragma omp parallel for if (!arts_omp_in_parallel()) firstprivate(wss)
1648 for (Index ilos = 0; ilos < nlos; ilos++) {
1650 if (failed)
continue;
1674 jacobian_quantities,
1681 if (geo_pos.nelem()) geo_pos_matrix(ilos, joker) = geo_pos;
1684 if (failed)
continue;
1687 out3 <<
" Not parallelizing los loop (" << nlos <<
" iterations, " << nf
1688 <<
" frequencies)\n";
1690 for (Index ilos = 0; ilos < nlos; ilos++) {
1692 if (failed)
continue;
1716 jacobian_quantities,
1723 if (geo_pos.nelem()) geo_pos_matrix(ilos, joker) = geo_pos;
1726 if (failed)
continue;
1731 "Run-time error in function: iyb_calc\n", fail_msg);
1735 const Index nq = iy_aux_array[0].nelem();
1738 for (Index q = 0; q < nq; q++) {
1739 iyb_aux[q].resize(niyb);
1741 for (Index ilos = 0; ilos < nlos; ilos++) {
1742 const Index row0 = ilos * nf * stokes_dim;
1743 for (Index iv = 0; iv < nf; iv++) {
1744 const Index row1 = row0 + iv * stokes_dim;
1745 const Index i1 =
min(iv, iy_aux_array[ilos][q].nrows() - 1);
1746 for (Index is = 0; is < stokes_dim; is++) {
1747 Index i2 =
min(is, iy_aux_array[ilos][q].ncols() - 1);
1748 iyb_aux[q][row1 + is] = iy_aux_array[ilos][q](i1, i2);
1756 const ConstTensor3View& iy_trans_old,
1757 const ConstTensor3View& iy_trans_new) {
1758 const Index nf = iy_trans_old.npages();
1759 const Index ns = iy_trans_old.ncols();
1766 iy_trans_total.resize(nf, ns, ns);
1768 for (Index iv = 0; iv < nf; iv++) {
1769 mult(iy_trans_total(iv, joker, joker),
1770 iy_trans_old(iv, joker, joker),
1771 iy_trans_new(iv, joker, joker));
1776 const ConstTensor3View& iy_trans,
1777 const ConstMatrixView& iy_old) {
1778 const Index nf = iy_trans.npages();
1779 const Index ns = iy_trans.ncols();
1785 iy_new.resize(nf, ns);
1787 for (Index iv = 0; iv < nf; iv++) {
1788 mult(iy_new(iv, joker), iy_trans(iv, joker, joker), iy_old(iv, joker));
1793 const ConstVectorView& los,
1794 const Index& atmosphere_dim) {
1795 los_mirrored.resize(2);
1797 if (atmosphere_dim == 1) {
1798 los_mirrored[0] = 180 - los[0];
1799 los_mirrored[1] = 180;
1800 }
else if (atmosphere_dim == 2) {
1801 los_mirrored[0] = 180 - fabs(los[0]);
1803 los_mirrored[1] = 180;
1805 los_mirrored[1] = 0;
1807 }
else if (atmosphere_dim == 3) {
1808 los_mirrored[0] = 180 - los[0];
1809 los_mirrored[1] = los[1] + 180;
1810 if (los_mirrored[1] > 180) {
1811 los_mirrored[1] -= 360;
1817 const Index& stokes_dim,
1818 const Numeric& rotangle) {
1829 if (stokes_dim > 2) {
1837 if (stokes_dim > 3) {
1846 const Index& stokes_dim) {
1850 Cs.resize(stokes_dim, stokes_dim);
1852 if (stokes_dim > 1 ) {
1853 Cs(0,1) = Cs(1,0) = 1;
1855 if (stokes_dim > 2 ) {
1856 Cs(0,2) = Cs(1,2) = Cs(2,0) = Cs(2,1) = 0;
1858 if (stokes_dim > 3 ) {
1859 Cs(0,3) = Cs(1,3) = Cs(2,3) = Cs(3,0) = Cs(3,1) = Cs(3,2) = 0;
1867 const Index& stokes_dim,
1868 const Numeric& rotangle) {
1872 L.resize(stokes_dim, stokes_dim);
1874 if (stokes_dim > 1 ) {
1876 const Numeric c2 = cos(alpha);
1877 L(0,1) = L(1,0) = 0;
1879 if (stokes_dim > 2 ) {
1880 const Numeric s2 = sin(alpha);
1881 L(0,2) = L(2,0) = 0;
1885 if (stokes_dim > 3 ) {
1886 L(0,3) = L(1,3) = L(2,3) = L(3,0) = L(3,1) = L(3,2) = 0;
1894 const Index& stokes_dim) {
1898 Cm.resize(stokes_dim, stokes_dim);
1900 if (stokes_dim > 1 ) {
1901 Cm(0,1) = Cm(1,0) = 0.5;
1903 if (stokes_dim > 2 ) {
1904 Cm(0,2) = Cm(1,2) = Cm(2,0) = Cm(2,1) = 0;
1906 if (stokes_dim > 3 ) {
1907 Cm(0,3) = Cm(1,3) = Cm(2,3) = Cm(3,0) = Cm(3,1) = Cm(3,2) = 0;
1916 const Index& atmosphere_dim,
1917 const ConstVectorView& lat_grid,
1918 const ConstVectorView& lat_true,
1919 const ConstVectorView& lon_true,
1920 const ConstVectorView& pos) {
1923 if (atmosphere_dim == 1) {
1931 else if (atmosphere_dim == 2) {
1932 ARTS_ASSERT(lat_true.nelem() == lat_grid.nelem());
1933 ARTS_ASSERT(lon_true.nelem() == lat_grid.nelem());
1936 gridpos(gp, lat_grid, pos[1]);
1938 lat =
interp(itw, lat_true, gp);
1939 lon =
interp(itw, lon_true, gp);
1950 ArrayOfTensor3& diy_dx,
1951 ArrayOfTensor3& diy_dpath,
1955 const Index& atmosphere_dim,
1957 const Vector& ppvar_p,
1958 const Vector& ppvar_t,
1959 const Matrix& ppvar_vmr,
1960 const Index& iy_agenda_call1,
1961 const Tensor3& iy_transmittance,
1962 const Agenda& water_p_eq_agenda,
1966 if (!iy_agenda_call1) {
1970 Y.resize(ns, diy_dpath[iq].npages());
1971 for (Index iv = 0; iv < nf; iv++) {
1972 X = transpose(diy_dpath[iq](joker, iv, joker));
1973 mult(Y, iy_transmittance(iv, joker, joker), X);
1974 diy_dpath[iq](joker, iv, joker) = transpose(Y);
1980 Tensor3 water_p_eq(0, 0, 0);
1983 for (Index iq = 0; iq < jacobian_quantities.
nelem(); iq++) {
1987 if (not(jacobian_quantities[iq] == Jacobian::Type::Sensor) and
1988 not(jacobian_quantities[iq] == Jacobian::Special::SurfaceString) and
1989 jac_species_i[iq] >= 0) {
1990 if (jacobian_quantities[iq].Mode() ==
"vmr") {
1993 else if (jacobian_quantities[iq].Mode() ==
"rel") {
1995 for (Index ip = 0; ip < np; ip++) {
1996 diy_dpath[iq](ip, joker, joker) *= ppvar_vmr(jac_species_i[iq], ip);
2000 else if (jacobian_quantities[iq].Mode() ==
"nd") {
2002 for (Index ip = 0; ip < np; ip++) {
2003 diy_dpath[iq](ip, joker, joker) /=
2008 else if (jacobian_quantities[iq].Mode() ==
"rh") {
2010 Tensor3 t_data(ppvar_t.nelem(), 1, 1);
2011 t_data(joker, 0, 0) = ppvar_t;
2013 for (Index ip = 0; ip < np; ip++) {
2014 diy_dpath[iq](ip, joker, joker) *= water_p_eq(ip, 0, 0) / ppvar_p[ip];
2018 else if (jacobian_quantities[iq].Mode() ==
"q") {
2020 diy_dpath[iq](joker, joker, joker) /= 0.622;
2030 for (Index iq = 0; iq < jacobian_quantities.
nelem(); iq++) {
2037 if (jacobian_quantities[iq] == Jacobian::Atm::Temperature) {
2039 for (Index ia = 0; ia < jacobian_quantities.
nelem(); ia++) {
2040 if (jac_species_i[ia] >= 0) {
2041 if (jacobian_quantities[ia].Mode() ==
"nd") {
2042 for (Index ip = 0; ip < np; ip++) {
2043 Matrix ddterm{diy_dpath[ia](ip, joker, joker)};
2044 ddterm *= ppvar_vmr(jac_species_i[ia], ip) *
2047 diy_dpath[iq](ip, joker, joker) += ddterm;
2049 }
else if (jacobian_quantities[ia].Mode() ==
"rh") {
2050 Tensor3 t_data(ppvar_t.nelem(), 1, 1);
2051 t_data(joker, 0, 0) = ppvar_t;
2053 if (water_p_eq.npages() == 0) {
2055 ws, water_p_eq, t_data, water_p_eq_agenda);
2058 Tensor3 water_p_eq1K;
2059 t_data(joker, 0, 0) += 1;
2061 ws, water_p_eq1K, t_data, water_p_eq_agenda);
2063 for (Index ip = 0; ip < np; ip++) {
2064 const Numeric p_eq = water_p_eq(ip, 0, 0);
2065 const Numeric p_eq1K = water_p_eq1K(ip, 0, 0);
2066 Matrix ddterm{diy_dpath[ia](ip, joker, joker)};
2067 ddterm *= ppvar_vmr(jac_species_i[ia], ip) *
2068 (ppvar_p[ip] / pow(p_eq, 2.0)) * (p_eq1K - p_eq);
2069 diy_dpath[iq](ip, joker, joker) += ddterm;
2079 jacobian_quantities[iq],
2088 ArrayOfTensor3& diy_dx,
2092 const Vector& f_grid,
2095 const Index& j_analytical_do,
2102 n = ppath.
nreal[np - 1];
2107 for (Index is = 0; is < ns; is++) {
2113 if (j_analytical_do) {
2122 for (Index ip = 0; ip < ppath.
np; ip++) {
2124 ppvar_iy(joker, joker, ip), iy_unit, f_grid, ppath.
nreal[ip], i_pol);
2130 ArrayOfArrayOfVector& iyb_aux_array,
2139 const Index& atmosphere_dim,
2141 const Index& cloudbox_on,
2142 const Index& stokes_dim,
2143 const Vector& f_grid,
2144 const Matrix& sensor_pos,
2145 const Matrix& sensor_los,
2146 const Matrix& transmitter_pos,
2147 const Matrix& mblock_dlos,
2148 const Sparse& sensor_response,
2149 const Vector& sensor_response_f,
2151 const Matrix& sensor_response_dlos,
2153 const Agenda& iy_main_agenda,
2154 const Agenda& jacobian_agenda,
2155 const Index& jacobian_do,
2160 const Index& mblock_index,
2162 const Index& j_analytical_do) {
2166 Vector iyb, iyb_error, yb(n1y);
2167 ArrayOfMatrix diyb_dx;
2168 Matrix geo_pos_matrix;
2172 iyb_aux_array[mblock_index],
2188 jacobian_quantities,
2196 const Index row0 = rowind.offset;
2198 mult(yb, sensor_response, iyb);
2204 for (Index i = 0; i < n1y; i++) {
2205 const Index ii = row0 + i;
2207 "One or several NaNs found in *y*.");
2208 y_f[ii] = sensor_response_f[i];
2209 y_pol[ii] = sensor_response_pol[i];
2210 y_pos(ii, joker) = sensor_pos(mblock_index, joker);
2211 y_los(ii, joker) = sensor_los(mblock_index, joker);
2212 y_los(ii, 0) += sensor_response_dlos(i, 0);
2213 if (sensor_response_dlos.ncols() > 1) {
2214 y_los(ii, 1) += sensor_response_dlos(i, 1);
2221 if (j_analytical_do) {
2223 mult(jacobian(rowind,
2224 Range(jacobian_indices[iq][0],
2225 jacobian_indices[iq][1] -
2226 jacobian_indices[iq][0] + 1)),
2235 ws, jacobian, mblock_index, iyb, yb, jacobian_agenda);
2239 if (!std::isnan(geo_pos_matrix(0, 0)))
2242 const Index nfs = f_grid.nelem() * stokes_dim;
2243 for (Index i = 0; i < n1y; i++) {
2245 Numeric rmax = -99e99;
2246 for (Index j = 0; j < sensor_response.ncols(); j++) {
2247 if (sensor_response(i, j) > rmax) {
2248 rmax = sensor_response(i, j);
2252 const auto jhit = Index(floor(jmax / nfs));
2253 y_geo(row0 + i, joker) = geo_pos_matrix(jhit, joker);
2258 catch (
const std::exception& e) {
2259#pragma omp critical(yCalc_fail)
2261 fail_msg = e.what();
2268 const Vector& f_grid,
2269 const Numeric& ze_tref,
2270 const Numeric& k2) {
2271 const Index nf = f_grid.nelem();
2276 Matrix complex_n(0, 0);
2284 for (Index iv = 0; iv < nf; iv++) {
2290 Complex n(complex_n(iv, 0), complex_n(iv, 1));
2292 Complex K = (n2 - Numeric(1.0)) / (n2 + Numeric(2.0));
2293 Numeric absK = abs(K);
2300 fac[iv] =
a * la * la * la * la / K2;
Array< Index > ArrayOfIndex
An array of Index.
Index TotalNumberOfElements(const Array< Array< base > > &aa)
Determine total number of elements in an ArrayOfArray.
base max(const Array< base > &x)
Max function.
base min(const Array< base > &x)
Min function.
Constants of physical expressions as constexpr.
Simple constexpr math that we can make use of in Arts.
int arts_omp_get_max_threads()
Wrapper for omp_get_max_threads.
void iy_cloudbox_agendaExecute(Workspace &ws, Matrix &iy, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Agenda &input_agenda)
void iy_surface_agendaExecute(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, const ArrayOfTensor4 &dsurface_rmatrix_dx, const ArrayOfMatrix &dsurface_emission_dx, const String &iy_unit, const Tensor3 &iy_transmittance, const Index iy_id, const Index cloudbox_on, const Index jacobian_do, const Agenda &iy_main_agenda, const Vector &f_grid, const EnergyLevelMap &nlte_field, const Vector &rtp_pos, const Vector &rtp_los, const Vector &rte_pos2, const Tensor3 &surface_props_data, const ArrayOfString &dsurface_names, const Agenda &input_agenda)
void iy_space_agendaExecute(Workspace &ws, Matrix &iy, const Vector &f_grid, const Vector &rtp_pos, const Vector &rtp_los, const Agenda &input_agenda)
void iy_main_agendaExecute(Workspace &ws, Matrix &iy, ArrayOfMatrix &iy_aux, Ppath &ppath, ArrayOfTensor3 &diy_dx, Vector &geo_pos, const Index iy_agenda_call1, const Tensor3 &iy_transmittance, const ArrayOfString &iy_aux_vars, const Index iy_id, const String &iy_unit, const Index cloudbox_on, const Index jacobian_do, const Vector &f_grid, const EnergyLevelMap &nlte_field, const Vector &rte_pos, const Vector &rte_los, const Vector &rte_pos2, const Agenda &input_agenda)
void jacobian_agendaExecute(Workspace &ws, Matrix &jacobian, const Index mblock_index, const Vector &iyb, const Vector &yb, const Agenda &input_agenda)
void propmat_clearsky_agendaExecute(Workspace &ws, PropagationMatrix &propmat_clearsky, StokesVector &nlte_source, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, ArrayOfStokesVector &dnlte_source_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfSpeciesTag &select_abs_species, const Vector &f_grid, const Vector &rtp_mag, const Vector &rtp_los, const Numeric rtp_pressure, const Numeric rtp_temperature, const EnergyLevelMap &rtp_nlte, const Vector &rtp_vmr, const Agenda &input_agenda)
void water_p_eq_agendaExecute(Workspace &ws, Tensor3 &water_p_eq_field, const Tensor3 &t_field, const Agenda &input_agenda)
This can be used to make arrays out of anything.
Index nelem() const ARTS_NOEXCEPT
bool is_gp_inside_cloudbox(const GridPos &gp_p, const GridPos &gp_lat, const GridPos &gp_lon, const ArrayOfIndex &cloudbox_limits, const bool &include_boundaries, const Index &atmosphere_dim)
void mult(MatrixView C, ConstMatrixView A, const Block &B)
#define ARTS_ASSERT(condition,...)
#define ARTS_USER_ERROR(...)
#define ARTS_USER_ERROR_IF(condition,...)
void distance3D(Numeric &l, const Numeric &r1, const Numeric &lat1, const Numeric &lon1, const Numeric &r2, const Numeric &lat2, const Numeric &lon2)
distance3D
void distance2D(Numeric &l, const Numeric &r1, const Numeric &lat1, const Numeric &r2, const Numeric &lat2)
distance2D
void cart2sph(Numeric &r, Numeric &lat, Numeric &lon, const Numeric &x, const Numeric &y, const Numeric &z, const Numeric &lat0, const Numeric &lon0, const Numeric &za0, const Numeric &aa0)
cart2sph
Numeric refell2r(ConstVectorView refellipsoid, const Numeric &lat)
refell2r
Numeric sphdist(const Numeric &lat1, const Numeric &lon1, const Numeric &lat2, const Numeric &lon2)
sphdist
void poslos2cart(Numeric &x, Numeric &z, Numeric &dx, Numeric &dz, const Numeric &r, const Numeric &lat, const Numeric &za)
poslos2cart
void cart2pol(Numeric &r, Numeric &lat, const Numeric &x, const Numeric &z, const Numeric &lat0, const Numeric &za0)
cart2pol
void gridpos(ArrayOfGridPos &gp, ConstVectorView old_grid, ConstVectorView new_grid, const Numeric &extpolfac)
Set up a grid position Array.
void interpweights(VectorView itw, const GridPos &tc)
Red 1D interpolation weights.
Numeric interp(ConstVectorView itw, ConstVectorView a, const GridPos &tc)
Red 1D Interpolate.
void gridpos_copy(GridPos &gp_new, const GridPos &gp_old)
gridpos_copy
void diy_from_path_to_rgrids(Tensor3View diy_dx, const RetrievalQuantity &jacobian_quantity, ConstTensor3View diy_dpath, const Index &atmosphere_dim, const Ppath &ppath, ConstVectorView ppath_p)
Maps jacobian data for points along the propagation path, to jacobian retrieval grid data.
#define FOR_ANALYTICAL_JACOBIANS_DO2(what_to_do)
#define FOR_ANALYTICAL_JACOBIANS_DO(what_to_do)
Array< RetrievalQuantity > ArrayOfRetrievalQuantity
Numeric fac(const Index n)
fac
Numeric AngIntegrate_trapezoid(ConstMatrixView Integrand, ConstVectorView za_grid, ConstVectorView aa_grid)
AngIntegrate_trapezoid.
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 Numeric c
Speed of light convenience name [m/s].
constexpr Numeric temperature_at_0c
Global constant, Temperature in Celsius of 0 Kelvin.
constexpr auto deg2rad(auto x) noexcept
Converts degrees to radians.
auto sind(auto x) noexcept
Returns the sine of the deg2rad of the input.
auto cosd(auto x) noexcept
Returns the cosine of the deg2rad of the input.
constexpr auto pow4(auto x) noexcept
power of four
void opt_prop_ScatSpecBulk(ArrayOfTensor5 &ext_mat, ArrayOfTensor4 &abs_vec, ArrayOfIndex &ptype, const ArrayOfArrayOfTensor5 &ext_mat_se, const ArrayOfArrayOfTensor4 &abs_vec_se, const ArrayOfArrayOfIndex &ptypes_se, ConstMatrixView pnds, ConstMatrixView t_ok)
Scattering species bulk extinction and absorption.
void pha_mat_1ScatElem(Tensor6View pha_mat, Index &ptype, VectorView t_ok, const SingleScatteringData &ssd, const Vector &T_array, const Matrix &pdir_array, const Matrix &idir_array, const Index &f_start, const Index &t_interp_order)
Preparing phase matrix from one scattering element.
void opt_prop_Bulk(Tensor5 &ext_mat, Tensor4 &abs_vec, Index &ptype, const ArrayOfTensor5 &ext_mat_ss, const ArrayOfTensor4 &abs_vec_ss, const ArrayOfIndex &ptypes_ss)
one-line descript
void opt_prop_NScatElems(ArrayOfArrayOfTensor5 &ext_mat, ArrayOfArrayOfTensor4 &abs_vec, ArrayOfArrayOfIndex &ptypes, Matrix &t_ok, const ArrayOfArrayOfSingleScatteringData &scat_data, const Index &stokes_dim, const Vector &T_array, const Matrix &dir_array, const Index &f_index, const Index &t_interp_order)
Extinction and absorption from all scattering elements.
Numeric dinvplanckdI(const Numeric &i, const Numeric &f)
dinvplanckdI
Numeric planck(const Numeric &f, const Numeric &t)
planck
Numeric invrayjean(const Numeric &i, const Numeric &f)
invrayjean
Numeric invplanck(const Numeric &i, const Numeric &f)
invplanck
Numeric dplanck_dt(const Numeric &f, const Numeric &t)
dplanck_dt
This file contains declerations of functions of physical character.
constexpr Numeric number_density(Numeric p, Numeric t) noexcept
number_density
void zaaa2cart(Numeric &dx, Numeric &dy, Numeric &dz, const Numeric &za, const Numeric &aa)
Converts zenith and azimuth angles to a cartesian unit vector.
Index ppath_what_background(const Ppath &ppath)
Returns the case number for the radiative background.
void add_za_aa(Numeric &za, Numeric &aa, const Numeric &za0, const Numeric &aa0, const Numeric &dza, const Numeric &daa)
Adds up zenith and azimuth angles.
void cart2zaaa(Numeric &za, Numeric &aa, const Numeric &dx, const Numeric &dy, const Numeric &dz)
Converts a cartesian directional vector to zenith and azimuth.
void find_tanpoint(Index &it, const Ppath &ppath)
Identifies the tangent point of a propagation path.
void ppath_calc(Workspace &ws, Ppath &ppath, const Agenda &ppath_step_agenda, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const Vector &f_grid, const Vector &refellipsoid, const Matrix &z_surface, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Vector &rte_pos, const Vector &rte_los, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const bool &ppath_inside_cloudbox_do, const Verbosity &verbosity)
This is the core for the WSM ppathStepByStep.
Propagation path structure and functions.
void complex_n_water_liebe93(Matrix &complex_n, const Vector &f_grid, const Numeric &t)
complex_n_water_liebe93
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 ze_cfac(Vector &fac, const Vector &f_grid, const Numeric &ze_tref, const Numeric &k2)
Calculates factor to convert back-scattering to Ze.
void iyb_calc_body(bool &failed, String &fail_msg, ArrayOfArrayOfMatrix &iy_aux_array, Workspace &ws, Ppath &ppath, Vector &iyb, ArrayOfMatrix &diyb_dx, Vector &geo_pos, const Index &mblock_index, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, const Matrix &sensor_pos, const Matrix &sensor_los, const Matrix &transmitter_pos, const Matrix &mblock_dlos, const String &iy_unit, const Agenda &iy_main_agenda, const Index &j_analytical_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfArrayOfIndex &jacobian_indices, const Vector &f_grid, const ArrayOfString &iy_aux_vars, const Index &ilos, const Index &nf)
void mueller_stokes2modif(Matrix &Cm, const Index &stokes_dim)
mueller_stokes2modif
void get_stepwise_scattersky_source(StokesVector &Sp, ArrayOfStokesVector &dSp_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, const ConstVectorView &ppath_1p_pnd, const ArrayOfMatrix &ppath_dpnd_dx, const Index ppath_1p_id, const ArrayOfArrayOfSingleScatteringData &scat_data, const ConstTensor7View &cloudbox_field, const ConstVectorView &za_grid, const ConstVectorView &aa_grid, const ConstMatrixView &ppath_line_of_sight, const GridPos &ppath_pressure, const Vector &temperature, const Index &atmosphere_dim, const bool &jacobian_do, const Index &t_interp_order)
Calculates the stepwise scattering source terms.
constexpr Numeric SPEED_OF_LIGHT
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 pos2true_latlon(Numeric &lat, Numeric &lon, const Index &atmosphere_dim, const ConstVectorView &lat_grid, const ConstVectorView &lat_true, const ConstVectorView &lon_true, const ConstVectorView &pos)
Determines the true alt and lon for an "ARTS position".
Numeric dotprod_with_los(const ConstVectorView &los, const Numeric &u, const Numeric &v, const Numeric &w, const Index &atmosphere_dim)
Calculates the dot product between a field and a LOS.
void apply_iy_unit2(Tensor3View J, const ConstMatrixView &iy, const String &iy_unit, const ConstVectorView &f_grid, const Numeric &n, const ArrayOfIndex &i_pol)
Largely as apply_iy_unit but operates on jacobian data.
void iy_transmittance_mult(Tensor3 &iy_trans_total, const ConstTensor3View &iy_trans_old, const ConstTensor3View &iy_trans_new)
Multiplicates iy_transmittance with transmissions.
void apply_iy_unit(MatrixView iy, const String &iy_unit, const ConstVectorView &f_grid, const Numeric &n, const ArrayOfIndex &i_pol)
Performs conversion from radiance to other units, as well as applies refractive index to fulfill the ...
void get_iy(Workspace &ws, Matrix &iy, const Index &cloudbox_on, const Vector &f_grid, const EnergyLevelMap &nlte_field, const Vector &rte_pos, const Vector &rte_los, const Vector &rte_pos2, const String &iy_unit, const Agenda &iy_main_agenda)
Basic call of iy_main_agenda.
void adjust_los(VectorView los, const Index &atmosphere_dim)
Ensures that the zenith and azimuth angles of a line-of-sight vector are inside defined ranges.
void defocusing_general_sub(Workspace &ws, Vector &pos, Vector &rte_los, Index &background, const Vector &rte_pos, const Numeric &lo0, const Agenda &ppath_step_agenda, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const Vector &f_grid, const Vector &refellipsoid, const Matrix &z_surface, const Verbosity &verbosity)
Just to avoid duplicatuion of code in defocusing_general.
void defocusing_sat2sat(Workspace &ws, Numeric &dlf, const Agenda &ppath_step_agenda, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const Vector &f_grid, const Vector &refellipsoid, const Matrix &z_surface, const Ppath &ppath, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Numeric &dza, const Verbosity &verbosity)
Calculates defocusing for limb measurements between two satellites.
void mueller_modif2stokes(Matrix &Cs, const Index &stokes_dim)
mueller_modif2stokes
void get_iy_of_background(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, const Tensor3 &iy_transmittance, const Index &iy_id, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Ppath &ppath, const Vector &rte_pos2, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, const Vector &f_grid, const String &iy_unit, const Tensor3 &surface_props_data, const Agenda &iy_main_agenda, const Agenda &iy_space_agenda, const Agenda &iy_surface_agenda, const Agenda &iy_cloudbox_agenda, const Index &iy_agenda_call1, const Verbosity &verbosity)
Determines iy of the "background" of a propgation path.
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 iyb_calc(Workspace &ws, Vector &iyb, ArrayOfVector &iyb_aux, ArrayOfMatrix &diyb_dx, Matrix &geo_pos_matrix, const Index &mblock_index, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, const Vector &f_grid, const Matrix &sensor_pos, const Matrix &sensor_los, const Matrix &transmitter_pos, const Matrix &mblock_dlos, const String &iy_unit, const Agenda &iy_main_agenda, const Index &j_analytical_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfArrayOfIndex &jacobian_indices, const ArrayOfString &iy_aux_vars, const Verbosity &verbosity)
Performs calculations for one measurement block, on iy-level.
void mueller_rotation(Matrix &L, const Index &stokes_dim, const Numeric &rotangle)
mueller_rotation
void rtmethods_unit_conversion(Matrix &iy, ArrayOfTensor3 &diy_dx, Tensor3 &ppvar_iy, const Index &ns, const Index &np, const Vector &f_grid, const Ppath &ppath, const ArrayOfRetrievalQuantity &jacobian_quantities, const Index &j_analytical_do, const String &iy_unit)
This function handles the unit conversion to be done at the end of some radiative transfer WSMs.
void bending_angle1d(Numeric &alpha, const Ppath &ppath)
Calculates the bending angle for a 1D atmosphere.
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 mirror_los(Vector &los_mirrored, const ConstVectorView &los, const Index &atmosphere_dim)
Determines the backward direction for a given line-of-sight.
void yCalc_mblock_loop_body(bool &failed, String &fail_msg, ArrayOfArrayOfVector &iyb_aux_array, Workspace &ws, Vector &y, Vector &y_f, ArrayOfIndex &y_pol, Matrix &y_pos, Matrix &y_los, Matrix &y_geo, Matrix &jacobian, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, const Vector &f_grid, const Matrix &sensor_pos, const Matrix &sensor_los, const Matrix &transmitter_pos, const Matrix &mblock_dlos, const Sparse &sensor_response, const Vector &sensor_response_f, const ArrayOfIndex &sensor_response_pol, const Matrix &sensor_response_dlos, const String &iy_unit, const Agenda &iy_main_agenda, const Agenda &jacobian_agenda, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfArrayOfIndex &jacobian_indices, const ArrayOfString &iy_aux_vars, const Verbosity &verbosity, const Index &mblock_index, const Index &n1y, const Index &j_analytical_do)
Performs calculations for one measurement block, on y-level.
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.
constexpr Numeric TEMP_0_C
void muellersparse_rotation(Sparse &H, const Index &stokes_dim, const Numeric &rotangle)
muellersparse_rotation
Vector get_stepwise_f_partials(const ConstVectorView &line_of_sight, const ConstVectorView &f_grid, const Jacobian::Atm wind_type, const Index &atmosphere_dim)
Computes the ratio that a partial derivative with regards to frequency relates to the wind of come co...
Range get_rowindex_for_mblock(const Sparse &sensor_response, const Index &mblock_index)
Returns the "range" of y corresponding to a measurement block.
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_stepwise_blackbody_radiation(VectorView B, VectorView dB_dT, const ConstVectorView &ppath_f_grid, const Numeric &ppath_temperature, const bool &do_temperature_derivative)
Get the blackbody radiation at propagation path point.
void defocusing_general(Workspace &ws, Numeric &dlf, const Agenda &ppath_step_agenda, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const Vector &f_grid, const Vector &refellipsoid, const Matrix &z_surface, const Ppath &ppath, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Numeric &dza, const Verbosity &verbosity)
Defocusing for arbitrary geometry (zenith angle part only)
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 interp_atmfield_gp2itw(Matrix &itw, const Index &atmosphere_dim, const ArrayOfGridPos &gp_p, const ArrayOfGridPos &gp_lat, const ArrayOfGridPos &gp_lon)
Converts atmospheric grid positions to weights for interpolation of an atmospheric field.
void itw2p(VectorView p_values, ConstVectorView p_grid, const ArrayOfGridPos &gp, ConstMatrixView itw)
Converts interpolation weights to pressures.
void interp_cloudfield_gp2itw(VectorView itw, GridPos &gp_p_out, GridPos &gp_lat_out, GridPos &gp_lon_out, const GridPos &gp_p_in, const GridPos &gp_lat_in, const GridPos &gp_lon_in, const Index &atmosphere_dim, const ArrayOfIndex &cloudbox_limits)
Converts atmospheric a grid position to weights for interpolation of a field defined ONLY inside the ...
void interp_atmfield_by_itw(VectorView x, const Index &atmosphere_dim, ConstTensor3View x_field, const ArrayOfGridPos &gp_p, const ArrayOfGridPos &gp_lat, const ArrayOfGridPos &gp_lon, ConstMatrixView itw)
Interpolates an atmospheric field with pre-calculated weights by interp_atmfield_gp2itw.
Header file for special_interp.cc.
EnergyLevelMap InterpToGridPos(Index atmosphere_dim, const ArrayOfGridPos &p, const ArrayOfGridPos &lat, const ArrayOfGridPos &lon) const
Structure to store a grid position.
The structure to describe a propagation path and releated quantities.
Matrix los
Line-of-sight at each ppath point.
ArrayOfGridPos gp_lon
Index position with respect to the longitude grid.
String background
Radiative background.
Index np
Number of points describing the ppath.
Matrix pos
The distance between start pos and the last position in pos.
ArrayOfGridPos gp_lat
Index position with respect to the latitude grid.
Numeric end_lstep
The distance between end pos and the first position in pos.
Vector start_pos
Start position.
Vector lstep
The length between ppath points.
Numeric start_lstep
Length between sensor and atmospheric boundary.
Numeric constant
The propagation path constant (only used for 1D)
Vector r
Radius of each ppath point.
ArrayOfGridPos gp_p
Index position with respect to the pressure grid.
Index dim
Atmospheric dimensionality.
Vector end_pos
End position.
Vector start_los
Start line-of-sight.
Vector nreal
The real part of the refractive index at each path position.
Vector end_los
End line-of-sight.