64 const Index& N_za_grid,
65 const Index& N_aa_grid,
66 const String& za_grid_type,
71 else if (N_aa_grid < 1) {
73 os <<
"N_aa_grid must be > 0 (even for 1D).";
74 throw std::runtime_error(os.str());
80 if (N_za_grid % 2 == 1) {
82 os <<
"N_za_grid must be even.";
83 throw runtime_error(os.str());
86 Index nph = N_za_grid / 2;
91 za_grid_weights.
resize(N_za_grid);
94 if (za_grid_type ==
"double_gauss") {
111 x[xtemp.
nelem() - 1] = acos((xtemp[0] + 1) / 2) /
DEG2RAD;
112 w[wtemp.
nelem() - 1] = wtemp[0] / 2;
114 for (
Index i = 0; i < xtemp.
nelem() - 1; i++) {
115 x[i] = acos((xtemp[xtemp.
nelem() - 1 - i] + 1) / 2.) /
DEG2RAD;
116 x[xtemp.
nelem() + i] = acos(1 - (xtemp[i + 1] + 1) / 2.) /
DEG2RAD;
118 w[i] = wtemp[wtemp.
nelem() - 1 - i] / 2;
119 w[wtemp.
nelem() + i] = wtemp[i + 1] / 2;
123 x[i] = acos((xtemp[xtemp.
nelem() - 1 - i] + 1) / 2.) /
DEG2RAD;
124 x[xtemp.
nelem() + i] = acos(1 - (xtemp[i] + 1) / 2.) /
DEG2RAD;
126 w[i] = wtemp[wtemp.
nelem() - 1 - i] / 2;
127 w[wtemp.
nelem() + i] = wtemp[i] / 2;
131 for (
Index i = 0; i < nph; i++) {
135 za_grid[za_grid.
nelem() - 1 - i] = 180 - x[i];
138 za_grid_weights[i] =
w[i];
139 za_grid_weights[za_grid_weights.
nelem() - 1 - i] =
w[i];
142 }
else if (za_grid_type ==
"linear") {
147 for (
Index i = 0; i < N_za_grid; i++) {
148 za_grid[i] = (x[i] + 1) * 90.;
152 za_grid_weights[i] =
w[i] * sin(za_grid[i] *
DEG2RAD);
154 }
else if (za_grid_type ==
"linear_mu") {
165 for (
Index i = 0; i < N_za_grid; i++) {
166 za_grid_temp[i] = acos(x[i]) /
DEG2RAD;
175 os <<
"The selected grid type is not implemented";
176 throw std::runtime_error(os.str());
182 if (za_grid[0] < 0) {
186 if (za_grid[za_grid.
nelem() - 1] > 180) {
187 za_grid[za_grid.
nelem() - 1] = 180.;
194 const Tensor4& irradiance_field,
195 const Tensor3& specific_heat_capacity,
200 irradiance_field.
npages(),
201 irradiance_field.
nrows());
214 for (
Index p = 0; p < irradiance_field.
npages(); p++) {
215 for (
Index r = 0; r < irradiance_field.
nrows(); r++) {
216 net_flux_b = (irradiance_field(
b - 1, p, r, 0) +
217 irradiance_field(
b - 1, p, r, 1));
218 net_flux_t = (irradiance_field(
b + 1, p, r, 0) +
219 irradiance_field(
b + 1, p, r, 1));
221 heating_rates(
b, p, r) = (net_flux_t - net_flux_b) /
222 (p_grid[
b + 1] - p_grid[
b - 1]) * g0 /
223 specific_heat_capacity(
b, p, r);
228 idx = irradiance_field.
nbooks();
231 for (
Index p = 0; p < irradiance_field.
npages(); p++) {
232 for (
Index r = 0; r < irradiance_field.
nrows(); r++) {
235 (irradiance_field(0, p, r, 0) + irradiance_field(0, p, r, 1));
237 (irradiance_field(1, p, r, 0) + irradiance_field(1, p, r, 1));
239 (irradiance_field(2, p, r, 0) + irradiance_field(0, p, r, 1));
241 heating_rates(0, p, r) = (-3 * net_flux_b + 4 * net_flux_c - net_flux_t) /
242 (p_grid[2] - p_grid[0]) * g0 /
243 specific_heat_capacity(0, p, r);
246 net_flux_t = (irradiance_field(idx - 1, p, r, 0) +
247 irradiance_field(idx - 1, p, r, 1));
248 net_flux_c = (irradiance_field(idx - 2, p, r, 0) +
249 irradiance_field(idx - 2, p, r, 1));
250 net_flux_b = (irradiance_field(idx - 3, p, r, 0) +
251 irradiance_field(idx - 3, p, r, 1));
253 heating_rates(idx - 1, p, r) =
254 -(-3 * net_flux_t + 4 * net_flux_c - net_flux_b) /
255 (p_grid[idx-1] - p_grid[idx-3]) * g0 / specific_heat_capacity(0, p, r);
265 const Vector& za_grid_weights,
271 Tensor4 radiance_field_aa_integrated;
276 radiance_field_aa_integrated =
278 radiance_field_aa_integrated *= 2 *
PI;
285 radiance_field.
nrows());
286 radiance_field_aa_integrated = 0.;
288 for (
Index b = 0;
b < radiance_field_aa_integrated.
nbooks();
b++) {
289 for (
Index p = 0; p < radiance_field_aa_integrated.
npages(); p++) {
290 for (
Index r = 0; r < radiance_field_aa_integrated.
nrows(); r++) {
291 for (
Index c = 0;
c < radiance_field_aa_integrated.
ncols();
c++) {
292 for (
Index i = 0; i < N_scat_aa - 1; i++) {
293 radiance_field_aa_integrated(
b, p, r,
c) +=
294 (radiance_field(
b, p, r,
c, i) +
295 radiance_field(
b, p, r,
c, i + 1)) /
296 2. *
abs(aa_grid[i + 1] - aa_grid[i]) *
DEG2RAD;
309 irradiance_field = 0;
314 for (
Index p = 0; p < irradiance_field.
npages(); p++) {
315 for (
Index r = 0; r < irradiance_field.
nrows(); r++) {
316 for (
Index i = 0; i < N_scat_za; i++) {
317 if (za_grid[i] <= 90.) {
318 irradiance_field(
b, p, r, 0) +=
319 radiance_field_aa_integrated(
b, p, r, i) *
320 cos(za_grid[i] *
DEG2RAD) * (-1.) * za_grid_weights[i];
322 irradiance_field(
b, p, r, 1) +=
323 radiance_field_aa_integrated(
b, p, r, i) *
324 cos(za_grid[i] *
DEG2RAD) * (-1.) * za_grid_weights[i];
335 const Tensor5& spectral_radiation_field,
339 "The length of f_grid does not match with\n"
340 " the first dimension of the spectral_radiation_field");
344 radiation_field.
resize(spectral_radiation_field.
nbooks(),
345 spectral_radiation_field.
npages(),
346 spectral_radiation_field.
nrows(),
347 spectral_radiation_field.
ncols());
351 for (
Index i = 0; i < spectral_radiation_field.
nshelves() - 1; i++) {
352 const Numeric df = f_grid[i + 1] - f_grid[i];
355 for (
Index p = 0; p < radiation_field.
npages(); p++) {
356 for (
Index r = 0; r < radiation_field.
nrows(); r++) {
358 radiation_field(
b, p, r,
c) +=
359 (spectral_radiation_field(i + 1,
b, p, r,
c) +
360 spectral_radiation_field(i,
b, p, r,
c)) /
372 const Tensor7& spectral_radiation_field,
376 "The length of f_grid does not match with\n"
377 " the first dimension of the spectral_radiation_field");
382 spectral_radiation_field.
nshelves(),
383 spectral_radiation_field.
nbooks(),
384 spectral_radiation_field.
npages(),
385 spectral_radiation_field.
nrows());
389 for (
Index i = 0; i < spectral_radiation_field.
nlibraries() - 1; i++) {
390 const Numeric df = f_grid[i + 1] - f_grid[i];
394 for (
Index p = 0; p < radiation_field.
npages(); p++) {
395 for (
Index r = 0; r < radiation_field.
nrows(); r++) {
397 radiation_field(s,
b, p, r,
c) +=
398 (spectral_radiation_field(i + 1, s,
b, p, r,
c, 0) +
399 spectral_radiation_field(i, s,
b, p, r,
c, 0)) /
411 Tensor5& spectral_irradiance_field,
412 const Tensor7& spectral_radiance_field,
415 const Vector& za_grid_weights,
418 const Index N_scat_za = spectral_radiance_field.
npages();
419 const Index N_scat_aa = spectral_radiance_field.
nrows();
421 Tensor5 iy_field_aa_integrated;
426 iy_field_aa_integrated =
428 iy_field_aa_integrated *= 2 *
PI;
435 spectral_radiance_field.
nbooks(),
436 spectral_radiance_field.
npages());
437 iy_field_aa_integrated = 0.;
439 for (
Index s = 0; s < iy_field_aa_integrated.
nshelves(); s++) {
441 for (
Index p = 0; p < iy_field_aa_integrated.
npages(); p++) {
442 for (
Index r = 0; r < iy_field_aa_integrated.
nrows(); r++) {
444 for (
Index i = 0; i < N_scat_aa - 1; i++) {
445 iy_field_aa_integrated(s,
b, p, r,
c) +=
446 (spectral_radiance_field(s,
b, p, r,
c, i, 0) +
447 spectral_radiance_field(s,
b, p, r,
c, i + 1, 0)) /
448 2. *
abs(aa_grid[i + 1] - aa_grid[i]) *
DEG2RAD;
461 spectral_radiance_field.
nbooks(),
463 spectral_irradiance_field = 0;
466 for (
Index s = 0; s < spectral_irradiance_field.
nshelves(); s++) {
468 for (
Index p = 0; p < spectral_irradiance_field.
npages(); p++) {
469 for (
Index r = 0; r < spectral_irradiance_field.
nrows(); r++) {
470 for (
Index i = 0; i < N_scat_za; i++) {
471 if (za_grid[i] <= 90.) {
472 spectral_irradiance_field(s,
b, p, r, 0) +=
473 iy_field_aa_integrated(s,
b, p, r, i) *
474 cos(za_grid[i] *
DEG2RAD) * (-1.) * za_grid_weights[i];
476 spectral_irradiance_field(s,
b, p, r, 1) +=
477 iy_field_aa_integrated(s,
b, p, r, i) *
478 cos(za_grid[i] *
DEG2RAD) * (-1.) * za_grid_weights[i];
490 Tensor7& spectral_radiance_field,
492 const Agenda& propmat_clearsky_agenda,
493 const Agenda& water_p_eq_agenda,
494 const Agenda& iy_space_agenda,
495 const Agenda& iy_surface_agenda,
496 const Agenda& iy_cloudbox_agenda,
497 const Index& stokes_dim,
499 const Index& atmosphere_dim,
515 const String& rt_integration_option,
516 const Tensor3& surface_props_data,
518 const Index& use_parallel_za [[maybe_unused]],
521 if (atmosphere_dim != 1)
522 throw runtime_error(
"This method only works for atmosphere_dim = 1.");
530 spectral_radiance_field.
resize(nf, nl, 1, 1, nza, 1, stokes_dim);
531 trans_field.
resize(nf, nl, nza);
534 const Index cloudbox_on = 0, ppath_inside_cloudbox_do = 0;
538 const String iy_unit =
"1";
541 const Index iy_agenda_call1 = 1;
542 const Tensor3 iy_transmittance(0, 0, 0);
543 const Index jacobian_do = 0;
546 const Numeric z_space = z_field(nl - 1, 0, 0) + 10;
564#pragma omp parallel for if (!arts_omp_in_parallel() && nza > 1 && \
565 use_parallel_za) firstprivate(wss)
566 for (
Index i = 0; i < nza; i++) {
567 if (failed)
continue;
572 Matrix iy, ppvar_vmr, ppvar_wind, ppvar_mag, ppvar_f;
575 Tensor4 ppvar_trans_cumulat, ppvar_trans_partial;
580 Vector rte_los(1, za_grid[i]);
581 Vector rte_pos(1, za_grid[i] < 90 ? z_surface(0, 0) : z_space);
589 ppath_inside_cloudbox_do,
595 ppath.
gp_p[ppath.
np - 1].idx == nl - 2);
633 propmat_clearsky_agenda,
635 rt_integration_option,
648 if (za_grid[i] < 90) {
649 spectral_radiance_field(
joker, i0, 0, 0, i, 0,
joker) =
651 spectral_radiance_field(
joker, nl - 1, 0, 0, i, 0,
joker) =
653 trans_field(
joker, 0, i) = ppvar_trans_partial(0,
joker, 0, 0);
654 trans_field(
joker, nl - 1, i) =
655 ppvar_trans_partial(ppath.
np - 1,
joker, 0, 0);
657 spectral_radiance_field(
joker, nl - 1, 0, 0, i, 0,
joker) =
659 spectral_radiance_field(
joker, i0, 0, 0, i, 0,
joker) =
661 trans_field(
joker, nl - 1, i) = ppvar_trans_partial(0,
joker, 0, 0);
662 trans_field(
joker, 0, i) =
663 ppvar_trans_partial(ppath.
np - 1,
joker, 0, 0);
667 for (
Index p = 1; p < ppath.
np - 1; p++) {
669 if (ppath.
gp_p[p].fd[0] < 1e-2) {
670 spectral_radiance_field(
673 trans_field(
joker, ppath.
gp_p[p].idx, i) =
674 ppvar_trans_partial(p,
joker, 0, 0);
681 for (
Index p = 0; p < i0; p++) {
682 spectral_radiance_field(
joker, p, 0, 0, i, 0,
joker) =
683 spectral_radiance_field(
joker, i0, 0, 0, i, 0,
joker);
684 trans_field(
joker, p, i) = trans_field(
joker, i0, i);
687 }
catch (
const std::exception& e) {
688#pragma omp critical(planep_setabort)
692 os <<
"Run-time error at nza #" << i <<
": \n" << e.what();
693#pragma omp critical(planep_push_fail_msg)
694 fail_msg.push_back(os.str());
699 if (fail_msg.
nelem()) {
701 for (
auto& msg : fail_msg) os << msg <<
'\n';
702 throw runtime_error(os.str());
708 Tensor7& spectral_radiance_field,
709 const Index& atmosphere_dim,
711 const Index& cloudbox_on,
715 if (atmosphere_dim > 1)
716 throw runtime_error(
"This method can only be used for 1D calculations.\n");
718 throw runtime_error(
"Cloudbox is off. This is not handled by this method.");
719 if (cloudbox_limits[0] || cloudbox_limits[1] != p_grid.
nelem() - 1)
721 "The cloudbox must cover all pressure levels "
722 "to use this method.");
725 spectral_radiance_field = cloudbox_field;
731 Tensor7& spectral_radiance_field,
732 const Agenda& propmat_clearsky_agenda,
733 const Agenda& water_p_eq_agenda,
734 const Agenda& iy_space_agenda,
735 const Agenda& iy_surface_agenda,
736 const Agenda& iy_cloudbox_agenda,
737 const Index& stokes_dim,
739 const Index& atmosphere_dim,
753 const Index& cloudbox_on,
758 const String& rt_integration_option,
759 const Tensor3& surface_props_data,
761 const Index& use_parallel_za [[maybe_unused]],
764 if (atmosphere_dim != 1)
765 throw runtime_error(
"This method only works for atmosphere_dim = 1.");
767 throw runtime_error(
"No ned to use this method with cloudbox=0.");
768 if (cloudbox_limits[0])
770 "The first element of *cloudbox_limits* must be zero "
771 "to use this method.");
779 spectral_radiance_field.
resize(nf, nl, 1, 1, nza, 1, stokes_dim);
782 spectral_radiance_field(
joker,
783 Range(0, cloudbox_limits[1] + 1),
788 joker) = cloudbox_field;
791 const Index ppath_inside_cloudbox_do = 0;
792 const String iy_unit =
"1";
795 const Index iy_agenda_call1 = 1;
796 const Tensor3 iy_transmittance(0, 0, 0);
797 const Index jacobian_do = 0;
800 const Numeric z_space = z_field(nl - 1, 0, 0) + 10;
810 const Index i0 = cloudbox_limits[1];
811 const Numeric z_top = z_field(i0 + 1, 0, 0);
818#pragma omp parallel for if (!arts_omp_in_parallel() && nza > 1 && \
819 use_parallel_za) firstprivate(wss)
820 for (
Index i = 0; i < nza; i++) {
821 if (failed)
continue;
826 Matrix iy, ppvar_vmr, ppvar_wind, ppvar_mag, ppvar_f;
829 Tensor4 ppvar_trans_cumulat, ppvar_trans_partial;
834 Vector rte_los(1, za_grid[i]);
835 Vector rte_pos(1, za_grid[i] < 90 ? z_top : z_space);
843 ppath_inside_cloudbox_do,
849 ppath.
gp_p[ppath.
np - 1].idx == nl - 2);
887 propmat_clearsky_agenda,
889 rt_integration_option,
903 if (za_grid[i] < 90) {
904 spectral_radiance_field(
joker, i0 + 1, 0, 0, i, 0,
joker) =
906 spectral_radiance_field(
joker, nl - 1, 0, 0, i, 0,
joker) =
909 spectral_radiance_field(
joker, nl - 1, 0, 0, i, 0,
joker) =
914 for (
Index p = 1; p < ppath.
np - 1; p++) {
916 if (ppath.
gp_p[p].fd[0] < 1e-2) {
917 spectral_radiance_field(
923 }
catch (
const std::exception& e) {
924#pragma omp critical(planep_setabort)
928 os <<
"Run-time error at nza #" << i <<
": \n" << e.what();
929#pragma omp critical(planep_push_fail_msg)
930 fail_msg.push_back(os.str());
935 if (fail_msg.
nelem()) {
937 for (
auto& msg : fail_msg) os << msg <<
'\n';
938 throw runtime_error(os.str());
Declarations required for the calculation of absorption coefficients.
Declarations for agendas.
Constants of physical expressions as constexpr.
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.
#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.
constexpr Numeric DEG2RAD
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
void abs(Sparse &A, const Sparse &B)
Absolute value of sparse matrix elements.
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.
Agenda get_iy_main_agenda(Workspace &ws, const String &option)
constexpr Numeric pi
The following mathematical constants are generated in python Decimal package by the code:
constexpr auto deg2rad(auto x) noexcept
Converts degrees to radians.
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.