61 const Index& N_za_grid,
62 const Index& N_aa_grid,
63 const String& za_grid_type,
68 else if (N_aa_grid < 1) {
70 os <<
"N_aa_grid must be > 0 (even for 1D).";
71 throw std::runtime_error(os.str());
77 if (N_za_grid % 2 == 1) {
79 os <<
"N_za_grid must be even.";
80 throw runtime_error(os.str());
83 Index nph = N_za_grid / 2;
88 za_grid_weights.
resize(N_za_grid);
91 if (za_grid_type ==
"double_gauss") {
108 x[xtemp.
nelem() - 1] = acos((xtemp[0] + 1) / 2) /
DEG2RAD;
109 w[wtemp.
nelem() - 1] = wtemp[0] / 2;
111 for (
Index i = 0; i < xtemp.
nelem() - 1; i++) {
112 x[i] = acos((xtemp[xtemp.
nelem() - 1 - i] + 1) / 2.) /
DEG2RAD;
113 x[xtemp.
nelem() + i] = acos(1 - (xtemp[i + 1] + 1) / 2.) /
DEG2RAD;
115 w[i] = wtemp[wtemp.
nelem() - 1 - i] / 2;
116 w[wtemp.
nelem() + i] = wtemp[i + 1] / 2;
120 x[i] = acos((xtemp[xtemp.
nelem() - 1 - i] + 1) / 2.) /
DEG2RAD;
121 x[xtemp.
nelem() + i] = acos(1 - (xtemp[i] + 1) / 2.) /
DEG2RAD;
123 w[i] = wtemp[wtemp.
nelem() - 1 - i] / 2;
124 w[wtemp.
nelem() + i] = wtemp[i] / 2;
128 for (
Index i = 0; i < nph; i++) {
132 za_grid[za_grid.
nelem() - 1 - i] = 180 - x[i];
135 za_grid_weights[i] =
w[i];
136 za_grid_weights[za_grid_weights.
nelem() - 1 - i] =
w[i];
139 }
else if (za_grid_type ==
"linear") {
144 for (
Index i = 0; i < N_za_grid; i++) {
145 za_grid[i] = (x[i] + 1) * 90.;
149 za_grid_weights[i] =
w[i] * sin(za_grid[i] *
DEG2RAD);
151 }
else if (za_grid_type ==
"linear_mu") {
162 for (
Index i = 0; i < N_za_grid; i++) {
163 za_grid_temp[i] = acos(x[i]) /
DEG2RAD;
172 os <<
"The selected grid type is not implemented";
173 throw std::runtime_error(os.str());
179 if (za_grid[0] < 0) {
183 if (za_grid[za_grid.
nelem() - 1] > 180) {
184 za_grid[za_grid.
nelem() - 1] = 180.;
191 const Tensor4& irradiance_field,
192 const Tensor3& specific_heat_capacity,
197 irradiance_field.
npages(),
198 irradiance_field.
nrows());
211 for (
Index p = 0; p < irradiance_field.
npages(); p++) {
212 for (
Index r = 0; r < irradiance_field.
nrows(); r++) {
213 net_flux_b = (irradiance_field(
b - 1, p, r, 0) +
214 irradiance_field(
b - 1, p, r, 1));
215 net_flux_t = (irradiance_field(
b + 1, p, r, 0) +
216 irradiance_field(
b + 1, p, r, 1));
218 heating_rates(
b, p, r) = (net_flux_t - net_flux_b) /
219 (p_grid[
b + 1] - p_grid[
b - 1]) * g0 /
220 specific_heat_capacity(
b, p, r);
225 idx = irradiance_field.
nbooks();
228 for (
Index p = 0; p < irradiance_field.
npages(); p++) {
229 for (
Index r = 0; r < irradiance_field.
nrows(); r++) {
232 (irradiance_field(0, p, r, 0) + irradiance_field(0, p, r, 1));
234 (irradiance_field(1, p, r, 0) + irradiance_field(1, p, r, 1));
236 (irradiance_field(2, p, r, 0) + irradiance_field(0, p, r, 1));
238 heating_rates(0, p, r) = (-3 * net_flux_b + 4 * net_flux_c - net_flux_t) /
239 (p_grid[2] - p_grid[0]) * g0 /
240 specific_heat_capacity(0, p, r);
243 net_flux_t = (irradiance_field(idx - 1, p, r, 0) +
244 irradiance_field(idx - 1, p, r, 1));
245 net_flux_c = (irradiance_field(idx - 2, p, r, 0) +
246 irradiance_field(idx - 2, p, r, 1));
247 net_flux_b = (irradiance_field(idx - 3, p, r, 0) +
248 irradiance_field(idx - 3, p, r, 1));
250 heating_rates(idx - 1, p, r) =
251 -(-3 * net_flux_t + 4 * net_flux_c - net_flux_b) /
252 (p_grid[2] - p_grid[0]) * g0 / specific_heat_capacity(0, p, r);
262 const Vector& za_grid_weights,
268 Tensor4 radiance_field_aa_integrated;
273 radiance_field_aa_integrated =
275 radiance_field_aa_integrated *= 2 *
PI;
282 radiance_field.
nrows());
283 radiance_field_aa_integrated = 0.;
285 for (
Index b = 0;
b < radiance_field_aa_integrated.
nbooks();
b++) {
286 for (
Index p = 0; p < radiance_field_aa_integrated.
npages(); p++) {
287 for (
Index r = 0; r < radiance_field_aa_integrated.
nrows(); r++) {
288 for (
Index c = 0;
c < radiance_field_aa_integrated.
ncols();
c++) {
289 for (
Index i = 0; i < N_scat_aa - 1; i++) {
290 radiance_field_aa_integrated(
b, p, r,
c) +=
291 (radiance_field(
b, p, r,
c, i) +
292 radiance_field(
b, p, r,
c, i + 1)) /
293 2. *
abs(aa_grid[i + 1] - aa_grid[i]) *
DEG2RAD;
306 irradiance_field = 0;
311 for (
Index p = 0; p < irradiance_field.
npages(); p++) {
312 for (
Index r = 0; r < irradiance_field.
nrows(); r++) {
313 for (
Index i = 0; i < N_scat_za; i++) {
314 if (za_grid[i] <= 90.) {
315 irradiance_field(
b, p, r, 0) +=
316 radiance_field_aa_integrated(
b, p, r, i) *
317 cos(za_grid[i] *
DEG2RAD) * (-1.) * za_grid_weights[i];
319 irradiance_field(
b, p, r, 1) +=
320 radiance_field_aa_integrated(
b, p, r, i) *
321 cos(za_grid[i] *
DEG2RAD) * (-1.) * za_grid_weights[i];
332 const Tensor5& spectral_radiation_field,
336 "The length of f_grid does not match with\n"
337 " the first dimension of the spectral_radiation_field");
341 radiation_field.
resize(spectral_radiation_field.
nbooks(),
342 spectral_radiation_field.
npages(),
343 spectral_radiation_field.
nrows(),
344 spectral_radiation_field.
ncols());
348 for (
Index i = 0; i < spectral_radiation_field.
nshelves() - 1; i++) {
349 const Numeric df = f_grid[i + 1] - f_grid[i];
352 for (
Index p = 0; p < radiation_field.
npages(); p++) {
353 for (
Index r = 0; r < radiation_field.
nrows(); r++) {
355 radiation_field(
b, p, r,
c) +=
356 (spectral_radiation_field(i + 1,
b, p, r,
c) +
357 spectral_radiation_field(i,
b, p, r,
c)) /
369 const Tensor7& spectral_radiation_field,
373 "The length of f_grid does not match with\n"
374 " the first dimension of the spectral_radiation_field");
379 spectral_radiation_field.
nshelves(),
380 spectral_radiation_field.
nbooks(),
381 spectral_radiation_field.
npages(),
382 spectral_radiation_field.
nrows());
386 for (
Index i = 0; i < spectral_radiation_field.
nlibraries() - 1; i++) {
387 const Numeric df = f_grid[i + 1] - f_grid[i];
391 for (
Index p = 0; p < radiation_field.
npages(); p++) {
392 for (
Index r = 0; r < radiation_field.
nrows(); r++) {
394 radiation_field(s,
b, p, r,
c) +=
395 (spectral_radiation_field(i + 1, s,
b, p, r,
c, 0) +
396 spectral_radiation_field(i, s,
b, p, r,
c, 0)) /
408 Tensor5& spectral_irradiance_field,
409 const Tensor7& spectral_radiance_field,
412 const Vector& za_grid_weights,
415 const Index N_scat_za = spectral_radiance_field.
npages();
416 const Index N_scat_aa = spectral_radiance_field.
nrows();
418 Tensor5 iy_field_aa_integrated;
423 iy_field_aa_integrated =
425 iy_field_aa_integrated *= 2 *
PI;
432 spectral_radiance_field.
nbooks(),
433 spectral_radiance_field.
npages());
434 iy_field_aa_integrated = 0.;
436 for (
Index s = 0; s < iy_field_aa_integrated.
nshelves(); s++) {
438 for (
Index p = 0; p < iy_field_aa_integrated.
npages(); p++) {
439 for (
Index r = 0; r < iy_field_aa_integrated.
nrows(); r++) {
441 for (
Index i = 0; i < N_scat_aa - 1; i++) {
442 iy_field_aa_integrated(s,
b, p, r,
c) +=
443 (spectral_radiance_field(s,
b, p, r,
c, i, 0) +
444 spectral_radiance_field(s,
b, p, r,
c, i + 1, 0)) /
445 2. *
abs(aa_grid[i + 1] - aa_grid[i]) *
DEG2RAD;
458 spectral_radiance_field.
nbooks(),
460 spectral_irradiance_field = 0;
463 for (
Index s = 0; s < spectral_irradiance_field.
nshelves(); s++) {
465 for (
Index p = 0; p < spectral_irradiance_field.
npages(); p++) {
466 for (
Index r = 0; r < spectral_irradiance_field.
nrows(); r++) {
467 for (
Index i = 0; i < N_scat_za; i++) {
468 if (za_grid[i] <= 90.) {
469 spectral_irradiance_field(s,
b, p, r, 0) +=
470 iy_field_aa_integrated(s,
b, p, r, i) *
471 cos(za_grid[i] *
DEG2RAD) * (-1.) * za_grid_weights[i];
473 spectral_irradiance_field(s,
b, p, r, 1) +=
474 iy_field_aa_integrated(s,
b, p, r, i) *
475 cos(za_grid[i] *
DEG2RAD) * (-1.) * za_grid_weights[i];
487 Tensor7& spectral_radiance_field,
489 const Agenda& propmat_clearsky_agenda,
490 const Agenda& water_p_eq_agenda,
491 const Agenda& iy_space_agenda,
492 const Agenda& iy_surface_agenda,
493 const Agenda& iy_cloudbox_agenda,
494 const Index& stokes_dim,
496 const Index& atmosphere_dim,
512 const String& rt_integration_option,
513 const Tensor3& surface_props_data,
515 const Index& use_parallel_za [[maybe_unused]],
518 if (atmosphere_dim != 1)
519 throw runtime_error(
"This method only works for atmosphere_dim = 1.");
527 spectral_radiance_field.
resize(nf, nl, 1, 1, nza, 1, stokes_dim);
528 trans_field.
resize(nf, nl, nza);
531 const Index cloudbox_on = 0, ppath_inside_cloudbox_do = 0;
535 const String iy_unit =
"1";
538 const Index iy_agenda_call1 = 1;
539 const Tensor3 iy_transmittance(0, 0, 0);
540 const Index jacobian_do = 0;
543 const Numeric z_space = z_field(nl - 1, 0, 0) + 10;
551 Agenda iy_main_agenda(ws);
559 iy_main_agenda.
set_name(
"iy_main_agenda");
560 iy_main_agenda.
check(ws, verbosity);
570#pragma omp parallel for if (!arts_omp_in_parallel() && nza > 1 && \
571 use_parallel_za) firstprivate(wss)
572 for (
Index i = 0; i < nza; i++) {
573 if (failed)
continue;
578 Matrix iy, ppvar_vmr, ppvar_wind, ppvar_mag, ppvar_f;
581 Tensor4 ppvar_trans_cumulat, ppvar_trans_partial;
586 Vector rte_los(1, za_grid[i]);
587 Vector rte_pos(1, za_grid[i] < 90 ? z_surface(0, 0) : z_space);
595 ppath_inside_cloudbox_do,
601 ppath.
gp_p[ppath.
np - 1].idx == nl - 2);
639 propmat_clearsky_agenda,
641 rt_integration_option,
654 if (za_grid[i] < 90) {
655 spectral_radiance_field(
joker, i0, 0, 0, i, 0,
joker) =
657 spectral_radiance_field(
joker, nl - 1, 0, 0, i, 0,
joker) =
659 trans_field(
joker, 0, i) = ppvar_trans_partial(0,
joker, 0, 0);
660 trans_field(
joker, nl - 1, i) =
661 ppvar_trans_partial(ppath.
np - 1,
joker, 0, 0);
663 spectral_radiance_field(
joker, nl - 1, 0, 0, i, 0,
joker) =
665 spectral_radiance_field(
joker, i0, 0, 0, i, 0,
joker) =
667 trans_field(
joker, nl - 1, i) = ppvar_trans_partial(0,
joker, 0, 0);
668 trans_field(
joker, 0, i) =
669 ppvar_trans_partial(ppath.
np - 1,
joker, 0, 0);
673 for (
Index p = 1; p < ppath.
np - 1; p++) {
675 if (ppath.
gp_p[p].fd[0] < 1e-2) {
676 spectral_radiance_field(
679 trans_field(
joker, ppath.
gp_p[p].idx, i) =
680 ppvar_trans_partial(p,
joker, 0, 0);
687 for (
Index p = 0; p < i0; p++) {
688 spectral_radiance_field(
joker, p, 0, 0, i, 0,
joker) =
689 spectral_radiance_field(
joker, i0, 0, 0, i, 0,
joker);
690 trans_field(
joker, p, i) = trans_field(
joker, i0, i);
693 }
catch (
const std::exception& e) {
694#pragma omp critical(planep_setabort)
698 os <<
"Run-time error at nza #" << i <<
": \n" << e.what();
699#pragma omp critical(planep_push_fail_msg)
700 fail_msg.push_back(os.str());
705 if (fail_msg.
nelem()) {
707 for (
auto& msg : fail_msg) os << msg <<
'\n';
708 throw runtime_error(os.str());
714 Tensor7& spectral_radiance_field,
715 const Index& atmosphere_dim,
717 const Index& cloudbox_on,
721 if (atmosphere_dim > 1)
722 throw runtime_error(
"This method can only be used for 1D calculations.\n");
724 throw runtime_error(
"Cloudbox is off. This is not handled by this method.");
725 if (cloudbox_limits[0] || cloudbox_limits[1] != p_grid.
nelem() - 1)
727 "The cloudbox must cover all pressure levels "
728 "to use this method.");
731 spectral_radiance_field = cloudbox_field;
737 Tensor7& spectral_radiance_field,
738 const Agenda& propmat_clearsky_agenda,
739 const Agenda& water_p_eq_agenda,
740 const Agenda& iy_space_agenda,
741 const Agenda& iy_surface_agenda,
742 const Agenda& iy_cloudbox_agenda,
743 const Index& stokes_dim,
745 const Index& atmosphere_dim,
759 const Index& cloudbox_on,
764 const String& rt_integration_option,
765 const Tensor3& surface_props_data,
767 const Index& use_parallel_za [[maybe_unused]],
770 if (atmosphere_dim != 1)
771 throw runtime_error(
"This method only works for atmosphere_dim = 1.");
773 throw runtime_error(
"No ned to use this method with cloudbox=0.");
774 if (cloudbox_limits[0])
776 "The first element of *cloudbox_limits* must be zero "
777 "to use this method.");
785 spectral_radiance_field.
resize(nf, nl, 1, 1, nza, 1, stokes_dim);
788 spectral_radiance_field(
joker,
789 Range(0, cloudbox_limits[1] + 1),
794 joker) = cloudbox_field;
797 const Index ppath_inside_cloudbox_do = 0;
798 const String iy_unit =
"1";
801 const Index iy_agenda_call1 = 1;
802 const Tensor3 iy_transmittance(0, 0, 0);
803 const Index jacobian_do = 0;
806 const Numeric z_space = z_field(nl - 1, 0, 0) + 10;
813 Agenda iy_main_agenda{ws};
814 iy_main_agenda.append(
"ppathPlaneParallel",
TokVal());
815 iy_main_agenda.append(
"iyEmissionStandard",
TokVal());
821 iy_main_agenda.set_name(
"iy_main_agenda");
822 iy_main_agenda.check(ws, verbosity);
825 const Index i0 = cloudbox_limits[1];
826 const Numeric z_top = z_field(i0 + 1, 0, 0);
833#pragma omp parallel for if (!arts_omp_in_parallel() && nza > 1 && \
834 use_parallel_za) firstprivate(wss)
835 for (
Index i = 0; i < nza; i++) {
836 if (failed)
continue;
841 Matrix iy, ppvar_vmr, ppvar_wind, ppvar_mag, ppvar_f;
844 Tensor4 ppvar_trans_cumulat, ppvar_trans_partial;
849 Vector rte_los(1, za_grid[i]);
850 Vector rte_pos(1, za_grid[i] < 90 ? z_top : z_space);
858 ppath_inside_cloudbox_do,
864 ppath.
gp_p[ppath.
np - 1].idx == nl - 2);
902 propmat_clearsky_agenda,
904 rt_integration_option,
918 if (za_grid[i] < 90) {
919 spectral_radiance_field(
joker, i0 + 1, 0, 0, i, 0,
joker) =
921 spectral_radiance_field(
joker, nl - 1, 0, 0, i, 0,
joker) =
924 spectral_radiance_field(
joker, nl - 1, 0, 0, i, 0,
joker) =
929 for (
Index p = 1; p < ppath.
np - 1; p++) {
931 if (ppath.
gp_p[p].fd[0] < 1e-2) {
932 spectral_radiance_field(
938 }
catch (
const std::exception& e) {
939#pragma omp critical(planep_setabort)
943 os <<
"Run-time error at nza #" << i <<
": \n" << e.what();
944#pragma omp critical(planep_push_fail_msg)
945 fail_msg.push_back(os.str());
950 if (fail_msg.
nelem()) {
952 for (
auto& msg : fail_msg) os << msg <<
'\n';
953 throw runtime_error(os.str());
Declarations required for the calculation of absorption coefficients.
Declarations for agendas.
void append(const String &methodname, const TokVal &keywordvalue)
Appends methods to an agenda.
void check(Workspace &ws_in, const Verbosity &verbosity)
Checks consistency of an agenda.
void push_back(const MRecord &n)
Append a new method to end of list.
void set_name(const String &nname)
Set agenda name.
Index nelem() const ARTS_NOEXCEPT
Index ncols() const noexcept
Index ncols() const noexcept
Index nrows() const noexcept
Index nbooks() const noexcept
Index npages() const noexcept
Index nrows() const noexcept
Index ncols() const noexcept
Index npages() const noexcept
Index nbooks() const noexcept
Index nshelves() const noexcept
Index npages() const noexcept
Index nrows() const noexcept
Index nlibraries() const noexcept
Index nvitrines() const noexcept
Index nshelves() const noexcept
Index nbooks() const noexcept
Index nelem() const noexcept
Returns the number of elements.
void resize(Index p, Index r, Index c)
Resize function.
void resize(Index b, Index p, Index r, Index c)
Resize function.
void resize(Index s, Index b, Index p, Index r, Index c)
Resize function.
void resize(Index l, Index v, Index s, Index b, Index p, Index r, Index c)
Resize function.
void resize(Index n)
Resize function.
std::shared_ptr< map< String, Index > > WsvMap_ptr
#define ARTS_ASSERT(condition,...)
bool gsl_integration_glfixed_table_alloc(Vector &x, Vector &w, Index n)
gsl_integration_glfixed_table_alloc
Contains the code to calculate Legendre polynomials.
void spectral_irradiance_fieldFromSpectralRadianceField(Tensor5 &spectral_irradiance_field, const Tensor7 &spectral_radiance_field, const Vector &za_grid, const Vector &aa_grid, const Vector &za_grid_weights, const Verbosity &)
WORKSPACE METHOD: spectral_irradiance_fieldFromSpectralRadianceField.
void heating_ratesFromIrradiance(Tensor3 &heating_rates, const Vector &p_grid, const Tensor4 &irradiance_field, const Tensor3 &specific_heat_capacity, const Numeric &g0, const Verbosity &)
WORKSPACE METHOD: heating_ratesFromIrradiance.
void irradiance_fieldFromRadiance(Tensor4 &irradiance_field, const Tensor5 &radiance_field, const Vector &za_grid, const Vector &aa_grid, const Vector &za_grid_weights, const Verbosity &)
WORKSPACE METHOD: irradiance_fieldFromRadiance.
void spectral_radiance_fieldClearskyPlaneParallel(Workspace &ws, Tensor7 &spectral_radiance_field, Tensor3 &trans_field, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &iy_space_agenda, const Agenda &iy_surface_agenda, const Agenda &iy_cloudbox_agenda, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Tensor3 &z_field, 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 Matrix &z_surface, const Numeric &ppath_lmax, const Numeric &rte_alonglos_v, const String &rt_integration_option, const Tensor3 &surface_props_data, const Vector &za_grid, const Index &use_parallel_za, const Verbosity &verbosity)
WORKSPACE METHOD: spectral_radiance_fieldClearskyPlaneParallel.
void spectral_radiance_fieldCopyCloudboxField(Tensor7 &spectral_radiance_field, const Index &atmosphere_dim, const Vector &p_grid, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Tensor7 &cloudbox_field, const Verbosity &)
WORKSPACE METHOD: spectral_radiance_fieldCopyCloudboxField.
void spectral_radiance_fieldExpandCloudboxField(Workspace &ws, Tensor7 &spectral_radiance_field, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const Agenda &iy_space_agenda, const Agenda &iy_surface_agenda, const Agenda &iy_cloudbox_agenda, const Index &stokes_dim, const Vector &f_grid, const Index &atmosphere_dim, const Vector &p_grid, const Tensor3 &z_field, 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 Matrix &z_surface, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Tensor7 &cloudbox_field, const Numeric &ppath_lmax, const Numeric &rte_alonglos_v, const String &rt_integration_option, const Tensor3 &surface_props_data, const Vector &za_grid, const Index &use_parallel_za, const Verbosity &verbosity)
WORKSPACE METHOD: spectral_radiance_fieldExpandCloudboxField.
void RadiationFieldSpectralIntegrate(Tensor4 &radiation_field, const Vector &f_grid, const Tensor5 &spectral_radiation_field, const Verbosity &)
WORKSPACE METHOD: RadiationFieldSpectralIntegrate.
void AngularGridsSetFluxCalc(Vector &za_grid, Vector &aa_grid, Vector &za_grid_weights, const Index &N_za_grid, const Index &N_aa_grid, const String &za_grid_type, const Verbosity &)
WORKSPACE METHOD: AngularGridsSetFluxCalc.
void ppathPlaneParallel(Ppath &ppath, const Index &atmosphere_dim, const Tensor3 &z_field, const Matrix &z_surface, const Index &cloudbox_on, const ArrayOfIndex &cloudbox_limits, const Index &ppath_inside_cloudbox_do, const Vector &rte_pos, const Vector &rte_los, const Numeric &ppath_lmax, const Verbosity &)
WORKSPACE METHOD: ppathPlaneParallel.
void iyEmissionStandard(Workspace &ws, Matrix &iy, ArrayOfMatrix &iy_aux, ArrayOfTensor3 &diy_dx, Vector &ppvar_p, Vector &ppvar_t, EnergyLevelMap &ppvar_nlte, Matrix &ppvar_vmr, Matrix &ppvar_wind, Matrix &ppvar_mag, Matrix &ppvar_f, Tensor3 &ppvar_iy, Tensor4 &ppvar_trans_cumulat, Tensor4 &ppvar_trans_partial, const Index &iy_id, 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 String &iy_unit, const ArrayOfString &iy_aux_vars, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Ppath &ppath, const Vector &rte_pos2, const Agenda &propmat_clearsky_agenda, const Agenda &water_p_eq_agenda, const String &rt_integration_option, const Agenda &iy_main_agenda, const Agenda &iy_space_agenda, const Agenda &iy_surface_agenda, const Agenda &iy_cloudbox_agenda, const Index &iy_agenda_call1, const Tensor3 &iy_transmittance, const Numeric &rte_alonglos_v, const Tensor3 &surface_props_data, const Verbosity &verbosity)
WORKSPACE METHOD: iyEmissionStandard.
void nlinspace(Vector &x, const Numeric start, const Numeric stop, const Index n)
nlinspace
void calculate_weights_linear(Vector &x, Vector &w, const Index nph)
calculate_weights_linear
NUMERIC Numeric
The type to use for all floating point numbers.
INDEX Index
The type to use for all integer numbers and indices.
Declarations having to do with the four output streams.
const map< String, Index > MdMap
The map associated with md_data.
Contains sorting routines.
The structure to describe a propagation path and releated quantities.
Index np
Number of points describing the ppath.
ArrayOfGridPos gp_p
Index position with respect to the pressure grid.
Index index_of_zsurface(const Numeric &z_surface, ConstVectorView z_profile)
Lccates the surface with respect to pressure levels.
This file contains the Workspace class.