ARTS 2.5.4 (git: 31ce4f0e)
m_abs.cc
Go to the documentation of this file.
1
2/* Copyright (C) 2000-2012
3 Stefan Buehler <sbuehler@ltu.se>
4 Patrick Eriksson <patrick.eriksson@chalmers.se>
5 Axel von Engeln <engeln@uni-bremen.de>
6 Thomas Kuhn <tkuhn@uni-bremen.de>
7
8 This program is free software; you can redistribute it and/or modify it
9 under the terms of the GNU General Public License as published by the
10 Free Software Foundation; either version 2, or (at your option) any
11 later version.
12
13 This program is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with this program; if not, write to the Free Software
20 Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
21 USA. */
22
23//
24
33#include <algorithm>
34#include <cmath>
35#include <cstddef>
36#include <memory>
37#include <utility>
38
39#include "absorption.h"
40#include "absorptionlines.h"
41#include "agenda_class.h"
42#include "agenda_set.h"
43#include "array.h"
44#include "arts.h"
45#include "arts_constants.h"
46#include "arts_omp.h"
47#include "artstime.h"
48#include "auto_md.h"
49#include "check_input.h"
50#include "debug.h"
51#include "depr.h"
52#include "file.h"
53#include "global_data.h"
54#include "hitran_species.h"
55#include "jacobian.h"
56#include "legacy_continua.h"
57#include "lineshape.h"
58#include "m_xml.h"
59#include "math_funcs.h"
60#include "matpackI.h"
61#include "messages.h"
62#include "methods.h"
63#include "montecarlo.h"
64#include "optproperties.h"
65#include "parameters.h"
66#include "physics_funcs.h"
67#include "rte.h"
68#include "species_tags.h"
69#include "xml_io.h"
70
71#ifdef ENABLE_NETCDF
72#include <netcdf.h>
73
74#include "nc_io.h"
75#endif
76
79inline constexpr Numeric PI=Constant::pi;
82
83/* Workspace method: Doxygen documentation will be auto-generated */
84void AbsInputFromRteScalars( // WS Output:
85 Vector& abs_p,
86 Vector& abs_t,
87 Matrix& abs_vmrs,
88 // WS Input:
89 const Numeric& rtp_pressure,
90 const Numeric& rtp_temperature,
91 const Vector& rtp_vmr,
92 const Verbosity&) {
93 // Prepare abs_p:
94 abs_p.resize(1);
95 abs_p = rtp_pressure;
96
97 // Prepare abs_t:
98 abs_t.resize(1);
99 abs_t = rtp_temperature;
100
101 // Prepare abs_vmrs:
102 abs_vmrs.resize(rtp_vmr.nelem(), 1);
103 abs_vmrs = rtp_vmr;
104}
105
106/* Workspace method: Doxygen documentation will be auto-generated */
108 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
109 // WS Input:
110 const ArrayOfAbsorptionLines& abs_lines,
111 const ArrayOfArrayOfSpeciesTag& abs_species,
112 const Verbosity&) {
113 // Size is set but inner size will now change from the original definition of species tags...
114 abs_lines_per_species.resize(abs_species.nelem());
115
116 // The inner arrays need to be emptied, because they may contain lines
117 // from a previous calculation
118 for (auto& lines : abs_lines_per_species) lines.resize(0);
119
120#pragma omp parallel for schedule(dynamic) if (!arts_omp_in_parallel())
121 for (Index ilines = 0; ilines < abs_lines.nelem(); ilines++) {
122 AbsorptionLines lines = abs_lines[ilines];
123
124 // Skip empty lines
125 if (lines.NumLines() == 0) continue;
126
127 // Loop all the tags
128 for (Index i = 0; i < abs_species.nelem() and lines.NumLines(); i++) {
129 for (auto& this_tag : abs_species[i]) {
130 // Test isotopologue, we have to hit the end of the list for no isotopologue or the exact value
131 if (not same_or_joker(this_tag.Isotopologue(), lines.Isotopologue()))
132 continue;
133
134 // If there is a frequency range, we have to check so that only selected lines are included
135 if (this_tag.lower_freq >= 0 or this_tag.upper_freq >= 0) {
136 const Numeric low = (this_tag.lower_freq >= 0)
137 ? this_tag.lower_freq
138 : std::numeric_limits<Numeric>::lowest();
139 const Numeric upp = (this_tag.upper_freq >= 0)
140 ? this_tag.upper_freq
142
143 // Fill up a copy of the line record to match with the wished frequency criteria
144 AbsorptionLines these_lines = lines;
145 these_lines.lines.resize(0);
146 for (Index k = lines.NumLines() - 1; k >= 0; k--)
147 if (low <= lines.lines[k].F0 and upp >= lines.lines[k].F0)
148 these_lines.AppendSingleLine(lines.PopLine(k));
149
150 // Append these lines after sorting them if there are any of them
151 if (these_lines.NumLines()) {
152 these_lines.ReverseLines();
153#pragma omp critical
154 abs_lines_per_species[i].push_back(these_lines);
155 }
156
157 // If this means we have deleted all lines, then we leave
158 if (lines.NumLines() == 0) goto leave_inner_loop;
159 } else {
160#pragma omp critical
161 abs_lines_per_species[i].push_back(lines);
162 goto leave_inner_loop;
163 }
164 }
165 }
166 leave_inner_loop : {}
167 }
168
169 abs_lines_per_species.shrink_to_fit();
170 for (auto& spec_band : abs_lines_per_species)
171 std::sort(spec_band.begin(), spec_band.end(), [](auto& a, auto& b) {
172 return a.lines.size() and b.lines.size() and
173 a.lines.front().F0 < b.lines.front().F0;
174 });
175}
176
177/* Workspace method: Doxygen documentation will be auto-generated */
180 Index& propmat_clearsky_agenda_checked,
181 // Control Parameters:
182 const String& basename,
183 const Verbosity& verbosity) {
185
186 // Invalidate agenda check flags
187 propmat_clearsky_agenda_checked = false;
188
189 // We want to make lists of included and excluded species:
190 ArrayOfString included(0), excluded(0);
191
192 tgs.resize(0);
193
194 for (Index i = 0; i < Index(Species::Species::FINAL); ++i) {
195 const String specname = Species::toShortName(Species::Species(i));
196
197 String filename = basename;
198 if (basename.length() && basename[basename.length() - 1] != '/')
199 filename += ".";
200 filename += specname;
201
202 try {
203 find_xml_file(filename, verbosity);
204 // Add to included list:
205 included.push_back(specname);
206
207 // Add this tag group to tgs:
208 tgs.emplace_back(ArrayOfSpeciesTag(specname));
209 } catch (const std::runtime_error& e) {
210 // The file for the species could not be found.
211 excluded.push_back(specname);
212 }
213 }
214
215 // Some nice output:
216 out2 << " Included Species (" << included.nelem() << "):\n";
217 for (Index i = 0; i < included.nelem(); ++i)
218 out2 << " " << included[i] << "\n";
219
220 out2 << " Excluded Species (" << excluded.nelem() << "):\n";
221 for (Index i = 0; i < excluded.nelem(); ++i)
222 out2 << " " << excluded[i] << "\n";
223}
224
225/* Workspace method: Doxygen documentation will be auto-generated */
226void abs_speciesDefineAll( // WS Output:
227 ArrayOfArrayOfSpeciesTag& abs_species,
228 Index& propmat_clearsky_agenda_checked,
229 // Control Parameters:
230 const Verbosity& verbosity) {
231 // Species lookup data:
232
233 // We want to make lists of all species
234 ArrayOfString specs(0);
235 for (Index i = 0; i < Index(Species::Species::FINAL); ++i) {
236 if (Species::Species(i) not_eq Species::Species::Bath) {
237 specs.emplace_back(Species::toShortName(Species::Species(i)));
238 }
239 }
240
241 // Set the values
242 abs_speciesSet(abs_species,
243 propmat_clearsky_agenda_checked,
244 specs,
245 verbosity);
246}
247
248/* Workspace method: Doxygen documentation will be auto-generated */
249void AbsInputFromAtmFields( // WS Output:
250 Vector& abs_p,
251 Vector& abs_t,
252 Matrix& abs_vmrs,
253 // WS Input:
254 const Index& atmosphere_dim,
255 const Vector& p_grid,
256 const Tensor3& t_field,
257 const Tensor4& vmr_field,
258 const Verbosity&) {
259 // First, make sure that we really have a 1D atmosphere:
260 ARTS_USER_ERROR_IF(1 != atmosphere_dim,
261 "Atmospheric dimension must be 1D, but atmosphere_dim is ",
262 atmosphere_dim,
263 ".")
264
265 abs_p = p_grid;
266 abs_t = t_field(joker, 0, 0);
267 abs_vmrs = vmr_field(joker, joker, 0, 0);
268}
269
271 const ArrayOfVector& abs_cont_parameters,
272 const ArrayOfString& abs_cont_models) {
273 std::ostringstream os;
274
275 for (Index i = 0; i < abs_cont_names.nelem(); ++i)
276 os << "abs_xsec_per_speciesAddConts: " << i
277 << " name : " << abs_cont_names[i] << "\n";
278
279 for (Index i = 0; i < abs_cont_parameters.nelem(); ++i)
280 os << "abs_xsec_per_speciesAddConts: " << i
281 << " param: " << abs_cont_parameters[i] << "\n";
282
283 for (Index i = 0; i < abs_cont_models.nelem(); ++i)
284 os << "abs_xsec_per_speciesAddConts: " << i
285 << " option: " << abs_cont_models[i] << "\n";
286
287 os << "The following variables must have the same dimension:\n"
288 << "abs_cont_names: " << abs_cont_names.nelem() << "\n"
289 << "abs_cont_parameters: " << abs_cont_parameters.nelem();
290
291 return os.str();
292}
293
294//======================================================================
295// Methods related to continua
296//======================================================================
297
298/* Workspace method: Doxygen documentation will be auto-generated */
299void abs_cont_descriptionInit( // WS Output:
300 ArrayOfString& abs_cont_names,
301 ArrayOfString& abs_cont_options,
302 ArrayOfVector& abs_cont_parameters,
303 const Verbosity& verbosity) {
305
306 abs_cont_names.resize(0);
307 abs_cont_options.resize(0);
308 abs_cont_parameters.resize(0);
309 out2 << " Initialized abs_cont_names \n"
310 " abs_cont_models\n"
311 " abs_cont_parameters.\n";
312}
313
314/* Workspace method: Doxygen documentation will be auto-generated */
315void abs_cont_descriptionAppend( // WS Output:
316 ArrayOfString& abs_cont_names,
317 ArrayOfString& abs_cont_models,
318 ArrayOfVector& abs_cont_parameters,
319 // Control Parameters:
320 const String& tagname,
321 const String& model,
322 const Vector& userparameters,
323 const Verbosity&) {
324 // First we have to check that name matches a continuum species tag.
325 check_continuum_model(tagname);
326
327 //cout << " + tagname: " << tagname << "\n";
328 //cout << " + model: " << model << "\n";
329 //cout << " + parameters: " << userparameters << "\n";
330
331 // Add name and parameters to the apropriate variables:
332 abs_cont_names.push_back(tagname);
333 abs_cont_models.push_back(model);
334 abs_cont_parameters.push_back(userparameters);
335}
336
337/* Workspace method: Doxygen documentation will be auto-generated */
339 StokesVector& nlte_source,
340 ArrayOfStokesVector& dnlte_source_dx,
341 // WS Input:
342 const ArrayOfMatrix& src_coef_per_species,
343 const ArrayOfMatrix& dsrc_coef_dx,
344 const ArrayOfRetrievalQuantity& jacobian_quantities,
345 const Vector& f_grid,
346 const Numeric& rtp_temperature,
347 const Verbosity&) {
348 // nlte_source has format
349 // [ abs_species, f_grid, stokes_dim ].
350 // src_coef_per_species has format ArrayOfMatrix (over species),
351 // where for each species the matrix has format [f_grid, abs_p].
352
353 Index n_species = src_coef_per_species.nelem(); // # species
354
355 ARTS_USER_ERROR_IF(not n_species, "Must have at least one species.")
356
357 Index n_f = src_coef_per_species[0].nrows(); // # frequencies
358
359 // # pressures must be 1:
360 ARTS_USER_ERROR_IF(1 not_eq src_coef_per_species[0].ncols(),
361 "Must have exactly one pressure.")
362
363 // Check frequency dimension of propmat_clearsky
364 ARTS_USER_ERROR_IF(nlte_source.NumberOfFrequencies() not_eq n_f,
365 "Frequency dimension of nlte_source does not\n"
366 "match abs_coef_per_species.")
367
368 const Vector B = planck(f_grid, rtp_temperature);
369
370 StokesVector sv(n_f, nlte_source.StokesDimensions());
371 for (Index si = 0; si < n_species; ++si) {
372 sv.Kjj() = src_coef_per_species[si](joker, 0);
373 sv *= B;
374 nlte_source.Kjj() += sv.Kjj();
375 }
376
377 // Jacobian
378 for (Index ii = 0; ii < jacobian_quantities.nelem(); ii++) {
379 const auto& deriv = jacobian_quantities[ii];
380
381 if (not deriv.propmattype()) continue;
382
383 if (deriv == Jacobian::Atm::Temperature) {
384 const Vector dB = dplanck_dt(f_grid, rtp_temperature);
385
386 for (Index si = 0; si < n_species; ++si) {
387 sv.Kjj() = src_coef_per_species[si](joker, 0);
388 sv *= dB;
389 dnlte_source_dx[ii].Kjj() += sv.Kjj();
390 }
391
392 sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
393 sv *= B;
394 dnlte_source_dx[ii].Kjj() += sv.Kjj();
395 } else if (is_frequency_parameter(deriv)) {
396 const Vector dB = dplanck_df(f_grid, rtp_temperature);
397
398 for (Index si = 0; si < n_species; ++si) {
399 sv.Kjj() = src_coef_per_species[si](joker, 0);
400 sv *= dB;
401 dnlte_source_dx[ii].Kjj() += sv.Kjj();
402 }
403
404 sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
405 sv *= B;
406 dnlte_source_dx[ii].Kjj() += sv.Kjj();
407 } else {
408 sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
409 sv *= B;
410 dnlte_source_dx[ii].Kjj() += sv.Kjj();
411 }
412 }
413}
414
415/* Workspace method: Doxygen documentation will be auto-generated */
416void propmat_clearskyInit( //WS Output
417 PropagationMatrix& propmat_clearsky,
418 StokesVector& nlte_source,
419 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
420 ArrayOfStokesVector& dnlte_source_dx,
421 //WS Input
422 const ArrayOfRetrievalQuantity& jacobian_quantities,
423 const Vector& f_grid,
424 const Index& stokes_dim,
425 const Index& propmat_clearsky_agenda_checked,
426 const Verbosity&) {
427 const Index nf = f_grid.nelem();
428 const Index nq = jacobian_quantities.nelem();
429
431 !propmat_clearsky_agenda_checked,
432 "You must call *propmat_clearsky_agenda_checkedCalc* before calling this method.")
433
434 ARTS_USER_ERROR_IF(not nf, "No frequencies");
435
436 ARTS_USER_ERROR_IF(stokes_dim < 1 or stokes_dim > 4,
437 "stokes_dim not in [1, 2, 3, 4]");
438
439 // Set size of propmat_clearsky or reset it's values
440 if (propmat_clearsky.StokesDimensions() == stokes_dim and
441 propmat_clearsky.NumberOfFrequencies() == nf) {
442 propmat_clearsky.SetZero();
443 } else {
444 propmat_clearsky = PropagationMatrix(nf, stokes_dim);
445 }
446
447 // Set size of dpropmat_clearsky_dx or reset it's values
448 if (dpropmat_clearsky_dx.nelem() not_eq nq) {
449 dpropmat_clearsky_dx =
450 ArrayOfPropagationMatrix(nq, PropagationMatrix(nf, stokes_dim));
451 } else {
452 for (auto& pm : dpropmat_clearsky_dx) {
453 if (pm.StokesDimensions() == stokes_dim and
454 pm.NumberOfFrequencies() == nf) {
455 pm.SetZero();
456 } else {
457 pm = PropagationMatrix(nf, stokes_dim);
458 }
459 }
460 }
461
462 // Set size of nlte_source or reset it's values
463 if (nlte_source.StokesDimensions() == stokes_dim and
464 nlte_source.NumberOfFrequencies() == nf) {
465 nlte_source.SetZero();
466 } else {
467 nlte_source = StokesVector(nf, stokes_dim);
468 }
469
470 // Set size of dnlte_source_dx or reset it's values
471 if (dnlte_source_dx.nelem() not_eq nq) {
472 dnlte_source_dx = ArrayOfStokesVector(nq, StokesVector(nf, stokes_dim));
473 } else {
474 for (auto& pm : dnlte_source_dx) {
475 if (pm.StokesDimensions() == stokes_dim and
476 pm.NumberOfFrequencies() == nf) {
477 pm.SetZero();
478 } else {
479 pm = StokesVector(nf, stokes_dim);
480 }
481 }
482 }
483}
484
485/* Workspace method: Doxygen documentation will be auto-generated */
487 PropagationMatrix& propmat_clearsky,
488 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
489 const Index& stokes_dim,
490 const Index& atmosphere_dim,
491 const Vector& f_grid,
492 const ArrayOfArrayOfSpeciesTag& abs_species,
493 const ArrayOfSpeciesTag& select_abs_species,
494 const ArrayOfRetrievalQuantity& jacobian_quantities,
495 const Vector& rtp_vmr,
496 const Vector& rtp_los,
497 const Vector& rtp_mag,
498 const Verbosity&) {
499 Index ife = -1;
500 for (Index sp = 0; sp < abs_species.nelem() && ife < 0; sp++) {
501 if (abs_species[sp].FreeElectrons()) {
502 ife = sp;
503 }
504 }
505
506 ARTS_USER_ERROR_IF(ife < 0,
507 "Free electrons not found in *abs_species* and "
508 "Faraday rotation can not be calculated.");
509
510 // Allow early exit for lookup table calculations
511 if (select_abs_species.nelem() and select_abs_species not_eq abs_species[ife]) return;
512
513 // All the physical constants joined into one static constant:
514 // (abs as e defined as negative)
515 static const Numeric FRconst =
519
520 const bool do_magn_jac = do_magnetic_jacobian(jacobian_quantities);
521 const Numeric dmag = magnetic_field_perturbation(jacobian_quantities);
522
524 stokes_dim < 3,
525 "To include Faraday rotation, stokes_dim >= 3 is required.")
527 atmosphere_dim == 1 && rtp_los.nelem() < 1,
528 "For applying propmat_clearskyAddFaraday, los needs to be specified\n"
529 "(at least zenith angle component for atmosphere_dim==1),\n"
530 "but it is not.\n")
532 atmosphere_dim > 1 && rtp_los.nelem() < 2,
533 "For applying propmat_clearskyAddFaraday, los needs to be specified\n"
534 "(both zenith and azimuth angle components for atmosphere_dim>1),\n"
535 "but it is not.\n")
536
537 const Numeric ne = rtp_vmr[ife];
538
539 if (ne != 0 && (rtp_mag[0] != 0 || rtp_mag[1] != 0 || rtp_mag[2] != 0)) {
540 // Include remaining terms, beside /f^2
541 const Numeric c1 =
542 2 * FRconst *
544 rtp_los, rtp_mag[0], rtp_mag[1], rtp_mag[2], atmosphere_dim);
545
546 Numeric dc1_u = 0.0, dc1_v = 0.0, dc1_w = 0.0;
547 if (do_magn_jac) {
548 dc1_u = (2 * FRconst *
549 dotprod_with_los(rtp_los,
550 rtp_mag[0] + dmag,
551 rtp_mag[1],
552 rtp_mag[2],
553 atmosphere_dim) -
554 c1) /
555 dmag;
556 dc1_v = (2 * FRconst *
557 dotprod_with_los(rtp_los,
558 rtp_mag[0],
559 rtp_mag[1] + dmag,
560 rtp_mag[2],
561 atmosphere_dim) -
562 c1) /
563 dmag;
564 dc1_w = (2 * FRconst *
565 dotprod_with_los(rtp_los,
566 rtp_mag[0],
567 rtp_mag[1],
568 rtp_mag[2] + dmag,
569 atmosphere_dim) -
570 c1) /
571 dmag;
572 }
573
574 for (Index iv = 0; iv < f_grid.nelem(); iv++) {
575 const Numeric f2 = f_grid[iv] * f_grid[iv];
576 const Numeric r = ne * c1 / f2;
577 propmat_clearsky.AddFaraday(r, iv);
578
579 // The Jacobian loop
580 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
581 if (is_frequency_parameter(jacobian_quantities[iq]))
582 dpropmat_clearsky_dx[iq].AddFaraday(-2.0 * ne * r / f_grid[iv], iv);
583 else if (jacobian_quantities[iq] == Jacobian::Atm::MagneticU)
584 dpropmat_clearsky_dx[iq].AddFaraday(ne * dc1_u / f2, iv);
585 else if (jacobian_quantities[iq] == Jacobian::Atm::MagneticV)
586 dpropmat_clearsky_dx[iq].AddFaraday(ne * dc1_v / f2, iv);
587 else if (jacobian_quantities[iq] == Jacobian::Atm::MagneticW)
588 dpropmat_clearsky_dx[iq].AddFaraday(ne * dc1_w / f2, iv);
589 else if (jacobian_quantities[iq] == Jacobian::Atm::Electrons)
590 dpropmat_clearsky_dx[iq].AddFaraday(r, iv);
591 else if (jacobian_quantities[iq] == abs_species[ife])
592 dpropmat_clearsky_dx[iq].AddFaraday(r, iv);
593 }
594 }
595 }
596}
597
598/* Workspace method: Doxygen documentation will be auto-generated */
600 // WS Output:
601 PropagationMatrix& propmat_clearsky,
602 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
603 // WS Input:
604 const Index& stokes_dim,
605 const Index& atmosphere_dim,
606 const Vector& f_grid,
607 const ArrayOfArrayOfSpeciesTag& abs_species,
608 const ArrayOfSpeciesTag& select_abs_species,
609 const ArrayOfRetrievalQuantity& jacobian_quantities,
610 const Vector& rtp_vmr,
611 const Vector& rtp_los,
612 const Numeric& rtp_temperature,
613 const ArrayOfArrayOfSingleScatteringData& scat_data,
614 const Index& scat_data_checked,
615 const Index& use_abs_as_ext,
616 // Verbosity object:
617 const Verbosity& verbosity) {
619
620 ARTS_USER_ERROR_IF(select_abs_species.nelem(), R"--(
621 We do not yet support select_abs_species for lookup table calculations
622 )--")
623
624 // (i)yCalc only checks scat_data_checked if cloudbox is on. It is off here,
625 // though, i.e. we need to check it here explicitly. (Also, cloudboxOff sets
626 // scat_data_checked=0 as it does not check it and as we ususally don't need
627 // scat_data for clearsky cases, hence don't want to check them by
628 // scat_data_checkedCalc in that case. This approach seems to be the more
629 // handy compared to cloudboxOff setting scat_data_checked=1 without checking
630 // it assuming we won't use it anyways.)
631 ARTS_USER_ERROR_IF(scat_data_checked != 1,
632 "The scat_data must be flagged to have "
633 "passed a consistency check (scat_data_checked=1).")
634
635 const Index ns = TotalNumberOfElements(scat_data);
636 Index np = 0;
637 for (Index sp = 0; sp < abs_species.nelem(); sp++) {
638 if (abs_species[sp].Particles()) {
639 np++;
640 }
641 }
642
644 np == 0,
645 "For applying propmat_clearskyAddParticles, *abs_species* needs to"
646 "contain species 'particles', but it does not.\n")
647
649 ns != np,
650 "Number of 'particles' entries in abs_species and of elements in\n"
651 "*scat_data* needs to be identical. But you have ",
652 np,
653 " 'particles' entries\n"
654 "and ",
655 ns,
656 " *scat_data* elements.\n")
657
659 atmosphere_dim == 1 && rtp_los.nelem() < 1,
660 "For applying *propmat_clearskyAddParticles*, *rtp_los* needs to be specified\n"
661 "(at least zenith angle component for atmosphere_dim==1),\n"
662 "but it is not.\n")
664 atmosphere_dim > 1 && rtp_los.nelem() < 2,
665 "For applying *propmat_clearskyAddParticles*, *rtp_los* needs to be specified\n"
666 "(both zenith and azimuth angle components for atmosphere_dim>1),\n"
667 "but it is not.\n")
668
669 // Use for rescaling vmr of particulates
670 Numeric rtp_vmr_sum = 0.0;
671
672 // Tests and setup partial derivatives
673 const bool do_jac_temperature = do_temperature_jacobian(jacobian_quantities);
674 const bool do_jac_frequencies = do_frequency_jacobian(jacobian_quantities);
675 const Numeric dT = temperature_perturbation(jacobian_quantities);
676
677 const Index na = abs_species.nelem();
678 Vector rtp_los_back;
679 mirror_los(rtp_los_back, rtp_los, atmosphere_dim);
680
681 // 170918 JM: along with transition to use of new-type (aka
682 // pre-f_grid-interpolated) scat_data, freq perturbation switched off. Typical
683 // clear-sky freq perturbations yield insignificant effects in particle
684 // properties. Hence, this feature is neglected here.
685 if (do_jac_frequencies) {
686 out1 << "WARNING:\n"
687 << "Frequency perturbation not available for absorbing particles.\n";
688 }
689
690 // creating temporary output containers
691 ArrayOfArrayOfTensor5 ext_mat_Nse;
692 ArrayOfArrayOfTensor4 abs_vec_Nse;
693 ArrayOfArrayOfIndex ptypes_Nse;
694 Matrix t_ok;
695
696 // preparing input in format needed
697 Vector T_array;
698 if (do_jac_temperature) {
699 T_array.resize(2);
700 T_array = rtp_temperature;
701 T_array[1] += dT;
702 } else {
703 T_array.resize(1);
704 T_array = rtp_temperature;
705 }
706 Matrix dir_array(1, 2);
707 dir_array(0, joker) = rtp_los_back;
708
709 // ext/abs per scat element for all freqs at once
710 opt_prop_NScatElems(ext_mat_Nse,
711 abs_vec_Nse,
712 ptypes_Nse,
713 t_ok,
714 scat_data,
715 stokes_dim,
716 T_array,
717 dir_array,
718 -1);
719
720 const Index nf = abs_vec_Nse[0][0].nbooks();
721 Tensor3 tmp(nf, stokes_dim, stokes_dim);
722
723 // Internal computations necessary since it relies on zero start
724 PropagationMatrix internal_propmat(propmat_clearsky.NumberOfFrequencies(),
725 propmat_clearsky.StokesDimensions());
726
727 // loop over the scat_data and link them with correct vmr_field entry according
728 // to the position of the particle type entries in abs_species.
729 Index sp = 0;
730 Index i_se_flat = 0;
731 for (Index i_ss = 0; i_ss < scat_data.nelem(); i_ss++) {
732 for (Index i_se = 0; i_se < scat_data[i_ss].nelem(); i_se++) {
733 // forward to next particle entry in abs_species
734 while (sp < na && not abs_species[sp].Particles()) sp++;
735 internal_propmat.SetZero();
736
737 // running beyond number of abs_species entries when looking for next
738 // particle entry. shouldn't happen, though.
739 ARTS_ASSERT(sp < na);
741 rtp_vmr[sp] < 0.,
742 "Negative absorbing particle 'vmr' (aka number density)"
743 " encountered:\n"
744 "scat species #",
745 i_ss,
746 ", scat elem #",
747 i_se,
748 " (vmr_field entry #",
749 sp,
750 ")\n")
751
752 if (rtp_vmr[sp] > 0.) {
753 ARTS_USER_ERROR_IF(t_ok(i_se_flat, 0) < 0.,
754 "Temperature interpolation error:\n"
755 "scat species #",
756 i_ss,
757 ", scat elem #",
758 i_se,
759 "\n")
760 if (use_abs_as_ext) {
761 if (nf > 1)
762 for (Index iv = 0; iv < f_grid.nelem(); iv++)
763 internal_propmat.AddAbsorptionVectorAtPosition(
764 abs_vec_Nse[i_ss][i_se](iv, 0, 0, joker), iv);
765 else
766 for (Index iv = 0; iv < f_grid.nelem(); iv++)
767 internal_propmat.AddAbsorptionVectorAtPosition(
768 abs_vec_Nse[i_ss][i_se](0, 0, 0, joker), iv);
769 } else {
770 if (nf > 1)
771 for (Index iv = 0; iv < f_grid.nelem(); iv++)
772 internal_propmat.SetAtPosition(
773 ext_mat_Nse[i_ss][i_se](iv, 0, 0, joker, joker), iv);
774 else
775 for (Index iv = 0; iv < f_grid.nelem(); iv++)
776 internal_propmat.SetAtPosition(
777 ext_mat_Nse[i_ss][i_se](0, 0, 0, joker, joker), iv);
778 }
779 propmat_clearsky += rtp_vmr[sp] * internal_propmat;
780 }
781
782 // For temperature derivatives (so we don't need to check it in jac loop)
783 if (do_jac_temperature) {
785 t_ok(i_se_flat, 1) < 0.,
786 "Temperature interpolation error (in perturbation):\n"
787 "scat species #",
788 i_ss,
789 ", scat elem #",
790 i_se,
791 "\n")
792 }
793
794 // For number density derivatives
795 if (jacobian_quantities.nelem()) rtp_vmr_sum += rtp_vmr[sp];
796
797 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
798 const auto& deriv = jacobian_quantities[iq];
799
800 if (not deriv.propmattype()) continue;
801
802 if (deriv == Jacobian::Atm::Temperature) {
803 if (use_abs_as_ext) {
804 tmp(joker, joker, 0) = abs_vec_Nse[i_ss][i_se](joker, 1, 0, joker);
805 tmp(joker, joker, 0) -= abs_vec_Nse[i_ss][i_se](joker, 0, 0, joker);
806 } else {
807 tmp = ext_mat_Nse[i_ss][i_se](joker, 1, 0, joker, joker);
808 tmp -= ext_mat_Nse[i_ss][i_se](joker, 0, 0, joker, joker);
809 }
810
811 tmp *= rtp_vmr[sp];
812 tmp /= dT;
813
814 if (nf > 1)
815 for (Index iv = 0; iv < f_grid.nelem(); iv++)
816 if (use_abs_as_ext)
817 dpropmat_clearsky_dx[iq].AddAbsorptionVectorAtPosition(
818 tmp(iv, joker, 0), iv);
819 else
820 dpropmat_clearsky_dx[iq].AddAtPosition(tmp(iv, joker, joker),
821 iv);
822 else
823 for (Index iv = 0; iv < f_grid.nelem(); iv++)
824 if (use_abs_as_ext)
825 dpropmat_clearsky_dx[iq].AddAbsorptionVectorAtPosition(
826 tmp(0, joker, 0), iv);
827 else
828 dpropmat_clearsky_dx[iq].AddAtPosition(tmp(0, joker, joker),
829 iv);
830 }
831
832 else if (deriv == Jacobian::Atm::Particulates) {
833 for (Index iv = 0; iv < f_grid.nelem(); iv++)
834 dpropmat_clearsky_dx[iq].AddAtPosition(internal_propmat, iv);
835 }
836
837 else if (deriv == abs_species[sp]) {
838 dpropmat_clearsky_dx[iq] += internal_propmat;
839 }
840 }
841 sp++;
842 i_se_flat++;
843 }
844 }
845
846 //checking that no further 'particle' entry left after all scat_data entries
847 //are processes. this is basically not necessary. but checking it anyway to
848 //really be safe. remove later, when more extensively tested.
849 while (sp < na) {
850 ARTS_ASSERT(abs_species[sp][0].Type() != Species::TagType::Particles);
851 sp++;
852 }
853
854 if (rtp_vmr_sum != 0.0) {
855 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
856 const auto& deriv = jacobian_quantities[iq];
857
858 if (not deriv.propmattype()) continue;
859
860 if (deriv == Jacobian::Atm::Particulates) {
861 dpropmat_clearsky_dx[iq] /= rtp_vmr_sum;
862 }
863 }
864 }
865}
866
868 const Vector& f_grid,
869 const Numeric& sparse_df,
870 const String& speedup_option,
871 // Verbosity object:
872 const Verbosity&) {
873 // Return empty for nothing
874 if (not f_grid.nelem()) {
875 sparse_f_grid.resize(0);
876 return;
877 };
878
879 switch (Options::toLblSpeedupOrThrow(speedup_option)) {
880 case Options::LblSpeedup::LinearIndependent:
881 sparse_f_grid = LineShape::linear_sparse_f_grid(f_grid, sparse_df);
883 break;
884 case Options::LblSpeedup::QuadraticIndependent:
885 sparse_f_grid = LineShape::triple_sparse_f_grid(f_grid, sparse_df);
886 break;
887 case Options::LblSpeedup::None:
888 sparse_f_grid.resize(0);
889 break;
890 case Options::LblSpeedup::FINAL: { /* Leave last */
891 }
892 }
893}
894
896 const Numeric& sparse_df,
897 const String& speedup_option,
898 // Verbosity object:
899 const Verbosity& verbosity) {
900 Vector sparse_f_grid;
902 sparse_f_grid, f_grid, sparse_df, speedup_option, verbosity);
903 return sparse_f_grid;
904}
905
906/* Workspace method: Doxygen documentation will be auto-generated */
907void propmat_clearskyAddConts( // Workspace reference:
908 // WS Output:
909 PropagationMatrix& propmat_clearsky,
910 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
911 // WS Input:
912 const ArrayOfArrayOfSpeciesTag& abs_species,
913 const ArrayOfSpeciesTag& select_abs_species,
914 const ArrayOfRetrievalQuantity& jacobian_quantities,
915 const Vector& f_grid,
916 const Numeric& rtp_pressure,
917 const Numeric& rtp_temperature,
918 const Vector& rtp_vmr,
919 const ArrayOfString& abs_cont_names,
920 const ArrayOfVector& abs_cont_parameters,
921 const ArrayOfString& abs_cont_models,
922 // Verbosity object:
923 const Verbosity& verbosity) {
925
926 // Size of problem
927 const Index nf = f_grid.nelem();
928 const Index nq = jacobian_quantities.nelem();
929 const Index ns = abs_species.nelem();
930
931 // Possible things that can go wrong in this code (excluding line parameters)
932 check_abs_species(abs_species);
933 // Check, that dimensions of abs_cont_names and
934 // abs_cont_parameters are consistent...
935 ARTS_USER_ERROR_IF(abs_cont_names.nelem() != abs_cont_parameters.nelem(),
937 abs_cont_names, abs_cont_parameters, abs_cont_models))
938 ARTS_USER_ERROR_IF(rtp_vmr.nelem() not_eq abs_species.nelem(),
939 "*rtp_vmr* must match *abs_species*")
940 ARTS_USER_ERROR_IF(propmat_clearsky.NumberOfFrequencies() not_eq nf,
941 "*f_grid* must match *propmat_clearsky*")
943 not nq and (nq not_eq dpropmat_clearsky_dx.nelem()),
944 "*dpropmat_clearsky_dx* must match derived form of *jacobian_quantities*")
946 not nq and bad_propmat(dpropmat_clearsky_dx, f_grid),
947 "*dpropmat_clearsky_dx* must have frequency dim same as *f_grid*")
949 "Negative frequency (at least one value).")
951 "Must be sorted and increasing.")
953 "Negative VMR (at least one value).")
954 ARTS_USER_ERROR_IF(rtp_temperature <= 0, "Non-positive temperature")
955 ARTS_USER_ERROR_IF(rtp_pressure <= 0, "Non-positive pressure")
956
957 // Jacobian overhead START
958 /* NOTE: if any of the functions inside continuum tags could
959 be made to give partial derivatives, then that would
960 speed things up. Also be aware that line specific
961 parameters cannot be retrieved while using these
962 models. */
963 const Numeric df = frequency_perturbation(jacobian_quantities);
964 const Numeric dt = temperature_perturbation(jacobian_quantities);
965
966 const bool do_jac = supports_continuum(
967 jacobian_quantities); // Throws runtime error if line parameters are wanted since we cannot know if the line is in the Continuum...
968 const bool do_wind_jac = do_wind_jacobian(jacobian_quantities);
969 const bool do_temp_jac = do_temperature_jacobian(jacobian_quantities);
970
971 Vector dfreq;
972 Vector dabs_t{rtp_temperature + dt};
973
974 if (do_wind_jac) {
975 dfreq.resize(f_grid.nelem());
976 for (Index iv = 0; iv < f_grid.nelem(); iv++) dfreq[iv] = f_grid[iv] + df;
977 }
978
979 Vector normal(f_grid.nelem());
980 Vector jacs_df, jacs_dt;
981
982 if (do_jac) {
983 if (do_wind_jac) {
985 !std::isnormal(df), "df must be >0 and not NaN or Inf: ", df)
986 jacs_df.resize(f_grid.nelem());
987 }
988 if (do_temp_jac) {
990 !std::isnormal(dt), "dt must be >0 and not NaN or Inf: ", dt)
991 jacs_dt.resize(f_grid.nelem());
992 }
993 }
994 // Jacobian overhead END
995
996 const Vector abs_p{rtp_pressure};
997 const Vector abs_t{rtp_temperature};
998
999 // We set abs_h2o, abs_n2, and abs_o2 later, because we only want to
1000 // do it if the parameters are really needed.
1001
1002 out3 << " Calculating continuum spectra.\n";
1003
1004 // Loop tag groups:
1005 for (Index ispecies = 0; ispecies < ns; ispecies++) {
1006 if (select_abs_species.nelem() and
1007 select_abs_species not_eq abs_species[ispecies])
1008 continue;
1009
1010 // Skip it if there are no species or there is Zeeman requested
1011 if (not abs_species[ispecies].nelem() or abs_species[ispecies].Zeeman())
1012 continue;
1013
1014 // Go through the tags in the current tag group to see if they
1015 // are continuum tags:
1016 for (Index s = 0; s < abs_species[ispecies].nelem(); ++s) {
1017 // Continuum tags in the sense that we talk about here
1018 // (including complete absorption models) are marked by a special type.
1019 if (abs_species[ispecies][s].Type() ==
1020 Species::TagType::PredefinedLegacy) {
1021 // We have identified a continuum tag!
1022
1023 // Get only the continuum name. The full tag name is something like:
1024 // H2O-HITRAN96Self-*-*. We want only the `H2O-HITRAN96Self' part:
1025 const String name = abs_species[ispecies][s].Isotopologue().FullName();
1026
1027 if (name == "O2-MPM2020") continue;
1028
1029 // Check, if we have parameters for this model. For
1030 // this, the model name must be listed in
1031 // abs_cont_names.
1032 const Index n =
1033 find(abs_cont_names.begin(), abs_cont_names.end(), name) -
1034 abs_cont_names.begin();
1035
1036 // n==abs_cont_names.nelem() indicates that
1037 // the name was not found.
1038 ARTS_USER_ERROR_IF(n == abs_cont_names.nelem(),
1039 "Cannot find model ",
1040 name,
1041 " in abs_cont_names.")
1042
1043 // Ok, the tag specifies a valid continuum model and
1044 // we have continuum parameters.
1045
1046 if (out3.sufficient_priority()) {
1047 ostringstream os;
1048 os << " Adding " << name << " to tag group " << ispecies << ".\n";
1049 out3 << os.str();
1050 }
1051
1052 // Set abs_h2o, abs_n2, and abs_o2 from the first matching species.
1053 const Index h2o =
1054 find_first_species(abs_species, Species::fromShortName("H2O"));
1055 const Vector h2o_vmr{h2o == -1 ? 0.0 : rtp_vmr[h2o]};
1056 const Index n2 =
1057 find_first_species(abs_species, Species::fromShortName("N2"));
1058 const Vector n2_vmr{n2 == -1 ? 0.0 : rtp_vmr[h2o]};
1059 const Index o2 =
1060 find_first_species(abs_species, Species::fromShortName("O2"));
1061 const Vector o2_vmr{o2 == -1 ? 0.0 : rtp_vmr[o2]};
1062
1063 // Add the continuum for this tag. The parameters in
1064 // this call should be clear. The vmr is in
1065 // abs_vmrs(i,Range(joker)). The other vmr variables,
1066 // abs_h2o, abs_n2, and abs_o2 contains the real vmr of H2O,
1067 // N2, nad O2, which are needed as additional information for
1068 // certain continua:
1069 // abs_h2o for
1070 // O2-PWR88, O2-PWR93, O2-PWR98,
1071 // O2-MPM85, O2-MPM87, O2-MPM89, O2-MPM92, O2-MPM93,
1072 // O2-TRE05,
1073 // O2-SelfContStandardType, O2-SelfContMPM93, O2-SelfContPWR93,
1074 // N2-SelfContMPM93, N2-DryContATM01,
1075 // N2-CIArotCKDMT252, N2-CIAfunCKDMT252
1076 // abs_n2 for
1077 // H2O-SelfContCKD24, H2O-ForeignContCKD24,
1078 // O2-v0v0CKDMT100,
1079 // CO2-ForeignContPWR93, CO2-ForeignContHo66
1080 // abs_o2 for
1081 // N2-CIArotCKDMT252, N2-CIAfunCKDMT252
1082 normal = 0.0;
1083 xsec_continuum_tag(normal,
1084 name,
1085 abs_cont_parameters[n],
1086 abs_cont_models[n],
1087 f_grid,
1088 abs_p,
1089 abs_t,
1090 n2_vmr,
1091 h2o_vmr,
1092 o2_vmr,
1093 rtp_vmr[ispecies],
1094 verbosity);
1095 if (!do_jac) {
1096 normal *=
1097 number_density(rtp_pressure, rtp_temperature) * rtp_vmr[ispecies];
1098 propmat_clearsky.Kjj() += normal;
1099 } else // The Jacobian block
1100 {
1101 // Needs a reseted block here...
1102 for (Index iv = 0; iv < f_grid.nelem(); iv++) {
1103 if (do_wind_jac) jacs_df[iv] = 0.0;
1104 if (do_temp_jac) jacs_dt[iv] = 0.0;
1105 }
1106
1107 // Frequency calculations
1108 if (do_wind_jac)
1109 xsec_continuum_tag(jacs_df,
1110 name,
1111 abs_cont_parameters[n],
1112 abs_cont_models[n],
1113 dfreq,
1114 abs_p,
1115 abs_t,
1116 n2_vmr,
1117 h2o_vmr,
1118 o2_vmr,
1119 rtp_vmr[ispecies],
1120 verbosity);
1121
1122 //Temperature calculations
1123 if (do_temp_jac)
1124 xsec_continuum_tag(jacs_dt,
1125 name,
1126 abs_cont_parameters[n],
1127 abs_cont_models[n],
1128 f_grid,
1129 abs_p,
1130 dabs_t,
1131 n2_vmr,
1132 h2o_vmr,
1133 o2_vmr,
1134 rtp_vmr[ispecies],
1135 verbosity);
1136
1137 Numeric nd = number_density(rtp_pressure, rtp_temperature);
1138 Numeric dnd_dt = dnumber_density_dt(rtp_pressure, rtp_temperature);
1139 for (Index iv = 0; iv < f_grid.nelem(); iv++) {
1140 propmat_clearsky.Kjj()[iv] += normal[iv] * nd * rtp_vmr[ispecies];
1141 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
1142 const auto& deriv = jacobian_quantities[iq];
1143
1144 if (not deriv.propmattype()) continue;
1145
1146 if (is_wind_parameter(deriv)) {
1147 dpropmat_clearsky_dx[iq].Kjj()[iv] +=
1148 (jacs_df[iv] - normal[iv]) / df * nd * rtp_vmr[ispecies];
1149 } else if (deriv == Jacobian::Atm::Temperature) {
1150 dpropmat_clearsky_dx[iq].Kjj()[iv] +=
1151 ((jacs_dt[iv] - normal[iv]) / dt * nd +
1152 normal[iv] * dnd_dt) *
1153 rtp_vmr[ispecies];
1154 } else if (deriv == abs_species[ispecies]) {
1155 dpropmat_clearsky_dx[iq].Kjj()[iv] +=
1156 normal[iv] * nd;
1157 }
1158 }
1159 }
1160 }
1161 }
1162 }
1163 }
1164}
1165
1166/* Workspace method: Doxygen documentation will be auto-generated */
1167void propmat_clearskyAddLines( // Workspace reference:
1168 // WS Output:
1169 PropagationMatrix& propmat_clearsky,
1170 StokesVector& nlte_source,
1171 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
1172 ArrayOfStokesVector& dnlte_source_dx,
1173 // WS Input:
1174 const Vector& f_grid,
1175 const ArrayOfArrayOfSpeciesTag& abs_species,
1176 const ArrayOfSpeciesTag& select_abs_species,
1177 const ArrayOfRetrievalQuantity& jacobian_quantities,
1178 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1180 const Numeric& rtp_pressure,
1181 const Numeric& rtp_temperature,
1182 const EnergyLevelMap& rtp_nlte,
1183 const Vector& rtp_vmr,
1184 const Index& nlte_do,
1185 const Index& lbl_checked,
1186 // WS User Generic inputs
1187 const Numeric& sparse_df,
1188 const Numeric& sparse_lim,
1189 const String& speedup_option,
1190 const Index& robust,
1191 // Verbosity object:
1192 const Verbosity& verbosity) {
1193 // Size of problem
1194 const Index nf = f_grid.nelem();
1195 const Index nq = jacobian_quantities.nelem();
1196 const Index ns = abs_species.nelem();
1197
1198 // Possible things that can go wrong in this code (excluding line parameters)
1199 ARTS_USER_ERROR_IF(not lbl_checked, "Must check LBL calculations")
1200 check_abs_species(abs_species);
1201 ARTS_USER_ERROR_IF(rtp_vmr.nelem() not_eq abs_species.nelem(),
1202 "*rtp_vmr* must match *abs_species*")
1203 ARTS_USER_ERROR_IF(propmat_clearsky.NumberOfFrequencies() not_eq nf,
1204 "*f_grid* must match *propmat_clearsky*")
1205 ARTS_USER_ERROR_IF(nlte_source.NumberOfFrequencies() not_eq nf,
1206 "*f_grid* must match *nlte_source*")
1208 not nq and (nq not_eq dpropmat_clearsky_dx.nelem()),
1209 "*dpropmat_clearsky_dx* must match derived form of *jacobian_quantities*")
1211 not nq and bad_propmat(dpropmat_clearsky_dx, f_grid),
1212 "*dpropmat_clearsky_dx* must have frequency dim same as *f_grid*")
1214 nlte_do and (nq not_eq dnlte_source_dx.nelem()),
1215 "*dnlte_source_dx* must match derived form of *jacobian_quantities* when non-LTE is on")
1217 nlte_do and bad_propmat(dnlte_source_dx, f_grid),
1218 "*dnlte_source_dx* must have frequency dim same as *f_grid* when non-LTE is on")
1220 "Negative frequency (at least one value).")
1221 ARTS_USER_ERROR_IF((any_cutoff(abs_lines_per_species) or speedup_option not_eq "None") and not is_increasing(f_grid),
1222 "Must be sorted and increasing if any cutoff or speedup is used.")
1224 "Negative VMR (at least one value).")
1226 "Negative NLTE (at least one value).")
1227 ARTS_USER_ERROR_IF(rtp_temperature <= 0, "Non-positive temperature")
1228 ARTS_USER_ERROR_IF(rtp_pressure <= 0, "Non-positive pressure")
1230 sparse_lim > 0 and sparse_df > sparse_lim,
1231 "If sparse grids are to be used, the limit must be larger than the grid-spacing.\n"
1232 "The limit is ",
1233 sparse_lim,
1234 " Hz and the grid_spacing is ",
1235 sparse_df,
1236 " Hz")
1237
1238 if (not nf) return;
1239
1240 // Deal with sparse computational grid
1241 const Vector f_grid_sparse = create_sparse_f_grid_internal(
1242 f_grid, sparse_df, speedup_option, verbosity);
1243 const Options::LblSpeedup speedup_type =
1244 f_grid_sparse.nelem() ? Options::toLblSpeedupOrThrow(speedup_option)
1245 : Options::LblSpeedup::None;
1247 sparse_lim <= 0 and speedup_type not_eq Options::LblSpeedup::None,
1248 "Must have a sparse limit if you set speedup_option")
1249
1250 // Calculations data
1251 LineShape::ComputeData com(f_grid, jacobian_quantities, nlte_do);
1252 LineShape::ComputeData sparse_com(
1253 f_grid_sparse, jacobian_quantities, nlte_do);
1254
1255 if (arts_omp_in_parallel()) {
1256 for (Index ispecies = 0; ispecies < ns; ispecies++) {
1257 if (select_abs_species.nelem() and
1258 select_abs_species not_eq abs_species[ispecies])
1259 continue;
1260
1261 // Skip it if there are no species or there is Zeeman requested
1262 if (not abs_species[ispecies].nelem() or abs_species[ispecies].Zeeman() or
1263 not abs_lines_per_species[ispecies].nelem())
1264 continue;
1265
1266 for (auto& band : abs_lines_per_species[ispecies]) {
1268 sparse_com,
1269 band,
1270 jacobian_quantities,
1271 rtp_nlte,
1272 band.BroadeningSpeciesVMR(rtp_vmr, abs_species),
1273 abs_species[ispecies],
1274 rtp_vmr[ispecies],
1275 isotopologue_ratios[band.Isotopologue()],
1276 rtp_pressure,
1277 rtp_temperature,
1278 0,
1279 sparse_lim,
1281 speedup_type,
1282 robust not_eq 0);
1283 }
1284 }
1285 } else { // In parallel
1286 const Index nbands = [](auto& lines) {
1287 Index n = 0;
1288 for (auto& abs_lines : lines) n += abs_lines.nelem();
1289 return n;
1290 }(abs_lines_per_species);
1291
1292 std::vector<LineShape::ComputeData> vcom(
1295 f_grid, jacobian_quantities, static_cast<bool>(nlte_do)});
1296 std::vector<LineShape::ComputeData> vsparse_com(
1299 f_grid_sparse, jacobian_quantities, static_cast<bool>(nlte_do)});
1300
1301#pragma omp parallel for schedule(dynamic)
1302 for (Index i = 0; i < nbands; i++) {
1303 const auto [ispecies, iband] =
1304 flat_index(i, abs_species, abs_lines_per_species);
1305
1306 if (select_abs_species.nelem() and
1307 select_abs_species not_eq abs_species[ispecies])
1308 continue;
1309
1310 // Skip it if there are no species or there is Zeeman requested
1311 if (not abs_species[ispecies].nelem() or abs_species[ispecies].Zeeman() or
1312 not abs_lines_per_species[ispecies].nelem())
1313 continue;
1314
1315 auto& band = abs_lines_per_species[ispecies][iband];
1317 vsparse_com[arts_omp_get_thread_num()],
1318 band,
1319 jacobian_quantities,
1320 rtp_nlte,
1321 band.BroadeningSpeciesVMR(rtp_vmr, abs_species),
1322 abs_species[ispecies],
1323 rtp_vmr[ispecies],
1324 isotopologue_ratios[band.Isotopologue()],
1325 rtp_pressure,
1326 rtp_temperature,
1327 0,
1328 sparse_lim,
1330 speedup_type,
1331 robust not_eq 0);
1332 }
1333
1334 for (auto& pcom: vcom) com += pcom;
1335 for (auto& pcom: vsparse_com) sparse_com += pcom;
1336 }
1337
1338 switch (speedup_type) {
1339 case Options::LblSpeedup::LinearIndependent:
1340 com.interp_add_even(sparse_com);
1341 break;
1342 case Options::LblSpeedup::QuadraticIndependent:
1343 com.interp_add_triplequad(sparse_com);
1344 break;
1345 case Options::LblSpeedup::None: /* Do nothing */
1346 break;
1347 case Options::LblSpeedup::FINAL: { /* Leave last */
1348 }
1349 }
1350
1351 // Sum up the propagation matrix
1352 propmat_clearsky.Kjj() += com.F.real();
1353
1354 // Sum up the Jacobian
1355 for (Index j = 0; j < nq; j++) {
1356 if (not jacobian_quantities[j].propmattype()) continue;
1357 dpropmat_clearsky_dx[j].Kjj() += com.dF.real()(joker, j);
1358 }
1359
1360 if (nlte_do) {
1361 // Sum up the source vector
1362 nlte_source.Kjj() += com.N.real();
1363
1364 // Sum up the Jacobian
1365 for (Index j = 0; j < nq; j++) {
1366 if (not jacobian_quantities[j].propmattype()) continue;
1367 dnlte_source_dx[j].Kjj() += com.dN.real()(joker, j);
1368 }
1369 }
1370}
1371
1372/* Workspace method: Doxygen documentation will be auto-generated */
1374 const Vector& f_grid,
1375 const Index& stokes_dim,
1376 const Verbosity&) {
1377 propmat_clearsky = PropagationMatrix(f_grid.nelem(), stokes_dim);
1378}
1379
1380/* Workspace method: Doxygen documentation will be auto-generated */
1382 const Verbosity&) {
1383 for (Index i = 0; i < propmat_clearsky.NumberOfFrequencies(); i++)
1384 if (propmat_clearsky.Kjj()[i] < 0.0) propmat_clearsky.SetAtPosition(0.0, i);
1385}
1386
1387/* Workspace method: Doxygen documentation will be auto-generated */
1391}
1392
1393/* Workspace method: Doxygen documentation will be auto-generated */
1397}
1398
1399#ifdef ENABLE_NETCDF
1400/* Workspace method: Doxygen documentation will be auto-generated */
1401/* Included by Claudia Emde, 20100707 */
1402void WriteMolTau( //WS Input
1403 const Vector& f_grid,
1404 const Tensor3& z_field,
1405 const Tensor7& propmat_clearsky_field,
1406 const Index& atmosphere_dim,
1407 //Keyword
1408 const String& filename,
1409 const Verbosity&) {
1410 int retval, ncid;
1411 int nlev_dimid, nlyr_dimid, nwvl_dimid, stokes_dimid, none_dimid;
1412 int dimids[4];
1413 int wvlmin_varid, wvlmax_varid, z_varid, wvl_varid, tau_varid;
1414
1415 ARTS_USER_ERROR_IF(atmosphere_dim != 1,
1416 "WriteMolTau can only be used for atmosphere_dim=1")
1417#pragma omp critical(netcdf__critical_region)
1418 {
1419 // Open file
1420 if ((retval = nc_create(filename.c_str(), NC_CLOBBER, &ncid)))
1421 nca_error(retval, "nc_create");
1422
1423 // Define dimensions
1424 if ((retval = nc_def_dim(ncid, "nlev", (int)z_field.npages(), &nlev_dimid)))
1425 nca_error(retval, "nc_def_dim");
1426
1427 if ((retval =
1428 nc_def_dim(ncid, "nlyr", (int)z_field.npages() - 1, &nlyr_dimid)))
1429 nca_error(retval, "nc_def_dim");
1430
1431 if ((retval = nc_def_dim(ncid, "nwvl", (int)f_grid.nelem(), &nwvl_dimid)))
1432 nca_error(retval, "nc_def_dim");
1433
1434 if ((retval = nc_def_dim(ncid, "none", 1, &none_dimid)))
1435 nca_error(retval, "nc_def_dim");
1436
1437 if ((retval = nc_def_dim(ncid,
1438 "nstk",
1439 (int)propmat_clearsky_field.nbooks(),
1440 &stokes_dimid)))
1441 nca_error(retval, "nc_def_dim");
1442
1443 // Define variables
1444 if ((retval = nc_def_var(
1445 ncid, "wvlmin", NC_DOUBLE, 1, &none_dimid, &wvlmin_varid)))
1446 nca_error(retval, "nc_def_var wvlmin");
1447
1448 if ((retval = nc_def_var(
1449 ncid, "wvlmax", NC_DOUBLE, 1, &none_dimid, &wvlmax_varid)))
1450 nca_error(retval, "nc_def_var wvlmax");
1451
1452 if ((retval = nc_def_var(ncid, "z", NC_DOUBLE, 1, &nlev_dimid, &z_varid)))
1453 nca_error(retval, "nc_def_var z");
1454
1455 if ((retval =
1456 nc_def_var(ncid, "wvl", NC_DOUBLE, 1, &nwvl_dimid, &wvl_varid)))
1457 nca_error(retval, "nc_def_var wvl");
1458
1459 dimids[0] = nlyr_dimid;
1460 dimids[1] = nwvl_dimid;
1461 dimids[2] = stokes_dimid;
1462 dimids[3] = stokes_dimid;
1463
1464 if ((retval =
1465 nc_def_var(ncid, "tau", NC_DOUBLE, 4, &dimids[0], &tau_varid)))
1466 nca_error(retval, "nc_def_var tau");
1467
1468 // Units
1469 if ((retval = nc_put_att_text(ncid, wvlmin_varid, "units", 2, "nm")))
1470 nca_error(retval, "nc_put_att_text");
1471
1472 if ((retval = nc_put_att_text(ncid, wvlmax_varid, "units", 2, "nm")))
1473 nca_error(retval, "nc_put_att_text");
1474
1475 if ((retval = nc_put_att_text(ncid, z_varid, "units", 2, "km")))
1476 nca_error(retval, "nc_put_att_text");
1477
1478 if ((retval = nc_put_att_text(ncid, wvl_varid, "units", 2, "nm")))
1479 nca_error(retval, "nc_put_att_text");
1480
1481 if ((retval = nc_put_att_text(ncid, tau_varid, "units", 1, "-")))
1482 nca_error(retval, "nc_put_att_text");
1483
1484 // End define mode. This tells netCDF we are done defining
1485 // metadata.
1486 if ((retval = nc_enddef(ncid))) nca_error(retval, "nc_enddef");
1487
1488 // Assign data
1489 double wvlmin[1];
1490 wvlmin[0] = SPEED_OF_LIGHT / f_grid[f_grid.nelem() - 1] * 1e9;
1491 if ((retval = nc_put_var_double(ncid, wvlmin_varid, &wvlmin[0])))
1492 nca_error(retval, "nc_put_var");
1493
1494 double wvlmax[1];
1495 wvlmax[0] = SPEED_OF_LIGHT / f_grid[0] * 1e9;
1496 if ((retval = nc_put_var_double(ncid, wvlmax_varid, &wvlmax[0])))
1497 nca_error(retval, "nc_put_var");
1498
1499 double z[z_field.npages()];
1500 for (int iz = 0; iz < z_field.npages(); iz++)
1501 z[iz] = z_field(z_field.npages() - 1 - iz, 0, 0) * 1e-3;
1502
1503 if ((retval = nc_put_var_double(ncid, z_varid, &z[0])))
1504 nca_error(retval, "nc_put_var");
1505
1506 double wvl[f_grid.nelem()];
1507 for (int iv = 0; iv < f_grid.nelem(); iv++)
1508 wvl[iv] = SPEED_OF_LIGHT / f_grid[f_grid.nelem() - 1 - iv] * 1e9;
1509
1510 if ((retval = nc_put_var_double(ncid, wvl_varid, &wvl[0])))
1511 nca_error(retval, "nc_put_var");
1512
1513 const Index zfnp = z_field.npages() - 1;
1514 const Index fgne = f_grid.nelem();
1515 const Index amfnb = propmat_clearsky_field.nbooks();
1516
1517 Tensor4 tau(zfnp, fgne, amfnb, amfnb, 0.);
1518
1519 // Calculate average tau for layers
1520 for (int is = 0; is < propmat_clearsky_field.nlibraries(); is++)
1521 for (int iz = 0; iz < zfnp; iz++)
1522 for (int iv = 0; iv < fgne; iv++)
1523 for (int is1 = 0; is1 < amfnb; is1++)
1524 for (int is2 = 0; is2 < amfnb; is2++)
1525 // sum up all species
1526 tau(iz, iv, is1, is2) +=
1527 0.5 *
1528 (propmat_clearsky_field(is,
1529 f_grid.nelem() - 1 - iv,
1530 is1,
1531 is2,
1532 z_field.npages() - 1 - iz,
1533 0,
1534 0) +
1535 propmat_clearsky_field(is,
1536 f_grid.nelem() - 1 - iv,
1537 is1,
1538 is2,
1539 z_field.npages() - 2 - iz,
1540 0,
1541 0)) *
1542 (z_field(z_field.npages() - 1 - iz, 0, 0) -
1543 z_field(z_field.npages() - 2 - iz, 0, 0));
1544
1545 if ((retval = nc_put_var_double(ncid, tau_varid, tau.get_c_array())))
1546 nca_error(retval, "nc_put_var");
1547
1548 // Close the file
1549 if ((retval = nc_close(ncid))) nca_error(retval, "nc_close");
1550 }
1551}
1552
1553#else
1554
1555void WriteMolTau( //WS Input
1556 const Vector& f_grid _U_,
1557 const Tensor3& z_field _U_,
1558 const Tensor7& propmat_clearsky_field _U_,
1559 const Index& atmosphere_dim _U_,
1560 //Keyword
1561 const String& filename _U_,
1562 const Verbosity&) {
1563 ARTS_USER_ERROR_IF(true,
1564 "The workspace method WriteMolTau is not available"
1565 "because ARTS was compiled without NetCDF support.");
1566}
1567
1568#endif /* ENABLE_NETCDF */
1569
1570void propmat_clearsky_agendaAuto(// Workspace reference:
1571 Workspace& ws,
1572 // WS Output:
1573 Agenda& propmat_clearsky_agenda,
1574 Index& propmat_clearsky_agenda_checked,
1575 // WS Input:
1576 const ArrayOfArrayOfSpeciesTag& abs_species,
1577 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1578 // WS Generic Input:
1579 const Numeric& H,
1580 const Numeric& T_extrapolfac,
1581 const Numeric& eta,
1582 const Numeric& extpolfac,
1583 const Numeric& force_p,
1584 const Numeric& force_t,
1585 const Index& ignore_errors,
1586 const Numeric& lines_sparse_df,
1587 const Numeric& lines_sparse_lim,
1588 const String& lines_speedup_option,
1589 const Index& manual_mag_field,
1590 const Index& no_negatives,
1591 const Numeric& theta,
1592 const Index& use_abs_as_ext,
1593 const Index& use_abs_lookup_ind,
1594 // Verbosity object:
1595 const Verbosity& verbosity) {
1596 using namespace AgendaManip;
1597
1598 propmat_clearsky_agenda_checked = 0; // In case of crash
1599
1600 AgendaCreator agenda(ws, "propmat_clearsky_agenda");
1601
1602 // Use bool because logic is easier
1603 const bool use_abs_lookup = static_cast<bool>(use_abs_lookup_ind);
1604
1605 const SpeciesTagTypeStatus any_species(abs_species);
1606 const AbsorptionTagTypesStatus any_lines(abs_lines_per_species);
1607
1608 // propmat_clearskyInit
1609 agenda.add("propmat_clearskyInit");
1610
1611 // propmat_clearskyAddFromLookup
1612 if (use_abs_lookup) {
1613 agenda.add("propmat_clearskyAddFromLookup",
1614 SetWsv{"extpolfac", extpolfac},
1615 SetWsv{"no_negatives", no_negatives});
1616 }
1617
1618 // propmat_clearskyAddLines
1619 if (not use_abs_lookup and any_species.Plain and
1620 (any_lines.population.LTE or any_lines.population.NLTE or
1621 any_lines.population.VibTemps)) {
1622 agenda.add("propmat_clearskyAddLines",
1623 SetWsv{"lines_sparse_df", lines_sparse_df},
1624 SetWsv{"lines_sparse_lim", lines_sparse_lim},
1625 SetWsv{"lines_speedup_option", lines_speedup_option},
1626 SetWsv{"no_negatives", no_negatives});
1627 }
1628
1629 // propmat_clearskyAddZeeman
1630 if (any_species.Zeeman and
1631 (any_lines.population.LTE or any_lines.population.NLTE or
1632 any_lines.population.VibTemps)) {
1633 agenda.add("propmat_clearskyAddZeeman",
1634 SetWsv{"manual_mag_field", manual_mag_field},
1635 SetWsv{"H", H},
1636 SetWsv{"theta", theta},
1637 SetWsv{"eta", eta});
1638 }
1639
1640 //propmat_clearskyAddHitranXsec
1641 if (not use_abs_lookup and any_species.XsecFit) {
1642 agenda.add("propmat_clearskyAddXsecFit",
1643 SetWsv{"force_p", force_p},
1644 SetWsv{"force_t", force_t});
1645 }
1646
1647 //propmat_clearskyAddOnTheFlyLineMixing
1648 if (not use_abs_lookup and any_species.Plain and
1649 (any_lines.population.ByMakarovFullRelmat or
1651 agenda.add("propmat_clearskyAddOnTheFlyLineMixing");
1652 }
1653
1654 //propmat_clearskyAddOnTheFlyLineMixingWithZeeman
1655 if (any_species.Zeeman and
1656 (any_lines.population.ByMakarovFullRelmat or
1658 agenda.add("propmat_clearskyAddOnTheFlyLineMixingWithZeeman");
1659 }
1660
1661 //propmat_clearskyAddCIA
1662 if (not use_abs_lookup and any_species.Cia) {
1663 agenda.add("propmat_clearskyAddCIA",
1664 SetWsv{"T_extrapolfac", T_extrapolfac},
1665 SetWsv{"ignore_errors", ignore_errors});
1666 }
1667
1668 //propmat_clearskyAddConts
1669 if (not use_abs_lookup and any_species.PredefinedLegacy) {
1670 agenda.add("propmat_clearskyAddConts");
1671 }
1672
1673 //propmat_clearskyAddPredefined
1674 if (not use_abs_lookup and any_species.PredefinedModern) {
1675 agenda.add("propmat_clearskyAddPredefined");
1676 }
1677
1678 //propmat_clearskyAddParticles
1679 if (any_species.Particles) {
1680 agenda.add("propmat_clearskyAddParticles",
1681 SetWsv{"use_abs_as_ext", use_abs_as_ext});
1682 }
1683
1684 //propmat_clearskyAddFaraday
1685 if (any_species.FreeElectrons) {
1686 agenda.add("propmat_clearskyAddFaraday");
1687 }
1688
1689 // propmat_clearskyAddHitranLineMixingLines
1690 if (not use_abs_lookup and any_species.Plain and
1691 (any_lines.population.ByHITRANFullRelmat or
1693 agenda.add("propmat_clearskyAddHitranLineMixingLines");
1694 }
1695
1696 // Extra check (should really never ever fail when species exist)
1697 propmat_clearsky_agenda = agenda.finalize();
1698 propmat_clearsky_agenda_checked = 1;
1699
1701 if (out3.sufficient_priority()) {
1702 out3 << "propmat_clearsky_agendaAuto sets propmat_clearsky_agenda to:\n\n"
1703 << propmat_clearsky_agenda << '\n';
1704 }
1705}
1706
1708 ArrayOfString& abs_cont_models,
1709 ArrayOfVector& abs_cont_parameters,
1710 const Verbosity& verbosity) {
1711 const Vector no_extras{};
1712
1713 //
1714 // initialize the continua tag structures
1715 //
1717 abs_cont_names, abs_cont_models, abs_cont_parameters, verbosity);
1718
1719 //
1720 // define the continua tags with the appropriate input
1721 //
1722
1723 //
1724 // H2O continuum
1725 //
1726 // Rosenkranz-type H2O-H2O continuum:
1727 abs_cont_descriptionAppend(abs_cont_names,
1728 abs_cont_models,
1729 abs_cont_parameters,
1730 "H2O-SelfContStandardType",
1731 "Rosenkranz",
1732 no_extras,
1733 verbosity);
1734
1735 // Rosenkranz-type H2O-dry air continuum
1736 abs_cont_descriptionAppend(abs_cont_names,
1737 abs_cont_models,
1738 abs_cont_parameters,
1739 "H2O-ForeignContStandardType",
1740 "Rosenkranz",
1741 no_extras,
1742 verbosity);
1743
1744 // Rosenkranz H2O full absorption model
1745 abs_cont_descriptionAppend(abs_cont_names,
1746 abs_cont_models,
1747 abs_cont_parameters,
1748 "H2O-PWR98",
1749 "Rosenkranz",
1750 no_extras,
1751 verbosity);
1752
1753 // Liebe 1993 H2O full absorption model
1754 abs_cont_descriptionAppend(abs_cont_names,
1755 abs_cont_models,
1756 abs_cont_parameters,
1757 "H2O-MPM93",
1758 "MPM93",
1759 no_extras,
1760 verbosity);
1761
1762 // Liebe 1989 H2O full absorption model
1763 abs_cont_descriptionAppend(abs_cont_names,
1764 abs_cont_models,
1765 abs_cont_parameters,
1766 "H2O-MPM89",
1767 "MPM89",
1768 no_extras,
1769 verbosity);
1770
1771 // Liebe 1987 H2O full absorption model
1772 abs_cont_descriptionAppend(abs_cont_names,
1773 abs_cont_models,
1774 abs_cont_parameters,
1775 "H2O-MPM87",
1776 "MPM87",
1777 no_extras,
1778 verbosity);
1779
1780 // Ma Tipping Continuum model
1781 abs_cont_descriptionAppend(abs_cont_names,
1782 abs_cont_models,
1783 abs_cont_parameters,
1784 "H2O-ForeignContMaTippingType",
1785 "MaTipping",
1786 no_extras,
1787 verbosity);
1788
1789 // MPM93-type H2O-air continuum
1790 abs_cont_descriptionAppend(abs_cont_names,
1791 abs_cont_models,
1792 abs_cont_parameters,
1793 "H2O-ContMPM93",
1794 "MPM93",
1795 no_extras,
1796 verbosity);
1797
1798 // CKD2.42 original parameters
1799 abs_cont_descriptionAppend(abs_cont_names,
1800 abs_cont_models,
1801 abs_cont_parameters,
1802 "H2O-SelfContCKD242",
1803 "CKD242",
1804 no_extras,
1805 verbosity);
1806
1807 // CKD2.42 original parameters
1808 abs_cont_descriptionAppend(abs_cont_names,
1809 abs_cont_models,
1810 abs_cont_parameters,
1811 "H2O-ForeignContCKD242",
1812 "CKD242",
1813 no_extras,
1814 verbosity);
1815
1816 // CKD2.4 original parameters
1817 abs_cont_descriptionAppend(abs_cont_names,
1818 abs_cont_models,
1819 abs_cont_parameters,
1820 "H2O-SelfContCKD24",
1821 "CKD24",
1822 no_extras,
1823 verbosity);
1824
1825 // CKD2.4 original parameters
1826 abs_cont_descriptionAppend(abs_cont_names,
1827 abs_cont_models,
1828 abs_cont_parameters,
1829 "H2O-ForeignContCKD24",
1830 "CKD24",
1831 no_extras,
1832 verbosity);
1833
1834 // CKD2.22 original parameters
1835 abs_cont_descriptionAppend(abs_cont_names,
1836 abs_cont_models,
1837 abs_cont_parameters,
1838 "H2O-SelfContCKD222",
1839 "CKD222",
1840 no_extras,
1841 verbosity);
1842
1843 // CKD2.22 original parameters
1844 abs_cont_descriptionAppend(abs_cont_names,
1845 abs_cont_models,
1846 abs_cont_parameters,
1847 "H2O-ForeignContCKD222",
1848 "CKD222",
1849 no_extras,
1850 verbosity);
1851
1852 //
1853 // CKD MT 1.00 self continuum original parameters:
1854 abs_cont_descriptionAppend(abs_cont_names,
1855 abs_cont_models,
1856 abs_cont_parameters,
1857 "H2O-SelfContCKDMT100",
1858 "CKDMT100",
1859 no_extras,
1860 verbosity);
1861 //
1862 // CKD MT 1.00 H2O foreign continuum original parameters:
1863 abs_cont_descriptionAppend(abs_cont_names,
1864 abs_cont_models,
1865 abs_cont_parameters,
1866 "H2O-ForeignContCKDMT100",
1867 "CKDMT100",
1868 no_extras,
1869 verbosity);
1870
1871 //
1872 // CKD MT 2.52 self continuum original parameters:
1873 abs_cont_descriptionAppend(abs_cont_names,
1874 abs_cont_models,
1875 abs_cont_parameters,
1876 "H2O-SelfContCKDMT252",
1877 "CKDMT252",
1878 no_extras,
1879 verbosity);
1880 //
1881 // CKD MT 2.52 H2O foreign continuum original parameters:
1882 abs_cont_descriptionAppend(abs_cont_names,
1883 abs_cont_models,
1884 abs_cont_parameters,
1885 "H2O-ForeignContCKDMT252",
1886 "CKDMT252",
1887 no_extras,
1888 verbosity);
1889
1890 //
1891 // CKD MT 3.20 self continuum original parameters:
1892 abs_cont_descriptionAppend(abs_cont_names,
1893 abs_cont_models,
1894 abs_cont_parameters,
1895 "H2O-SelfContCKDMT320",
1896 "CKDMT320",
1897 no_extras,
1898 verbosity);
1899 //
1900 // CKD MT 3.20 H2O foreign continuum original parameters:
1901 abs_cont_descriptionAppend(abs_cont_names,
1902 abs_cont_models,
1903 abs_cont_parameters,
1904 "H2O-ForeignContCKDMT320",
1905 "CKDMT320",
1906 no_extras,
1907 verbosity);
1908 //
1909 // Cruz-Pol et al. 98 Continuum:
1910 abs_cont_descriptionAppend(abs_cont_names,
1911 abs_cont_models,
1912 abs_cont_parameters,
1913 "H2O-CP98",
1914 "CruzPol",
1915 no_extras,
1916 verbosity);
1917
1918 //
1919 // Pardo et al. model 2001 Continuum:
1920 abs_cont_descriptionAppend(abs_cont_names,
1921 abs_cont_models,
1922 abs_cont_parameters,
1923 "H2O-ForeignContATM01",
1924 "ATM",
1925 no_extras,
1926 verbosity);
1927
1928 //
1929 // O2 continuum
1930 //
1931 // Standard O2-air continuum
1932 abs_cont_descriptionAppend(abs_cont_names,
1933 abs_cont_models,
1934 abs_cont_parameters,
1935 "O2-SelfContStandardType",
1936 "Rosenkranz",
1937 no_extras,
1938 verbosity);
1939
1940 //
1941 // Liebe 1993 O2-air continuum
1942 abs_cont_descriptionAppend(abs_cont_names,
1943 abs_cont_models,
1944 abs_cont_parameters,
1945 "O2-SelfContMPM93",
1946 "MPM93",
1947 no_extras,
1948 verbosity);
1949
1950 // Rosenkranz 1998 full oxygen continuum
1951 abs_cont_descriptionAppend(abs_cont_names,
1952 abs_cont_models,
1953 abs_cont_parameters,
1954 "O2-PWR98",
1955 "Rosenkranz",
1956 no_extras,
1957 verbosity);
1958
1959 // Rosenkranz 1993 full oxygen continuum
1960 abs_cont_descriptionAppend(abs_cont_names,
1961 abs_cont_models,
1962 abs_cont_parameters,
1963 "O2-PWR93",
1964 "Rosenkranz",
1965 no_extras,
1966 verbosity);
1967
1968 // Rosenkranz 1988 full oxygen continuum
1969 abs_cont_descriptionAppend(abs_cont_names,
1970 abs_cont_models,
1971 abs_cont_parameters,
1972 "O2-PWR88",
1973 "Rosenkranz",
1974 no_extras,
1975 verbosity);
1976
1977 // Liebe 1993 O2 full absorption model
1978 abs_cont_descriptionAppend(abs_cont_names,
1979 abs_cont_models,
1980 abs_cont_parameters,
1981 "O2-MPM93",
1982 "MPM93",
1983 no_extras,
1984 verbosity);
1985
1986 // Liebe 1992 O2 full absorption model
1987 abs_cont_descriptionAppend(abs_cont_names,
1988 abs_cont_models,
1989 abs_cont_parameters,
1990 "O2-MPM92",
1991 "MPM92",
1992 no_extras,
1993 verbosity);
1994
1995 // Liebe 1989 O2 full absorption model
1996 abs_cont_descriptionAppend(abs_cont_names,
1997 abs_cont_models,
1998 abs_cont_parameters,
1999 "O2-MPM89",
2000 "MPM89",
2001 no_extras,
2002 verbosity);
2003
2004 // Liebe 1987 O2 full absorption model
2005 abs_cont_descriptionAppend(abs_cont_names,
2006 abs_cont_models,
2007 abs_cont_parameters,
2008 "O2-MPM87",
2009 "MPM87",
2010 no_extras,
2011 verbosity);
2012
2013 // Liebe 1985 O2 full absorption model
2014 abs_cont_descriptionAppend(abs_cont_names,
2015 abs_cont_models,
2016 abs_cont_parameters,
2017 "O2-MPM85",
2018 "MPM85",
2019 no_extras,
2020 verbosity);
2021
2022 // Tretyakov etal. 2005 modernization of Liebe 1993 O2 full absorption model
2023 abs_cont_descriptionAppend(abs_cont_names,
2024 abs_cont_models,
2025 abs_cont_parameters,
2026 "O2-TRE05",
2027 "TRE05",
2028 no_extras,
2029 verbosity);
2030
2031 //
2032 // Thibault et al. from CKD F77 program:
2033 abs_cont_descriptionAppend(abs_cont_names,
2034 abs_cont_models,
2035 abs_cont_parameters,
2036 "O2-CIAfunCKDMT100",
2037 "CKDMT100",
2038 no_extras,
2039 verbosity);
2040
2041 // Mate et al. from CKD F77 program:
2042 abs_cont_descriptionAppend(abs_cont_names,
2043 abs_cont_models,
2044 abs_cont_parameters,
2045 "O2-v0v0CKDMT100",
2046 "CKDMT100",
2047 no_extras,
2048 verbosity);
2049
2050 // Mlawer et al. from CKD F77 program:
2051 abs_cont_descriptionAppend(abs_cont_names,
2052 abs_cont_models,
2053 abs_cont_parameters,
2054 "O2-v1v0CKDMT100",
2055 "CKDMT100",
2056 no_extras,
2057 verbosity);
2058
2059 // Mlawer et al. from CKD F77 program:
2060 abs_cont_descriptionAppend(abs_cont_names,
2061 abs_cont_models,
2062 abs_cont_parameters,
2063 "O2-v1v0CKDMT100",
2064 "CKDMT100",
2065 no_extras,
2066 verbosity);
2067
2068 // Greenblatt et al. from CKD F77 program:
2069 abs_cont_descriptionAppend(abs_cont_names,
2070 abs_cont_models,
2071 abs_cont_parameters,
2072 "O2-visCKDMT252",
2073 "CKDMT252",
2074 no_extras,
2075 verbosity);
2076
2077 //
2078
2079 //
2080 // N2 continuum
2081 //
2082 // Rosenkranz N2-N2 continuum (only N2-N2 broadening;
2083 abs_cont_descriptionAppend(abs_cont_names,
2084 abs_cont_models,
2085 abs_cont_parameters,
2086 "N2-SelfContStandardType",
2087 "Rosenkranz",
2088 no_extras,
2089 verbosity);
2090
2091 // MPM93 N2-N2 continuum (only N2-N2 broadening);
2092 abs_cont_descriptionAppend(abs_cont_names,
2093 abs_cont_models,
2094 abs_cont_parameters,
2095 "N2-SelfContMPM93",
2096 "MPM93",
2097 no_extras,
2098 verbosity);
2099
2100 // Rosenkranz N2-N2 continuum (only N2-N2 broadening);
2101 abs_cont_descriptionAppend(abs_cont_names,
2102 abs_cont_models,
2103 abs_cont_parameters,
2104 "N2-SelfContPWR93",
2105 "Rosenkranz",
2106 no_extras,
2107 verbosity);
2108
2109 // Borysow-Frommhold original parameters
2110 abs_cont_descriptionAppend(abs_cont_names,
2111 abs_cont_models,
2112 abs_cont_parameters,
2113 "N2-SelfContBorysow",
2114 "BF86",
2115 no_extras,
2116 verbosity);
2117
2118 // Borysow-Frommhold from CKD F77 program:
2119 abs_cont_descriptionAppend(abs_cont_names,
2120 abs_cont_models,
2121 abs_cont_parameters,
2122 "N2-CIArotCKDMT100",
2123 "CKDMT100",
2124 no_extras,
2125 verbosity);
2126 //
2127 //
2128 // Lafferty et al. from CKD F77 program:
2129 abs_cont_descriptionAppend(abs_cont_names,
2130 abs_cont_models,
2131 abs_cont_parameters,
2132 "N2-CIAfunCKDMT100",
2133 "CKDMT100",
2134 no_extras,
2135 verbosity);
2136
2137 // Borysow-Frommhold from CKD F77 program:
2138 abs_cont_descriptionAppend(abs_cont_names,
2139 abs_cont_models,
2140 abs_cont_parameters,
2141 "N2-CIArotCKDMT252",
2142 "CKDMT252",
2143 no_extras,
2144 verbosity);
2145 //
2146 //
2147 // Lafferty et al. from CKD F77 program:
2148 abs_cont_descriptionAppend(abs_cont_names,
2149 abs_cont_models,
2150 abs_cont_parameters,
2151 "N2-CIAfunCKDMT252",
2152 "CKDMT252",
2153 no_extras,
2154 verbosity);
2155
2156 //
2157 // Pardo et al. model 2001 Continuum:
2158 abs_cont_descriptionAppend(abs_cont_names,
2159 abs_cont_models,
2160 abs_cont_parameters,
2161 "N2-DryContATM01",
2162 "ATM",
2163 no_extras,
2164 verbosity);
2165
2166 //
2167
2168 //
2169 // CO2 continuum
2170 //
2171 // Rosenkranz CO2-CO2 continuum:
2172 abs_cont_descriptionAppend(abs_cont_names,
2173 abs_cont_models,
2174 abs_cont_parameters,
2175 "CO2-SelfContPWR93",
2176 "Rosenkranz",
2177 no_extras,
2178 verbosity);
2179
2180 // Rosenkranz CO2-N2 continuum:
2181 abs_cont_descriptionAppend(abs_cont_names,
2182 abs_cont_models,
2183 abs_cont_parameters,
2184 "CO2-ForeignContPWR93",
2185 "Rosenkranz",
2186 no_extras,
2187 verbosity);
2188
2189 // Ho et al 1966 CO2-CO2 continuum:
2190 abs_cont_descriptionAppend(abs_cont_names,
2191 abs_cont_models,
2192 abs_cont_parameters,
2193 "CO2-SelfContHo66",
2194 "Ho66",
2195 no_extras,
2196 verbosity);
2197
2198 // Ho et al 1966 CO2-N2 continuum:
2199 abs_cont_descriptionAppend(abs_cont_names,
2200 abs_cont_models,
2201 abs_cont_parameters,
2202 "CO2-ForeignContHo66",
2203 "Ho66",
2204 no_extras,
2205 verbosity);
2206
2207 // far line wing continuum from CKD F77 program:
2208 abs_cont_descriptionAppend(abs_cont_names,
2209 abs_cont_models,
2210 abs_cont_parameters,
2211 "CO2-CKDMT100",
2212 "CKDMT100",
2213 no_extras,
2214 verbosity);
2215
2216 // far line wing continuum from CKD F77 program:
2217 abs_cont_descriptionAppend(abs_cont_names,
2218 abs_cont_models,
2219 abs_cont_parameters,
2220 "CO2-CKDMT252",
2221 "CKDMT252",
2222 no_extras,
2223 verbosity);
2224
2225 // far line wing continuum from CKD F77 program:
2226 abs_cont_descriptionAppend(abs_cont_names,
2227 abs_cont_models,
2228 abs_cont_parameters,
2229 "CO2-CKD241",
2230 "CKD241",
2231 no_extras,
2232 verbosity);
2233
2234 //
2235 // Hydrometeors
2236 //
2237
2238 // Rosenkranz liquid cloud
2239 abs_cont_descriptionAppend(abs_cont_names,
2240 abs_cont_models,
2241 abs_cont_parameters,
2242 "liquidcloud-MPM93",
2243 "MPM93",
2244 no_extras,
2245 verbosity);
2246
2247 // Ellison liquid cloud (2007, no_extras, verbosity );
2248 abs_cont_descriptionAppend(abs_cont_names,
2249 abs_cont_models,
2250 abs_cont_parameters,
2251 "liquidcloud-ELL07",
2252 "ELL07",
2253 no_extras,
2254 verbosity);
2255
2256 // Rosenkranz ice cloud
2257 abs_cont_descriptionAppend(abs_cont_names,
2258 abs_cont_models,
2259 abs_cont_parameters,
2260 "icecloud-MPM93",
2261 "MPM93",
2262 no_extras,
2263 verbosity);
2264
2265 // Rosenkranz rain
2266 abs_cont_descriptionAppend(abs_cont_names,
2267 abs_cont_models,
2268 abs_cont_parameters,
2269 "rain-MPM93",
2270 "MPM93",
2271 no_extras,
2272 verbosity);
2273}
Declarations required for the calculation of absorption coefficients.
AbsorptionSpeciesBandIndex flat_index(Index i, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species)
Get a flat index pair for species and band.
Contains the absorption namespace.
Declarations for agendas.
This file contains the definition of Array.
Index TotalNumberOfElements(const Array< Array< base > > &aa)
Determine total number of elements in an ArrayOfArray.
Definition: array.h:241
The global header file for ARTS.
Constants of physical expressions as constexpr.
int arts_omp_get_max_threads()
Wrapper for omp_get_max_threads.
Definition: arts_omp.cc:46
int arts_omp_get_thread_num()
Wrapper for omp_get_thread_num.
Definition: arts_omp.cc:76
bool arts_omp_in_parallel()
Wrapper for omp_in_parallel.
Definition: arts_omp.cc:62
Header file for helper functions for OpenMP.
Stuff related to time in ARTS.
The Agenda class.
Definition: agenda_class.h:69
This can be used to make arrays out of anything.
Definition: array.h:48
Index nelem() const ARTS_NOEXCEPT
Definition: array.h:92
MatrixView real()
Get a view of the real part of the matrix.
VectorView real()
Get a view of the real part of the vector.
Index npages() const
Returns the number of pages.
Definition: matpackIII.h:140
Index nlibraries() const noexcept
Definition: matpackVII.h:153
Index nbooks() const noexcept
Definition: matpackVII.h:156
Index nelem() const noexcept
Returns the number of elements.
Definition: matpackI.h:536
The Matrix class.
Definition: matpackI.h:1261
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1010
void AddAbsorptionVectorAtPosition(const ConstVectorView &x, const Index iv=0, const Index iz=0, const Index ia=0)
Adds as a Stokes vector at position.
Index NumberOfFrequencies() const
The number of frequencies of the propagation matrix.
VectorView Kjj(const Index iz=0, const Index ia=0)
Vector view to diagonal elements.
void SetZero()
Sets all data to zero.
Index StokesDimensions() const
The stokes dimension of the propagation matrix.
void AddFaraday(const Numeric &rot, const Index iv=0, const Index iz=0, const Index ia=0)
Adds the Faraday rotation to the PropagationMatrix at required position.
void SetAtPosition(const PropagationMatrix &x, const Index iv=0, const Index iz=0, const Index ia=0)
Set the At Position object.
Stokes vector is as Propagation matrix but only has 4 possible values.
The Tensor3 class.
Definition: matpackIII.h:346
const Numeric * get_c_array() const ARTS_NOEXCEPT
Conversion to plain C-array.
Definition: matpackIV.cc:342
The Tensor4 class.
Definition: matpackIV.h:429
The Tensor7 class.
Definition: matpackVII.h:2399
The Vector class.
Definition: matpackI.h:899
void resize(Index n)
Resize function.
Definition: matpackI.cc:390
Workspace class.
Definition: workspace_ng.h:53
#define _U_
Definition: config.h:180
Helper macros for debugging.
#define ARTS_ASSERT(condition,...)
Definition: debug.h:83
#define ARTS_USER_ERROR_IF(condition,...)
Definition: debug.h:134
void find_xml_file(String &filename, const Verbosity &verbosity)
Find an xml file.
Definition: file.cc:373
This file contains basic functions to handle ASCII files.
bool is_wind_parameter(const RetrievalQuantity &t) noexcept
Returns if the Retrieval quantity is a wind parameter in propagation matrix calculations.
Definition: jacobian.cc:1002
bool do_frequency_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants a frequency derivative.
Definition: jacobian.cc:1175
bool supports_continuum(const ArrayOfRetrievalQuantity &js)
Returns if the array supports continuum derivatives.
Definition: jacobian.cc:1097
bool do_magnetic_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants a magnetic derivative.
Definition: jacobian.cc:1179
bool do_temperature_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants the temperature derivative.
Definition: jacobian.cc:1150
Numeric magnetic_field_perturbation(const ArrayOfRetrievalQuantity &js) noexcept
Returns the magnetic field perturbation if it exists.
Definition: jacobian.cc:1199
bool do_wind_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants a wind-based frequency derivative derivative.
Definition: jacobian.cc:1171
bool is_frequency_parameter(const RetrievalQuantity &t) noexcept
Returns if the Retrieval quantity is a frequency parameter in propagation matrix calculations.
Definition: jacobian.cc:1006
Numeric temperature_perturbation(const ArrayOfRetrievalQuantity &js) noexcept
Returns the temperature perturbation if it exists.
Definition: jacobian.cc:1183
Numeric frequency_perturbation(const ArrayOfRetrievalQuantity &js) noexcept
Returns the frequency perturbation if it exists.
Definition: jacobian.cc:1191
Routines for setting up the jacobian.
#define abs(x)
void xsec_continuum_tag(MatrixView xsec, const String &name, ConstVectorView parameters, const String &model, ConstVectorView f_grid, ConstVectorView abs_p, ConstVectorView abs_t, ConstVectorView abs_n2, ConstVectorView abs_h2o, ConstVectorView abs_o2, ConstVectorView vmr, const Verbosity &verbosity)
Calculates model absorption for one continuum or full model tag.
#define ns
void check_continuum_model(const String &name)
An auxiliary functions that checks if a given continuum model is listed in Species::Isotopologues.
#define max(a, b)
bool is_increasing(ConstVectorView x)
Checks if a vector is sorted and strictly increasing.
Definition: logic.cc:218
void sparse_f_gridFromFrequencyGrid(Vector &sparse_f_grid, const Vector &f_grid, const Numeric &sparse_df, const String &speedup_option, const Verbosity &)
WORKSPACE METHOD: sparse_f_gridFromFrequencyGrid.
Definition: m_abs.cc:867
constexpr Numeric ELECTRON_MASS
Definition: m_abs.cc:78
void isotopologue_ratiosInitFromHitran(SpeciesIsotopologueRatios &isotopologue_ratios, const Verbosity &)
WORKSPACE METHOD: isotopologue_ratiosInitFromHitran.
Definition: m_abs.cc:1394
constexpr Numeric SPEED_OF_LIGHT
Definition: m_abs.cc:80
Vector create_sparse_f_grid_internal(const Vector &f_grid, const Numeric &sparse_df, const String &speedup_option, const Verbosity &verbosity)
Definition: m_abs.cc:895
void LegacyContinuaInit(ArrayOfString &abs_cont_names, ArrayOfString &abs_cont_models, ArrayOfVector &abs_cont_parameters, const Verbosity &verbosity)
WORKSPACE METHOD: LegacyContinuaInit.
Definition: m_abs.cc:1707
void propmat_clearskyInit(PropagationMatrix &propmat_clearsky, StokesVector &nlte_source, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, ArrayOfStokesVector &dnlte_source_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, const Vector &f_grid, const Index &stokes_dim, const Index &propmat_clearsky_agenda_checked, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyInit.
Definition: m_abs.cc:416
void nlte_sourceFromTemperatureAndSrcCoefPerSpecies(StokesVector &nlte_source, ArrayOfStokesVector &dnlte_source_dx, const ArrayOfMatrix &src_coef_per_species, const ArrayOfMatrix &dsrc_coef_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, const Vector &f_grid, const Numeric &rtp_temperature, const Verbosity &)
Definition: m_abs.cc:338
void abs_speciesDefineAllInScenario(ArrayOfArrayOfSpeciesTag &tgs, Index &propmat_clearsky_agenda_checked, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: abs_speciesDefineAllInScenario.
Definition: m_abs.cc:178
void abs_cont_descriptionAppend(ArrayOfString &abs_cont_names, ArrayOfString &abs_cont_models, ArrayOfVector &abs_cont_parameters, const String &tagname, const String &model, const Vector &userparameters, const Verbosity &)
WORKSPACE METHOD: abs_cont_descriptionAppend.
Definition: m_abs.cc:315
constexpr Numeric VACUUM_PERMITTIVITY
Definition: m_abs.cc:81
void abs_speciesDefineAll(ArrayOfArrayOfSpeciesTag &abs_species, Index &propmat_clearsky_agenda_checked, const Verbosity &verbosity)
WORKSPACE METHOD: abs_speciesDefineAll.
Definition: m_abs.cc:226
String continua_model_error_message(const ArrayOfString &abs_cont_names, const ArrayOfVector &abs_cont_parameters, const ArrayOfString &abs_cont_models)
Definition: m_abs.cc:270
void abs_cont_descriptionInit(ArrayOfString &abs_cont_names, ArrayOfString &abs_cont_options, ArrayOfVector &abs_cont_parameters, const Verbosity &verbosity)
WORKSPACE METHOD: abs_cont_descriptionInit.
Definition: m_abs.cc:299
void AbsInputFromRteScalars(Vector &abs_p, Vector &abs_t, Matrix &abs_vmrs, const Numeric &rtp_pressure, const Numeric &rtp_temperature, const Vector &rtp_vmr, const Verbosity &)
Definition: m_abs.cc:84
void propmat_clearskyForceNegativeToZero(PropagationMatrix &propmat_clearsky, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyForceNegativeToZero.
Definition: m_abs.cc:1381
void propmat_clearskyAddLines(PropagationMatrix &propmat_clearsky, StokesVector &nlte_source, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, ArrayOfStokesVector &dnlte_source_dx, const Vector &f_grid, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfSpeciesTag &select_abs_species, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const SpeciesIsotopologueRatios &isotopologue_ratios, const Numeric &rtp_pressure, const Numeric &rtp_temperature, const EnergyLevelMap &rtp_nlte, const Vector &rtp_vmr, const Index &nlte_do, const Index &lbl_checked, const Numeric &sparse_df, const Numeric &sparse_lim, const String &speedup_option, const Index &robust, const Verbosity &verbosity)
WORKSPACE METHOD: propmat_clearskyAddLines.
Definition: m_abs.cc:1167
void propmat_clearskyAddConts(PropagationMatrix &propmat_clearsky, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfSpeciesTag &select_abs_species, const ArrayOfRetrievalQuantity &jacobian_quantities, const Vector &f_grid, const Numeric &rtp_pressure, const Numeric &rtp_temperature, const Vector &rtp_vmr, const ArrayOfString &abs_cont_names, const ArrayOfVector &abs_cont_parameters, const ArrayOfString &abs_cont_models, const Verbosity &verbosity)
WORKSPACE METHOD: propmat_clearskyAddConts.
Definition: m_abs.cc:907
void propmat_clearskyZero(PropagationMatrix &propmat_clearsky, const Vector &f_grid, const Index &stokes_dim, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyZero.
Definition: m_abs.cc:1373
void propmat_clearskyAddParticles(PropagationMatrix &propmat_clearsky, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &f_grid, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfSpeciesTag &select_abs_species, const ArrayOfRetrievalQuantity &jacobian_quantities, const Vector &rtp_vmr, const Vector &rtp_los, const Numeric &rtp_temperature, const ArrayOfArrayOfSingleScatteringData &scat_data, const Index &scat_data_checked, const Index &use_abs_as_ext, const Verbosity &verbosity)
WORKSPACE METHOD: propmat_clearskyAddParticles.
Definition: m_abs.cc:599
constexpr Numeric PI
Definition: m_abs.cc:79
void propmat_clearskyAddFaraday(PropagationMatrix &propmat_clearsky, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, const Index &stokes_dim, const Index &atmosphere_dim, const Vector &f_grid, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfSpeciesTag &select_abs_species, const ArrayOfRetrievalQuantity &jacobian_quantities, const Vector &rtp_vmr, const Vector &rtp_los, const Vector &rtp_mag, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyAddFaraday.
Definition: m_abs.cc:486
void WriteMolTau(const Vector &f_grid, const Tensor3 &z_field, const Tensor7 &propmat_clearsky_field, const Index &atmosphere_dim, const String &filename, const Verbosity &)
WORKSPACE METHOD: WriteMolTau.
Definition: m_abs.cc:1402
void AbsInputFromAtmFields(Vector &abs_p, Vector &abs_t, Matrix &abs_vmrs, const Index &atmosphere_dim, const Vector &p_grid, const Tensor3 &t_field, const Tensor4 &vmr_field, const Verbosity &)
WORKSPACE METHOD: AbsInputFromAtmFields.
Definition: m_abs.cc:249
void abs_lines_per_speciesCreateFromLines(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfAbsorptionLines &abs_lines, const ArrayOfArrayOfSpeciesTag &abs_species, const Verbosity &)
WORKSPACE METHOD: abs_lines_per_speciesCreateFromLines.
Definition: m_abs.cc:107
void isotopologue_ratiosInitFromBuiltin(SpeciesIsotopologueRatios &isotopologue_ratios, const Verbosity &)
WORKSPACE METHOD: isotopologue_ratiosInitFromBuiltin.
Definition: m_abs.cc:1388
constexpr Numeric ELECTRON_CHARGE
Definition: m_abs.cc:77
void propmat_clearsky_agendaAuto(Workspace &ws, Agenda &propmat_clearsky_agenda, Index &propmat_clearsky_agenda_checked, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Numeric &H, const Numeric &T_extrapolfac, const Numeric &eta, const Numeric &extpolfac, const Numeric &force_p, const Numeric &force_t, const Index &ignore_errors, const Numeric &lines_sparse_df, const Numeric &lines_sparse_lim, const String &lines_speedup_option, const Index &manual_mag_field, const Index &no_negatives, const Numeric &theta, const Index &use_abs_as_ext, const Index &use_abs_lookup_ind, const Verbosity &verbosity)
WORKSPACE METHOD: propmat_clearsky_agendaAuto.
Definition: m_abs.cc:1570
void abs_speciesSet(ArrayOfArrayOfSpeciesTag &abs_species, Index &propmat_clearsky_agenda_checked, const ArrayOfString &names, const Verbosity &verbosity)
WORKSPACE METHOD: abs_speciesSet.
Workspace methods and template functions for supergeneric XML IO.
constexpr bool any_negative(const MatpackType &var) noexcept
Checks for negative values.
Definition: math_funcs.h:147
Implementation of Matrix, Vector, and such stuff.
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
const Joker joker
Declarations having to do with the four output streams.
#define CREATE_OUT1
Definition: messages.h:204
#define CREATE_OUT3
Definition: messages.h:206
#define CREATE_OUT2
Definition: messages.h:205
Declaration of the class MdRecord.
Index nelem(const Lines &l)
Number of lines.
bool any_cutoff(const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species)
constexpr Numeric pi
The following mathematical constants are generated in python Decimal package by the code:
constexpr Numeric vacuum_permittivity
Vacuum permittivity [F/m].
constexpr Numeric electron_mass
Mass of resting electron [kg].
constexpr Numeric k
Boltzmann constant convenience name [J/K].
constexpr Numeric speed_of_light
Speed of light [m/s] From: https://en.wikipedia.org/wiki/2019_redefinition_of_SI_base_units Date: 201...
constexpr Numeric e
Elementary charge convenience name [C].
constexpr Numeric elementary_charge
Elementary charge [C] From: https://en.wikipedia.org/wiki/2019_redefinition_of_SI_base_units Date: 20...
SpeciesIsotopologueRatios isotopologue_ratios()
Temperature
Definition: jacobian.h:59
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 Zeeman::Polarization zeeman_polarization, const Options::LblSpeedup speedup_type, const bool robust) ARTS_NOEXCEPT
Compute the absorption of an absorption band.
Definition: lineshape.cc:3534
bool good_linear_sparse_f_grid(const Vector &f_grid_dense, const Vector &f_grid_sparse) noexcept
Definition: lineshape.cc:3637
Vector triple_sparse_f_grid(const Vector &f_grid, const Numeric &sparse_df) noexcept
Definition: lineshape.cc:3651
Vector linear_sparse_f_grid(const Vector &f_grid, const Numeric &sparse_df) ARTS_NOEXCEPT
Definition: lineshape.cc:3614
Online Transmission UseSurfaceRtprop TransmitterReceiverPath FreeElectrons
Definition: arts_options.h:114
constexpr bool same_or_joker(const IsotopeRecord &ir1, const IsotopeRecord &ir2) noexcept
constexpr IsotopologueRatios isotopologue_ratiosInitFromBuiltin()
Implements Zeeman modeling.
Definition: zeemandata.cc:333
void nca_error(const int e, const String s)
Throws a runtime error for the given NetCDF error code.
Definition: nc_io.cc:622
This file contains basic functions to handle NetCDF data files.
void opt_prop_NScatElems(ArrayOfArrayOfTensor5 &ext_mat, ArrayOfArrayOfTensor4 &abs_vec, ArrayOfArrayOfIndex &ptypes, Matrix &t_ok, const ArrayOfArrayOfSingleScatteringData &scat_data, const Index &stokes_dim, const Vector &T_array, const Matrix &dir_array, const Index &f_index, const Index &t_interp_order)
Extinction and absorption from all scattering elements.
Scattering database structure and functions.
This file contains header information for the dealing with command line parameters.
Numeric dplanck_df(const Numeric &f, const Numeric &t)
dplanck_df
Numeric planck(const Numeric &f, const Numeric &t)
planck
Numeric dplanck_dt(const Numeric &f, const Numeric &t)
dplanck_dt
This file contains declerations of functions of physical character.
constexpr Numeric dnumber_density_dt(Numeric p, Numeric t) noexcept
dnumber_density_dT
Definition: physics_funcs.h:84
constexpr Numeric number_density(Numeric p, Numeric t) noexcept
number_density
Definition: physics_funcs.h:70
Array< StokesVector > ArrayOfStokesVector
Array< PropagationMatrix > ArrayOfPropagationMatrix
bool bad_propmat(const Array< T > &main, const Vector &f_grid, const Index sd=1) noexcept
Checks if a Propagation Matrix or something similar has good grids.
Numeric dotprod_with_los(const ConstVectorView &los, const Numeric &u, const Numeric &v, const Numeric &w, const Index &atmosphere_dim)
Calculates the dot product between a field and a LOS.
Definition: rte.cc:671
void mirror_los(Vector &los_mirrored, const ConstVectorView &los, const Index &atmosphere_dim)
Determines the backward direction for a given line-of-sight.
Definition: rte.cc:1812
Declaration of functions in rte.cc.
Index find_first_species(const ArrayOfArrayOfSpeciesTag &specs, Species::Species spec) noexcept
void check_abs_species(const ArrayOfArrayOfSpeciesTag &abs_species)
AbsorptionPopulationTagTypeStatus population
SingleLine PopLine(Index) noexcept
Pops a single line.
Array< SingleLine > lines
A list of individual lines.
Index NumLines() const noexcept
Number of lines.
Species::IsotopeRecord Isotopologue() const noexcept
Isotopologue Index.
void ReverseLines() noexcept
Reverses the order of the internal lines.
void AppendSingleLine(SingleLine &&sl)
Appends a single line to the absorption lines.
Helper class to create an agenda.
Definition: agenda_set.h:144
void add(const std::string_view method, Input &&... input)
Add a method with as many inputs as you want. These inputs must be of type Wsv.
Definition: agenda_set.h:155
Agenda finalize()
Check the agenda and return a copy of it (ignoring/touching all agenda input/output not dealt with)
Definition: agenda_set.cc:178
Main computational data for the line shape and strength calculations.
Definition: lineshape.h:793
void interp_add_even(const ComputeData &sparse) ARTS_NOEXCEPT
Add a sparse grid to this grid via linear interpolation.
Definition: lineshape.cc:3677
ComplexVector N
Definition: lineshape.h:794
void interp_add_triplequad(const ComputeData &sparse) ARTS_NOEXCEPT
Add a sparse grid to this grid via square interpolation.
Definition: lineshape.cc:3722
ComplexVector F
Definition: lineshape.h:794
ComplexMatrix dF
Definition: lineshape.h:795
ComplexMatrix dN
Definition: lineshape.h:795
Struct to test of an ArrayOfArrayOfSpeciesTag contains a tagtype.
Definition: species_tags.h:173
#define a
#define b
This file contains basic functions to handle XML data files.