69 const Index& stokes_dim,
71 const Index& atmosphere_dim,
83 const Index& cloudbox_on,
89 const Index& scat_data_checked,
91 const Index& jacobian_do,
94 const Matrix& iy_transmitter,
95 const Agenda& propmat_clearsky_agenda,
96 const Agenda& water_p_eq_agenda,
98 const Index& trans_in_jacobian,
100 const Index& t_interp_order,
101 const Verbosity& verbosity [[maybe_unused]]) {
103 const Index j_analytical_do = jacobian_do ?
104 do_analytical_jacobian<1>(jacobian_quantities) : 0;
108 const Index ns = stokes_dim;
111 const Index nq = j_analytical_do ? jacobian_quantities.
nelem() : 0;
120 "ppath.background is invalid. Check your "
121 "calculation of *ppath*?");
123 "The propagation path ends inside or at boundary of "
124 "the cloudbox.\nFor this method, *ppath* must be "
125 "calculated in this way:\n ppathCalc( cloudbox_on = 0 ).");
128 "The scattering data must be flagged to have\n"
129 "passed a consistency check (scat_data_checked=1).");
131 "*pnd_field* and *scat_data* inconsistent regarding total number of"
132 " scattering elements.");
137 "revision before safe to use.");
140 "*dpnd_field_dx* not properly initialized:\n"
141 "Number of elements in dpnd_field_dx must be equal number of jacobian"
142 " quantities.\n(Note: jacobians have to be defined BEFORE *pnd_field*"
143 " is calculated/set.");
150 "The size of *iy* returned from *iy_transmitter_agenda* is\n"
152 " expected size = [", nf,
",", stokes_dim,
"]\n"
153 " size of iy = [", iy_transmitter.
nrows(),
",", iy_transmitter.
ncols(),
"]\n")
154 for (
Index iv = 0; iv < nf; iv++)
156 "The *iy* returned from *iy_transmitter_agenda* "
157 "must have the value 1 in the first column.");
171 diy_dx = j_analytical_do ?
184 Index auxBackScat = -1;
185 Index auxAbSpAtte = -1;
186 Index auxPartAtte = -1;
188 for (
Index i = 0; i < naux; i++) {
189 iy_aux[i].resize(nf * np, ns);
191 if (iy_aux_vars[i] ==
"Radiative background") {
194 }
else if( iy_aux_vars[i] ==
"Backscattering" ) {
197 }
else if( iy_aux_vars[i] ==
"Abs species extinction" ) {
200 }
else if( iy_aux_vars[i] ==
"Particle extinction" ) {
205 "The only allowed strings in *iy_aux_vars* are:\n"
206 " \"Radiative background\"\n"
207 " \"Backscattering\"\n"
210 "but you have selected: \"", iy_aux_vars[i],
"\"")
228 Tensor5 Pe(ne, np, nf, ns, ns, 0);
234 const Index nf_ssd = scat_data[0][0].pha_mat_data.nlibraries();
235 const Index duplicate_freqs = ((nf == nf_ssd) ? 0 : 1);
236 Tensor6 pha_mat_1se(nf_ssd, 1, 1, 1, ns, ns);
237 Vector t_ok(1), t_array(1);
238 Matrix mlos_sca(1, 2), mlos_inc(1, 2);
242 if (np == 1 && rbi == 1) {
272 ppvar_f, ppath, f_grid, atmosphere_dim, rte_alonglos_v, ppvar_wind);
285 clear2cloudy.resize(np);
286 for (
Index ip = 0; ip < np; ip++)
287 clear2cloudy[ip] = -1;
301 Index temperature_derivative_position = -1;
304 if (trans_in_jacobian && j_analytical_do) {
305 dK_this_dx.resize(nq);
306 dK_past_dx.resize(nq);
317 if (jacobian_quantities[iq] == Jacobian::Atm::Temperature) {
318 temperature_derivative_position = iq;
319 do_hse = jacobian_quantities[iq].Subtag() ==
"HSE on";
324 for (
Index ip = 0; ip < np; ip++) {
331 propmat_clearsky_agenda,
334 ppvar_mag(
joker, ip),
337 ppvar_vmr(
joker, ip),
340 trans_in_jacobian && j_analytical_do);
342 if (trans_in_jacobian && j_analytical_do)
351 trans_in_jacobian && j_analytical_do);
353 if (clear2cloudy[ip] + 1) {
364 ppvar_t[
Range(ip, 1)],
366 trans_in_jacobian && jacobian_do);
368 if (
abs(pext_scaling - 1) > 1e-6) {
370 if (trans_in_jacobian && j_analytical_do) {
376 if (auxAbSpAtte >= 0) {
377 for (
Index iv = 0; iv < nf; iv++) {
378 iy_aux[auxAbSpAtte](iv*np+ip,
joker) = K_this.Kjj()[iv];
381 if (auxPartAtte >= 0) {
382 for (
Index iv = 0; iv < nf; iv++) {
383 iy_aux[auxPartAtte](iv*np+ip,
joker) = Kp.
Kjj()[iv];
389 if (trans_in_jacobian && j_analytical_do)
397 mlos_sca(0,
joker) = los_sca;
401 if (atmosphere_dim == 3) {
406 mlos_inc(0,
joker) = los_inc;
409 for (
Index i_ss = 0; i_ss < scat_data.
nelem(); i_ss++)
410 for (
Index i_se = 0; i_se < scat_data[i_ss].
nelem(); i_se++) {
414 if (ppvar_pnd(i_se_flat, ip) != 0)
416 else if (j_analytical_do)
417 for (
Index iq = 0; iq < nq && !val_pnd; iq++)
418 if (jac_scat_i[iq] >= 0)
419 if (ppvar_dpnd_dx[iq](i_se_flat, ip) != 0) {
427 scat_data[i_ss][i_se],
428 ppvar_t[
Range(ip, 1)],
433 if (t_ok[0] not_eq 0)
435 for (
Index iv = 0; iv < nf; iv++)
443 "Interpolation error for (flat-array) scattering"
444 " element #", i_se_flat,
"\n"
445 "at location/temperature point #", ip,
"\n")
456 do_hse ? ppath.
lstep[ip - 1] / (2.0 * ppvar_t[ip - 1]) : 0;
458 do_hse ? ppath.
lstep[ip - 1] / (2.0 * ppvar_t[ip]) : 0;
469 temperature_derivative_position);
472 swap(K_past, K_this);
473 swap(dK_past_dx, dK_this_dx);
486 lvl_rad[0] = iy_transmitter;
488 rad_inc = iy_transmitter;
506 for (
Index ip = 0; ip < np; ip++) {
507 for (
Index iv = 0; iv < nf; iv++) {
508 for (
Index is = 0; is < stokes_dim; is++) {
509 iy(iv*np+ip, is) = lvl_rad[ip](iv, is);
510 if (j_analytical_do) {
512 diy_dpath[iq](ip, iv*np+ip2, is) =
513 dlvl_rad[ip][ip2][iq](iv, is););
520 if (j_analytical_do) {
521 const Index iy_agenda_call1 = 1;
522 const Tensor3 iy_transmittance(0, 0, 0);
543 if (auxBackScat >= 0) {
544 for (
Index ip = 0; ip < np; ip++) {
545 for (
Index iv = 0; iv < nf; iv++) {
565 const Index& atmgeom_checked,
566 const Index& atmfields_checked,
567 const String& iy_unit_radar,
569 const Index& stokes_dim,
571 const Index& cloudbox_on,
572 const Index& cloudbox_checked,
575 const Index& sensor_checked,
576 const Index& jacobian_do,
578 const Agenda& iy_radar_agenda,
605 "The atmospheric fields must be flagged to have "
606 "passed a consistency check (atmfields_checked=1).");
608 "The atmospheric geometry must be flagged to have "
609 "passed a consistency check (atmgeom_checked=1).");
611 "The cloudbox must be flagged to have "
612 "passed a consistency check (cloudbox_checked=1).");
614 "The sensor variables must be flagged to have "
615 "passed a consistency check (sensor_checked=1).");
618 bool is_z =
max(range_bins) > 1;
620 "The vector *range_bins* must contain strictly "
621 "increasing values.");
623 "The vector *range_bins* is not allowed to contain "
626 "The main length of *instrument_pol_array* must match "
627 "the number of frequencies.");
632 const Numeric jfac = 10 * log10(exp(1.0));
633 if (iy_unit_radar ==
"1") {
634 }
else if (iy_unit_radar ==
"Ze") {
635 ze_cfac(cfac, f_grid, ze_tref, k2);
636 }
else if (iy_unit_radar ==
"dBZe") {
637 ze_cfac(cfac, f_grid, ze_tref, k2);
638 ze_min =
pow(10.0, dbze_min / 10);
641 "For this method, *iy_unit_radar* must be set to \"1\", "
642 "\"Ze\" or \"dBZe\".");
653 for (
Index i = 0; i < nf; i++) {
654 const Index ni = instrument_pol_array[i].
nelem();
655 npolcum[i + 1] = npolcum[i] + ni;
657 for (
Index j = 0; j < ni; j++) {
658 W[i][j].resize(stokes_dim);
659 stokes2pol(W[i][j], stokes_dim, instrument_pol_array[i][j], 0.5);
665 const Index ntot = npos * npolcum[nf] * nbins;
677 for (
Index i = 0; i < naux; i++) {
678 y_aux[i].resize(ntot);
684 Index j_analytical_do = 0;
693 jacobian_indices[jacobian_indices.
nelem() - 1][1] + 1);
695 njq = jacobian_quantities.
nelem();
707 for (
Index p = 0; p < npos; p++) {
726 sensor_pos(p,
joker),
727 sensor_los(p,
joker),
733 "A path consisting of a single point found. "
734 "This is not allowed.");
737 "The size of *iy* returned from *iy_radar_agenda* "
738 "is not correct (for this method).");
746 for (
Index i = 1; i < np; i++) {
747 range[i] = range[i - 1] + ppath.
lstep[i - 1] *
752 const Numeric range_end1 =
min(range[0], range[np - 1]);
753 const Numeric range_end2 =
max(range[0], range[np - 1]);
757 if (!(range_bins[
b] >= range_end2 ||
758 range_bins[
b + 1] <= range_end1))
762 Numeric blim2 =
min(range_bins[
b + 1], range_end2);
769 hbin /= (blim2 - blim1);
771 for (
Index iv = 0; iv < nf; iv++) {
775 if (j_analytical_do) {
787 if (j_analytical_do) {
792 auxvar[
a].resize(np);
795 for (
Index ip = 0; ip < instrument_pol_array[iv].
nelem(); ip++) {
797 mult(refl, I, W[iv][ip]);
798 if (j_analytical_do) {
800 k < drefl[iq].nrows();
806 if (iy_aux_vars[
a] ==
"Backscattering") {
807 mult(auxvar[
a], A[
a], W[iv][ip]);
809 for (
Index j = 0; j < np; j++) {
810 auxvar[
a][j] = A[
a](j, 0);
816 Index iout = nbins * (p * npolcum[nf] + npolcum[iv] + ip) +
b;
817 y[iout] = cfac[iv] * (hbin * refl);
819 if (j_analytical_do) {
821 for (
Index k = 0; k < drefl[iq].nrows(); k++) {
822 jacobian(iout, jacobian_indices[iq][0] + k) =
823 cfac[iv] * (hbin * drefl[iq](k,
joker));
824 if (iy_unit_radar ==
"dBZe") {
825 jacobian(iout, jacobian_indices[iq][0] + k) *=
826 jfac /
max(y[iout], ze_min);
831 if (iy_unit_radar ==
"dBZe") {
832 y[iout] = y[iout] <= ze_min ? dbze_min : 10 * log10(y[iout]);
837 if (iy_aux_vars[
a] ==
"Backscattering") {
838 y_aux[
a][iout] = cfac[iv] * (hbin * auxvar[
a]);
839 if (iy_unit_radar ==
"dBZe") {
840 y_aux[
a][iout] = y_aux[
a][iout] <= ze_min
842 : 10 * log10(y_aux[
a][iout]);
845 y_aux[
a][iout] = hbin * auxvar[
a];
850 if (geo_pos.
nelem()) y_geo(iout,
joker) = geo_pos;
859 for (
Index iv = 0; iv < nf; iv++) {
860 for (
Index ip = 0; ip < instrument_pol_array[iv].
nelem(); ip++) {
861 const Index iout = nbins * (p * npolcum[nf] + npolcum[iv] + ip) +
b;
862 y_f[iout] = f_grid[iv];
863 y_pol[iout] = instrument_pol_array[iv][ip];
877 Tensor4& particle_bulkprop_field,
879 const Index& atmosphere_dim,
887 const Index& atmfields_checked,
888 const Index& atmgeom_checked,
890 const Agenda& propmat_clearsky_agenda,
897 const Index& fill_clutter,
901 const Index& do_atten_abs,
902 const Index& do_atten_hyd,
903 const Numeric& atten_hyd_scaling,
913 "The atmospheric fields must be flagged to have\n"
914 "passed a consistency check (atmfields_checked=1).");
916 "The atmospheric geometry must be flagged to have\n"
917 "passed a consistency check (atmgeom_checked=1).");
921 "*dbze_noise* not covered by invtable[0]." );
923 "*dbze_noise* not covered by invtable[1]." );
924 chk_atm_field(
"GIN reflectivities", dBZe, atmosphere_dim, p_grid,
933 if (h_clutter.
nrows() > 1 || h_clutter.
ncols() > 1) {
934 chk_atm_surface(
"GIN h_clutter", h_clutter, atmosphere_dim, lat_grid, lon_grid);
935 hclutterm = h_clutter;
937 hclutterm.
resize(nlat, nlon);
938 hclutterm = h_clutter(0,0);
942 particle_bulkprop_field.
resize(2, np, nlat, nlon);
943 particle_bulkprop_field = 0;
944 particle_bulkprop_names = scat_species;
949 const Numeric extrap_fac = 100;
956#pragma omp parallel for if (!arts_omp_in_parallel() && nlat + nlon > 2) \
957 firstprivate(wss) collapse(2)
958 for (
Index ilat = 0; ilat < nlat; ilat++) {
959 for (
Index ilon = 0; ilon < nlon; ilon++) {
960 if (fail_msg.
nelem() != 0)
continue;
963 if (
max(dBZe(
joker, ilat, ilon)) > dbze_noise) {
966 Numeric k_part_above = 0, k_abs_above = 0;
970 for (
Index ip = np - 1; ip >= 0; ip--) {
972 if (z_field(ip, ilat, ilon) >= z_surface(ilat, ilon) +
973 hclutterm(ilat, ilon)) {
976 dBZe(ip, ilat, ilon) + dbze_corr_abs + dbze_corr_hyd;
979 const Numeric t = t_field(ip, ilat, ilon);
982 Index phase = -1, it = -1;
983 if (dbze > dbze_noise) {
985 phase = t >= t_phase ? 0 : 1;
993 else if (t >= invtable[phase].get_numeric_grid(
1000 it = gp.
fd[0] < 0.5 ? gp.
idx : gp.
idx + 1;
1007 if (dBZe(ip, ilat, ilon) > dbze_noise && do_atten_hyd &&
1008 dbze_corr_hyd < atten_hyd_max) {
1012 "Max 20 dB extrapolation allowed. Found a value\n"
1013 "exceeding this limit.");
1022 interp(itw, invtable[phase].data(1,
joker, it), gp);
1025 0.5 * hfac * (k_part_above + k_this) *
1026 (z_field(ip + 1, ilat, ilon) - z_field(ip, ilat, ilon));
1031 dbze_corr_hyd =
min(dbze_corr_hyd, atten_hyd_max);
1033 k_part_above = k_this;
1059 t_field(ip, ilat, ilon),
1060 rtp_nlte_local_dummy,
1061 vmr_field(
joker, ip, ilat, ilon),
1062 propmat_clearsky_agenda);
1063 k_this = propmat.
Kjj()[0];
1066 0.5 * hfac * (k_abs_above + k_this) *
1067 (z_field(ip + 1, ilat, ilon) - z_field(ip, ilat, ilon));
1071 k_abs_above = k_this;
1076 if (dBZe(ip, ilat, ilon) > dbze_noise) {
1078 dbze = dBZe(ip, ilat, ilon) + dbze_corr_abs + dbze_corr_hyd;
1083 "Max 20 dB extrapolation allowed. Found a value\n"
1084 "exceeding this limit.");
1094 (
interp(itw, invtable[phase].data(0,
joker, it), gp)));
1098 else if (wc > wc_clip)
1100 particle_bulkprop_field(phase, ip, ilat, ilon) = wc;
1106 particle_bulkprop_field(0, ip, ilat, ilon) =
1107 particle_bulkprop_field(0, ip + 1, ilat, ilon);
1108 particle_bulkprop_field(1, ip, ilat, ilon) =
1109 particle_bulkprop_field(1, ip + 1, ilat, ilon);
1115 }
catch (
const std::exception& e) {
1117 os <<
"Run-time error at ilat: " << ilat <<
" ilon: " << ilon <<
":\n"
1119#pragma omp critical(particle_bulkpropRadarOnionPeeling_latlon)
1120 fail_msg.push_back(os.str());
1125 if (fail_msg.
nelem()) {
1127 for (
const auto& msg : fail_msg) os << msg <<
'\n';
1142 const Index& i_species,
1155 const Index iss = i_species;
1159 "First element in T_grid of scattering species is "
1160 "below 220 K.\nThat seems to be too low for liquid.\n"
1161 "Please note that the onion peeling function expects "
1162 "rain\nas the first scattering element.");
1164 "This method requires that *f_grid* has length 1.");
1166 "*i_species* must either be 0 or 1.");
1168 "*scat_data* must contain data for exactly two "
1169 "scattering species.");
1171 "*scat_data* and *scat_species* are inconsistent in size.");
1173 "*scat_data* and *scat_meta* are inconsistent in size.");
1175 "*scat_data* and *scat_meta* have inconsistent sizes.");
1177 "*scat_data* should just contain one frequency.");
1179 "The PSD applied must be of 1-moment type.");
1182 if (invtable.empty())
1184 invtable[iss].set_name(
"Radar inversion table created by *RadarOnionPeelingTableCalc*");
1185 invtable[iss].resize(2, ndb, nt);
1186 invtable[iss].data = 0;
1187 invtable[iss].set_grid_name(0,
"Radiative properties");
1188 invtable[iss].set_grid(0, {
"Log of water content",
"Extinction"});
1189 invtable[iss].set_grid_name(1,
"Radar reflectivity");
1190 invtable[iss].set_grid(1, dbze_grid);
1191 invtable[iss].set_grid_name(2,
"Temperature");
1192 invtable[iss].set_grid(2, t_grid);
1200 for (
Index i=0; i<nse; i++) {
1202 "So far only TRO is handled by this method.");
1203 const Index ib = scat_data[iss][i].za_grid.
nelem() - 1;
1205 "All za_grid in scat_data must end with 180.\n"
1206 "This is not the case for scat_data[", iss,
"][",
1208 gridpos(gp, scat_data[iss][i].T_grid, t_grid, 1);
1210 interp(e(
joker,i), itw, scat_data[iss][i].ext_mat_data(0,
joker,0,0,0), gp);
1211 interp(
b(
joker,i), itw, scat_data[iss][i].pha_mat_data(0,
joker,ib,0,0,0,0), gp);
1224 const Vector& pnd_agenda_input_t = t_grid;
1225 Matrix pnd_agenda_input(nt, 1);
1229 ze_cfac(cfac, f_grid, ze_tref, k2);
1233 pnd_agenda_input = wc_grid[
w];
1242 pnd_agenda_array_input_names[iss],
1247 for (
Index t=0; t<nt; t++) {
1248 for (
Index i=0; i<nse; i++) {
1250 "A negative back-scattering found for scat_species ", iss,
1251 ",\ntemperature ", t_grid[t],
"K and element ", i)
1253 "A negative extinction found for scat_species ", iss,
1254 ",\ntemperature ", t_grid[t],
"K and element ", i)
1256 "A negative PSD value found for scat_species ", iss,
1257 ",\ntemperature ", t_grid[t],
"K and ", wc_grid[
w],
1259 D(0,
w,t) += pnd_data(t,i) *
b(t,i);
1260 D(1,
w,t) += pnd_data(t,i) * e(t,i);
1263 D(0,
w,t) = 10 * log10(cfac[0] * D(0,
w,t));
1273 for (
Index t=0; t<nt; t++) {
1276 cout << wc_grid[
w] <<
" " << D(0,
w,t) << endl;
1279 "A case found of non-increasing dBZe.\n"
1280 "Found for scat_species ", iss,
" and ", t_grid[t],
"K.")
1282 if (D(0,0,t) > dbze_grid[0]) {
1284 cout << wc_grid[
w] <<
" " << D(0,
w,t) << endl;
1287 "A case found where start of dbze_grid not covered.\n"
1288 "Found for scat_species ", iss,
" and ", t_grid[t],
"K.")
1290 if (D(0,nwc-1,t) < dbze_grid[ndb-1]) {
1292 cout << wc_grid[
w] <<
" " << D(0,
w,t) << endl;
1295 "A case found where end of dbze_grid not covered.\n"
1296 "Found for scat_species ", iss,
" and ", t_grid[t],
"K.")
1301 interp(invtable[iss].data(0,
joker,t), itw, wc_log, gp);
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.
The global header file for ARTS.
Constants of physical expressions as constexpr.
Header file for helper functions for OpenMP.
void pnd_agenda_arrayExecute(Workspace &ws, Matrix &pnd_data, Tensor3 &dpnd_data_dx, const Index agenda_array_index, const Vector &pnd_agenda_input_t, const Matrix &pnd_agenda_input, const ArrayOfString &pnd_agenda_input_names, const ArrayOfString &dpnd_data_dx_names, const ArrayOfAgenda &input_agenda_array)
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 iy_radar_agendaExecute(Workspace &ws, Matrix &iy, ArrayOfMatrix &iy_aux, Ppath &ppath, ArrayOfTensor3 &diy_dx, Vector &geo_pos, const ArrayOfString &iy_aux_vars, const Index iy_id, const Index cloudbox_on, const Index jacobian_do, const Vector &rte_pos, const Vector &rte_los, const Agenda &input_agenda)
This can be used to make arrays out of anything.
Index nelem() const ARTS_NOEXCEPT
Index nrows() const noexcept
Index ncols() const noexcept
Index npages() const
Returns the number of pages.
Index nrows() const
Returns the number of rows.
Index ncols() const
Returns the number of columns.
Index nbooks() const noexcept
Index nelem() const noexcept
Returns the number of elements.
bool empty() const noexcept
Returns true if variable size is zero.
void resize(Index r, Index c)
Resize function.
VectorView Kjj(const Index iz=0, const Index ia=0)
Vector view to diagonal elements.
Stokes vector is as Propagation matrix but only has 4 possible values.
void resize(Index b, Index p, Index r, Index c)
Resize function.
void resize(Index n)
Resize function.
void mult(MatrixView C, ConstMatrixView A, const Block &B)
#define ARTS_USER_ERROR(...)
#define ARTS_USER_ERROR_IF(condition,...)
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.
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.
void jac_ranges_indices(ArrayOfArrayOfIndex &jis, bool &any_affine, const ArrayOfRetrievalQuantity &jqs, const bool &before_affine)
Determines the index range inside x and the Jacobian for each retrieval quantity.
#define FOR_ANALYTICAL_JACOBIANS_DO(what_to_do)
Array< RetrievalQuantity > ArrayOfRetrievalQuantity
bool is_increasing(ConstVectorView x)
Checks if a vector is sorted and strictly increasing.
Header file for logic.cc.
void VectorLogSpace(Vector &x, const Numeric &start, const Numeric &stop, const Numeric &step, const Verbosity &verbosity)
WORKSPACE METHOD: VectorLogSpace.
void yRadar(Workspace &ws, Vector &y, Vector &y_f, ArrayOfIndex &y_pol, Matrix &y_pos, Matrix &y_los, ArrayOfVector &y_aux, Matrix &y_geo, Matrix &jacobian, const Index &atmgeom_checked, const Index &atmfields_checked, const String &iy_unit_radar, const ArrayOfString &iy_aux_vars, const Index &stokes_dim, const Vector &f_grid, const Index &cloudbox_on, const Index &cloudbox_checked, const Matrix &sensor_pos, const Matrix &sensor_los, const Index &sensor_checked, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Agenda &iy_radar_agenda, const ArrayOfArrayOfIndex &instrument_pol_array, const Vector &range_bins, const Numeric &ze_tref, const Numeric &k2, const Numeric &dbze_min, const Verbosity &)
WORKSPACE METHOD: yRadar.
const Index GFIELD3_T_GRID
const Index GFIELD3_DB_GRID
constexpr Numeric SPEED_OF_LIGHT
void RadarOnionPeelingTableCalc(Workspace &ws, ArrayOfGriddedField3 &invtable, const Vector &f_grid, const ArrayOfString &scat_species, const ArrayOfArrayOfSingleScatteringData &scat_data, const ArrayOfArrayOfScatteringMetaData &scat_meta, const ArrayOfAgenda &pnd_agenda_array, const ArrayOfArrayOfString &pnd_agenda_array_input_names, const Index &i_species, const Vector &dbze_grid, const Vector &t_grid, const Numeric &wc_min, const Numeric &wc_max, const Numeric &ze_tref, const Numeric &k2, const Verbosity &verbosity)
WORKSPACE METHOD: RadarOnionPeelingTableCalc.
void particle_bulkpropRadarOnionPeeling(Workspace &ws, Tensor4 &particle_bulkprop_field, ArrayOfString &particle_bulkprop_names, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &t_field, const Tensor3 &z_field, const Tensor4 &vmr_field, const Matrix &z_surface, const Index &atmfields_checked, const Index &atmgeom_checked, const Vector &f_grid, const Agenda &propmat_clearsky_agenda, const ArrayOfString &scat_species, const ArrayOfGriddedField3 &invtable, const Matrix &incangles, const Tensor3 &dBZe, const Numeric &dbze_noise, const Matrix &h_clutter, const Index &fill_clutter, const Numeric &t_phase, const Numeric &wc_max, const Numeric &wc_clip, const Index &do_atten_abs, const Index &do_atten_hyd, const Numeric &atten_hyd_scaling, const Numeric &atten_hyd_max, const Verbosity &)
WORKSPACE METHOD: particle_bulkpropRadarOnionPeeling.
constexpr Numeric LOG10_EULER_NUMBER
void iyRadarSingleScat(Workspace &ws, Matrix &iy, ArrayOfMatrix &iy_aux, ArrayOfTensor3 &diy_dx, Vector &ppvar_p, Vector &ppvar_t, Matrix &ppvar_vmr, Matrix &ppvar_wind, Matrix &ppvar_mag, Matrix &ppvar_pnd, Matrix &ppvar_f, 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 Tensor4 &pnd_field, const ArrayOfTensor4 &dpnd_field_dx, const ArrayOfString &scat_species, const ArrayOfArrayOfSingleScatteringData &scat_data, const Index &scat_data_checked, 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 Numeric &rte_alonglos_v, const Index &trans_in_jacobian, const Numeric &pext_scaling, const Index &t_interp_order, const Verbosity &verbosity)
WORKSPACE METHOD: iyRadarSingleScat.
Numeric last(ConstVectorView x)
last
Array< Tensor3 > ArrayOfTensor3
An array of Tensor3.
void abs(Sparse &A, const Sparse &B)
Absolute value of sparse matrix elements.
void transform(VectorView y, double(&my_func)(double), ConstVectorView x)
A generic transform function for vectors, which can be used to implement mathematical functions opera...
NUMERIC Numeric
The type to use for all floating point numbers.
INDEX Index
The type to use for all integer numbers and indices.
void swap(ComplexVector &v1, ComplexVector &v2) noexcept
Swaps two objects.
Declarations having to do with the four output streams.
constexpr Numeric log10_euler
Ten's logarithm of Euler's number.
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...
auto cosd(auto x) noexcept
Returns the cosine of the deg2rad of the input.
auto row_vec(matpack::vector auto &&x)
Map the input to a non-owning const-correct Eigen Map representing a column vector.
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 error_if_limb_ppath(const Ppath &ppath)
Throws an error if ppath altitudes not are strictly increasing or decreasing.
Index ppath_what_background(const Ppath &ppath)
Returns the case number for the radiative background.
Stuff related to the propagation matrix.
Numeric pow(const Rational base, Numeric exp)
Power of.
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 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 ConstVectorView &ppath_f_grid, const ConstVectorView &ppath_magnetic_field, const ConstVectorView &ppath_line_of_sight, const EnergyLevelMap &ppath_nlte, const ConstVectorView &ppath_vmrs, const Numeric &ppath_temperature, const Numeric &ppath_pressure, const bool &jacobian_do)
Gets the clearsky propgation matrix and NLTE contributions.
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 mirror_los(Vector &los_mirrored, const ConstVectorView &los, const Index &atmosphere_dim)
Determines the backward direction for a given line-of-sight.
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 integration_bin_by_vecmult(VectorView h, ConstVectorView x_g_in, const Numeric &limit1, const Numeric &limit2)
integration_bin_by_vecmult
void stokes2pol(VectorView w, const Index &stokes_dim, const Index &ipol_1based, const Numeric nv)
stokes2pol
Sensor modelling functions.
Structure to store a grid position.
std::array< Numeric, 2 > fd
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.
Matrix pos
The distance between start pos and the last position in pos.
Numeric end_lstep
The distance between end pos and the first position in pos.
Vector lstep
The length between ppath points.
Vector ngroup
The group index of refraction.
Radiation Vector for Stokes dimension 1-4.
Eigen::VectorXd Vec(size_t i) const
Return Vector at position by copy.
Class to keep track of Transmission Matrices for Stokes Dim 1-4.
ArrayOfTransmissionMatrix bulk_backscatter(const ConstTensor5View &Pe, const ConstMatrixView &pnd)
Bulk back-scattering.
ArrayOfArrayOfTransmissionMatrix bulk_backscatter_derivative(const ConstTensor5View &Pe, const ArrayOfMatrix &dpnd_dx)
Derivatives of bulk back-scattering
ArrayOfTransmissionMatrix cumulative_transmission(const ArrayOfTransmissionMatrix &T, const CumulativeTransmission type)
Accumulate the transmission matrix over all layers.
void set_backscatter_radiation_vector(ArrayOfRadiationVector &I, ArrayOfArrayOfArrayOfRadiationVector &dI, const RadiationVector &I_incoming, const ArrayOfTransmissionMatrix &T, const ArrayOfTransmissionMatrix &PiTf, const ArrayOfTransmissionMatrix &PiTr, const ArrayOfTransmissionMatrix &Z, const ArrayOfArrayOfTransmissionMatrix &dT1, const ArrayOfArrayOfTransmissionMatrix &dT2, const ArrayOfArrayOfTransmissionMatrix &dZ, const BackscatterSolver solver)
Set the backscatter 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.
@ CommutativeTransmission
Array< ArrayOfRadiationVector > ArrayOfArrayOfRadiationVector
Array< RadiationVector > ArrayOfRadiationVector
Array< TransmissionMatrix > ArrayOfTransmissionMatrix