35#include "matpack_data.h"
53 Index& doit_za_grid_size,
57 const Index& N_za_grid,
58 const Index& N_aa_grid,
59 const String& za_grid_opt_file,
67 "N_aa_grid must be > 0 (even for 1D / DISORT cases).")
75 "N_za_grid must be >= 0.")
76 doit_za_grid_size = N_za_grid;
78 if (za_grid_opt_file ==
"")
83 "N_za_grid must be >1 or =0 (the latter only allowed for RT4).")
92 Index& doit_conv_flag,
93 Index& doit_iteration_counter,
94 Tensor6& cloudbox_field_mono,
96 const Tensor6& cloudbox_field_mono_old,
98 const Vector& epsilon,
99 const Index& max_iterations,
100 const Index& throw_nonconv_error,
107 "Convergence flag is non-zero, which means that this\n"
108 "WSM is not used correctly. *doit_conv_flagAbs* should\n"
109 "be used only in *doit_conv_test_agenda*\n");
111 const Index N_p = cloudbox_field_mono.nvitrines();
112 const Index N_lat = cloudbox_field_mono.nshelves();
113 const Index N_lon = cloudbox_field_mono.nbooks();
114 const Index N_za = cloudbox_field_mono.npages();
115 const Index N_aa = cloudbox_field_mono.nrows();
116 const Index stokes_dim = cloudbox_field_mono.ncols();
120 "You have to specify limiting values for the "
121 "convergence test for each Stokes component "
122 "separately. That means that *epsilon* must "
123 "have *stokes_dim* elements!");
127 cloudbox_field_mono_old, N_p, N_lat, N_lon, N_za, N_aa, stokes_dim),
128 "The fields (Tensor6) *cloudbox_field* and \n"
129 "*cloudbox_field_old* which are compared in the \n"
130 "convergence test do not have the same size.\n");
134 doit_iteration_counter += 1;
135 out2 <<
" Number of DOIT iteration: " << doit_iteration_counter <<
"\n";
137 if (doit_iteration_counter > max_iterations) {
139 out <<
"Method does not converge (number of iterations \n"
140 <<
"is > " << max_iterations <<
"). Either the cloud "
141 <<
"particle number density \n"
142 <<
"is too large or the numerical setup for the DOIT \n"
143 <<
"calculation is not correct. In case of limb \n"
144 <<
"simulations please make sure that you use an \n"
145 <<
"optimized zenith angle grid. \n"
146 <<
"*cloudbox_field* might be wrong.\n";
147 if (throw_nonconv_error != 0) {
153 out1 <<
"Warning in DOIT calculation (output set to NaN):\n" << out.str();
154 cloudbox_field_mono = NAN;
157 out1 <<
"Warning in DOIT calculation (output equals current status):\n"
162 for (Index p_index = 0; p_index < N_p; p_index++) {
163 for (Index lat_index = 0; lat_index < N_lat; lat_index++) {
164 for (Index lon_index = 0; lon_index < N_lon; lon_index++) {
165 for (Index za_index = 0; za_index < N_za; za_index++) {
166 for (Index aa_index = 0; aa_index < N_aa; aa_index++) {
167 for (Index stokes_index = 0; stokes_index < stokes_dim;
169 Numeric diff = (cloudbox_field_mono(p_index,
175 cloudbox_field_mono_old(p_index,
186 if (abs(diff) > epsilon[stokes_index]) {
187 out1 <<
"difference: " << diff <<
"\n";
205 Index& doit_conv_flag,
206 Index& doit_iteration_counter,
207 Tensor6& cloudbox_field_mono,
209 const Tensor6& cloudbox_field_mono_old,
210 const Vector& f_grid,
211 const Index& f_index,
213 const Vector& epsilon,
214 const Index& max_iterations,
215 const Index& throw_nonconv_error,
223 "Convergence flag is non-zero, which means that this \n"
224 "WSM is not used correctly. *doit_conv_flagAbs* should\n"
225 "be used only in *doit_conv_test_agenda*\n");
227 const Index N_p = cloudbox_field_mono.nvitrines();
228 const Index N_lat = cloudbox_field_mono.nshelves();
229 const Index N_lon = cloudbox_field_mono.nbooks();
230 const Index N_za = cloudbox_field_mono.npages();
231 const Index N_aa = cloudbox_field_mono.nrows();
232 const Index stokes_dim = cloudbox_field_mono.ncols();
236 "You have to specify limiting values for the "
237 "convergence test for each Stokes component "
238 "separately. That means that *epsilon* must "
239 "have *stokes_dim* elements!");
243 cloudbox_field_mono_old, N_p, N_lat, N_lon, N_za, N_aa, stokes_dim),
244 "The fields (Tensor6) *cloudbox_field* and \n"
245 "*cloudbox_field_old* which are compared in the \n"
246 "convergence test do not have the same size.\n");
255 "*f_index* is greater than number of elements in the\n"
256 "frequency grid.\n");
260 doit_iteration_counter += 1;
261 out2 <<
" Number of DOIT iteration: " << doit_iteration_counter <<
"\n";
264 if (doit_iteration_counter > max_iterations) {
266 out <<
"At frequency " << f_grid[f_index] <<
" GHz \n"
267 <<
"method does not converge (number of iterations \n"
268 <<
"is > " << max_iterations <<
"). Either the particle"
269 <<
" number density \n"
270 <<
"is too large or the numerical setup for the DOIT \n"
271 <<
"calculation is not correct. In case of limb \n"
272 <<
"simulations please make sure that you use an \n"
273 <<
"optimized zenith angle grid. \n"
274 <<
"*cloudbox_field* might be wrong.\n";
275 if (throw_nonconv_error != 0) {
281 out1 <<
"Warning in DOIT calculation (output set to NaN):\n" << out.str();
282 cloudbox_field_mono = NAN;
285 out1 <<
"Warning in DOIT calculation (output equals current status):\n"
290 for (Index p_index = 0; p_index < N_p; p_index++) {
291 for (Index lat_index = 0; lat_index < N_lat; lat_index++) {
292 for (Index lon_index = 0; lon_index < N_lon; lon_index++) {
293 for (Index za_index = 0; za_index < N_za; za_index++) {
294 for (Index aa_index = 0; aa_index < N_aa; aa_index++) {
295 for (Index stokes_index = 0; stokes_index < stokes_dim;
297 Numeric diff = cloudbox_field_mono(p_index,
303 cloudbox_field_mono_old(p_index,
313 Numeric diff_bt =
invrayjean(diff, f_grid[f_index]);
316 if (abs(diff_bt) > epsilon[stokes_index]) {
317 out1 <<
"BT difference: " << diff_bt <<
" in stokes dim "
318 << stokes_index <<
"\n";
340 Index& doit_conv_flag,
341 Index& doit_iteration_counter,
342 Tensor6& cloudbox_field_mono,
344 const Tensor6& cloudbox_field_mono_old,
345 const Vector& f_grid,
346 const Index& f_index,
348 const Vector& epsilon,
349 const Index& max_iterations,
350 const Index& throw_nonconv_error,
358 "Convergence flag is non-zero, which means that this \n"
359 "WSM is not used correctly. *doit_conv_flagAbs* should\n"
360 "be used only in *doit_conv_test_agenda*\n");
362 const Index N_p = cloudbox_field_mono.nvitrines();
363 const Index N_lat = cloudbox_field_mono.nshelves();
364 const Index N_lon = cloudbox_field_mono.nbooks();
365 const Index N_za = cloudbox_field_mono.npages();
366 const Index N_aa = cloudbox_field_mono.nrows();
367 const Index stokes_dim = cloudbox_field_mono.ncols();
371 "You have to specify limiting values for the "
372 "convergence test for each Stokes component "
373 "separately. That means that *epsilon* must "
374 "have *stokes_dim* elements!");
378 cloudbox_field_mono_old, N_p, N_lat, N_lon, N_za, N_aa, stokes_dim),
379 "The fields (Tensor6) *cloudbox_field* and \n"
380 "*cloudbox_field_old* which are compared in the \n"
381 "convergence test do not have the same size.\n");
390 "*f_index* is greater than number of elements in the\n"
391 "frequency grid.\n");
395 doit_iteration_counter += 1;
396 out2 <<
" Number of DOIT iteration: " << doit_iteration_counter <<
"\n";
398 if (doit_iteration_counter > max_iterations) {
400 out <<
"Method does not converge (number of iterations \n"
401 <<
"is > " << max_iterations <<
"). Either the"
402 <<
" particle number density \n"
403 <<
"is too large or the numerical setup for the DOIT \n"
404 <<
"calculation is not correct. In case of limb \n"
405 <<
"simulations please make sure that you use an \n"
406 <<
"optimized zenith angle grid. \n";
407 if (throw_nonconv_error != 0) {
413 out1 <<
"Warning in DOIT calculation (output set to NaN):\n" << out.str();
414 cloudbox_field_mono = NAN;
417 out1 <<
"Warning in DOIT calculation (output equals current status):\n"
426 for (Index i = 0; i < epsilon.nelem(); i++) {
427 for (Index p_index = 0; p_index < N_p; p_index++) {
428 for (Index lat_index = 0; lat_index < N_lat; lat_index++) {
429 for (Index lon_index = 0; lon_index < N_lon; lon_index++) {
430 for (Index za_index = 0; za_index < N_za; za_index++) {
431 for (Index aa_index = 0; aa_index < N_aa; aa_index++) {
434 p_index, lat_index, lon_index, za_index, aa_index, i) -
435 cloudbox_field_mono_old(p_index,
448 lqs[i] = sqrt(lqs[i]);
449 lqs[i] /= (Numeric)(N_p * N_lat * N_lon * N_za * N_aa);
454 if (lqs[i] >= epsilon[i]) doit_conv_flag = 0;
457 out1 <<
"lqs [I]: " << lqs[0] <<
"\n";
464 Tensor6& cloudbox_field_mono,
467 const Agenda& doit_scat_field_agenda,
468 const Agenda& doit_rte_agenda,
469 const Agenda& doit_conv_test_agenda,
470 const Index& accelerated,
477 chk_not_empty(
"doit_scat_field_agenda", doit_scat_field_agenda);
479 chk_not_empty(
"doit_conv_test_agenda", doit_conv_test_agenda);
481 for (Index
v = 0;
v < cloudbox_field_mono.nvitrines();
v++)
482 for (Index s = 0; s < cloudbox_field_mono.nshelves(); s++)
483 for (Index
b = 0;
b < cloudbox_field_mono.nbooks();
b++)
484 for (Index p = 0; p < cloudbox_field_mono.npages(); p++)
485 for (Index r = 0; r < cloudbox_field_mono.nrows(); r++)
486 for (Index
c = 0;
c < cloudbox_field_mono.ncols();
c++)
488 "*cloudbox_field_mono* contains at least one NaN value.\n"
489 "This indicates an improper initialization of *cloudbox_field*.");
496 Tensor6 cloudbox_field_mono_old_local;
497 Index doit_conv_flag_local;
498 Index doit_iteration_counter_local;
502 Tensor6 doit_scat_field_local(cloudbox_field_mono.nvitrines(),
503 cloudbox_field_mono.nshelves(),
504 cloudbox_field_mono.nbooks(),
505 cloudbox_field_mono.npages(),
506 cloudbox_field_mono.nrows(),
507 cloudbox_field_mono.ncols(),
510 doit_conv_flag_local = 0;
511 doit_iteration_counter_local = 0;
513 ArrayOfTensor6 acceleration_input;
515 acceleration_input.resize(4);
517 while (doit_conv_flag_local == 0) {
519 cloudbox_field_mono_old_local = cloudbox_field_mono;
524 out2 <<
" Execute doit_scat_field_agenda. \n";
526 ws, doit_scat_field_local, cloudbox_field_mono, doit_scat_field_agenda);
529 out2 <<
" Execute doit_rte_agenda. \n";
531 ws, cloudbox_field_mono, doit_scat_field_local, doit_rte_agenda);
535 doit_conv_flag_local,
536 doit_iteration_counter_local,
538 cloudbox_field_mono_old_local,
539 doit_conv_test_agenda);
542 if (accelerated > 0 && doit_conv_flag_local == 0) {
543 acceleration_input[(doit_iteration_counter_local - 1) % 4] =
546 if (doit_iteration_counter_local % 4 == 0) {
548 cloudbox_field_mono, acceleration_input, accelerated, verbosity);
558 Tensor6& cloudbox_field_mono,
560 const Tensor6& doit_scat_field,
563 const Agenda& propmat_clearsky_agenda,
564 const Tensor4& vmr_field,
566 const Agenda& spt_calc_agenda,
567 const Vector& za_grid,
568 const Tensor4& pnd_field,
570 const Agenda& ppath_step_agenda,
571 const Numeric& ppath_lmax,
572 const Numeric& ppath_lraytrace,
573 const Vector& p_grid,
574 const Tensor3& z_field,
575 const Vector& refellipsoid,
577 const Tensor3& t_field,
578 const Vector& f_grid,
579 const Index& f_index,
580 const Agenda& surface_rtprop_agenda,
581 const Index& doit_za_interp,
587 <<
" cloudbox_fieldUpdate1D: Radiative transfer calculation in cloudbox\n";
588 out2 <<
" ------------------------------------------------------------- \n";
597 "The cloudbox dimension is not 1D! \n"
598 "Do you really want to do a 1D calculation? \n"
599 "If not, use *cloudbox_fieldUpdateSeq3D*.\n");
602 const Index N_scat_za = za_grid.nelem();
605 "The range of *za_grid* must [0 180].");
608 "The length of *p_grid* must be >= 2.");
611 chk_size(
"z_field", z_field, p_grid.nelem(), 1, 1);
612 chk_size(
"t_field", t_field, p_grid.nelem(), 1, 1);
621 "*f_index* is greater than number of elements in the\n"
622 "frequency grid.\n");
625 "Interpolation method is not defined. Use \n"
626 "*doit_za_interpSet*.\n");
628 const Index stokes_dim = doit_scat_field.ncols();
633 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
641 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
655 out3 <<
"Calculate optical properties of individual scattering elements\n";
665 Tensor5 ext_mat_field(cloudbox_limits[1] - cloudbox_limits[0] + 1,
671 Tensor4 abs_vec_field(
672 cloudbox_limits[1] - cloudbox_limits[0] + 1, 1, 1, stokes_dim, 0.);
674 Tensor6 cloudbox_field_old(cloudbox_field_mono);
677 Index aa_index_local = 0;
680 for (Index za_index_local = 0; za_index_local < N_scat_za; za_index_local++) {
699 for (Index p_index = cloudbox_limits[0]; p_index <= cloudbox_limits[1];
701 if ((p_index != 0) || (za_grid[za_index_local] <= 90.)) {
710 propmat_clearsky_agenda,
723 surface_rtprop_agenda,
735 Tensor6& cloudbox_field_mono,
736 Tensor6& doit_scat_field,
740 const Agenda& propmat_clearsky_agenda,
741 const Tensor4& vmr_field,
743 const Agenda& spt_calc_agenda,
744 const Vector& za_grid,
745 const Vector& aa_grid,
746 const Tensor4& pnd_field,
748 const Agenda& ppath_step_agenda,
749 const Numeric& ppath_lmax,
750 const Numeric& ppath_lraytrace,
751 const Vector& p_grid,
752 const Tensor3& z_field,
753 const Vector& refellipsoid,
755 const Tensor3& t_field,
756 const Vector& f_grid,
757 const Index& f_index,
758 const Agenda& surface_rtprop_agenda,
759 const Index& doit_za_interp,
760 const Index& normalize,
761 const Numeric& norm_error_threshold,
762 const Index& norm_debug,
768 <<
" cloudbox_fieldUpdateSeq1D: Radiative transfer calculation in cloudbox\n";
769 out2 <<
" ------------------------------------------------------------- \n";
774 chk_not_empty(
"propmat_clearsky_agenda", propmat_clearsky_agenda);
779 "The cloudbox dimension is not 1D! \n"
780 "Do you really want to do a 1D calculation? \n"
781 "For 3D use *cloudbox_fieldUpdateSeq3D*.\n");
784 const Index N_scat_za = za_grid.nelem();
787 "The range of *za_grid* must [0 180].");
790 "The length of *p_grid* must be >= 2.");
793 chk_size(
"z_field", z_field, p_grid.nelem(), 1, 1);
794 chk_size(
"t_field", t_field, p_grid.nelem(), 1, 1);
803 "*f_index* is greater than number of elements in the\n"
804 "frequency grid.\n");
807 "Interpolation method is not defined. Use \n"
808 "*doit_za_interpSet*.\n");
810 const Index stokes_dim = doit_scat_field.ncols();
815 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
823 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
837 out3 <<
"Calculate optical properties of individual scattering elements\n";
847 Tensor5 ext_mat_field(cloudbox_limits[1] - cloudbox_limits[0] + 1,
853 Tensor4 abs_vec_field(
854 cloudbox_limits[1] - cloudbox_limits[0] + 1, 1, 1, stokes_dim, 0.);
859 180. - asin((refellipsoid[0] + z_field(cloudbox_limits[0], 0, 0)) /
860 (refellipsoid[0] + z_field(cloudbox_limits[1], 0, 0))) *
870 Matrix cloudbox_field_limb;
873 Index aa_index_local = 0;
876 Tensor4 si, sei, si_corr;
887 norm_error_threshold,
893 for (Index za_index_local = 0; za_index_local < N_scat_za; za_index_local++) {
913 if (za_grid[za_index_local] <= 90.) {
919 for (Index p_index = cloudbox_limits[1] - 1;
920 p_index >= cloudbox_limits[0];
929 propmat_clearsky_agenda,
942 surface_rtprop_agenda,
946 }
else if (za_grid[za_index_local] >= theta_lim) {
950 for (Index p_index = cloudbox_limits[0] + 1;
951 p_index <= cloudbox_limits[1];
960 propmat_clearsky_agenda,
973 surface_rtprop_agenda,
987 bool conv_flag =
false;
989 while (!conv_flag && limb_it < 10) {
991 cloudbox_field_limb =
992 cloudbox_field_mono(joker, 0, 0, za_index_local, 0, joker);
993 for (Index p_index = cloudbox_limits[0]; p_index <= cloudbox_limits[1];
1001 cloudbox_field_mono,
1007 propmat_clearsky_agenda,
1020 surface_rtprop_agenda,
1027 for (Index p_index = 0;
1028 conv_flag && p_index < cloudbox_field_mono.nvitrines();
1030 for (Index stokes_index = 0; conv_flag && stokes_index < stokes_dim;
1032 Numeric diff = cloudbox_field_mono(
1033 p_index, 0, 0, za_index_local, 0, stokes_index) -
1034 cloudbox_field_limb(p_index, stokes_index);
1039 Numeric diff_bt =
invrayjean(diff, f_grid[f_index]);
1040 if (abs(diff_bt) > epsilon[stokes_index]) {
1041 out2 <<
"Limb BT difference: " << diff_bt <<
" in stokes dim "
1042 << stokes_index <<
"\n";
1048 out2 <<
"Limb iterations: " << limb_it <<
"\n";
1057 Tensor6& cloudbox_field_mono,
1059 const Tensor6& doit_scat_field,
1062 const Agenda& propmat_clearsky_agenda,
1063 const Tensor4& vmr_field,
1065 const Agenda& spt_calc_agenda,
1066 const Vector& za_grid,
1067 const Vector& aa_grid,
1068 const Tensor4& pnd_field,
1070 const Agenda& ppath_step_agenda,
1071 const Numeric& ppath_lmax,
1072 const Numeric& ppath_lraytrace,
1073 const Vector& p_grid,
1074 const Vector& lat_grid,
1075 const Vector& lon_grid,
1076 const Tensor3& z_field,
1077 const Vector& refellipsoid,
1079 const Tensor3& t_field,
1080 const Vector& f_grid,
1081 const Index& f_index,
1082 const Index& doit_za_interp,
1088 <<
" cloudbox_fieldUpdateSeq3D: Radiative transfer calculatiuon in cloudbox.\n";
1089 out2 <<
" ------------------------------------------------------------- \n";
1094 chk_not_empty(
"propmat_clearsky_agenda", propmat_clearsky_agenda);
1099 "The cloudbox dimension is not 3D! \n"
1100 "Do you really want to do a 3D calculation? \n"
1101 "For 1D use *cloudbox_fieldUpdateSeq1D*.\n");
1104 const Index N_scat_za = za_grid.nelem();
1107 "The range of *za_grid* must [0 180].");
1110 const Index N_scat_aa = aa_grid.nelem();
1113 "The range of *za_grid* must [0 360].");
1120 "z_field", z_field, p_grid.nelem(), lat_grid.nelem(), lon_grid.nelem());
1122 "t_field", t_field, p_grid.nelem(), lat_grid.nelem(), lon_grid.nelem());
1131 "*f_index* is greater than number of elements in the\n"
1132 "frequency grid.\n");
1135 "Interpolation method is not defined. Use \n"
1136 "*doit_za_interpSet*.\n");
1138 const Index stokes_dim = doit_scat_field.ncols();
1143 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
1144 (cloudbox_limits[3] - cloudbox_limits[2]) + 1,
1145 (cloudbox_limits[5] - cloudbox_limits[4]) + 1,
1151 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
1152 (cloudbox_limits[3] - cloudbox_limits[2]) + 1,
1153 (cloudbox_limits[5] - cloudbox_limits[4]) + 1,
1165 out3 <<
"Calculate optical properties of individual scattering elements\n";
1175 const Index p_low = cloudbox_limits[0];
1176 const Index p_up = cloudbox_limits[1];
1177 const Index lat_low = cloudbox_limits[2];
1178 const Index lat_up = cloudbox_limits[3];
1179 const Index lon_low = cloudbox_limits[4];
1180 const Index lon_up = cloudbox_limits[5];
1184 Tensor5 ext_mat_field(p_up - p_low + 1,
1185 lat_up - lat_low + 1,
1186 lon_up - lon_low + 1,
1190 Tensor4 abs_vec_field(p_up - p_low + 1,
1191 lat_up - lat_low + 1,
1192 lon_up - lon_low + 1,
1197 for (Index za_index = 0; za_index < N_scat_za; za_index++) {
1200 for (Index aa_index = 1; aa_index < N_scat_aa; aa_index++) {
1219 Vector stokes_vec(stokes_dim, 0.);
1221 Numeric theta_lim = 180. - asin((refellipsoid[0] + z_field(p_low, 0, 0)) /
1222 (refellipsoid[0] + z_field(p_up, 0, 0))) *
1226 if (za_grid[za_index] <= 90.) {
1232 for (Index p_index = p_up - 1; p_index >= p_low; p_index--) {
1233 for (Index lat_index = lat_low; lat_index <= lat_up; lat_index++) {
1234 for (Index lon_index = lon_low; lon_index <= lon_up; lon_index++) {
1236 cloudbox_field_mono,
1246 propmat_clearsky_agenda,
1267 else if (za_grid[za_index] > theta_lim) {
1271 for (Index p_index = p_low + 1; p_index <= p_up; p_index++) {
1272 for (Index lat_index = lat_low; lat_index <= lat_up; lat_index++) {
1273 for (Index lon_index = lon_low; lon_index <= lon_up; lon_index++) {
1275 cloudbox_field_mono,
1285 propmat_clearsky_agenda,
1314 else if (za_grid[za_index] > 90. && za_grid[za_index] < theta_lim) {
1315 for (Index p_index = p_low; p_index <= p_up; p_index++) {
1320 if (!(p_index == 0 && za_grid[za_index] > 90.)) {
1321 for (Index lat_index = lat_low; lat_index <= lat_up; lat_index++) {
1322 for (Index lon_index = lon_low; lon_index <= lon_up;
1325 cloudbox_field_mono,
1335 propmat_clearsky_agenda,
1360 cloudbox_field_mono(joker, joker, joker, joker, 0, joker) =
1361 cloudbox_field_mono(joker, joker, joker, joker, N_scat_aa - 1, joker);
1368 Tensor6& cloudbox_field_mono,
1372 const Tensor6& doit_scat_field,
1375 const Agenda& propmat_clearsky_agenda,
1376 const Tensor4& vmr_field,
1378 const Agenda& spt_calc_agenda,
1379 const Vector& za_grid,
1380 const Tensor4& pnd_field,
1382 const Vector& p_grid,
1383 const Tensor3& z_field,
1385 const Tensor3& t_field,
1386 const Vector& f_grid,
1387 const Index& f_index,
1393 <<
" cloudbox_fieldUpdateSeq1DPP: Radiative transfer calculation in cloudbox.\n";
1395 <<
" --------------------------------------------------------------------- \n";
1397 const Index stokes_dim = doit_scat_field.ncols();
1403 "The dimension of stokes vector must be"
1407 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
1415 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
1428 const Index N_scat_za = za_grid.nelem();
1433 out3 <<
"Calculate optical properties of individual scattering elements\n";
1444 Tensor5 ext_mat_field(cloudbox_limits[1] - cloudbox_limits[0] + 1,
1450 Tensor4 abs_vec_field(
1451 cloudbox_limits[1] - cloudbox_limits[0] + 1, 1, 1, stokes_dim, 0.);
1454 for (za_index = 0; za_index < N_scat_za; za_index++) {
1473 Vector stokes_vec(stokes_dim, 0.);
1475 if (za_grid[za_index] <= 90) {
1484 for (Index p_index = cloudbox_limits[1] - 1;
1485 p_index >= cloudbox_limits[0];
1488 cloudbox_field_mono,
1494 propmat_clearsky_agenda,
1505 }
else if (za_grid[za_index] > 90) {
1509 for (Index p_index = cloudbox_limits[0] + 1;
1510 p_index <= cloudbox_limits[1];
1513 cloudbox_field_mono,
1519 propmat_clearsky_agenda,
1537 Tensor6& doit_scat_field,
1538 Tensor7& cloudbox_field,
1539 Index& doit_is_initialized,
1541 const Index& stokes_dim,
1542 const Index& atmosphere_dim,
1543 const Vector& f_grid,
1544 const Vector& za_grid,
1545 const Vector& aa_grid,
1546 const Index& doit_za_grid_size,
1547 const Index& cloudbox_on,
1552 doit_is_initialized = 0;
1553 out0 <<
" Cloudbox is off, DOIT calculation will be skipped.\n";
1562 "The dimension of stokes vector must be"
1568 const Index N_scat_za = za_grid.nelem();
1574 "For accurate results, za_grid must have "
1575 "more than 15 elements.");
1576 if (N_scat_za > 100) {
1578 out1 <<
"Warning: za_grid is very large, which means that the \n"
1579 <<
"calculation will be very slow.\n";
1583 "The range of *za_grid* must [0 180].");
1586 "*za_grid* must be increasing.");
1589 const Index N_scat_aa = aa_grid.nelem();
1592 "For accurate results, aa_grid must have "
1593 "more than 5 elements.");
1594 if (N_scat_aa > 100) {
1596 out1 <<
"Warning: aa_grid is very large which means that the \n"
1597 <<
"calculation will be very slow.\n";
1601 "The range of *aa_grid* must [0 360].");
1604 "*doit_za_grid_size* must be greater than 15 for accurate results");
1605 if (doit_za_grid_size > 100) {
1607 out1 <<
"Warning: doit_za_grid_size is very large which means that the \n"
1608 <<
"calculation will be very slow.\n";
1612 "*cloudbox_limits* is a vector which contains the"
1613 "upper and lower limit of the cloud for all "
1614 "atmospheric dimensions. So its dimension must"
1615 "be 2 x *atmosphere_dim*");
1619 const Index Nf = f_grid.nelem();
1620 const Index Np_cloud = cloudbox_limits[1] - cloudbox_limits[0] + 1;
1621 const Index Nza = za_grid.
nelem();
1622 const Index Ns = stokes_dim;
1625 if (atmosphere_dim == 1) {
1626 cloudbox_field.resize(Nf, Np_cloud, 1, 1, Nza, 1, Ns);
1627 doit_scat_field.resize(Np_cloud, 1, 1, Nza, 1, Ns);
1628 }
else if (atmosphere_dim == 3) {
1629 const Index Nlat_cloud = cloudbox_limits[3] - cloudbox_limits[2] + 1;
1630 const Index Nlon_cloud = cloudbox_limits[5] - cloudbox_limits[4] + 1;
1631 const Index Naa = aa_grid.
nelem();
1633 cloudbox_field.resize(Nf, Np_cloud, Nlat_cloud, Nlon_cloud, Nza, Naa, Ns);
1634 doit_scat_field.resize(Np_cloud, Nlat_cloud, Nlon_cloud, Nza, Naa, Ns);
1637 "Scattering calculations are not possible for a 2D"
1638 "atmosphere. If you want to do scattering calculations"
1639 "*atmosphere_dim* has to be either 1 or 3");
1642 cloudbox_field = NAN;
1643 doit_scat_field = NAN;
1644 doit_is_initialized = 1;
1649 Tensor6& cloudbox_field_mono,
1650 const Vector& p_grid_orig,
1651 const Vector& p_grid,
1654 Tensor6 cloudbox_field_mono_opt(
1655 p_grid_orig.nelem(), 1, 1, cloudbox_field_mono.npages(), 1, 1);
1657 Matrix itw_z(Z_gp.
nelem(), 2);
1660 Vector p_grid_cloudbox {p_grid[Range(
1661 cloudbox_limits[0], cloudbox_limits[1] - cloudbox_limits[0] + 1)]};
1663 os <<
"There is a problem with the pressure grid interpolation";
1667 gridpos(Z_gp, p_grid_cloudbox, p_grid_orig);
1671 for (Index i = 0; i < cloudbox_field_mono.npages(); i++) {
1672 interp(cloudbox_field_mono_opt(joker, 0, 0, i, 0, 0),
1674 cloudbox_field_mono(joker, 0, 0, i, 0, 0),
1677 cloudbox_field_mono = cloudbox_field_mono_opt;
1690 Tensor6& cloudbox_field_mono,
1691 Tensor7& pha_mat_doit,
1693 Vector& p_grid_orig,
1694 const Vector& f_grid,
1695 const Index& f_index,
1696 const Agenda& propmat_clearsky_agenda,
1697 const Numeric& tau_scat_max,
1698 const Numeric& sgl_alb_max,
1699 const Index& cloudbox_size_max,
1705 " This method currently only works for unpolarized radiation "
1706 "( stokes_dim = 1)");
1710 pnd_field.nbooks() != cloudbox_limits[1] - cloudbox_limits[0] + 1,
1711 " ScatSpeciesMerge has to be applied in order to use this method");
1713 bool was_too_much =
false;
1714 Numeric tau_scat_max_internal = tau_scat_max;
1716 p_grid_orig = p_grid;
1717 vector<Numeric> z_grid_new;
1718 Vector ext_mat(cloudbox_limits[1] - cloudbox_limits[0] + 1);
1719 Vector abs_vec(cloudbox_limits[1] - cloudbox_limits[0] + 1);
1720 Vector scat_vec(cloudbox_limits[1] - cloudbox_limits[0] + 1);
1722 z_grid_new.reserve(1000);
1723 for (Index i = 0; i < cloudbox_limits[0]; i++)
1724 z_grid_new.push_back(z_field(i, 0, 0));
1730 const Vector rtp_mag_dummy(3, 0);
1731 const Vector ppath_los_dummy;
1732 StokesVector nlte_dummy;
1733 ArrayOfPropagationMatrix partial_dummy;
1734 ArrayOfStokesVector partial_nlte_dummy;
1736 PropagationMatrix cur_propmat_clearsky;
1738 Index scat_data_insert_offset = 0;
1739 Vector single_scattering_albedo(cloudbox_limits[1] - cloudbox_limits[0] + 1,
1741 for (Index k = cloudbox_limits[0]; k < cloudbox_limits[1] + 1; ++k) {
1742 Index cloudbox_index = k - cloudbox_limits[0];
1743 Numeric abs_coeff = 0;
1745 ext_mat[cloudbox_index] =
1746 scat_data_local[cloudbox_index].ext_mat_data(0, 0, 0, 0, 0);
1747 abs_vec[cloudbox_index] =
1748 scat_data_local[cloudbox_index].abs_vec_data(0, 0, 0, 0, 0);
1749 scat_vec[cloudbox_index] =
1750 ext_mat[cloudbox_index] - abs_vec[cloudbox_index];
1753 cur_propmat_clearsky,
1759 Vector{f_grid[Range(f_index, 1)]},
1765 Vector{vmr_field(joker, k, 0, 0)},
1766 propmat_clearsky_agenda);
1767 abs_coeff += cur_propmat_clearsky.Kjj()[0];
1768 abs_coeff /= (Numeric)vmr_field.nbooks();
1769 single_scattering_albedo[cloudbox_index] =
1770 scat_vec[cloudbox_index] / (ext_mat[cloudbox_index] + abs_coeff);
1776 scat_data_insert_offset = 0;
1777 for (Index k = cloudbox_limits[0]; k < cloudbox_limits[1]; ++k) {
1778 Index cloudbox_index = k - cloudbox_limits[0];
1779 const Numeric sgl_alb = (single_scattering_albedo[cloudbox_index] +
1780 single_scattering_albedo[cloudbox_index + 1]) /
1782 const Numeric scat_opt_thk =
1783 (z_field(k + 1, 0, 0) - z_field(k, 0, 0)) *
1784 ((scat_vec[cloudbox_index] + scat_vec[cloudbox_index + 1]) / 2);
1785 if (scat_opt_thk > tau_scat_max_internal && sgl_alb > sgl_alb_max) {
1786 Index factor = (Index)ceil(scat_opt_thk / tau_scat_max_internal);
1787 for (Index j = 1; j < factor; j++) {
1788 scat_data_insert_offset++;
1793 if (scat_data_insert_offset + cloudbox_limits[1] - cloudbox_limits[0] + 1 >
1794 cloudbox_size_max) {
1795 tau_scat_max_internal += 0.01;
1796 was_too_much =
true;
1798 }
while (scat_data_insert_offset + cloudbox_limits[1] - cloudbox_limits[0] +
1801 scat_data_insert_offset = 0;
1805 <<
"Warning: Pressure grid optimization with the thresholds tau_scat_max = "
1806 << tau_scat_max <<
" and sgl_slb_max = " << sgl_alb_max
1807 <<
" has lead to an enhancement larger than the value of"
1808 <<
" cloudbox_size_max = " << cloudbox_size_max
1809 <<
". This is why I changed tau_scat_max to " << tau_scat_max_internal
1810 <<
". This may lead to errors of a too coarse grid! \n";
1816 for (Index k = cloudbox_limits[0]; k < cloudbox_limits[1]; ++k) {
1817 Index cloudbox_index = k - cloudbox_limits[0];
1818 const Numeric sgl_alb = (single_scattering_albedo[cloudbox_index] +
1819 single_scattering_albedo[cloudbox_index + 1]) /
1821 const Numeric scat_opt_thk =
1822 (z_field(k + 1, 0, 0) - z_field(k, 0, 0)) *
1823 ((scat_vec[cloudbox_index] + scat_vec[cloudbox_index + 1]) / 2);
1824 z_grid_new.push_back(z_field(k, 0, 0));
1826 if (scat_opt_thk > tau_scat_max_internal && sgl_alb > sgl_alb_max) {
1827 Index factor = (Index)ceil(scat_opt_thk / tau_scat_max_internal);
1829 (z_field(k + 1, 0, 0) - z_field(k, 0, 0)) / (Numeric)factor;
1831 scat_data_local[cloudbox_index + scat_data_insert_offset + 1];
1833 scat_data_local[cloudbox_index + scat_data_insert_offset];
1835 for (Index j = 1; j < factor; j++) {
1836 z_grid_new.push_back(z_field(k, 0, 0) + (Numeric)j * step);
1838 const Numeric weight = (Numeric)j / (Numeric)factor;
1840 Tensor7 weightednextLayerPhamat = nextLayer.
pha_mat_data;
1841 Tensor5 weightednextLayerExtmat = nextLayer.
ext_mat_data;
1842 Tensor5 weightednextLayerAbsvec = nextLayer.
abs_vec_data;
1844 weightednextLayerPhamat *= weight;
1845 weightednextLayerExtmat *= weight;
1846 weightednextLayerAbsvec *= weight;
1857 scat_data_local.insert(scat_data_local.begin() + cloudbox_index +
1858 scat_data_insert_offset + 1,
1859 std::move(newLayer));
1861 scat_data_insert_offset++;
1866 cloudbox_limits_opt[0] = cloudbox_limits[0];
1867 cloudbox_limits_opt[1] = scat_data_insert_offset + cloudbox_limits[1];
1868 const Index cloudbox_opt_size =
1869 cloudbox_limits_opt[1] - cloudbox_limits_opt[0] + 1;
1870 out3 <<
"Frequency: " << f_grid[f_index]
1871 <<
" old: " << cloudbox_limits[1] - cloudbox_limits[0] + 1
1872 <<
" new: " << cloudbox_opt_size <<
" Factor: "
1873 << (Numeric)cloudbox_opt_size /
1874 (Numeric)(cloudbox_limits[1] - cloudbox_limits[0] + 1)
1877 for (Index i = cloudbox_limits[1]; i < z_field.npages(); i++)
1878 z_grid_new.push_back(z_field(i, 0, 0));
1880 Vector z_grid(z_grid_new.size());
1881 for (Index i = 0; i < z_grid.nelem(); i++) z_grid[i] = z_grid_new[i];
1882 p_grid_orig = p_grid[Range(cloudbox_limits[0],
1883 cloudbox_limits[1] - cloudbox_limits[0] + 1)];
1888 Tensor3 t_field_new(z_grid.nelem(), 1, 1);
1889 Vector p_grid_opt(z_grid.nelem());
1890 Tensor6 cloudbox_field_mono_opt(
1891 cloudbox_opt_size, 1, 1, cloudbox_field_mono.npages(), 1, 1);
1892 Tensor7 pha_mat_doit_opt(cloudbox_opt_size,
1893 pha_mat_doit.nvitrines(),
1895 pha_mat_doit.nbooks(),
1896 pha_mat_doit.npages(),
1900 Matrix itw_z(Z_gp.
nelem(), 2);
1902 os <<
"At the current frequency " << f_grid[f_index]
1903 <<
" there was an error while interpolating the fields to the new z_field";
1907 gridpos(Z_gp, z_field(joker, 0, 0), z_grid);
1912 interp(t_field_new(joker, 0, 0), itw_z, t_field(joker, 0, 0), Z_gp);
1914 for (Index k = cloudbox_limits_opt[0]; k < cloudbox_limits_opt[1]; k++) {
1915 Index i = k - cloudbox_limits[0];
1916 scat_data_local[i].T_grid = t_field_new(i, 0, 0);
1920 interp(p_grid_opt, itw_z, p_grid, Z_gp);
1923 Tensor4 vmr_field_opt(vmr_field.nbooks(), p_grid_opt.nelem(), 1, 1);
1924 for (Index i = 0; i < vmr_field.nbooks(); i++)
1926 vmr_field_opt(i, joker, 0, 0), itw_z, vmr_field(i, joker, 0, 0), Z_gp);
1930 Matrix itw_z_2(Z_gp_2.
nelem(), 2);
1932 Range(cloudbox_limits[0], cloudbox_limits[1] - cloudbox_limits[0] + 1);
1933 Range r2 = Range(cloudbox_limits_opt[0], cloudbox_opt_size);
1936 z_field(Range(cloudbox_limits[0],
1937 cloudbox_limits[1] - cloudbox_limits[0] + 1),
1940 z_grid[Range(cloudbox_limits_opt[0], cloudbox_opt_size)]);
1943 for (Index i = 0; i < cloudbox_field_mono.npages(); i++) {
1944 interp(cloudbox_field_mono_opt(joker, 0, 0, i, 0, 0),
1946 cloudbox_field_mono(joker, 0, 0, i, 0, 0),
1949 for (Index i = 0; i < pha_mat_doit.nvitrines(); i++) {
1950 for (Index j = 0; j < pha_mat_doit.nbooks(); j++) {
1951 for (Index k = 0; k < pha_mat_doit.npages(); k++) {
1952 interp(pha_mat_doit_opt(joker, i, 0, j, k, 0, 0),
1954 pha_mat_doit(joker, i, 0, j, k, 0, 0),
1961 pnd_field.resize(cloudbox_opt_size, cloudbox_opt_size, 1, 1);
1963 for (Index i = 0; i < cloudbox_opt_size; i++) pnd_field(i, i, 0, 0) = 1.;
1966 p_grid = p_grid_opt;
1967 t_field = t_field_new;
1968 cloudbox_limits = cloudbox_limits_opt;
1969 cloudbox_field_mono = cloudbox_field_mono_opt;
1970 pha_mat_doit = pha_mat_doit_opt;
1971 z_field.resize(z_grid.nelem(), 1, 1);
1972 z_field(joker, 0, 0) = z_grid;
1973 vmr_field = vmr_field_opt;
1978 const Index& doit_iteration_counter,
1979 const Tensor6& cloudbox_field_mono,
1980 const Index& f_index,
1985 if (!frequencies.
nelem() || !iterations.
nelem())
return;
1991 os <<
"doit_iteration_f" << f_index <<
"_i" << doit_iteration_counter;
1994 if (frequencies[0] == -1 && iterations[0] == -1) {
1996 os.str() +
".xml", cloudbox_field_mono,
FILE_TYPE_ASCII, 0, verbosity);
1999 for (Index j = 0; j < frequencies.
nelem(); j++) {
2000 if (f_index == frequencies[j] || (!j && frequencies[j] == -1)) {
2002 if (iterations[0] == -1) {
2004 cloudbox_field_mono,
2012 for (Index i = 0; i < iterations.
nelem(); i++) {
2013 if (doit_iteration_counter == iterations[i])
2015 cloudbox_field_mono,
2028 Tensor6& doit_scat_field,
2030 const Agenda& pha_mat_spt_agenda,
2031 const Tensor6& cloudbox_field_mono,
2032 const Tensor4& pnd_field,
2033 const Tensor3& t_field,
2034 const Index& atmosphere_dim,
2036 const Vector& za_grid,
2037 const Vector& aa_grid,
2038 const Index& doit_za_grid_size,
2039 const Tensor7& pha_mat_doit,
2052 const Index Nza = za_grid.nelem();
2055 "The range of *za_grid* must [0 180].");
2058 const Index Naa = aa_grid.nelem();
2061 "The range of *aa_grid* must [0 360].");
2064 const Index stokes_dim = doit_scat_field.ncols();
2075 if (atmosphere_dim == 1) {
2077 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
2084 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
2090 }
else if (atmosphere_dim == 3) {
2092 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
2093 (cloudbox_limits[3] - cloudbox_limits[2]) + 1,
2094 (cloudbox_limits[5] - cloudbox_limits[4]) + 1,
2099 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
2100 (cloudbox_limits[3] - cloudbox_limits[2]) + 1,
2101 (cloudbox_limits[5] - cloudbox_limits[4]) + 1,
2107 "The atmospheric dimension must be 1D or 3D \n"
2108 "for scattering calculations using the DOIT\n"
2109 "module, but it is not. The value of *atmosphere_dim*\n"
2110 "is ", atmosphere_dim,
".")
2114 "*cloudbox_limits* is a vector which contains the"
2115 "upper and lower limit of the cloud for all "
2116 "atmospheric dimensions. So its dimension must"
2117 "be 2 x *atmosphere_dim*");
2122 "The zenith angle grids for the computation of\n"
2123 "the scattering integral and the RT part must \n"
2124 "be equal. Check definitions in \n"
2125 "*DOAngularGridsSet*. The keyword \n"
2126 "'za_grid_opt_file' should be empty. \n");
2131 Tensor4 pha_mat_local(
2132 doit_za_grid_size, aa_grid.nelem(), stokes_dim, stokes_dim, 0.);
2134 Tensor5 pha_mat_spt_local(pnd_field.nbooks(),
2142 Vector grid_stepsize(2);
2143 grid_stepsize[0] = 180. / (Numeric)(doit_za_grid_size - 1);
2145 grid_stepsize[1] = 360. / (Numeric)(Naa - 1);
2148 Tensor3 product_field(Nza, Naa, stokes_dim, 0);
2150 out2 <<
" Calculate the scattered field\n";
2152 if (atmosphere_dim == 1) {
2155 for (Index p_index = 0; p_index <= cloudbox_limits[1] - cloudbox_limits[0];
2158 for (Index za_index_local = 0; za_index_local < Nza; za_index_local++) {
2160 out3 <<
"Multiplication of phase matrix with incoming"
2161 <<
" intensities \n";
2167 for (Index za_in = 0; za_in < Nza; ++za_in) {
2168 for (Index aa_in = 0; aa_in < Naa; ++aa_in) {
2172 for (Index i = 0; i < stokes_dim; i++) {
2173 for (Index j = 0; j < stokes_dim; j++) {
2174 product_field(za_in, aa_in, i) +=
2176 p_index, za_index_local, 0, za_in, aa_in, i, j) *
2177 cloudbox_field_mono(p_index, 0, 0, za_in, 0, j);
2186 for (Index i = 0; i < stokes_dim; i++) {
2187 doit_scat_field(p_index, 0, 0, za_index_local, 0, i) =
2192 for (Index i = 0; i < stokes_dim; i++) {
2193 doit_scat_field(p_index, 0, 0, za_index_local, 0, i) =
2205 else if (atmosphere_dim == 3) {
2210 for (Index p_index = 0; p_index <= cloudbox_limits[1] - cloudbox_limits[0];
2212 for (Index lat_index = 0;
2213 lat_index <= cloudbox_limits[3] - cloudbox_limits[2];
2215 for (Index lon_index = 0;
2216 lon_index <= cloudbox_limits[5] - cloudbox_limits[4];
2218 Numeric rtp_temperature_local =
2219 t_field(p_index + cloudbox_limits[0],
2220 lat_index + cloudbox_limits[2],
2221 lon_index + cloudbox_limits[4]);
2223 for (Index aa_index_local = 1; aa_index_local < Naa;
2225 for (Index za_index_local = 0; za_index_local < Nza;
2227 out3 <<
"Calculate phase matrix \n";
2235 rtp_temperature_local,
2236 pha_mat_spt_agenda);
2251 for (Index za_in = 0; za_in < Nza; ++za_in) {
2252 for (Index aa_in = 0; aa_in < Naa; ++aa_in) {
2255 for (Index i = 0; i < stokes_dim; i++) {
2256 for (Index j = 0; j < stokes_dim; j++) {
2257 product_field(za_in, aa_in, i) +=
2258 pha_mat_local(za_in, aa_in, i, j) *
2259 cloudbox_field_mono(p_index,
2273 for (Index i = 0; i < stokes_dim; i++) {
2274 doit_scat_field(p_index,
2291 doit_scat_field(joker, joker, joker, joker, 0, joker) =
2292 doit_scat_field(joker, joker, joker, joker, Naa - 1, joker);
2299 Tensor6& doit_scat_field,
2301 const Agenda& pha_mat_spt_agenda,
2302 const Tensor6& cloudbox_field_mono,
2303 const Tensor4& pnd_field,
2304 const Tensor3& t_field,
2305 const Index& atmosphere_dim,
2307 const Vector& za_grid,
2308 const Vector& aa_grid,
2309 const Index& doit_za_grid_size,
2310 const Index& doit_za_interp,
2311 const Tensor7& pha_mat_doit,
2322 const Index Nza = za_grid.nelem();
2325 "The range of *za_grid* must [0 180].");
2328 const Index Naa = aa_grid.nelem();
2331 "The range of *aa_grid* must [0 360].");
2334 const Index stokes_dim = doit_scat_field.ncols();
2345 if (atmosphere_dim == 1) {
2347 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
2354 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
2360 }
else if (atmosphere_dim == 3) {
2362 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
2363 (cloudbox_limits[3] - cloudbox_limits[2]) + 1,
2364 (cloudbox_limits[5] - cloudbox_limits[4]) + 1,
2369 (cloudbox_limits[1] - cloudbox_limits[0]) + 1,
2370 (cloudbox_limits[3] - cloudbox_limits[2]) + 1,
2371 (cloudbox_limits[5] - cloudbox_limits[4]) + 1,
2377 "The atmospheric dimension must be 1D or 3D \n"
2378 "for scattering calculations using the DOIT\n"
2379 "module, but it is not. The value of *atmosphere_dim*\n"
2380 "is ", atmosphere_dim,
".")
2384 "Interpolation method is not defined. Use \n"
2385 "*doit_za_interpSet*.\n");
2388 "*cloudbox_limits* is a vector which contains the"
2389 "upper and lower limit of the cloud for all "
2390 "atmospheric dimensions. So its dimension must"
2391 "be 2 x *atmosphere_dim*");
2394 "*doit_za_grid_size* must be greater than 15 for"
2395 "accurate results");
2396 if (doit_za_grid_size > 100) {
2398 out1 <<
"Warning: doit_za_grid_size is very large which means that the \n"
2399 <<
"calculation will be very slow. The recommended value is 19.\n";
2405 Tensor4 pha_mat_local(
2406 doit_za_grid_size, aa_grid.nelem(), stokes_dim, stokes_dim, 0.);
2408 Tensor5 pha_mat_spt_local(pnd_field.nbooks(),
2417 nlinspace(za_g, 0, 180, doit_za_grid_size);
2422 gridpos(gp_za_i, za_grid, za_g);
2424 Matrix itw_za_i(doit_za_grid_size, 2);
2428 Matrix cloudbox_field_int(doit_za_grid_size, stokes_dim, 0);
2433 gridpos(gp_za, za_g, za_grid);
2435 Matrix itw_za(Nza, 2);
2439 Matrix doit_scat_field_org(doit_za_grid_size, stokes_dim, 0);
2443 Vector grid_stepsize(2);
2444 grid_stepsize[0] = 180. / (Numeric)(doit_za_grid_size - 1);
2446 grid_stepsize[1] = 360. / (Numeric)(Naa - 1);
2449 Tensor3 product_field(doit_za_grid_size, Naa, stokes_dim, 0);
2451 if (atmosphere_dim == 1) {
2454 for (Index p_index = 0; p_index <= cloudbox_limits[1] - cloudbox_limits[0];
2457 for (Index i = 0; i < stokes_dim; i++) {
2458 if (doit_za_interp == 0) {
2459 interp(cloudbox_field_int(joker, i),
2461 cloudbox_field_mono(p_index, 0, 0, joker, 0, i),
2463 }
else if (doit_za_interp == 1) {
2465 for (Index za = 0; za < za_g.nelem(); za++) {
2466 cloudbox_field_int(za, i) =
2468 cloudbox_field_mono(p_index, 0, 0, joker, 0, i),
2479 for (Index za_index_local = 0; za_index_local < doit_za_grid_size;
2481 out3 <<
"Multiplication of phase matrix with incoming"
2482 <<
" intensities \n";
2488 for (Index za_in = 0; za_in < doit_za_grid_size; za_in++) {
2489 for (Index aa_in = 0; aa_in < Naa; ++aa_in) {
2493 for (Index i = 0; i < stokes_dim; i++) {
2494 for (Index j = 0; j < stokes_dim; j++) {
2495 product_field(za_in, aa_in, i) +=
2497 p_index, za_index_local, 0, za_in, aa_in, i, j) *
2498 cloudbox_field_int(za_in, j);
2505 out3 <<
"Compute integral. \n";
2507 for (Index i = 0; i < stokes_dim; i++) {
2508 doit_scat_field_org(za_index_local, i) =
2513 for (Index i = 0; i < stokes_dim; i++) {
2514 doit_scat_field_org(za_index_local, i) =
2526 for (Index i = 0; i < stokes_dim; i++) {
2527 if (doit_za_interp == 0)
2529 interp(doit_scat_field(p_index, 0, 0, joker, 0, i),
2531 doit_scat_field_org(joker, i),
2535 for (Index za = 0; za < za_grid.nelem(); za++) {
2536 doit_scat_field(p_index, 0, 0, za, 0, i) =
interp_poly(
2537 za_g, doit_scat_field_org(joker, i), za_grid[za], gp_za[za]);
2545 else if (atmosphere_dim == 3) {
2547 for (Index p_index = 0; p_index <= cloudbox_limits[1] - cloudbox_limits[0];
2549 for (Index lat_index = 0;
2550 lat_index <= cloudbox_limits[3] - cloudbox_limits[2];
2552 for (Index lon_index = 0;
2553 lon_index <= cloudbox_limits[5] - cloudbox_limits[4];
2555 Numeric rtp_temperature_local =
2556 t_field(p_index + cloudbox_limits[0],
2557 lat_index + cloudbox_limits[2],
2558 lon_index + cloudbox_limits[4]);
2561 for (Index aa_index_local = 1; aa_index_local < Naa;
2564 for (Index i = 0; i < stokes_dim; i++) {
2566 cloudbox_field_int(joker, i),
2568 cloudbox_field_mono(
2569 p_index, lat_index, lon_index, joker, aa_index_local, i),
2573 for (Index za_index_local = 0; za_index_local < doit_za_grid_size;
2575 out3 <<
"Calculate phase matrix \n";
2583 rtp_temperature_local,
2584 pha_mat_spt_agenda);
2599 out3 <<
"Multiplication of phase matrix with"
2600 <<
"incoming intensity \n";
2602 for (Index za_in = 0; za_in < doit_za_grid_size; za_in++) {
2603 for (Index aa_in = 0; aa_in < Naa; ++aa_in) {
2606 for (Index i = 0; i < stokes_dim; i++) {
2607 for (Index j = 0; j < stokes_dim; j++) {
2608 product_field(za_in, aa_in, i) +=
2609 pha_mat_local(za_in, aa_in, i, j) *
2610 cloudbox_field_int(za_in, j);
2616 out3 <<
"Compute the integral \n";
2618 for (Index i = 0; i < stokes_dim; i++) {
2619 doit_scat_field_org(za_index_local, i) =
2628 for (Index i = 0; i < stokes_dim; i++) {
2631 p_index, lat_index, lon_index, joker, aa_index_local, i),
2633 doit_scat_field_org(joker, i),
2640 doit_scat_field(joker, joker, joker, joker, 0, joker) =
2641 doit_scat_field(joker, joker, joker, joker, Naa - 1, joker);
2643 out2 <<
" Finished scattered field.\n";
2648 Vector& doit_za_grid_opt,
2650 const Tensor6& cloudbox_field_mono,
2651 const Vector& za_grid,
2652 const Index& doit_za_interp,
2663 cloudbox_field_mono,
2664 cloudbox_field_mono.nvitrines(),
2669 cloudbox_field_mono.ncols());
2671 ARTS_USER_ERROR_IF (cloudbox_field_mono.ncols() < 1 || cloudbox_field_mono.ncols() > 4,
2672 "The last dimension of *cloudbox_field* corresponds\n"
2673 "to the Stokes dimension, therefore the number of\n"
2674 "columns in *cloudbox_field* must be a number between\n"
2675 "1 and 4, but it is not!");
2678 "Interpolation method is not defined. Use \n"
2679 "*doit_za_interpSet*.\n");
2681 if (za_grid.nelem() < 500) {
2682 out1 <<
"Warning: The fine grid (*za_grid*) has less than\n"
2683 <<
"500 grid points which is likely not sufficient for\n"
2684 <<
"grid_optimization\n";
2693 Matrix cloudbox_field_opt_mat;
2694 cloudbox_field_opt_mat = 0.;
2698 cloudbox_field_opt_mat,
2700 cloudbox_field_mono,
2707 const Index& atmosphere_dim,
2714 "Polynomial interpolation is only implemented for\n"
2715 "1D DOIT calculations as \n"
2716 "in 3D there can be numerical problems.\n"
2717 "Please use 'linear' interpolation method.");
2719 if (method ==
"linear")
2721 else if (method ==
"polynomial")
2725 "Possible interpolation methods are 'linear' "
2726 "and 'polynomial'.\n");
2732 Tensor7& cloudbox_field,
2733 const Index& atmfields_checked,
2734 const Index& atmgeom_checked,
2735 const Index& cloudbox_checked,
2736 const Index& scat_data_checked,
2737 const Index& cloudbox_on,
2738 const Vector& f_grid,
2739 const Agenda& doit_mono_agenda,
2740 const Index& doit_is_initialized,
2748 out0 <<
" Cloudbox is off, DOIT calculation will be skipped.\n";
2757 "The atmospheric fields must be flagged to have "
2758 "passed a consistency check (atmfields_checked=1).");
2760 "The atmospheric geometry must be flagged to have "
2761 "passed a consistency check (atmgeom_checked=1).");
2763 "The cloudbox must be flagged to have "
2764 "passed a consistency check (cloudbox_checked=1).");
2767 if (!cloudbox_on)
return;
2770 "The scattering data must be flagged to have "
2771 "passed a consistency check (scat_data_checked=1).");
2782 "Initialization method *DoitInit* has to be called before "
2788 const Index nf = f_grid.nelem();
2792 bool failed =
false;
2796#pragma omp parallel for if (!arts_omp_in_parallel() && nf > 1) firstprivate(wss)
2797 for (Index f_index = 0; f_index < nf; f_index++) {
2799 cloudbox_field(f_index, joker, joker, joker, joker, joker, joker) = NAN;
2805 os <<
"Frequency: " << f_grid[f_index] / 1e9 <<
" GHz \n";
2808 Tensor6 cloudbox_field_mono_local{
2809 cloudbox_field(f_index, joker, joker, joker, joker, joker, joker)};
2811 cloudbox_field_mono_local,
2815 cloudbox_field(f_index, joker, joker, joker, joker, joker, joker) =
2816 cloudbox_field_mono_local;
2817 }
catch (
const std::exception& e) {
2818 cloudbox_field(f_index, joker, joker, joker, joker, joker, joker) = NAN;
2820 os <<
"Error for f_index = " << f_index <<
" (" << f_grid[f_index]
2823#pragma omp critical(DoitCalc_fail)
2826 fail_msg = os.str();
2838 Tensor7& cloudbox_field,
2839 const Index& atmfields_checked,
2840 const Index& atmgeom_checked,
2841 const Index& cloudbox_checked,
2842 const Index& doit_is_initialized,
2843 const Agenda& iy_main_agenda,
2844 const Index& atmosphere_dim,
2845 const Vector& lat_grid,
2846 const Vector& lon_grid,
2847 const Tensor3& z_field,
2849 const Index& cloudbox_on,
2851 const Vector& f_grid,
2852 const Index& stokes_dim,
2853 const Vector& za_grid,
2854 const Vector& aa_grid,
2855 const Index& rigorous,
2856 const Numeric& maxratio,
2860 "The atmospheric fields must be flagged to have "
2861 "passed a consistency check (atmfields_checked=1).");
2863 "The atmospheric geometry must be flagged to have "
2864 "passed a consistency check (atmgeom_checked=1).");
2866 "The cloudbox must be flagged to have "
2867 "passed a consistency check (cloudbox_checked=1).");
2870 if (!cloudbox_on)
return;
2874 "Initialization method *DoitInit* has to be called before "
2875 "*DoitGetIncoming*.")
2878 const String iy_unit =
"1";
2880 Index Nf = f_grid.
nelem();
2881 Index Np_cloud = cloudbox_limits[1] - cloudbox_limits[0] + 1;
2882 Index Nza = za_grid.
nelem();
2887 ArrayOfMatrix iy_aux;
2888 ArrayOfTensor3 diy_dx;
2890 Tensor3 iy_transmittance(0, 0, 0);
2897 "The atmospheric dimensionality must be 1 or 3.");
2899 "*za_grid* must include 0 and 180 degrees as endpoints.");
2902 if (atmosphere_dim == 1) {
2904 Vector los(1), pos(1);
2908 for (Index boundary = 0; boundary <= 1; boundary++) {
2909 const Index boundary_index =
2910 boundary ? cloudbox_field.nvitrines() - 1 : 0;
2911 pos[0] = z_field(cloudbox_limits[boundary], 0, 0);
2915 los[0] = za_grid[0];
2937 cloudbox_field(joker, boundary_index, 0, 0, 0, 0, joker) = iy;
2939 for (Index za_index = 1; za_index < Nza; za_index++) {
2940 los[0] = za_grid[za_index];
2962 cloudbox_field(joker, boundary_index, 0, 0, za_index, 0, joker) = iy;
2965 for (Index fi = 0; fi < Nf; fi++) {
2968 fi, boundary_index, 0, 0, za_index, 0, 0) >
2970 cloudbox_field(fi, boundary_index, 0, 0, za_index - 1, 0, 0) /
2972 fi, boundary_index, 0, 0, za_index, 0, 0) <
2974 "ERROR: Radiance difference between "
2975 "interpolation points is too large (factor ", maxratio,
2977 "to safely interpolate. This might be due to "
2978 "za_grid being too coarse or the radiance field "
2979 "being a step-like function.\n"
2980 "Happens at boundary ", boundary_index,
2981 " between zenith angles ", za_grid[za_index - 1],
2982 " and ", za_grid[za_index],
"deg\n"
2983 "for frequency #", fi,
", where radiances are ",
2984 cloudbox_field(fi, boundary_index, 0, 0, za_index - 1, 0, 0),
2986 cloudbox_field(fi, boundary_index, 0, 0, za_index, 0, 0),
2996 Index Naa = aa_grid.nelem();
2999 "*aa_grid* must include 0 and 360 degrees as endpoints.");
3001 Index Nlat_cloud = cloudbox_limits[3] - cloudbox_limits[2] + 1;
3002 Index Nlon_cloud = cloudbox_limits[5] - cloudbox_limits[4] + 1;
3008 for (Index i = 0; i < Naa; i++) aa_g[i] = aa_grid[i] - 180;
3011 Vector los(2), pos(3);
3015 for (Index boundary = 0; boundary <= 1; boundary++) {
3016 const Index boundary_index =
3017 boundary ? cloudbox_field.nvitrines() - 1 : 0;
3018 for (Index lat_index = 0; lat_index < Nlat_cloud; lat_index++) {
3019 for (Index lon_index = 0; lon_index < Nlon_cloud; lon_index++) {
3020 pos[2] = lon_grid[lon_index + cloudbox_limits[4]];
3021 pos[1] = lat_grid[lat_index + cloudbox_limits[2]];
3022 pos[0] = z_field(cloudbox_limits[boundary],
3023 lat_index + cloudbox_limits[2],
3024 lon_index + cloudbox_limits[4]);
3026 for (Index za_index = 0; za_index < Nza; za_index++) {
3027 for (Index aa_index = 0; aa_index < Naa; aa_index++) {
3028 los[0] = za_grid[za_index];
3029 los[1] = aa_g[aa_index];
3034 if ((za_index != 0 && za_index != (Nza - 1)) || aa_index == 0) {
3056 cloudbox_field(joker,
3070 for (Index boundary = 0; boundary <= 1; boundary++) {
3071 const Index boundary_index = boundary ? cloudbox_field.nshelves() - 1 : 0;
3072 for (Index p_index = 0; p_index < Np_cloud; p_index++) {
3073 for (Index lon_index = 0; lon_index < Nlon_cloud; lon_index++) {
3074 pos[2] = lon_grid[lon_index + cloudbox_limits[4]];
3075 pos[1] = lat_grid[cloudbox_limits[boundary + 2]];
3076 pos[0] = z_field(p_index + cloudbox_limits[0],
3077 cloudbox_limits[boundary + 2],
3078 lon_index + cloudbox_limits[4]);
3080 for (Index za_index = 0; za_index < Nza; za_index++) {
3081 for (Index aa_index = 0; aa_index < Naa; aa_index++) {
3082 los[0] = za_grid[za_index];
3083 los[1] = aa_g[aa_index];
3087 if ((za_index != 0 && za_index != (Nza - 1)) || aa_index == 0) {
3109 cloudbox_field(joker,
3123 for (Index boundary = 0; boundary <= 1; boundary++) {
3124 const Index boundary_index = boundary ? cloudbox_field.nbooks() - 1 : 0;
3125 for (Index p_index = 0; p_index < Np_cloud; p_index++) {
3126 for (Index lat_index = 0; lat_index < Nlat_cloud; lat_index++) {
3127 pos[2] = lon_grid[cloudbox_limits[boundary + 4]];
3128 pos[1] = lat_grid[lat_index + cloudbox_limits[2]];
3129 pos[0] = z_field(p_index + cloudbox_limits[0],
3130 lat_index + cloudbox_limits[2],
3131 cloudbox_limits[boundary + 4]);
3133 for (Index za_index = 0; za_index < Nza; za_index++) {
3134 for (Index aa_index = 0; aa_index < Naa; aa_index++) {
3135 los[0] = za_grid[za_index];
3136 los[1] = aa_g[aa_index];
3140 if ((za_index != 0 && za_index != (Nza - 1)) || aa_index == 0) {
3162 cloudbox_field(joker,
3179 Tensor7& cloudbox_field,
3181 const Index& atmfields_checked,
3182 const Index& atmgeom_checked,
3183 const Index& cloudbox_checked,
3184 const Index& doit_is_initialized,
3185 const Agenda& iy_main_agenda,
3186 const Index& atmosphere_dim,
3187 const Vector& lat_grid,
3188 const Vector& lon_grid,
3189 const Tensor3& z_field,
3192 const Vector& f_grid,
3193 const Index& stokes_dim,
3194 const Vector& za_grid,
3195 const Vector& aa_grid,
3199 "The atmospheric fields must be flagged to have "
3200 "passed a consistency check (atmfields_checked=1).");
3202 "The atmospheric geometry must be flagged to have "
3203 "passed a consistency check (atmgeom_checked=1).");
3205 "The cloudbox must be flagged to have "
3206 "passed a consistency check (cloudbox_checked=1).");
3209 if (!cloudbox_on)
return;
3213 "Initialization method *DoitInit* has to be called before "
3214 "*DoitGetIncoming1DAtm*.")
3217 const String iy_unit =
"1";
3219 Index Np_cloud = cloudbox_limits[1] - cloudbox_limits[0] + 1;
3220 Index Nlat_cloud = cloudbox_limits[3] - cloudbox_limits[2] + 1;
3221 Index Nlon_cloud = cloudbox_limits[5] - cloudbox_limits[4] + 1;
3222 Index Nza = za_grid.
nelem();
3223 Index Naa = aa_grid.nelem();
3228 ArrayOfMatrix iy_aux;
3229 ArrayOfTensor3 diy_dx;
3231 Tensor3 iy_transmittance(0, 0, 0);
3237 "The atmospheric dimensionality must be 3.");
3239 "*za_grid* must include 0 and 180 degrees as endpoints.");
3241 "*aa_grid* must include 0 and 360 degrees as endpoints.");
3253 for (Index i = 0; i < Naa; i++) aa_g[i] = aa_grid[i] - 180;
3260 Vector los(2), pos(3);
3263 pos[1] = lat_grid[cloudbox_limits[2]];
3264 pos[2] = lon_grid[cloudbox_limits[4]];
3265 los[1] = aa_g[aa_index];
3268 for (Index p_index = 0; p_index < Np_cloud; p_index++) {
3270 cloudbox_limits[0] + p_index, cloudbox_limits[2], cloudbox_limits[4]);
3272 for (Index za_index = 0; za_index < Nza; za_index++) {
3273 los[0] = za_grid[za_index];
3295 for (Index aa = 0; aa < Naa; aa++) {
3298 for (Index lat = 0; lat < Nlat_cloud; lat++) {
3299 for (Index lon = 0; lon < Nlon_cloud; lon++) {
3300 cloudbox_field(joker, 0, lat, lon, za_index, aa, joker) = iy;
3305 else if (p_index == Np_cloud - 1)
3306 for (Index lat = 0; lat < Nlat_cloud; lat++) {
3307 for (Index lon = 0; lon < Nlon_cloud; lon++) {
3308 cloudbox_field(joker,
3309 cloudbox_field.nvitrines(),
3319 for (Index lat = 0; lat < 2; lat++) {
3320 for (Index lon = 0; lon < Nlon_cloud; lon++) {
3321 const Index boundary_index =
3322 lat ? cloudbox_field.nshelves() - 1 : 0;
3324 joker, p_index, boundary_index, lon, za_index, aa, joker) = iy;
3329 for (Index lat = 0; lat < Nlat_cloud; lat++) {
3330 for (Index lon = 0; lon < 2; lon++) {
3331 const Index boundary_index = lon ? cloudbox_field.nbooks() - 1 : 0;
3333 joker, p_index, lat, boundary_index, za_index, aa, joker) = iy;
3344 const Vector& za_grid,
3345 const Vector& f_grid,
3346 const Index& atmosphere_dim,
3347 const Index& stokes_dim,
3349 const Index& doit_is_initialized,
3350 const Tensor7& cloudbox_field_precalc,
3355 "This method is currently only implemented for 1D atmospheres!\n")
3359 "Initialization method *DoitInit* has to be called before "
3360 "*cloudbox_fieldSetFromPrecalc*")
3363 Index nf = f_grid.nelem();
3364 Index nza = za_grid.nelem();
3365 Index np = cloudbox_limits[1] - cloudbox_limits[0] + 1;
3368 "cloudbox_field_precalc has wrong size in frequency dimension.\n",
3369 nf,
" frequency points are expected, but cloudbox_field_precalc "
3370 "contains ", cloudbox_field_precalc.nlibraries(),
3371 "frequency points.\n")
3373 "cloudbox_field_precalc has wrong size in pressure level dimension.\n",
3374 np,
" pressure levels expected, but cloudbox_field_precalc "
3375 "contains ", cloudbox_field_precalc.nvitrines(),
3376 "pressure levels.\n")
3378 "cloudbox_field_precalc has wrong size in polar angle dimension.\n",
3379 nza,
" angles expected, but cloudbox_field_precalc "
3380 "contains ", cloudbox_field_precalc.npages(),
"angles.\n")
3382 "cloudbox_field_precalc has wrong stokes dimension.\n"
3383 "Dimension ", stokes_dim,
3384 " expected, but cloudbox_field_precalc is dimesnion ",
3385 cloudbox_field_precalc.ncols(),
".\n")
3393 Index first_upwell = 0;
3396 while (za_grid[first_upwell] < 90.) first_upwell++;
3398 Range downwell(0, first_upwell);
3399 Range upwell(first_upwell, za_grid.nelem() - first_upwell);
3402 cloudbox_field(joker, Range(1, np - 2), 0, 0, joker, 0, joker) =
3403 cloudbox_field_precalc(joker, Range(1, np - 2), 0, 0, joker, 0, joker);
3410 cloudbox_field(joker, np - 1, 0, 0, upwell, 0, joker) =
3411 cloudbox_field_precalc(joker, np - 1, 0, 0, upwell, 0, joker);
3413 if (cloudbox_limits[0] != 0)
3416 cloudbox_field(joker, 0, 0, 0, downwell, 0, joker) =
3417 cloudbox_field_precalc(joker, 0, 0, 0, downwell, 0, joker);
3423 cloudbox_field(joker, 0, 0, 0, joker, 0, joker) =
3424 cloudbox_field_precalc(joker, 0, 0, 0, joker, 0, joker);
3429 const Vector& p_grid,
3430 const Vector& lat_grid,
3431 const Vector& lon_grid,
3433 const Index& atmosphere_dim,
3434 const Index& cloudbox_on,
3435 const Index& doit_is_initialized,
3436 const Index& all_frequencies,
3441 if (!cloudbox_on)
return;
3445 "Initialization method *DoitInit* has to be called before "
3446 "*cloudbox_fieldSetClearsky*.")
3449 <<
" Interpolate boundary clearsky field to obtain the initial field.\n";
3454 if (atmosphere_dim == 1) {
3455 const Index nf = all_frequencies ? cloudbox_field.nlibraries() : 1;
3457 for (Index f_index = 0; f_index < nf; f_index++) {
3458 Index N_za = cloudbox_field.npages();
3459 Index N_aa = cloudbox_field.nrows();
3460 Index N_i = cloudbox_field.ncols();
3468 ArrayOfGridPos p_gp((cloudbox_limits[1] - cloudbox_limits[0]) + 1);
3471 p_grid[Range(cloudbox_limits[0],
3473 (cloudbox_limits[1] - cloudbox_limits[0]))],
3474 p_grid[Range(cloudbox_limits[0],
3475 (cloudbox_limits[1] - cloudbox_limits[0]) + 1)]);
3477 Matrix itw((cloudbox_limits[1] - cloudbox_limits[0]) + 1, 2);
3480 Tensor6 scat_i_p(2, 1, 1, N_za, 1, N_i);
3481 scat_i_p(0, joker, joker, joker, joker, joker) =
3482 cloudbox_field(f_index, 0, joker, joker, joker, joker, joker);
3483 scat_i_p(1, joker, joker, joker, joker, joker) =
3484 cloudbox_field(f_index,
3485 cloudbox_field.nvitrines() - 1,
3492 for (Index za_index = 0; za_index < N_za; ++za_index) {
3493 for (Index aa_index = 0; aa_index < N_aa; ++aa_index) {
3494 for (Index i = 0; i < N_i; ++i) {
3495 VectorView target_field = cloudbox_field(
3496 f_index, Range(joker), 0, 0, za_index, aa_index, i);
3498 ConstVectorView source_field =
3499 scat_i_p(Range(joker), 0, 0, za_index, aa_index, i);
3501 interp(target_field, itw, source_field, p_gp);
3506 }
else if (atmosphere_dim == 3) {
3508 "Error in cloudbox_fieldSetClearsky: For 3D "
3509 "all_frequencies option is not implemented \n");
3511 for (Index f_index = 0; f_index < cloudbox_field.nvitrines(); f_index++) {
3512 Index N_p = cloudbox_field.nvitrines();
3513 Index N_lat = cloudbox_field.nshelves();
3514 Index N_lon = cloudbox_field.nbooks();
3515 Index N_za = cloudbox_field.npages();
3516 Index N_aa = cloudbox_field.nrows();
3517 Index N_i = cloudbox_field.ncols();
3519 Tensor6 scat_i_p(2, N_lat, N_lon, N_za, N_aa, N_i);
3520 scat_i_p(0, joker, joker, joker, joker, joker) =
3521 cloudbox_field(f_index, 0, joker, joker, joker, joker, joker);
3522 scat_i_p(1, joker, joker, joker, joker, joker) =
3523 cloudbox_field(f_index, N_p - 1, joker, joker, joker, joker, joker);
3525 Tensor6 scat_i_lat(N_p, 2, N_lon, N_za, N_aa, N_i);
3526 scat_i_lat(joker, 0, joker, joker, joker, joker) =
3527 cloudbox_field(f_index, joker, 0, joker, joker, joker, joker);
3528 scat_i_lat(joker, 1, joker, joker, joker, joker) =
3529 cloudbox_field(f_index, joker, N_lat - 1, joker, joker, joker, joker);
3531 Tensor6 scat_i_lon(N_p, N_lat, 2, N_za, N_aa, N_i);
3532 scat_i_lon(joker, joker, 0, joker, joker, joker) =
3533 cloudbox_field(f_index, joker, joker, 0, joker, joker, joker);
3534 scat_i_lon(joker, joker, 1, joker, joker, joker) =
3535 cloudbox_field(f_index, joker, joker, N_lon - 1, joker, joker, joker);
3539 ArrayOfGridPos p_gp((cloudbox_limits[1] - cloudbox_limits[0]) + 1);
3540 ArrayOfGridPos lat_gp((cloudbox_limits[3] - cloudbox_limits[2]) + 1);
3541 ArrayOfGridPos lon_gp((cloudbox_limits[5] - cloudbox_limits[4]) + 1);
3548 p_grid[Range(cloudbox_limits[0],
3550 (cloudbox_limits[1] - cloudbox_limits[0]))],
3551 p_grid[Range(cloudbox_limits[0],
3552 (cloudbox_limits[1] - cloudbox_limits[0]) + 1)]);
3554 lat_grid[Range(cloudbox_limits[2],
3556 (cloudbox_limits[3] - cloudbox_limits[2]))],
3557 lat_grid[Range(cloudbox_limits[2],
3558 (cloudbox_limits[3] - cloudbox_limits[2]) + 1)]);
3560 lon_grid[Range(cloudbox_limits[4],
3562 (cloudbox_limits[5] - cloudbox_limits[4]))],
3563 lon_grid[Range(cloudbox_limits[4],
3564 (cloudbox_limits[5] - cloudbox_limits[4]) + 1)]);
3569 Matrix itw_p((cloudbox_limits[1] - cloudbox_limits[0]) + 1, 2);
3570 Matrix itw_lat((cloudbox_limits[3] - cloudbox_limits[2]) + 1, 2);
3571 Matrix itw_lon((cloudbox_limits[5] - cloudbox_limits[4]) + 1, 2);
3578 for (Index lat_index = 0;
3579 lat_index <= (cloudbox_limits[3] - cloudbox_limits[2]);
3581 for (Index lon_index = 0;
3582 lon_index <= (cloudbox_limits[5] - cloudbox_limits[4]);
3584 for (Index za_index = 0; za_index < N_za; ++za_index) {
3585 for (Index aa_index = 0; aa_index < N_aa; ++aa_index) {
3586 for (Index i = 0; i < N_i; ++i) {
3587 VectorView target_field = cloudbox_field(f_index,
3595 ConstVectorView source_field = scat_i_p(
3596 Range(joker), lat_index, lon_index, za_index, aa_index, i);
3598 interp(target_field, itw_p, source_field, p_gp);
3605 for (Index p_index = 0;
3606 p_index <= (cloudbox_limits[1] - cloudbox_limits[0]);
3608 for (Index lon_index = 0;
3609 lon_index <= (cloudbox_limits[5] - cloudbox_limits[4]);
3611 for (Index za_index = 0; za_index < N_za; ++za_index) {
3612 for (Index aa_index = 0; aa_index < N_aa; ++aa_index) {
3613 for (Index i = 0; i < N_i; ++i) {
3614 VectorView target_field = cloudbox_field(f_index,
3622 ConstVectorView source_field = scat_i_lat(
3623 p_index, Range(joker), lon_index, za_index, aa_index, i);
3625 interp(target_field, itw_lat, source_field, lat_gp);
3632 for (Index p_index = 0;
3633 p_index <= (cloudbox_limits[1] - cloudbox_limits[0]);
3635 for (Index lat_index = 0;
3636 lat_index <= (cloudbox_limits[3] - cloudbox_limits[2]);
3638 for (Index za_index = 0; za_index < N_za; ++za_index) {
3639 for (Index aa_index = 0; aa_index < N_aa; ++aa_index) {
3640 for (Index i = 0; i < N_i; ++i) {
3641 VectorView target_field = cloudbox_field(f_index,
3649 ConstVectorView source_field = scat_i_lon(
3650 p_index, lat_index, Range(joker), za_index, aa_index, i);
3652 interp(target_field, itw_lon, source_field, lon_gp);
3664 Tensor7& cloudbox_field,
3666 const Vector& p_grid,
3667 const Vector& lat_grid,
3668 const Vector& lon_grid,
3670 const Index& atmosphere_dim,
3671 const Index& stokes_dim,
3673 const Vector& cloudbox_field_values,
3677 Tensor6 cloudbox_field_mono(cloudbox_field.nvitrines(),
3678 cloudbox_field.nshelves(),
3679 cloudbox_field.nbooks(),
3680 cloudbox_field.npages(),
3681 cloudbox_field.nrows(),
3682 cloudbox_field.ncols());
3684 for (Index f_index = 0; f_index < cloudbox_field.nlibraries(); f_index++) {
3685 cloudbox_field_mono =
3686 cloudbox_field(f_index, joker, joker, joker, joker, joker, joker);
3695 cloudbox_field_values,
3698 cloudbox_field(f_index, joker, joker, joker, joker, joker, joker) =
3699 cloudbox_field_mono;
3705 Tensor7& cloudbox_field,
3707 const Vector& p_grid,
3708 const Vector& lat_grid,
3709 const Vector& lon_grid,
3711 const Index& atmosphere_dim,
3712 const Index& stokes_dim,
3714 const Matrix& cloudbox_field_values,
3719 "Number of rows in *cloudbox_field_values* has to be equal"
3720 " the frequency dimension of *cloudbox_field*.");
3722 Tensor6 cloudbox_field_mono(cloudbox_field.nvitrines(),
3723 cloudbox_field.nshelves(),
3724 cloudbox_field.nbooks(),
3725 cloudbox_field.npages(),
3726 cloudbox_field.nrows(),
3727 cloudbox_field.ncols());
3729 for (Index f_index = 0; f_index < cloudbox_field.nlibraries(); f_index++) {
3730 cloudbox_field_mono =
3731 cloudbox_field(f_index, joker, joker, joker, joker, joker, joker);
3740 Vector{cloudbox_field_values(f_index, joker)},
3743 cloudbox_field(f_index, joker, joker, joker, joker, joker, joker) =
3744 cloudbox_field_mono;
3750 Tensor6& cloudbox_field_mono,
3752 const Vector& p_grid,
3753 const Vector& lat_grid,
3754 const Vector& lon_grid,
3756 const Index& atmosphere_dim,
3757 const Index& stokes_dim,
3759 const Vector& cloudbox_field_values,
3764 out2 <<
" Set initial field to constant values: " << cloudbox_field_values
3774 "The dimension of stokes vector must be"
3778 "Length of *cloudbox_field_values* has to be equal"
3781 "*cloudbox_limits* is a vector which contains the"
3782 "upper and lower limit of the cloud for all "
3783 "atmospheric dimensions. So its dimension must"
3784 "be 2 x *atmosphere_dim*.");
3786 for (Index i = 0; i < stokes_dim; i++) {
3787 cloudbox_field_mono(joker, joker, joker, joker, joker, i) =
3788 cloudbox_field_values[i];
Declarations for agendas.
This file contains the definition of Array.
The global header file for ARTS.
Constants of physical expressions as constexpr.
void doit_rte_agendaExecute(Workspace &ws, Tensor6 &cloudbox_field_mono, const Tensor6 &doit_scat_field, const Agenda &input_agenda)
void doit_scat_field_agendaExecute(Workspace &ws, Tensor6 &doit_scat_field, const Tensor6 &cloudbox_field_mono, const Agenda &input_agenda)
void doit_mono_agendaExecute(Workspace &ws, Tensor6 &cloudbox_field_mono, const Vector &f_grid, const Index f_index, 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 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 doit_conv_test_agendaExecute(Workspace &ws, Index &doit_conv_flag, Index &doit_iteration_counter, const Tensor6 &cloudbox_field_mono, const Tensor6 &cloudbox_field_mono_old, const Agenda &input_agenda)
void pha_mat_spt_agendaExecute(Workspace &ws, Tensor5 &pha_mat_spt, const Index za_index, const Index scat_lat_index, const Index scat_lon_index, const Index scat_p_index, const Index aa_index, const Numeric rtp_temperature, const Agenda &input_agenda)
Index nelem() const ARTS_NOEXCEPT
#define ARTS_ASSERT(condition,...)
#define ARTS_USER_ERROR(...)
#define ARTS_USER_ERROR_IF(condition,...)
void za_gridOpt(Vector &za_grid_opt, Matrix &cloudbox_field_opt, ConstVectorView za_grid_fine, ConstTensor6View cloudbox_field_mono, const Numeric &acc, const Index &scat_za_interp)
Optimize the zenith angle grid.
void doit_scat_fieldNormalize(Workspace &ws, Tensor6 &doit_scat_field, const Tensor6 &cloudbox_field_mono, const ArrayOfIndex &cloudbox_limits, const Agenda &spt_calc_agenda, const Index &atmosphere_dim, const Vector &za_grid, const Vector &aa_grid, const Tensor4 &pnd_field, const Tensor3 &t_field, const Numeric &norm_error_threshold, const Index &norm_debug, const Verbosity &verbosity)
Normalization of scattered field.
void cloud_fieldsCalc(Workspace &ws, Tensor5View ext_mat_field, Tensor4View abs_vec_field, const Agenda &spt_calc_agenda, const Index &za_index, const Index &aa_index, const ArrayOfIndex &cloudbox_limits, ConstTensor3View t_field, ConstTensor4View pnd_field, const Verbosity &verbosity)
Calculate ext_mat, abs_vec for all points inside the cloudbox for one.
void cloud_ppath_update1D_noseq(Workspace &ws, Tensor6View cloudbox_field_mono, const Index &p_index, const Index &za_index, ConstVectorView za_grid, const ArrayOfIndex &cloudbox_limits, ConstTensor6View cloudbox_field_mono_old, ConstTensor6View doit_scat_field, const Agenda &propmat_clearsky_agenda, ConstTensor4View vmr_field, const Agenda &ppath_step_agenda, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, ConstVectorView p_grid, ConstTensor3View z_field, ConstVectorView refellipsoid, ConstTensor3View t_field, ConstVectorView f_grid, const Index &f_index, ConstTensor5View ext_mat_field, ConstTensor4View abs_vec_field, const Agenda &surface_rtprop_agenda, const Index &scat_za_interp, const Verbosity &verbosity)
Calculation of radiation field along a propagation path step for specified zenith direction and press...
void cloudbox_field_ngAcceleration(Tensor6 &cloudbox_field_mono, const ArrayOfTensor6 &acceleration_input, const Index &accelerated, const Verbosity &)
Convergence acceleration.
void cloud_ppath_update3D(Workspace &ws, Tensor6View cloudbox_field_mono, const Index &p_index, const Index &lat_index, const Index &lon_index, const Index &za_index, const Index &aa_index, ConstVectorView za_grid, ConstVectorView aa_grid, const ArrayOfIndex &cloudbox_limits, ConstTensor6View doit_scat_field, const Agenda &propmat_clearsky_agenda, ConstTensor4View vmr_field, const Agenda &ppath_step_agenda, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, ConstVectorView p_grid, ConstVectorView lat_grid, ConstVectorView lon_grid, ConstTensor3View z_field, ConstVectorView refellipsoid, ConstTensor3View t_field, ConstVectorView f_grid, const Index &f_index, ConstTensor5View ext_mat_field, ConstTensor4View abs_vec_field, const Index &, const Verbosity &verbosity)
Radiative transfer calculation along a path inside the cloudbox (3D).
void cloud_ppath_update1D_planeparallel(Workspace &ws, Tensor6View cloudbox_field_mono, const Index &p_index, const Index &za_index, ConstVectorView za_grid, const ArrayOfIndex &cloudbox_limits, ConstTensor6View doit_scat_field, const Agenda &propmat_clearsky_agenda, ConstTensor4View vmr_field, ConstVectorView p_grid, ConstTensor3View z_field, ConstTensor3View t_field, ConstVectorView f_grid, const Index &f_index, ConstTensor5View ext_mat_field, ConstTensor4View abs_vec_field, const Verbosity &verbosity)
Radiative transfer calculation inside cloudbox for planeparallel case.
void cloud_ppath_update1D(Workspace &ws, Tensor6View cloudbox_field_mono, const Index &p_index, const Index &za_index, ConstVectorView za_grid, const ArrayOfIndex &cloudbox_limits, ConstTensor6View doit_scat_field, const Agenda &propmat_clearsky_agenda, ConstTensor4View vmr_field, const Agenda &ppath_step_agenda, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, ConstVectorView p_grid, ConstTensor3View z_field, ConstVectorView refellipsoid, ConstTensor3View t_field, ConstVectorView f_grid, const Index &f_index, ConstTensor5View ext_mat_field, ConstTensor4View abs_vec_field, const Agenda &surface_rtprop_agenda, const Index &scat_za_interp, const Verbosity &verbosity)
Calculates radiation field along a propagation path step for specified zenith direction and pressure ...
Radiative transfer in cloudbox.
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_poly(ConstVectorView x, ConstVectorView y, const Numeric &x_i, const GridPos &gp)
Polynomial interpolation.
Numeric interp(ConstVectorView itw, ConstVectorView a, const GridPos &tc)
Red 1D Interpolate.
Array< RetrievalQuantity > ArrayOfRetrievalQuantity
void DoitCalc(Workspace &ws, Tensor7 &cloudbox_field, const Index &atmfields_checked, const Index &atmgeom_checked, const Index &cloudbox_checked, const Index &scat_data_checked, const Index &cloudbox_on, const Vector &f_grid, const Agenda &doit_mono_agenda, const Index &doit_is_initialized, const Verbosity &verbosity)
WORKSPACE METHOD: DoitCalc.
void cloudbox_fieldUpdateSeq3D(Workspace &ws, Tensor6 &cloudbox_field_mono, const Tensor6 &doit_scat_field, const ArrayOfIndex &cloudbox_limits, const Agenda &propmat_clearsky_agenda, const Tensor4 &vmr_field, const Agenda &spt_calc_agenda, const Vector &za_grid, const Vector &aa_grid, const Tensor4 &pnd_field, const Agenda &ppath_step_agenda, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const Vector &refellipsoid, const Tensor3 &t_field, const Vector &f_grid, const Index &f_index, const Index &doit_za_interp, const Verbosity &verbosity)
WORKSPACE METHOD: cloudbox_fieldUpdateSeq3D.
void doit_za_grid_optCalc(Vector &doit_za_grid_opt, const Tensor6 &cloudbox_field_mono, const Vector &za_grid, const Index &doit_za_interp, const Numeric &acc, const Verbosity &verbosity)
WORKSPACE METHOD: doit_za_grid_optCalc.
void DoitGetIncoming1DAtm(Workspace &ws, Tensor7 &cloudbox_field, Index &cloudbox_on, const Index &atmfields_checked, const Index &atmgeom_checked, const Index &cloudbox_checked, const Index &doit_is_initialized, const Agenda &iy_main_agenda, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const EnergyLevelMap &nlte_field, const ArrayOfIndex &cloudbox_limits, const Vector &f_grid, const Index &stokes_dim, const Vector &za_grid, const Vector &aa_grid, const Verbosity &)
WORKSPACE METHOD: DoitGetIncoming1DAtm.
void DoitWriteIterationFields(const Index &doit_iteration_counter, const Tensor6 &cloudbox_field_mono, const Index &f_index, const ArrayOfIndex &iterations, const ArrayOfIndex &frequencies, const Verbosity &verbosity)
WORKSPACE METHOD: DoitWriteIterationFields.
void doit_conv_flagLsq(Index &doit_conv_flag, Index &doit_iteration_counter, Tensor6 &cloudbox_field_mono, const Tensor6 &cloudbox_field_mono_old, const Vector &f_grid, const Index &f_index, const Vector &epsilon, const Index &max_iterations, const Index &throw_nonconv_error, const Verbosity &verbosity)
WORKSPACE METHOD: doit_conv_flagLsq.
void cloudbox_fieldUpdate1D(Workspace &ws, Tensor6 &cloudbox_field_mono, const Tensor6 &doit_scat_field, const ArrayOfIndex &cloudbox_limits, const Agenda &propmat_clearsky_agenda, const Tensor4 &vmr_field, const Agenda &spt_calc_agenda, const Vector &za_grid, const Tensor4 &pnd_field, const Agenda &ppath_step_agenda, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Vector &p_grid, const Tensor3 &z_field, const Vector &refellipsoid, const Tensor3 &t_field, const Vector &f_grid, const Index &f_index, const Agenda &surface_rtprop_agenda, const Index &doit_za_interp, const Verbosity &verbosity)
WORKSPACE METHOD: cloudbox_fieldUpdate1D.
void doit_conv_flagAbs(Index &doit_conv_flag, Index &doit_iteration_counter, Tensor6 &cloudbox_field_mono, const Tensor6 &cloudbox_field_mono_old, const Vector &epsilon, const Index &max_iterations, const Index &throw_nonconv_error, const Verbosity &verbosity)
WORKSPACE METHOD: doit_conv_flagAbs.
void cloudbox_field_monoOptimizeReverse(Tensor6 &cloudbox_field_mono, const Vector &p_grid_orig, const Vector &p_grid, const ArrayOfIndex &cloudbox_limits, const Verbosity &)
WORKSPACE METHOD: cloudbox_field_monoOptimizeReverse.
void doit_za_interpSet(Index &doit_za_interp, const Index &atmosphere_dim, const String &method, const Verbosity &)
WORKSPACE METHOD: doit_za_interpSet.
constexpr Numeric RAD2DEG
void doit_scat_fieldCalcLimb(Workspace &ws, Tensor6 &doit_scat_field, const Agenda &pha_mat_spt_agenda, const Tensor6 &cloudbox_field_mono, const Tensor4 &pnd_field, const Tensor3 &t_field, const Index &atmosphere_dim, const ArrayOfIndex &cloudbox_limits, const Vector &za_grid, const Vector &aa_grid, const Index &doit_za_grid_size, const Index &doit_za_interp, const Tensor7 &pha_mat_doit, const Verbosity &verbosity)
WORKSPACE METHOD: doit_scat_fieldCalcLimb.
void cloudbox_fieldSetConstPerFreq(Tensor7 &cloudbox_field, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const ArrayOfIndex &cloudbox_limits, const Index &atmosphere_dim, const Index &stokes_dim, const Matrix &cloudbox_field_values, const Verbosity &verbosity)
WORKSPACE METHOD: cloudbox_fieldSetConstPerFreq.
void cloudbox_field_monoIterate(Workspace &ws, Tensor6 &cloudbox_field_mono, const Agenda &doit_scat_field_agenda, const Agenda &doit_rte_agenda, const Agenda &doit_conv_test_agenda, const Index &accelerated, const Verbosity &verbosity)
WORKSPACE METHOD: cloudbox_field_monoIterate.
void cloudbox_fieldSetFromPrecalc(Tensor7 &cloudbox_field, const Vector &za_grid, const Vector &f_grid, const Index &atmosphere_dim, const Index &stokes_dim, const ArrayOfIndex &cloudbox_limits, const Index &doit_is_initialized, const Tensor7 &cloudbox_field_precalc, const Verbosity &)
WORKSPACE METHOD: cloudbox_fieldSetFromPrecalc.
void DoitGetIncoming(Workspace &ws, Tensor7 &cloudbox_field, const Index &atmfields_checked, const Index &atmgeom_checked, const Index &cloudbox_checked, const Index &doit_is_initialized, const Agenda &iy_main_agenda, const Index &atmosphere_dim, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Vector &f_grid, const Index &stokes_dim, const Vector &za_grid, const Vector &aa_grid, const Index &rigorous, const Numeric &maxratio, const Verbosity &)
WORKSPACE METHOD: DoitGetIncoming.
void DOAngularGridsSet(Index &doit_za_grid_size, Vector &aa_grid, Vector &za_grid, const Index &N_za_grid, const Index &N_aa_grid, const String &za_grid_opt_file, const Verbosity &verbosity)
WORKSPACE METHOD: DOAngularGridsSet.
void cloudbox_field_monoSetConst(Tensor6 &cloudbox_field_mono, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const ArrayOfIndex &cloudbox_limits, const Index &atmosphere_dim, const Index &stokes_dim, const Vector &cloudbox_field_values, const Verbosity &verbosity)
WORKSPACE METHOD: cloudbox_field_monoSetConst.
void cloudbox_fieldSetConst(Tensor7 &cloudbox_field, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const ArrayOfIndex &cloudbox_limits, const Index &atmosphere_dim, const Index &stokes_dim, const Vector &cloudbox_field_values, const Verbosity &verbosity)
WORKSPACE METHOD: cloudbox_fieldSetConst.
void cloudbox_fieldSetClearsky(Tensor7 &cloudbox_field, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const ArrayOfIndex &cloudbox_limits, const Index &atmosphere_dim, const Index &cloudbox_on, const Index &doit_is_initialized, const Index &all_frequencies, const Verbosity &verbosity)
WORKSPACE METHOD: cloudbox_fieldSetClearsky.
void OptimizeDoitPressureGrid(Workspace &ws, Vector &p_grid, Tensor4 &pnd_field, Tensor3 &t_field, ArrayOfArrayOfSingleScatteringData &scat_data_mono, Tensor3 &z_field, ArrayOfIndex &cloudbox_limits, Tensor6 &cloudbox_field_mono, Tensor7 &pha_mat_doit, Tensor4 &vmr_field, Vector &p_grid_orig, const Vector &f_grid, const Index &f_index, const Agenda &propmat_clearsky_agenda, const Numeric &tau_scat_max, const Numeric &sgl_alb_max, const Index &cloudbox_size_max, const Verbosity &verbosity)
WORKSPACE METHOD: OptimizeDoitPressureGrid.
void doit_conv_flagAbsBT(Index &doit_conv_flag, Index &doit_iteration_counter, Tensor6 &cloudbox_field_mono, const Tensor6 &cloudbox_field_mono_old, const Vector &f_grid, const Index &f_index, const Vector &epsilon, const Index &max_iterations, const Index &throw_nonconv_error, const Verbosity &verbosity)
WORKSPACE METHOD: doit_conv_flagAbsBT.
void cloudbox_fieldUpdateSeq1D(Workspace &ws, Tensor6 &cloudbox_field_mono, Tensor6 &doit_scat_field, const ArrayOfIndex &cloudbox_limits, const Agenda &propmat_clearsky_agenda, const Tensor4 &vmr_field, const Agenda &spt_calc_agenda, const Vector &za_grid, const Vector &aa_grid, const Tensor4 &pnd_field, const Agenda &ppath_step_agenda, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Vector &p_grid, const Tensor3 &z_field, const Vector &refellipsoid, const Tensor3 &t_field, const Vector &f_grid, const Index &f_index, const Agenda &surface_rtprop_agenda, const Index &doit_za_interp, const Index &normalize, const Numeric &norm_error_threshold, const Index &norm_debug, const Verbosity &verbosity)
WORKSPACE METHOD: cloudbox_fieldUpdateSeq1D.
void cloudbox_fieldUpdateSeq1DPP(Workspace &ws, Tensor6 &cloudbox_field_mono, Index &za_index, const Tensor6 &doit_scat_field, const ArrayOfIndex &cloudbox_limits, const Agenda &propmat_clearsky_agenda, const Tensor4 &vmr_field, const Agenda &spt_calc_agenda, const Vector &za_grid, const Tensor4 &pnd_field, const Vector &p_grid, const Tensor3 &z_field, const Tensor3 &t_field, const Vector &f_grid, const Index &f_index, const Verbosity &verbosity)
WORKSPACE METHOD: cloudbox_fieldUpdateSeq1DPP.
void DoitInit(Tensor6 &doit_scat_field, Tensor7 &cloudbox_field, Index &doit_is_initialized, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &f_grid, const Vector &za_grid, const Vector &aa_grid, const Index &doit_za_grid_size, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Verbosity &verbosity)
WORKSPACE METHOD: DoitInit.
void doit_scat_fieldCalc(Workspace &ws, Tensor6 &doit_scat_field, const Agenda &pha_mat_spt_agenda, const Tensor6 &cloudbox_field_mono, const Tensor4 &pnd_field, const Tensor3 &t_field, const Index &atmosphere_dim, const ArrayOfIndex &cloudbox_limits, const Vector &za_grid, const Vector &aa_grid, const Index &doit_za_grid_size, const Tensor7 &pha_mat_doit, const Verbosity &verbosity)
WORKSPACE METHOD: doit_scat_fieldCalc.
Template functions for general supergeneric ws methods.
void pha_matCalc(Tensor4 &pha_mat, const Tensor5 &pha_mat_spt, const Tensor4 &pnd_field, const Index &atmosphere_dim, const Index &scat_p_index, const Index &scat_lat_index, const Index &scat_lon_index, const Verbosity &)
WORKSPACE METHOD: pha_matCalc.
void nlinspace(Vector &x, const Numeric start, const Numeric stop, const Index n)
nlinspace
Numeric AngIntegrate_trapezoid(ConstMatrixView Integrand, ConstVectorView za_grid, ConstVectorView aa_grid)
AngIntegrate_trapezoid.
Numeric AngIntegrate_trapezoid_opti(ConstMatrixView Integrand, ConstVectorView za_grid, ConstVectorView aa_grid, ConstVectorView grid_stepsize)
AngIntegrate_trapezoid_opti.
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 auto rad2deg(auto x) noexcept
Converts radians to degrees.
Numeric invrayjean(const Numeric &i, const Numeric &f)
invrayjean
This file contains declerations of functions of physical character.
Propagation path structure and functions.
Declaration of functions in rte.cc.
void p2gridpos(ArrayOfGridPos &gp, ConstVectorView old_pgrid, ConstVectorView new_pgrid, const Numeric &extpolfac)
Calculates grid positions for pressure values.
Header file for special_interp.cc.
The structure to describe a propagation path and releated quantities.
Auxiliary header stuff related to workspace variable groups.
This file contains basic functions to handle XML data files.
void xml_write_to_file(const String &filename, const T &type, const FileType ftype, const Index no_clobber, const Verbosity &verbosity)
Write data to XML file.
void xml_read_from_file(const String &filename, T &type, const Verbosity &verbosity)
Reads data from XML file.