ARTS 2.5.11 (git: 725533f0)
m_radiation_field.cc
Go to the documentation of this file.
1
9#include "absorption.h"
10#include "arts.h"
11#include "arts_omp.h"
12#include "auto_md.h"
13#include "debug.h"
14#include "lineshape.h"
15#include "logic.h"
16#include "physics_funcs.h"
17#include "ppath.h"
18#include "propmat_field.h"
19#include "radiation_field.h"
20
22 Workspace& ws,
23 Matrix& line_irradiance,
24 Tensor3& line_transmission,
25 const ArrayOfArrayOfSpeciesTag& abs_species,
26 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
27 const EnergyLevelMap& nlte_field,
28 const Tensor4& vmr_field,
29 const Tensor3& t_field,
30 const Tensor3& z_field,
31 const Vector& p_grid,
32 const Vector& refellipsoid,
33 const Tensor3& surface_props_data,
34 const Agenda& ppath_agenda,
35 const Agenda& iy_main_agenda,
36 const Agenda& iy_space_agenda,
37 const Agenda& iy_surface_agenda,
38 const Agenda& iy_cloudbox_agenda,
39 const Agenda& propmat_clearsky_agenda,
40 const Numeric& df,
41 const Index& nz,
42 const Index& nf,
43 const Numeric& r,
44 const Verbosity& verbosity)
45{
46 ARTS_USER_ERROR_IF (abs_lines_per_species.nelem() not_eq 1,
47 "Only for one species...");
48 ARTS_USER_ERROR_IF (nf % 2 not_eq 1,
49 "Must hit line center, nf % 2 must be 1.");
50 const Index nl = nelem(abs_lines_per_species);
51 const Index np = p_grid.nelem();
52
53 // Compute variables
54 ArrayOfTensor3 diy_dx;
55 FieldOfPropagationMatrix propmat_field;
56 FieldOfStokesVector absorption_field, additional_source_field;
57 ArrayOfPpath ppath_field;
58
59 // Check that the lines and nf is correct
60 Vector f_grid(nf * nl);
61 Index il=0;
62 for (auto& lines: abs_lines_per_species) {
63 for (auto& band: lines) {
64 for (Index k=0; k<band.NumLines(); k++) {
65 nlinspace(f_grid[Range(il * nf, nf)], band.lines[k].F0 * (1 - df), band.lines[k].F0 * (1 + df), nf);
66 il++;
67 }
68 }
69 }
70
71 ARTS_USER_ERROR_IF(not is_increasing(f_grid), "Frequency grid is not increasing, abs_lines_per_species must be sorted and no overlap is allowed");
72
74 ppath_field,
75 ppath_agenda,
76 -1,
77 1e99,
78 1,
79 z_field,
80 f_grid,
81 0,
82 1,
83 0,
84 Vector(1, 0),
85 Vector(1, 0),
86 Vector(0),
87 refellipsoid,
88 1,
89 nz,
90 verbosity);
91
92 ArrayOfArrayOfIndex sorted_index;
93 ArrayOfVector cos_zenith_angles;
94 sorted_index_of_ppath_field(sorted_index, cos_zenith_angles, ppath_field);
95
96 for (Index ip = 0; ip < np; ip++)
98 "Your lineshape integration does normalize. Increase nf and decrease df until it does.",
99 test_integrate_zenith(cos_zenith_angles[ip], sorted_index[ip]));
100
102 propmat_field,
103 absorption_field,
104 additional_source_field,
105 1,
106 f_grid,
107 p_grid,
108 z_field,
109 t_field,
110 nlte_field,
111 vmr_field,
113 propmat_clearsky_agenda);
114
116 nl, Array<Eigen::VectorXcd>(np, Eigen::VectorXcd(nf * nl)));
117 for (Index ip=0; ip<np; ip++) {
118 il=0;
119 for (auto& lines: abs_lines_per_species) {
120 for (auto& band: lines) {
121 const Numeric DC = band.DopplerConstant(t_field(ip, 0, 0));
122 const Vector vmrs = band.BroadeningSpeciesVMR(vmr_field(joker, ip, 0, 0), abs_species);
123 for (Index k=0; k<band.NumLines(); k++) {
124 const auto X = band.ShapeParameters(k, t_field(ip, 0, 0), p_grid[ip], vmrs);
125 LineShape::Calculator ls(band.lineshapetype, band.lines[k].F0, X, DC, 0, false);
126 for (Index iv=0; iv<f_grid.nelem(); iv++) {
127 lineshapes[il][ip][iv] = ls(f_grid[iv]);
128 }
129 il++;
130 }
131 }
132 }
133 }
134 for (auto& aols : lineshapes)
135 for (auto& ls : aols)
137 "Your lineshape integration does not normalize. Increase nf and decrease df until it does.",
138 test_integrate_convolved(ls, f_grid));
139
140 // Counting the path index so we can make the big loop parallel
141 ArrayOfArrayOfIndex counted_path_index(ppath_field.nelem());
142 ArrayOfIndex counter(np, 0);
143 for (Index i = 0; i < ppath_field.nelem(); i++) {
144 const Ppath& path = ppath_field[i];
145 counted_path_index[i].resize(path.np);
146 for (Index ip_path = 0; ip_path < path.np; ip_path++) {
147 const Index ip_grid = grid_index_from_gp(path.gp_p[ip_path]);
148 counted_path_index[i][ip_path] = counter[ip_grid];
149 counter[ip_grid]++;
150 }
151 }
152
153 line_irradiance = Matrix(nl, np, 0.0);
154 line_transmission = Tensor3(1, nl, np, 0.0);
155
156 ArrayOfMatrix line_radiance(np);
157 for (Index i = 0; i < np; i++)
158 line_radiance[i].resize(sorted_index[i].nelem(), nl);
159
161
162#pragma omp parallel for if (not arts_omp_in_parallel()) \
163 schedule(guided) default(shared) firstprivate(wss, il)
164 for (Index i = 0; i < ppath_field.nelem(); i++) {
165 const Ppath& path = ppath_field[i];
166
167 thread_local ArrayOfRadiationVector lvl_rad;
168 thread_local ArrayOfRadiationVector src_rad;
169 thread_local ArrayOfTransmissionMatrix lyr_tra;
170 thread_local ArrayOfTransmissionMatrix tot_tra;
171
173 lvl_rad,
174 src_rad,
175 lyr_tra,
176 tot_tra,
177 propmat_field,
178 absorption_field,
179 additional_source_field,
180 f_grid,
181 t_field,
182 nlte_field,
183 path,
184 iy_main_agenda,
185 iy_space_agenda,
186 iy_surface_agenda,
187 iy_cloudbox_agenda,
188 surface_props_data,
189 verbosity);
190
191 for (Index ip_path = 0; ip_path < path.np; ip_path++) {
192 const Index ip_grid = grid_index_from_gp(path.gp_p[ip_path]);
193 for (il = 0; il < nl; il++)
194 line_radiance[ip_grid](counted_path_index[i][ip_path], il) =
196 lvl_rad[ip_path], lineshapes[il][ip_grid], f_grid);
197 }
198 }
199
200 for (Index ip = 0; ip < np; ip++) {
201 for (il = 0; il < nl; il++) {
202 line_irradiance(il, ip) = integrate_zenith(line_radiance[ip](joker, il),
203 cos_zenith_angles[ip],
204 sorted_index[ip]);
205 }
206 }
207
208 if (r > 0) {
209 const FieldOfTransmissionMatrix transmat_field =
211 for (Index ip = 0; ip < np; ip++)
212 for (il = 0; il < nl; il++)
213 line_transmission(0, il, ip) = integrate_convolved(
214 transmat_field(ip, 0, 0), lineshapes[il][ip], f_grid);
215 }
216}
217
Declarations required for the calculation of absorption coefficients.
The global header file for ARTS.
Header file for helper functions for OpenMP.
The Agenda class.
Index nelem() const ARTS_NOEXCEPT
Definition array.h:75
Creates a 3D field of a base unit.
Definition field.h:16
Line shape calculator.
Definition lineshape.h:661
Workspace class.
Helper macros for debugging.
#define ARTS_USER_ERROR_IF(condition,...)
Definition debug.h:137
Array< RetrievalQuantity > ArrayOfRetrievalQuantity
Definition jacobian.h:529
void ppath_fieldFromDownUpLimbGeoms(Workspace &ws, ArrayOfPpath &ppath_field, const Agenda &ppath_agenda, const Numeric &ppath_lmax, const Numeric &ppath_lraytrace, const Index &atmgeom_checked, const Tensor3 &z_field, const Vector &f_grid, const Index &cloudbox_on, const Index &cloudbox_checked, const Index &ppath_inside_cloudbox_do, const Vector &rte_pos, const Vector &rte_los, const Vector &rte_pos2, const Vector &refellipsoid, const Index &atmosphere_dim, const Index &zenith_angles_per_position, const Verbosity &verbosity)
WORKSPACE METHOD: ppath_fieldFromDownUpLimbGeoms.
Definition m_ppath.cc:1673
void line_irradianceCalcForSingleSpeciesNonOverlappingLinesPseudo2D(Workspace &ws, Matrix &line_irradiance, Tensor3 &line_transmission, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const EnergyLevelMap &nlte_field, const Tensor4 &vmr_field, const Tensor3 &t_field, const Tensor3 &z_field, const Vector &p_grid, const Vector &refellipsoid, const Tensor3 &surface_props_data, const Agenda &ppath_agenda, const Agenda &iy_main_agenda, const Agenda &iy_space_agenda, const Agenda &iy_surface_agenda, const Agenda &iy_cloudbox_agenda, const Agenda &propmat_clearsky_agenda, const Numeric &df, const Index &nz, const Index &nf, const Numeric &r, const Verbosity &verbosity)
WORKSPACE METHOD: line_irradianceCalcForSingleSpeciesNonOverlappingLinesPseudo2D.
void nlinspace(Vector &x, const Numeric start, const Numeric stop, const Index n)
nlinspace
This file contains declerations of functions of physical character.
Propagation path structure and functions.
void field_of_propagation(Workspace &ws, FieldOfPropagationMatrix &propmat_field, FieldOfStokesVector &absorption_field, FieldOfStokesVector &additional_source_field, const Index &stokes_dim, const Vector &f_grid, const Vector &p_grid, const Tensor3 &z_field, const Tensor3 &t_field, const EnergyLevelMap &nlte_field, const Tensor4 &vmr_field, const ArrayOfRetrievalQuantity &jacobian_quantities, const Agenda &propmat_clearsky_agenda)
Creates a field of propagation matrices, absorption vectors, and source vectors.
FieldOfTransmissionMatrix transmat_field_calc_from_propmat_field(const FieldOfPropagationMatrix &propmat_field, const Numeric &r)
Get a field of transmission matrices from the propagation matrix field.
void emission_from_propmat_field(Workspace &ws, ArrayOfRadiationVector &lvl_rad, ArrayOfRadiationVector &src_rad, ArrayOfTransmissionMatrix &lyr_tra, ArrayOfTransmissionMatrix &tot_tra, const FieldOfPropagationMatrix &propmat_field, const FieldOfStokesVector &absorption_field, const FieldOfStokesVector &additional_source_field, const Vector &f_grid, const Tensor3 &t_field, const EnergyLevelMap &nlte_field, const Ppath &ppath, const Agenda &iy_main_agenda, const Agenda &iy_space_agenda, const Agenda &iy_surface_agenda, const Agenda &iy_cloudbox_agenda, const Tensor3 &surface_props_data, const Verbosity &verbosity)
Computes the radiation and transmission from fields of atmospheric propagation.
Implements a propagation matrix field.
void sorted_index_of_ppath_field(ArrayOfArrayOfIndex &sorted_index, ArrayOfVector &cosza, const ArrayOfPpath &ppath_field)
Get sorting of zenith angles in field of ppath.
Numeric test_integrate_convolved(const Eigen::Ref< Eigen::VectorXcd > &F, const Vector &f)
Integrate the line shape.
void error_in_integrate(const String &error_msg, const Numeric &value_that_should_be_unity)
Throws an error if integration values are bad.
Numeric integrate_zenith(const ConstVectorView &j, const Vector &cosza, const Array< Index > &sorted_index)
Convolve source function with 1D sphere and integrate.
Index grid_index_from_gp(const GridPos &gp)
Get a discrete position from grid pos.
Numeric test_integrate_zenith(const Vector &cosza, const Array< Index > &sorted_index)
Integrate cos(za) over the angles.
Numeric integrate_convolved(const RadiationVector &I, const Eigen::VectorXcd &F, const Vector &f)
Convolve intensity and line shape and integrate.
Radiation field calculations.
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.