ARTS 2.5.11 (git: 6827797f)
propmat_field.cc
Go to the documentation of this file.
1
13#include "propmat_field.h"
14#include "matpack_data.h"
15#include "rte.h"
16#include "special_interp.h"
17#include "transmissionmatrix.h"
18
20 FieldOfPropagationMatrix& propmat_field,
21 FieldOfStokesVector& absorption_field,
22 FieldOfStokesVector& additional_source_field,
23 const Index& stokes_dim,
24 const Vector& f_grid,
25 const Vector& p_grid,
26 const Tensor3& z_field,
27 const Tensor3& t_field,
28 const EnergyLevelMap& nlte_field,
29 const Tensor4& vmr_field,
30 const ArrayOfRetrievalQuantity& jacobian_quantities,
31 const Agenda& propmat_clearsky_agenda)
32{
33 const Index nalt = z_field.npages();
34 const Index nlat = z_field.nrows();
35 const Index nlon = z_field.ncols();
36 const Index nq = jacobian_quantities.nelem();
37 const Index nf = f_grid.nelem();
38
40 "Does not support Jacobian calculations at this time");
41 ARTS_USER_ERROR_IF (stokes_dim not_eq 1,
42 "Only for stokes_dim 1 at this time.");
43
44 // Compute variables
45 const Vector mag_field = Vector(3, 0);
46 const Vector los = Vector(2, 0);
47 const Vector tmp(0);
48 ArrayOfStokesVector dS_dx(nq);
49 ArrayOfPropagationMatrix dK_dx(nq);
50
51 propmat_field = FieldOfPropagationMatrix(
52 nalt, nlat, nlon, PropagationMatrix(nf, stokes_dim));
53 absorption_field =
54 FieldOfStokesVector(nalt, nlat, nlon, StokesVector(nf, stokes_dim));
55 additional_source_field =
56 FieldOfStokesVector(nalt, nlat, nlon, StokesVector(nf, stokes_dim));
57
59
60#pragma omp parallel for if (not arts_omp_in_parallel()) collapse(3) \
61 firstprivate(wss)
62 for (Index i = 0; i < nalt; i++) {
63 for (Index j = 0; j < nlat; j++) {
64 for (Index k = 0; k < nlon; k++) {
65 thread_local Index itmp;
67 wss,
68 propmat_field(i, j, k),
69 additional_source_field(i, j, k),
70 itmp,
71 dK_dx,
72 dS_dx,
73 propmat_clearsky_agenda,
74 jacobian_quantities,
75 f_grid,
76 mag_field,
77 los,
78 nlte_field(i, j, k),
79 Vector{vmr_field(joker, i, j, k)},
80 t_field(i, j, k),
81 p_grid[i],
82 false);
83 absorption_field(i, j, k) = propmat_field(i, j, k);
84 }
85 }
86 }
87}
88
90 const FieldOfPropagationMatrix& propmat_field, const Numeric& r)
91{
92 FieldOfTransmissionMatrix transmat_field(
93 propmat_field.npages(), propmat_field.nrows(), propmat_field.ncols());
94 for (size_t ip = 0; ip < propmat_field.npages(); ip++)
95 for (size_t ir = 0; ir < propmat_field.nrows(); ir++)
96 for (size_t ic = 0; ic < propmat_field.ncols(); ic++)
97 transmat_field(ip, ir, ic) =
98 TransmissionMatrix(propmat_field(ip, ir, ic), r);
99 return transmat_field;
100}
101
103 Workspace& ws,
104 ArrayOfRadiationVector& lvl_rad,
105 ArrayOfRadiationVector& src_rad,
108 const FieldOfPropagationMatrix& propmat_field,
109 const FieldOfStokesVector& absorption_field,
110 const FieldOfStokesVector& additional_source_field,
111 const Vector& f_grid,
112 const Tensor3& t_field,
113 const EnergyLevelMap& nlte_field,
114 const Ppath& ppath,
115 const Agenda& iy_main_agenda,
116 const Agenda& iy_space_agenda,
117 const Agenda& iy_surface_agenda,
118 const Agenda& iy_cloudbox_agenda,
119 const Tensor3& surface_props_data,
120 const Verbosity& verbosity)
121{
122 // Size of problem
123 const Index nf = f_grid.nelem();
124 const Index ns = propmat_field(0, 0, 0).StokesDimensions();
125 const Index np = ppath.np;
126
127 // Current limitations
128 ARTS_USER_ERROR_IF (ns not_eq 1, "Only for stokes_dim 1");
129 ARTS_USER_ERROR_IF (ppath.dim not_eq 1, "Only for atmosphere_dim 1");
130
131 // Size of compute variables
132 lvl_rad = ArrayOfRadiationVector(np, RadiationVector(nf, ns));
133 src_rad = ArrayOfRadiationVector(np, RadiationVector(nf, ns));
134 lyr_tra = ArrayOfTransmissionMatrix(np, TransmissionMatrix(nf, ns));
135
136 RadiationVector J_add_dummy;
137 ArrayOfRadiationVector dJ_add_dummy;
138
139 // Size radiative variables always used
140 Vector B(nf);
141 PropagationMatrix K_this(nf, ns), K_past(nf, ns);
142
143 // Temporary empty variables to fit available function handles
144 Vector vtmp(0);
145 ArrayOfTensor3 t3tmp;
147 ArrayOfRadiationVector rvtmp(0);
148 ArrayOfPropagationMatrix pmtmp(0);
149 ArrayOfStokesVector svtmp(0);
151
152 // Loop ppath points and determine radiative properties
153 for (Index ip = 0; ip < np; ip++) {
155 B, vtmp, f_grid, interp_atmfield_by_gp(1, t_field, ppath.gp_p[ip]), false);
156 K_this = propmat_field(ppath.gp_p[ip]);
157 const StokesVector S(additional_source_field(ppath.gp_p[ip]));
158 const StokesVector a(absorption_field(ppath.gp_p[ip]));
159
160 if (ip)
161 stepwise_transmission(lyr_tra[ip],
162 tmtmp,
163 tmtmp,
164 K_past,
165 K_this,
166 pmtmp,
167 pmtmp,
168 ppath.lstep[ip - 1],
169 0,
170 0,
171 -1);
172
173 stepwise_source(src_rad[ip],
174 rvtmp,
175 J_add_dummy,
176 K_this,
177 a,
178 S,
179 pmtmp,
180 svtmp,
181 svtmp,
182 B,
183 vtmp,
184 rqtmp,
185 false);
186
187 swap(K_past, K_this);
188 }
189
190 // In case of backwards RT necessary
192 const Tensor3 iy_trans_new = tot_tra[np - 1];
193
194 // Radiative background
195 Matrix iy;
197 iy,
198 t3tmp,
199 iy_trans_new,
200 0,
201 0,
202 rqtmp,
203 ppath,
204 Vector{0},
205 1,
206 nlte_field,
207 0,
208 1,
209 f_grid,
210 "1",
211 surface_props_data,
212 iy_main_agenda,
213 iy_space_agenda,
214 iy_surface_agenda,
215 iy_cloudbox_agenda,
216 1,
217 verbosity);
218 lvl_rad[np - 1] = iy;
219
220 // Radiative transfer calculations
221 for (Index ip = np - 2; ip >= 0; ip--)
222 update_radiation_vector(lvl_rad[ip] = lvl_rad[ip + 1],
223 rvtmp,
224 rvtmp,
225 src_rad[ip],
226 src_rad[ip + 1],
227 rvtmp,
228 rvtmp,
229 lyr_tra[ip + 1],
230 tot_tra[ip],
231 tmtmp,
232 tmtmp,
233 PropagationMatrix(),
234 PropagationMatrix(),
235 ArrayOfPropagationMatrix(),
236 ArrayOfPropagationMatrix(),
237 Numeric(),
238 Vector(),
239 Vector(),
240 0,
241 0,
243}
The Agenda class.
Definition: agenda_class.h:52
This can be used to make arrays out of anything.
Definition: array.h:31
Index nelem() const ARTS_NOEXCEPT
Definition: array.h:75
Creates a 3D field of a base unit.
Definition: field.h:16
size_t ncols() const
Number of columns.
Definition: field.h:159
size_t npages() const
Number of pages.
Definition: field.h:153
size_t nrows() const
Number of rows.
Definition: field.h:156
Workspace class.
Definition: workspace_ng.h:36
#define ARTS_USER_ERROR_IF(condition,...)
Definition: debug.h:135
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.
Field3D< StokesVector > FieldOfStokesVector
Definition: propmat_field.h:26
Field3D< PropagationMatrix > FieldOfPropagationMatrix
Definition: propmat_field.h:25
void get_iy_of_background(Workspace &ws, Matrix &iy, ArrayOfTensor3 &diy_dx, const Tensor3 &iy_transmittance, const Index &iy_id, const Index &jacobian_do, const ArrayOfRetrievalQuantity &jacobian_quantities, const Ppath &ppath, const Vector &rte_pos2, const Index &atmosphere_dim, const EnergyLevelMap &nlte_field, const Index &cloudbox_on, const Index &stokes_dim, const Vector &f_grid, const String &iy_unit, const Tensor3 &surface_props_data, 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 Verbosity &verbosity)
Determines iy of the "background" of a propgation path.
Definition: rte.cc:713
void get_stepwise_clearsky_propmat(Workspace &ws, PropagationMatrix &K, StokesVector &S, Index &lte, ArrayOfPropagationMatrix &dK_dx, ArrayOfStokesVector &dS_dx, const Agenda &propmat_clearsky_agenda, const ArrayOfRetrievalQuantity &jacobian_quantities, const Vector &ppath_f_grid, const Vector &ppath_magnetic_field, const Vector &ppath_line_of_sight, const EnergyLevelMap &ppath_nlte, const Vector &ppath_vmrs, const Numeric &ppath_temperature, const Numeric &ppath_pressure, const bool &jacobian_do)
Gets the clearsky propgation matrix and NLTE contributions.
Definition: rte.cc:1116
void get_stepwise_blackbody_radiation(VectorView B, VectorView dB_dT, const ConstVectorView &ppath_f_grid, const Numeric &ppath_temperature, const bool &do_temperature_derivative)
Get the blackbody radiation at propagation path point.
Definition: rte.cc:1101
Declaration of functions in rte.cc.
void interp_atmfield_by_gp(VectorView x, const Index &atmosphere_dim, ConstTensor3View x_field, const ArrayOfGridPos &gp_p, const ArrayOfGridPos &gp_lat, const ArrayOfGridPos &gp_lon)
Interpolates an atmospheric field given the grid positions.
Header file for special_interp.cc.
The structure to describe a propagation path and releated quantities.
Definition: ppath_struct.h:17
Index np
Number of points describing the ppath.
Definition: ppath_struct.h:21
Vector lstep
The length between ppath points.
Definition: ppath_struct.h:39
ArrayOfGridPos gp_p
Index position with respect to the pressure grid.
Definition: ppath_struct.h:51
Index dim
Atmospheric dimensionality.
Definition: ppath_struct.h:19
Radiation Vector for Stokes dimension 1-4.
Class to keep track of Transmission Matrices for Stokes Dim 1-4.
ArrayOfTransmissionMatrix cumulative_transmission(const ArrayOfTransmissionMatrix &T, const CumulativeTransmission type)
Accumulate the transmission matrix over all layers.
void stepwise_source(RadiationVector &J, ArrayOfRadiationVector &dJ, RadiationVector &J_add, const PropagationMatrix &K, const StokesVector &a, const StokesVector &S, const ArrayOfPropagationMatrix &dK, const ArrayOfStokesVector &da, const ArrayOfStokesVector &dS, const ConstVectorView &B, const ConstVectorView &dB_dT, const ArrayOfRetrievalQuantity &jacobian_quantities, const bool &jacobian_do)
Set the stepwise source.
void update_radiation_vector(RadiationVector &I, ArrayOfRadiationVector &dI1, ArrayOfRadiationVector &dI2, const RadiationVector &J1, const RadiationVector &J2, const ArrayOfRadiationVector &dJ1, const ArrayOfRadiationVector &dJ2, const TransmissionMatrix &T, const TransmissionMatrix &PiT, const ArrayOfTransmissionMatrix &dT1, const ArrayOfTransmissionMatrix &dT2, const PropagationMatrix &K1, const PropagationMatrix &K2, const ArrayOfPropagationMatrix &dK1, const ArrayOfPropagationMatrix &dK2, const Numeric r, const Vector &dr1, const Vector &dr2, const Index ia, const Index iz, const RadiativeTransferSolver solver)
Update the Radiation Vector.
void stepwise_transmission(TransmissionMatrix &T, ArrayOfTransmissionMatrix &dT1, ArrayOfTransmissionMatrix &dT2, const PropagationMatrix &K1, const PropagationMatrix &K2, const ArrayOfPropagationMatrix &dK1, const ArrayOfPropagationMatrix &dK2, const Numeric &r, const Numeric &dr_dtemp1, const Numeric &dr_dtemp2, const Index temp_deriv_pos)
Set the stepwise transmission matrix.
Stuff related to the transmission matrix.
Array< RadiationVector > ArrayOfRadiationVector
#define a
Array< TransmissionMatrix > ArrayOfTransmissionMatrix