ARTS 2.5.0 (git: 9ee3ac6c)
test_lineshape.cc
Go to the documentation of this file.
1#include "absorption.h"
2#include "artstime.h"
3#include "gui/plot.h"
4#include "lineshape.h"
5#include "species.h"
6
7namespace {
8 constexpr Index N = 10001;
9 constexpr Index M = 25;
10 constexpr Numeric dG0 = 10;
11 constexpr Numeric dD0 = 1;
12 constexpr Numeric dG2 = 1;
13 constexpr Numeric dD2 = 1;
14 constexpr Numeric dETA = 1e-4;
15 constexpr Numeric dFVC = 1e-3;
16 constexpr Numeric dY = 1e-10;
17 constexpr Numeric dG = 1e-14;
18 constexpr Numeric dDV = 1e-1;
19 constexpr Numeric dVMR = 1e-4;
20 constexpr Numeric dT = 0.1;
21 constexpr Numeric dH = 1e-8;
22 constexpr Numeric df = 1e1;
23 constexpr Numeric dF0 = 1e1;
24 constexpr Numeric dI0 = 1e1;
25 static_assert (N > 1);
26} // namespace
27
28ENUMCLASS (ModificationInternal, char,
29 None,
30 T,
31 f,
32 H,
33 F0,
34 I0,
35 SelfG0,
36 SelfD0,
37 SelfG2,
38 SelfD2,
39 SelfETA,
40 SelfFVC,
41 SelfY,
42 SelfG,
43 SelfDV,
44 AirG0,
45 AirD0,
46 AirG2,
47 AirD2,
48 AirETA,
49 AirFVC,
50 AirY,
51 AirG,
52 AirDV,
53 SelfVMR,
54 AirVMR
55 )
56
57template <ModificationInternal mod>
58LineShape::Model lineshape_model() {
59 LineShape::Model model(2);
60 model[0].G0() = LineShape::ModelParameters(LineShape::TemperatureModel::T1, 20e3, 0.8);
61 model[0].D0() = LineShape::ModelParameters(LineShape::TemperatureModel::T5, 200, 0.8);
62 model[0].G2() = LineShape::ModelParameters(LineShape::TemperatureModel::T3, 200, 20.);
63 model[0].D2() = LineShape::ModelParameters(LineShape::TemperatureModel::T0, 20.);
64 model[0].G2() = LineShape::ModelParameters(LineShape::TemperatureModel::T3, 200., 20.);
65 model[0].ETA() = LineShape::ModelParameters(LineShape::TemperatureModel::T0, 0.1);
66 model[0].FVC() = LineShape::ModelParameters(LineShape::TemperatureModel::T2, 1e2, .2, 1.);
67 model[0].Y() = LineShape::ModelParameters(LineShape::TemperatureModel::LM_AER, 1e-6, 1e-7, 1e-8, 1e-9);
68 model[0].G() = LineShape::ModelParameters(LineShape::TemperatureModel::LM_AER, 1e-12, 1e-11, 1e-10, 1e-9);
69 model[0].DV() = LineShape::ModelParameters(LineShape::TemperatureModel::DPL, 1e2, 0.8, -1e1, 0.3);
70
71 model[1].G0() = LineShape::ModelParameters(LineShape::TemperatureModel::T1, 1.2*20e3, 0.8);
72 model[1].D0() = LineShape::ModelParameters(LineShape::TemperatureModel::T5, 1.2*200, 0.8);
73 model[1].G2() = LineShape::ModelParameters(LineShape::TemperatureModel::T3, 1.2*200, 1.2*20);
74 model[1].D2() = LineShape::ModelParameters(LineShape::TemperatureModel::T0, 1.2*20);
75 model[1].G2() = LineShape::ModelParameters(LineShape::TemperatureModel::T3, 1.2*200, 1.2*20);
76 model[1].ETA() = LineShape::ModelParameters(LineShape::TemperatureModel::T0, 1.2*0.1);
77 model[1].FVC() = LineShape::ModelParameters(LineShape::TemperatureModel::T2, 1.2*1e2, 1.2*.2, 1.2*1);
78 model[1].Y() = LineShape::ModelParameters(LineShape::TemperatureModel::LM_AER, 1.2*1e-6, 1.2*1e-7, 1.2*1e-8, 1.2*1e-9);
79 model[1].G() = LineShape::ModelParameters(LineShape::TemperatureModel::LM_AER, 1.2*1e-12, 1.2*1e-11, 1.2*1e-10, 1.2*1e-9);
80 model[1].DV() = LineShape::ModelParameters(LineShape::TemperatureModel::DPL, 1.2*1e2, 0.8, 1.2*-1e1, 0.3);
81
82 // Modify
83 if constexpr (mod == ModificationInternal::SelfG0) model[0].G0().X0 += ::dG0;
84 if constexpr (mod == ModificationInternal::SelfD0) model[0].D0().X0 += ::dD0;
85 if constexpr (mod == ModificationInternal::SelfG2) model[0].G2().X0 += ::dG2;
86 if constexpr (mod == ModificationInternal::SelfD2) model[0].D2().X0 += ::dD2;
87 if constexpr (mod == ModificationInternal::SelfETA) model[0].ETA().X0 += ::dETA;
88 if constexpr (mod == ModificationInternal::SelfFVC) model[0].FVC().X0 += ::dFVC;
89 if constexpr (mod == ModificationInternal::SelfY) model[0].Y().X1 += ::dY;
90 if constexpr (mod == ModificationInternal::SelfG) model[0].G().X1 += ::dG;
91 if constexpr (mod == ModificationInternal::SelfDV) model[0].DV().X0 += ::dDV;
92 if constexpr (mod == ModificationInternal::AirG0) model[1].G0().X0 += ::dG0;
93 if constexpr (mod == ModificationInternal::AirD0) model[1].D0().X0 += ::dD0;
94 if constexpr (mod == ModificationInternal::AirG2) model[1].G2().X0 += ::dG2;
95 if constexpr (mod == ModificationInternal::AirD2) model[1].D2().X0 += ::dD2;
96 if constexpr (mod == ModificationInternal::AirETA) model[1].ETA().X0 += ::dETA;
97 if constexpr (mod == ModificationInternal::AirFVC) model[1].FVC().X0 += ::dFVC;
98 if constexpr (mod == ModificationInternal::AirY) model[1].Y().X1 += ::dY;
99 if constexpr (mod == ModificationInternal::AirG) model[1].G().X1 += ::dG;
100 if constexpr (mod == ModificationInternal::AirDV) model[1].DV().X0 += ::dDV;
101
102 return model;
103}
104
113};
114
115template <ModificationInternal mod, LineShape::Type ls_type>
117 QuantumNumbers upp, low;
118 upp[QuantumNumberType::J] = 3;
119 low[QuantumNumberType::J] = 2;
120 QuantumIdentifier qid = QuantumIdentifier(Species::Tag("H2O-161").Isotopologue(), upp, low);
121
122 Absorption::SingleLine line(1e9, 1e10, 1e-20, 1., 3., 1e-14,
123 Zeeman::Model(2, 1.5), lineshape_model<mod>(),
124 {low[QuantumNumberType::J]}, {upp[QuantumNumberType::J]});
125
126 if constexpr (mod == ModificationInternal::F0) line.F0() += ::dF0;
127 if constexpr (mod == ModificationInternal::I0) line.I0() += ::dI0;
128
129 AbsorptionLines band(true, true,
131 Absorption::PopulationType::LTE, Absorption::NormalizationType::None,
132 ls_type, 296, -1, -1, qid, {QuantumNumberType::J},
133 {Species::Species::Water, Species::Species::Bath}, {line});
134
135 const Numeric f0 = ((mod == ModificationInternal::f) ? 500e6 + ::df : 500e6) + (ls_type == LineShape::Type::DP ? 498e6 : 0);
136 const Numeric df_ = ls_type == LineShape::Type::DP ? 4e6/Numeric(::N - 1) : 1e9/Numeric(::N - 1);
137
138 const Vector f_grid(f0, ::N, df_);
139 Vector vmr({0.3, 0.7});
140 if constexpr (mod == ModificationInternal::SelfVMR) vmr[0] += ::dVMR;
141 else if constexpr (mod == ModificationInternal::AirVMR) vmr[1] += ::dVMR;
142
143 const Numeric P = 1e3;
144 const Numeric T = (mod == ModificationInternal::T) ? 255 + ::dT : 255;
145 const Numeric H = (mod == ModificationInternal::H) ? 50e-6 + ::dH : 50e-6;
146
147 return {qid, band, f_grid, vmr, P, T, H};
148}
149
150template <LineShape::Type ls_type>
151void test_ls() {
152 std::cout << "Line Shape: " << ls_type << "\n";
153 Index j=0;
154
155 Vector f_grid(::N);
156 ComplexVector F(::N), N(::N);
157 ComplexMatrix dF(::N, ::M, 0), dN(::N, ::M, 0);
158 ArrayOfComplexVector dF_mod(::M, ComplexVector(::N, 0)), dN_mod(::M, ComplexVector(::N, 0));
159 ComplexMatrix dFnull(0, 0), dNnull(0, 0);
160 ArrayOfRetrievalQuantity jacobian_quantities(::M);
161 ArrayOfRetrievalQuantity jacobian_quantities_null(0);
162 EnergyLevelMap nlte;
163
164 // Ignore sparsity here
165 const Vector f_grid_sparse(0);
166 LineShape::ComputeData sparse_com(f_grid_sparse, {}, false);
167
168 constexpr Index TN = 1;
169
170 // Compute
171 {
172 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::None, ls_type>();
173 f_grid = f_g; f_grid -= 1e9; f_grid /= 1e6;
174 jacobian_quantities[0].Target() = JacobianTarget(Jacobian::Atm::Temperature);
175 jacobian_quantities[1].Target() = JacobianTarget(Jacobian::Atm::WindMagnitude);
176 jacobian_quantities[2].Target() = JacobianTarget(Jacobian::Atm::MagneticMagnitude);
177 jacobian_quantities[3].Target() = JacobianTarget(Jacobian::Line::Center, qid, qid.Species());
178 jacobian_quantities[4].Target() = JacobianTarget(Jacobian::Line::Strength, qid, qid.Species());
179 jacobian_quantities[5].Target() = JacobianTarget(Jacobian::Line::ShapeG0X0, qid, qid.Species());
180 jacobian_quantities[6].Target() = JacobianTarget(Jacobian::Line::ShapeD0X0, qid, qid.Species());
181 jacobian_quantities[7].Target() = JacobianTarget(Jacobian::Line::ShapeG2X0, qid, qid.Species());
182 jacobian_quantities[8].Target() = JacobianTarget(Jacobian::Line::ShapeD2X0, qid, qid.Species());
183 jacobian_quantities[9].Target() = JacobianTarget(Jacobian::Line::ShapeETAX0, qid, qid.Species());
184 jacobian_quantities[10].Target() = JacobianTarget(Jacobian::Line::ShapeFVCX0, qid, qid.Species());
185 jacobian_quantities[11].Target() = JacobianTarget(Jacobian::Line::ShapeYX1, qid, qid.Species());
186 jacobian_quantities[12].Target() = JacobianTarget(Jacobian::Line::ShapeGX1, qid, qid.Species());
187 jacobian_quantities[13].Target() = JacobianTarget(Jacobian::Line::ShapeDVX0, qid, qid.Species());
188 jacobian_quantities[14].Target() = JacobianTarget(Jacobian::Line::ShapeG0X0, qid, Species::Species::Bath);
189 jacobian_quantities[15].Target() = JacobianTarget(Jacobian::Line::ShapeD0X0, qid, Species::Species::Bath);
190 jacobian_quantities[16].Target() = JacobianTarget(Jacobian::Line::ShapeG2X0, qid, Species::Species::Bath);
191 jacobian_quantities[17].Target() = JacobianTarget(Jacobian::Line::ShapeD2X0, qid, Species::Species::Bath);
192 jacobian_quantities[18].Target() = JacobianTarget(Jacobian::Line::ShapeETAX0, qid, Species::Species::Bath);
193 jacobian_quantities[19].Target() = JacobianTarget(Jacobian::Line::ShapeFVCX0, qid, Species::Species::Bath);
194 jacobian_quantities[20].Target() = JacobianTarget(Jacobian::Line::ShapeYX1, qid, Species::Species::Bath);
195 jacobian_quantities[21].Target() = JacobianTarget(Jacobian::Line::ShapeGX1, qid, Species::Species::Bath);
196 jacobian_quantities[22].Target() = JacobianTarget(Jacobian::Line::ShapeDVX0, qid, Species::Species::Bath);
197 qid = QuantumIdentifier(qid.Isotopologue());
198 jacobian_quantities[23].Target() = JacobianTarget(Jacobian::Line::VMR, qid, qid.Species()); // Self is increasing
199 jacobian_quantities[24].Target() = JacobianTarget(Jacobian::Line::VMR, QuantumIdentifier(Species::Isotopologues[0]), Species::Species::Bath); // Another species... main point: it's air
200 LineShape::ComputeData com(f_g, jacobian_quantities, false);
201 ArrayOfTimeStep dt(TN);
202 std::cout << "all derivs...\n";
203 for (Index i=0; i<TN; i++) {
204 Time start;
205 LineShape::compute(com, sparse_com, band, jacobian_quantities, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
206 dt[i] = Time() - start;
207 }
208 std::sort(dt.begin(), dt.end());
209 std::cerr << "Slow: " << dt[TN-1] << '\n' << "Mean: " << mean(dt) << '\n' << "Med.: " << median(dt) << '\n' << "Fast: " << dt[0] << '\n';
210
211 LineShape::ComputeData com2(f_g, jacobian_quantities, false);
212 LineShape::compute(com2, sparse_com, band, jacobian_quantities, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
213 F = com2.F;
214 N = com2.N;
215 dF = com2.dF;
216 dN = com2.dN;
217 }
218 {
219 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::T, ls_type>();
220
221 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
222 ArrayOfTimeStep dt(TN);
223 std::cout << "\njust forward calcs...\n";
224 for (Index i=0; i<TN; i++) {
225 Time start;
226 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
227 dt[i] = Time() - start;
228 }
229 std::sort(dt.begin(), dt.end());
230 std::cerr << "Slow: " << dt[TN-1] << '\n' << "Mean: " << mean(dt) << '\n' << "Med.: " << median(dt) << '\n' << "Fast: " << dt[0] << '\n';
231
232 LineShape::ComputeData com2(f_g, jacobian_quantities_null, false);
233 LineShape::compute(com2, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
234 dF_mod[j] = com2.F;
235 dN_mod[j] = com2.N;
236
237 dF_mod[j] -= F;
238 dF_mod[j] /= ::dT;
239 j++;
240 }
241 {
242 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::f, ls_type>();
243 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
244 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
245 dF_mod[j] = com.F;
246 dN_mod[j] = com.N;
247 dF_mod[j] -= F;
248 dF_mod[j] /= ::df;
249 j++;
250 }
251 {
252 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::H, ls_type>();
253 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
254 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
255 dF_mod[j] = com.F;
256 dN_mod[j] = com.N;
257 dF_mod[j] -= F;
258 dF_mod[j] /= ::dH;
259 j++;
260 }
261 {
262 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::F0, ls_type>();
263 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
264 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
265 dF_mod[j] = com.F;
266 dN_mod[j] = com.N;
267 dF_mod[j] -= F;
268 dF_mod[j] /= ::dF0;
269 j++;
270 }
271 {
272 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::I0, ls_type>();
273 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
274 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
275 dF_mod[j] = com.F;
276 dN_mod[j] = com.N;
277 dF_mod[j] -= F;
278 dF_mod[j] /= ::dI0;
279 j++;
280 }
281 {
282 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::SelfG0, ls_type>();
283 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
284 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
285 dF_mod[j] = com.F;
286 dN_mod[j] = com.N;
287 dF_mod[j] -= F;
288 dF_mod[j] /= ::dG0;
289 j++;
290 }
291 {
292 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::SelfD0, ls_type>();
293 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
294 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
295 dF_mod[j] = com.F;
296 dN_mod[j] = com.N;
297 dF_mod[j] -= F;
298 dF_mod[j] /= ::dD0;
299 j++;
300 }
301 {
302 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::SelfG2, ls_type>();
303 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
304 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
305 dF_mod[j] = com.F;
306 dN_mod[j] = com.N;
307 dF_mod[j] -= F;
308 dF_mod[j] /= ::dG2;
309 j++;
310 }
311 {
312 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::SelfD2, ls_type>();
313 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
314 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
315 dF_mod[j] = com.F;
316 dN_mod[j] = com.N;
317 dF_mod[j] -= F;
318 dF_mod[j] /= ::dD2;
319 j++;
320 }
321 {
322 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::SelfETA, ls_type>();
323 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
324 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
325 dF_mod[j] = com.F;
326 dN_mod[j] = com.N;
327 dF_mod[j] -= F;
328 dF_mod[j] /= ::dETA;
329 j++;
330 }
331 {
332 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::SelfFVC, ls_type>();
333 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
334 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
335 dF_mod[j] = com.F;
336 dN_mod[j] = com.N;
337 dF_mod[j] -= F;
338 dF_mod[j] /= ::dFVC;
339 j++;
340 }
341 {
342 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::SelfY, ls_type>();
343 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
344 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
345 dF_mod[j] = com.F;
346 dN_mod[j] = com.N;
347 dF_mod[j] -= F;
348 dF_mod[j] /= ::dY;
349 j++;
350 }
351 {
352 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::SelfG, ls_type>();
353 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
354 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
355 dF_mod[j] = com.F;
356 dN_mod[j] = com.N;
357 dF_mod[j] -= F;
358 dF_mod[j] /= ::dG;
359 j++;
360 }
361 {
362 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::SelfDV, ls_type>();
363 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
364 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
365 dF_mod[j] = com.F;
366 dN_mod[j] = com.N;
367 dF_mod[j] -= F;
368 dF_mod[j] /= ::dDV;
369 j++;
370 }
371 {
372 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::AirG0, ls_type>();
373 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
374 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
375 dF_mod[j] = com.F;
376 dN_mod[j] = com.N;
377 dF_mod[j] -= F;
378 dF_mod[j] /= ::dG0;
379 j++;
380 }
381 {
382 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::AirD0, ls_type>();
383 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
384 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
385 dF_mod[j] = com.F;
386 dN_mod[j] = com.N;
387 dF_mod[j] -= F;
388 dF_mod[j] /= ::dD0;
389 j++;
390 }
391 {
392 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::AirG2, ls_type>();
393 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
394 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
395 dF_mod[j] = com.F;
396 dN_mod[j] = com.N;
397 dF_mod[j] -= F;
398 dF_mod[j] /= ::dG2;
399 j++;
400 }
401 {
402 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::AirD2, ls_type>();
403 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
404 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
405 dF_mod[j] = com.F;
406 dN_mod[j] = com.N;
407 dF_mod[j] -= F;
408 dF_mod[j] /= ::dD2;
409 j++;
410 }
411 {
412 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::AirETA, ls_type>();
413 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
414 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
415 dF_mod[j] = com.F;
416 dN_mod[j] = com.N;
417 dF_mod[j] -= F;
418 dF_mod[j] /= ::dETA;
419 j++;
420 }
421 {
422 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::AirFVC, ls_type>();
423 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
424 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
425 dF_mod[j] = com.F;
426 dN_mod[j] = com.N;
427 dF_mod[j] -= F;
428 dF_mod[j] /= ::dFVC;
429 j++;
430 }
431 {
432 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::AirY, ls_type>();
433 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
434 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
435 dF_mod[j] = com.F;
436 dN_mod[j] = com.N;
437 dF_mod[j] -= F;
438 dF_mod[j] /= ::dY;
439 j++;
440 }
441 {
442 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::AirG, ls_type>();
443 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
444 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
445 dF_mod[j] = com.F;
446 dN_mod[j] = com.N;
447 dF_mod[j] -= F;
448 dF_mod[j] /= ::dG;
449 j++;
450 }
451 {
452 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::AirDV, ls_type>();
453 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
454 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
455 dF_mod[j] = com.F;
456 dN_mod[j] = com.N;
457 dF_mod[j] -= F;
458 dF_mod[j] /= ::dDV;
459 j++;
460 }
461 {
462 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::SelfVMR, ls_type>();
463 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
464 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
465 dF_mod[j] = com.F;
466 dN_mod[j] = com.N;
467 dF_mod[j] -= F;
468 dF_mod[j] /= ::dVMR;
469 j++;
470 }
471 {
472 auto [qid, band, f_g, vmr, P, T, H] = line_model<ModificationInternal::AirVMR, ls_type>();
473 LineShape::ComputeData com(f_g, jacobian_quantities_null, false);
474 LineShape::compute(com, sparse_com, band, jacobian_quantities_null, nlte, vmr, {}, 1, 1, P, T, H, 0, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
475 dF_mod[j] = com.F;
476 dN_mod[j] = com.N;
477 dF_mod[j] -= F;
478 dF_mod[j] /= ::dVMR;
479 j++;
480 }
481
482 ARTSGUI::plot(f_grid, F.real(), f_grid, F.imag());
483 for (Index i=0; i<::M; i++) {
484 bool all_zero = true;
485 for (Index iv=0; iv<::N; iv++) all_zero = all_zero and dF(iv, i) == Complex(0, 0) and dF_mod[i][iv] == Complex(0, 0);
486 if (all_zero) continue; // plot only if some are non-zero
487
488 std::cout << jacobian_quantities[i].Target() << '\n';
489 ARTSGUI::plot(f_grid, dF(joker, i).real(), f_grid, dF(joker, i).imag(), f_grid, dF_mod[i].real(), f_grid, dF_mod[i].imag());
490 }
491 std::cout << '\n';
492}
493
494void test_norm() {
495
496 [[maybe_unused]]
497 auto [qid, band, f_grid2, vmr, P, T, H] = line_model<ModificationInternal::None, LineShape::Type::VP>();
498 Vector f_grid(::N);
499 Numeric f0 = f_grid2[0];
500 for (Index i=0; i<::N; i++) {
501 f_grid[i] = f0;
502 f0 += 3 * (f_grid2[1] - f_grid2[0]);
503 }
504 Vector N_nonorm(::N);
505 Vector N_vvh(::N);
506 Vector N_vvw(::N);
507 Vector N_rq(::N);
508 Vector N_sfs(::N);
509
510 ArrayOfVector dN_nonorm(3, Vector(::N));
511 ArrayOfVector dN_vvh(3, Vector(::N));
512 ArrayOfVector dN_vvw(3, Vector(::N));
513 ArrayOfVector dN_rq(3, Vector(::N));
514 ArrayOfVector dN_sfs(3, Vector(::N));
515
516 ArrayOfVector dN_nonorm_mod(3, Vector(::N));
517 ArrayOfVector dN_vvh_mod(3, Vector(::N));
518 ArrayOfVector dN_vvw_mod(3, Vector(::N));
519 ArrayOfVector dN_rq_mod(3, Vector(::N));
520 ArrayOfVector dN_sfs_mod(3, Vector(::N));
521
522 LineShape::Normalizer lsn_nonorm(Absorption::NormalizationType::None, band.Line(0).F0(), T);
523 LineShape::Normalizer lsn_vvh(Absorption::NormalizationType::VVH, band.Line(0).F0(), T);
524 LineShape::Normalizer lsn_vvw(Absorption::NormalizationType::VVW, band.Line(0).F0(), T);
525 LineShape::Normalizer lsn_rq(Absorption::NormalizationType::RQ, band.Line(0).F0(), T);
526 LineShape::Normalizer lsn_sfs(Absorption::NormalizationType::SFS, band.Line(0).F0(), T);
527 for (Index i=0; i<::N; i++) {
528 const Numeric& f = f_grid[i];
529
530 N_nonorm[i] = lsn_nonorm(f);
531 N_vvh[i] = lsn_vvh(f);
532 N_vvw[i] = lsn_vvw(f);
533 N_rq[i] = lsn_rq(f);
534 N_sfs[i] = lsn_sfs(f);
535
536 dN_nonorm[0][i] = lsn_nonorm.dNdT(T, f);
537 dN_vvh[0][i] = lsn_vvh.dNdT(T, f);
538 dN_vvw[0][i] = lsn_vvw.dNdT(T, f);
539 dN_rq[0][i] = lsn_rq.dNdT(T, f);
540 dN_sfs[0][i] = lsn_sfs.dNdT(T, f);
541
542 dN_nonorm[1][i] = lsn_nonorm.dNdf(f);
543 dN_vvh[1][i] = lsn_vvh.dNdf(f);
544 dN_vvw[1][i] = lsn_vvw.dNdf(f);
545 dN_rq[1][i] = lsn_rq.dNdf(f);
546 dN_sfs[1][i] = lsn_sfs.dNdf(f);
547
548 dN_nonorm[2][i] = lsn_nonorm.dNdF0();
549 dN_vvh[2][i] = lsn_vvh.dNdF0();
550 dN_vvw[2][i] = lsn_vvw.dNdF0();
551 dN_rq[2][i] = lsn_rq.dNdF0();
552 dN_sfs[2][i] = lsn_sfs.dNdF0();
553
554 // Do derivatives last...
555 dN_nonorm_mod[1][i] = lsn_nonorm(f+::df);
556 dN_vvh_mod[1][i] = lsn_vvh(f+::df);
557 dN_vvw_mod[1][i] = lsn_vvw(f+::df);
558 dN_rq_mod[1][i] = lsn_rq(f+::df);
559 dN_sfs_mod[1][i] = lsn_sfs(f+::df);
560 }
561
562 LineShape::Normalizer lsn_nonormdT(Absorption::NormalizationType::None, band.Line(0).F0(), T+::dT);
563 LineShape::Normalizer lsn_vvhdT(Absorption::NormalizationType::VVH, band.Line(0).F0(), T+::dT);
564 LineShape::Normalizer lsn_vvwdT(Absorption::NormalizationType::VVW, band.Line(0).F0(), T+::dT);
565 LineShape::Normalizer lsn_rqdT(Absorption::NormalizationType::RQ, band.Line(0).F0(), T+::dT);
566 LineShape::Normalizer lsn_sfsdT(Absorption::NormalizationType::SFS, band.Line(0).F0(), T+::dT);
567 for (Index i=0; i<::N; i++) {
568 const Numeric& f = f_grid[i];
569
570 dN_nonorm_mod[0][i] = lsn_nonormdT(f);
571 dN_vvh_mod[0][i] = lsn_vvhdT(f);
572 dN_vvw_mod[0][i] = lsn_vvwdT(f);
573 dN_rq_mod[0][i] = lsn_rqdT(f);
574 dN_sfs_mod[0][i] = lsn_sfsdT(f);
575 }
576
577 LineShape::Normalizer lsn_nonormdF0(Absorption::NormalizationType::None, band.Line(0).F0()+::dF0, T);
578 LineShape::Normalizer lsn_vvhdF0(Absorption::NormalizationType::VVH, band.Line(0).F0()+::dF0, T);
579 LineShape::Normalizer lsn_vvwdF0(Absorption::NormalizationType::VVW, band.Line(0).F0()+::dF0, T);
580 LineShape::Normalizer lsn_rqdF0(Absorption::NormalizationType::RQ, band.Line(0).F0()+::dF0, T);
581 LineShape::Normalizer lsn_sfsdF0(Absorption::NormalizationType::SFS, band.Line(0).F0()+::dF0, T);
582 for (Index i=0; i<::N; i++) {
583 const Numeric& f = f_grid[i];
584
585 dN_nonorm_mod[2][i] = lsn_nonormdF0(f);
586 dN_vvh_mod[2][i] = lsn_vvhdF0(f);
587 dN_vvw_mod[2][i] = lsn_vvwdF0(f);
588 dN_rq_mod[2][i] = lsn_rqdF0(f);
589 dN_sfs_mod[2][i] = lsn_sfsdF0(f);
590 }
591
592 dN_nonorm_mod[0] -= N_nonorm; dN_nonorm_mod[0] /= ::dT;
593 dN_nonorm_mod[1] -= N_nonorm; dN_nonorm_mod[1] /= ::df;
594 dN_nonorm_mod[2] -= N_nonorm; dN_nonorm_mod[2] /= ::dF0;
595 dN_vvh_mod[0] -= N_vvh; dN_vvh_mod[0] /= ::dT;
596 dN_vvh_mod[1] -= N_vvh; dN_vvh_mod[1] /= ::df;
597 dN_vvh_mod[2] -= N_vvh; dN_vvh_mod[2] /= ::dF0;
598 dN_vvw_mod[0] -= N_vvw; dN_vvw_mod[0] /= ::dT;
599 dN_vvw_mod[1] -= N_vvw; dN_vvw_mod[1] /= ::df;
600 dN_vvw_mod[2] -= N_vvw; dN_vvw_mod[2] /= ::dF0;
601 dN_rq_mod[0] -= N_rq; dN_rq_mod[0] /= ::dT;
602 dN_rq_mod[1] -= N_rq; dN_rq_mod[1] /= ::df;
603 dN_rq_mod[2] -= N_rq; dN_rq_mod[2] /= ::dF0;
604 dN_sfs_mod[0] -= N_sfs; dN_sfs_mod[0] /= ::dT;
605 dN_sfs_mod[1] -= N_sfs; dN_sfs_mod[1] /= ::df;
606 dN_sfs_mod[2] -= N_sfs; dN_sfs_mod[2] /= ::dF0;
607
609 f_grid /= band.F0(0);
610
611 ARTSGUI::PlotConfig::X = "Frequency [line center]";
612 ARTSGUI::PlotConfig::Y = "Norm";
613 ARTSGUI::PlotConfig::Legend = {"None", "VVH", "VVW", "RQ", "SFS"};
614 ARTSGUI::PlotConfig::Title = "Normalizers";
615 ARTSGUI::plot(f_grid, N_nonorm, f_grid, N_vvh, f_grid, N_vvw, f_grid, N_rq, f_grid, N_sfs);
616
617 ARTSGUI::PlotConfig::Legend = {"Analytical", "Perturbed"};
618 for (Index i=0; i<3; i++) {
619 if (i == 0) {
620 ARTSGUI::PlotConfig::Y = "Temperature Derivative";
621 } else if (i == 1) {
622 ARTSGUI::PlotConfig::Y = "Frequency Derivative";
623 } else if (i == 2) {
624 ARTSGUI::PlotConfig::Y = "Line Center Derivative";
625 } else {
626 ARTSGUI::PlotConfig::Y = "Unknown Derivative";
627 }
628
629 if (max(dN_nonorm[i]) not_eq 0 or min(dN_nonorm[i]) not_eq 0 or max(dN_nonorm_mod[i]) not_eq 0 or min(dN_nonorm_mod[i]) not_eq 0) {
630 ARTSGUI::PlotConfig::Title = "No Normalizer";
631 ARTSGUI::plot(f_grid, dN_nonorm[i], f_grid, dN_nonorm_mod[i]);
632 }
633
634 if (max(dN_vvh[i]) not_eq 0 or min(dN_vvh[i]) not_eq 0 or max(dN_vvh_mod[i]) not_eq 0 or min(dN_vvh_mod[i]) not_eq 0) {
635 ARTSGUI::PlotConfig::Title = "VVH Normalizer";
636 ARTSGUI::plot(f_grid, dN_vvh[i], f_grid, dN_vvh_mod[i]);
637 }
638
639 if (max(dN_vvw[i]) not_eq 0 or min(dN_vvw[i]) not_eq 0 or max(dN_vvw_mod[i]) not_eq 0 or min(dN_vvw_mod[i]) not_eq 0) {
640 ARTSGUI::PlotConfig::Title = "VVW Normalizer";
641 ARTSGUI::plot(f_grid, dN_vvw[i], f_grid, dN_vvw_mod[i]);
642 }
643
644 if (max(dN_rq[i]) not_eq 0 or min(dN_rq[i]) not_eq 0 or max(dN_rq_mod[i]) not_eq 0 or min(dN_rq_mod[i]) not_eq 0) {
645 ARTSGUI::PlotConfig::Title = "RQ Normalizer";
646 ARTSGUI::plot(f_grid, dN_rq[i], f_grid, dN_rq_mod[i]);
647 }
648
649 if (max(dN_sfs[i]) not_eq 0 or min(dN_sfs[i]) not_eq 0 or max(dN_sfs_mod[i]) not_eq 0 or min(dN_sfs_mod[i]) not_eq 0) {
650 ARTSGUI::PlotConfig::Title = "SFS Normalizer";
651 ARTSGUI::plot(f_grid, dN_sfs[i], f_grid, dN_sfs_mod[i]);
652 }
653 }
654}
655
657 Numeric I0=1e10, T0=296, T=255, F0=1e9, E0=1e-20, r=1.0;
659 single_partition_function(T, Species::Tag("H2O-161").Isotopologue()),
660 single_partition_function(T0, Species::Tag("H2O-161").Isotopologue()),
661 dsingle_partition_function_dT(T, Species::Tag("H2O-161").Isotopologue()), r, 0, 0);
662 LineShape::LocalThermodynamicEquilibrium dxdI0(I0+1e4, T0, T, F0, E0,
663 single_partition_function(T, Species::Tag("H2O-161").Isotopologue()),
664 single_partition_function(T0, Species::Tag("H2O-161").Isotopologue()),
665 dsingle_partition_function_dT(T, Species::Tag("H2O-161").Isotopologue()), r, 0, 0);
666 LineShape::LocalThermodynamicEquilibrium dxdT(I0, T0, T+0.01, F0, E0,
667 single_partition_function(T+0.01, Species::Tag("H2O-161").Isotopologue()),
668 single_partition_function(T0, Species::Tag("H2O-161").Isotopologue()),
669 dsingle_partition_function_dT(T+0.01, Species::Tag("H2O-161").Isotopologue()), r, 0, 0);
670 LineShape::LocalThermodynamicEquilibrium dxdF0(I0, T0, T, F0+100e3, E0,
671 single_partition_function(T, Species::Tag("H2O-161").Isotopologue()),
672 single_partition_function(T0, Species::Tag("H2O-161").Isotopologue()),
673 dsingle_partition_function_dT(T, Species::Tag("H2O-161").Isotopologue()), r, 0, 0);
674
675 std::cout << x.S << '\n';
676 std::cout << x.dSdI0() << ' ' << 'v' << 's' << ' ' << (dxdI0.S - x.S) / 1e4 << '\n';
677 std::cout << x.dSdT() << ' ' << 'v' << 's' << ' ' << (dxdT.S - x.S) / 0.01 << '\n';
678 std::cout << x.dSdF0() << ' ' << 'v' << 's' << ' ' << (dxdF0.S - x.S) / 100e3 << '\n';
679}
680
682 auto [qid, band, f_gp, vmr, P, T, H] = line_model<ModificationInternal::None, LineShape::Type::HTP>();
683 P *= 1000;
684 band.F0(0) = 1500e9;
685 Vector f_g(500e9, 20001, 10e7);
686 auto band2 = band;
687 band2.Cutoff(Absorption::CutoffType::ByLine);
688 band2.CutoffFreqValue(750e9);
689 auto band3 = band2;
690 band3.CutoffFreqValue(333.33e9);
691
692 EnergyLevelMap nlte;
693
694 const Numeric dfreq = 45e9;
695
696 LineShape::ComputeData com_full(f_g, {}, false);
697 LineShape::ComputeData sparse_com_full(Vector(0), {}, false);
698 LineShape::compute(com_full, sparse_com_full, band, {}, nlte,
699 vmr, {}, 1, 1, P, T, H, dfreq, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
700 LineShape::compute(com_full, sparse_com_full, band2, {}, nlte,
701 vmr, {}, 1, 1, P, T, H, dfreq, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
702 LineShape::compute(com_full, sparse_com_full, band3, {}, nlte,
703 vmr, {}, 1, 1, P, T, H, dfreq, true, Zeeman::Polarization::Pi, Options::LblSpeedup::None);
704
705 Vector sparse3_f_grid = LineShape::triple_sparse_f_grid(f_g, 30e9);
706 LineShape::ComputeData com3(f_g, {}, false);
707 LineShape::ComputeData sparse_com3(sparse3_f_grid, {}, false);
708 LineShape::compute(com3, sparse_com3, band, {}, nlte,
709 vmr, {}, 1, 1, P, T, H, dfreq, true, Zeeman::Polarization::Pi, Options::LblSpeedup::QuadraticIndependent);
710 LineShape::compute(com3, sparse_com3, band2, {}, nlte,
711 vmr, {}, 1, 1, P, T, H, dfreq, true, Zeeman::Polarization::Pi, Options::LblSpeedup::QuadraticIndependent);
712 LineShape::compute(com3, sparse_com3, band3, {}, nlte,
713 vmr, {}, 1, 1, P, T, H, dfreq, true, Zeeman::Polarization::Pi, Options::LblSpeedup::QuadraticIndependent);
714
715 Vector sparsel_f_grid = LineShape::linear_sparse_f_grid(f_g, 10e9);
716 std::cout << sparsel_f_grid <<'\n';
717 LineShape::ComputeData coml(f_g, {}, false);
718 LineShape::ComputeData sparse_coml(sparsel_f_grid, {}, false);
719 LineShape::compute(coml, sparse_coml, band, {}, nlte,
720 vmr, {}, 1, 1, P, T, H, dfreq, true, Zeeman::Polarization::Pi, Options::LblSpeedup::LinearIndependent);
721 LineShape::compute(coml, sparse_coml, band2, {}, nlte,
722 vmr, {}, 1, 1, P, T, H, dfreq, true, Zeeman::Polarization::Pi, Options::LblSpeedup::LinearIndependent);
723 LineShape::compute(coml, sparse_coml, band3, {}, nlte,
724 vmr, {}, 1, 1, P, T, H, dfreq, true, Zeeman::Polarization::Pi, Options::LblSpeedup::LinearIndependent);
725 Vector f = f_g;
726 f /= 1e9;
727
728 std::cout << f[0] << ' ' << f[f_g.nelem() - 1] << '\n';
729 std::cout << com_full.f_grid.nelem()<< ' ' << sparse_com3.f_grid.nelem()<< ' ' << sparse_coml.f_grid.nelem() << '\n';
730
731 ARTSGUI::plot(f, com_full.F.real(), f, com3.F.real(), f, coml.F.real());
732
733 com3.interp_add_triplequad(sparse_com3);
734 coml.interp_add_even(sparse_coml);
735
736 ARTSGUI::plot(f, com_full.F.real(),
737 f, com3.F.real(),
738 f, coml.F.real());
739}
740
741int main() try {
742 make_wigner_ready(20, 30, 6);
743 std::cout << std::setprecision(15);
744
745// test_ls<LineShape::Type::DP>();
746// test_ls<LineShape::Type::LP>();
747 test_ls<LineShape::Type::VP>();
748// test_ls<LineShape::Type::SDVP>();
749// test_ls<LineShape::Type::HTP>();
750//
751 test_norm();
752//
754//
755 test_sparse();
756//
757 return EXIT_SUCCESS;
758} catch (std::exception& e) {
759 std::cerr << e.what() << '\n';
760 return EXIT_FAILURE;
761}
Declarations required for the calculation of absorption coefficients.
Numeric E0
TimeStep mean(const ArrayOfTimeStep &dt)
Returns the mean time step.
Definition: artstime.cc:190
TimeStep median(ArrayOfTimeStep dt)
Returns the median time step.
Definition: artstime.cc:179
Stuff related to time in ARTS.
Computations and data for a single absorption line.
Numeric F0() const noexcept
Central frequency.
This can be used to make arrays out of anything.
Definition: array.h:107
The ComplexMatrix class.
The ComplexVector class.
Index nelem() const ARTS_NOEXCEPT
Returns the number of elements.
Definition: matpackI.cc:48
Main line shape model class.
constexpr Numeric dNdF0() const noexcept
Definition: lineshape.h:604
constexpr Numeric dNdf(Numeric f) const noexcept
Definition: lineshape.h:603
constexpr Numeric dNdT(Numeric T, Numeric f) const noexcept
Definition: lineshape.h:602
Container class for Quantum Numbers.
Definition: quantum.h:112
Class to handle time in ARTS.
Definition: artstime.h:41
The Vector class.
Definition: matpackI.h:876
Main Zeeman Model.
Definition: zeemandata.h:357
Jacobian::Target JacobianTarget
Definition: jacobian.h:340
#define min(a, b)
#define max(a, b)
Numeric dsingle_partition_function_dT(const Numeric &T, const Species::IsotopeRecord &ir)
Computes the partition function temperature derivative.
Definition: linescaling.cc:28
Numeric single_partition_function(const Numeric &T, const Species::IsotopeRecord &ir)
Computes the partition function at one temperature.
Definition: linescaling.cc:23
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
constexpr Numeric imag(Complex c) noexcept
imag
const Joker joker
std::complex< Numeric > Complex
constexpr Numeric real(Complex c) noexcept
real
Temperature
Definition: jacobian.h:57
MagneticMagnitude
Definition: jacobian.h:59
WindMagnitude
Definition: jacobian.h:58
Vector triple_sparse_f_grid(const Vector &f_grid, const Numeric &sparse_df) noexcept
Definition: lineshape.cc:3209
Vector linear_sparse_f_grid(const Vector &f_grid, const Numeric &sparse_df) ARTS_NOEXCEPT
Definition: lineshape.cc:3174
void compute(ComputeData &com, ComputeData &sparse_com, const AbsorptionLines &band, const ArrayOfRetrievalQuantity &jacobian_quantities, const EnergyLevelMap &nlte, const Vector &vmrs, const ArrayOfSpeciesTag &self_tag, const Numeric &self_vmr, const Numeric &isot_ratio, const Numeric &P, const Numeric &T, const Numeric &H, const Numeric &sparse_lim, const bool do_zeeman, const Zeeman::Polarization zeeman_polarization, const Options::LblSpeedup speedup_type) ARTS_NOEXCEPT
Compute the line shape in its entirety.
Definition: lineshape.cc:3116
X3 LineCenter QuadraticIndependent
Definition: constants.h:577
X3 LineCenter None
Definition: constants.h:576
constexpr std::array Isotopologues
Definition: isotopologues.h:50
constexpr Rational start(Rational Ju, Rational Jl, Polarization type) noexcept
Gives the lowest M for a polarization type of this transition.
Definition: zeemandata.h:78
invlib::Vector< ArtsVector > Vector
invlib wrapper type for ARTS vectors.
Definition: oem.h:31
Quantum::Identifier QuantumIdentifier
Definition: quantum.h:471
#define N
Definition: rng.cc:164
#define M
Definition: rng.cc:165
QuantumIdentifier qid
AbsorptionLines band
ComplexVector N
Definition: lineshape.h:687
ComplexVector F
Definition: lineshape.h:687
ComplexMatrix dF
Definition: lineshape.h:688
ComplexMatrix dN
Definition: lineshape.h:688
constexpr Numeric dSdF0() const noexcept
Definition: lineshape.h:378
constexpr Numeric dSdT() const noexcept
Definition: lineshape.h:376
constexpr Numeric dSdI0() const noexcept
Definition: lineshape.h:377
Coefficients and temperature model for SingleSpeciesModel.
InternalData line_model()
void test_sparse()
void test_norm()
void test_lte_strength()
ENUMCLASS(ModificationInternal, char, None, T, f, H, F0, I0, SelfG0, SelfD0, SelfG2, SelfD2, SelfETA, SelfFVC, SelfY, SelfG, SelfDV, AirG0, AirD0, AirG2, AirD2, AirETA, AirFVC, AirY, AirG, AirDV, SelfVMR, AirVMR) template< ModificationInternal mod > LineShape
void test_ls()
int main()
Index make_wigner_ready(int largest, int fastest, int size)
Ready Wigner.