ARTS 2.5.4 (git: 4c0d3b4d)
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 "absorption.h"
34#include "array.h"
35#include "arts.h"
36#include "auto_md.h"
37#include "check_input.h"
38#include "depr.h"
39#include "file.h"
40#include "global_data.h"
41#include "hitran_species.h"
42#include "jacobian.h"
43#include "legacy_continua.h"
44#include "lineshape.h"
45#include "m_xml.h"
46#include "math_funcs.h"
47#include "matpackI.h"
48#include "messages.h"
49#include "montecarlo.h"
50#include "optproperties.h"
51#include "parameters.h"
52#include "physics_funcs.h"
53#include "rte.h"
54#include "xml_io.h"
55#include <algorithm>
56#include <cmath>
57
58#ifdef ENABLE_NETCDF
59#include <netcdf.h>
60#include "nc_io.h"
61#endif
62
63extern const Numeric ELECTRON_CHARGE;
64extern const Numeric ELECTRON_MASS;
65extern const Numeric PI;
66extern const Numeric SPEED_OF_LIGHT;
67extern const Numeric VACUUM_PERMITTIVITY;
68
69/* Workspace method: Doxygen documentation will be auto-generated */
70void AbsInputFromRteScalars( // WS Output:
71 Vector& abs_p,
72 Vector& abs_t,
73 Matrix& abs_vmrs,
74 // WS Input:
75 const Numeric& rtp_pressure,
76 const Numeric& rtp_temperature,
77 const Vector& rtp_vmr,
78 const Verbosity&) {
79 // Prepare abs_p:
80 abs_p.resize(1);
81 abs_p = rtp_pressure;
82
83 // Prepare abs_t:
84 abs_t.resize(1);
85 abs_t = rtp_temperature;
86
87 // Prepare abs_vmrs:
88 abs_vmrs.resize(rtp_vmr.nelem(), 1);
89 abs_vmrs = rtp_vmr;
90}
91
92/* Workspace method: Doxygen documentation will be auto-generated */
94 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
95 // WS Input:
96 const ArrayOfAbsorptionLines& abs_lines,
97 const ArrayOfArrayOfSpeciesTag& tgs,
98 const Verbosity&) {
99 // Size is set but inner size will now change from the original definition of species tags...
100 abs_lines_per_species.resize(tgs.nelem());
101
102 // The inner arrays need to be emptied, because they may contain lines
103 // from a previous calculation
104 for (auto &tg : abs_lines_per_species)
105 tg.resize(0);
106
107 // Take copies because we have to support frequency ranges, so might have to delete
108 for (AbsorptionLines lines: abs_lines) {
109
110 // Skip empty lines
111 if (lines.NumLines() == 0) continue;
112
113 // Loop all the tags
114 for (Index i=0; i<tgs.nelem() and lines.NumLines(); i++) {
115 for (auto& this_tag: tgs[i]) {
116 // Test isotopologue, we have to hit the end of the list for no isotopologue or the exact value
117 if (not same_or_joker(this_tag.Isotopologue(), lines.Isotopologue()))
118 continue;
119
120 // If there is a frequency range, we have to check so that only selected lines are included
121 if (this_tag.lower_freq >= 0 or this_tag.upper_freq >= 0) {
122 const Numeric low = (this_tag.lower_freq >= 0) ? this_tag.lower_freq : std::numeric_limits<Numeric>::lowest();
123 const Numeric upp = (this_tag.upper_freq >= 0) ? this_tag.upper_freq : std::numeric_limits<Numeric>::max();
124
125 // Fill up a copy of the line record to match with the wished frequency criteria
126 AbsorptionLines these_lines = lines;
127 these_lines.lines.resize(0);
128 for (Index k=lines.NumLines()-1; k>=0; k--)
129 if (low <= lines.lines[k].F0 and upp >= lines.lines[k].F0)
130 these_lines.AppendSingleLine(lines.PopLine(k));
131
132 // Append these lines after sorting them if there are any of them
133 if (these_lines.NumLines()) {
134 these_lines.ReverseLines();
135 abs_lines_per_species[i].push_back(these_lines);
136 }
137
138 // If this means we have deleted all lines, then we leave
139 if (lines.NumLines() == 0)
140 goto leave_inner_loop;
141 }
142 else {
143 abs_lines_per_species[i].push_back(lines);
144 goto leave_inner_loop;
145 }
146 }
147 }
148 leave_inner_loop: {}
149 }
150}
151
152/* Workspace method: Doxygen documentation will be auto-generated */
155 Index& propmat_clearsky_agenda_checked,
156 Index& abs_xsec_agenda_checked,
157 // Control Parameters:
158 const String& basename,
159 const Verbosity& verbosity) {
161
162 // Invalidate agenda check flags
163 propmat_clearsky_agenda_checked = false;
164 abs_xsec_agenda_checked = false;
165
166 // We want to make lists of included and excluded species:
167 ArrayOfString included(0), excluded(0);
168
169 tgs.resize(0);
170
171 for (Index i = 0; i < Index(Species::Species::FINAL); ++i) {
172 const String specname = Species::toShortName(Species::Species(i));
173
174 String filename = basename;
175 if (basename.length() && basename[basename.length() - 1] != '/')
176 filename += ".";
177 filename += specname;
178
179 try {
180 find_xml_file(filename, verbosity);
181 // Add to included list:
182 included.push_back(specname);
183
184 // Add this tag group to tgs:
185 tgs.emplace_back(ArrayOfSpeciesTag(specname));
186 } catch (const std::runtime_error& e) {
187 // The file for the species could not be found.
188 excluded.push_back(specname);
189 }
190 }
191
192 // Some nice output:
193 out2 << " Included Species (" << included.nelem() << "):\n";
194 for (Index i = 0; i < included.nelem(); ++i)
195 out2 << " " << included[i] << "\n";
196
197 out2 << " Excluded Species (" << excluded.nelem() << "):\n";
198 for (Index i = 0; i < excluded.nelem(); ++i)
199 out2 << " " << excluded[i] << "\n";
200}
201
202/* Workspace method: Doxygen documentation will be auto-generated */
203void abs_speciesDefineAll( // WS Output:
204 ArrayOfArrayOfSpeciesTag& abs_species,
205 Index& propmat_clearsky_agenda_checked,
206 Index& abs_xsec_agenda_checked,
207 // Control Parameters:
208 const Verbosity& verbosity) {
209 // Species lookup data:
210
211 // We want to make lists of all species
212 ArrayOfString specs(0);
213 for (Index i = 0; i < Index(Species::Species::FINAL); ++i) {
214 if (Species::Species(i) not_eq Species::Species::Bath) {
215 specs.emplace_back(Species::toShortName(Species::Species(i)));
216 }
217 }
218
219 // Set the values
220 abs_speciesSet(abs_species, abs_xsec_agenda_checked, propmat_clearsky_agenda_checked, specs, verbosity);
221}
222
223/* Workspace method: Doxygen documentation will be auto-generated */
224void AbsInputFromAtmFields( // WS Output:
225 Vector& abs_p,
226 Vector& abs_t,
227 Matrix& abs_vmrs,
228 // WS Input:
229 const Index& atmosphere_dim,
230 const Vector& p_grid,
231 const Tensor3& t_field,
232 const Tensor4& vmr_field,
233 const Verbosity&) {
234 // First, make sure that we really have a 1D atmosphere:
235 ARTS_USER_ERROR_IF (1 != atmosphere_dim,
236 "Atmospheric dimension must be 1D, but atmosphere_dim is ",
237 atmosphere_dim, ".")
238
239 abs_p = p_grid;
240 abs_t = t_field(joker, 0, 0);
241 abs_vmrs = vmr_field(joker, joker, 0, 0);
242}
243
244/* Workspace method: Doxygen documentation will be auto-generated */
245void abs_coefCalcFromXsec( // WS Output:
246 Matrix& abs_coef,
247 ArrayOfMatrix& dabs_coef_dx,
248 ArrayOfMatrix& abs_coef_per_species,
249 // WS Input:
250 const ArrayOfMatrix& abs_xsec_per_species,
251 const ArrayOfArrayOfMatrix& dabs_xsec_per_species_dx,
252 const ArrayOfArrayOfSpeciesTag& abs_species,
253 const ArrayOfRetrievalQuantity& jacobian_quantities,
254 const Matrix& abs_vmrs,
255 const Vector& abs_p,
256 const Vector& abs_t,
257 const Verbosity& verbosity) {
259
260 // Check that abs_vmrs and abs_xsec_per_species really have compatible
261 // dimensions. In abs_vmrs there should be one row for each tg:
262 ARTS_USER_ERROR_IF(abs_vmrs.nrows() != abs_xsec_per_species.nelem(),
263 "Variable abs_vmrs must have compatible dimension to abs_xsec_per_species.\n"
264 "abs_vmrs.nrows() = ", abs_vmrs.nrows(), "\n"
265 "abs_xsec_per_species.nelem() = ", abs_xsec_per_species.nelem())
266
267 // Check that number of altitudes are compatible. We only check the
268 // first element, this is possilble because within arts all elements
269 // are on the same altitude grid.
270 ARTS_USER_ERROR_IF(abs_vmrs.ncols() != abs_xsec_per_species[0].ncols(),
271 "Variable abs_vmrs must have same numbers of altitudes as abs_xsec_per_species.\n"
272 "abs_vmrs.ncols() = ", abs_vmrs.ncols(), "\n"
273 "abs_xsec_per_species[0].ncols() = ", abs_xsec_per_species[0].ncols());
274
275 // Check dimensions of abs_p and abs_t:
276 chk_size("abs_p", abs_p, abs_vmrs.ncols());
277 chk_size("abs_t", abs_t, abs_vmrs.ncols());
278
279 // Initialize abs_coef and abs_coef_per_species. The array dimension of abs_coef_per_species
280 // is the same as that of abs_xsec_per_species. The dimension of abs_coef should
281 // be equal to one of the abs_xsec_per_species enries.
282
283 abs_coef.resize(abs_xsec_per_species[0].nrows(),
284 abs_xsec_per_species[0].ncols());
285 abs_coef = 0;
286
287 dabs_coef_dx.resize(jacobian_quantities.nelem());
288
289 for (Index ii = 0; ii < jacobian_quantities.nelem(); ii++) {
290 const auto& deriv = jacobian_quantities[ii];
291
292 if (not deriv.propmattype()) continue;
293
294 dabs_coef_dx[ii].resize(abs_xsec_per_species[0].nrows(),
295 abs_xsec_per_species[0].ncols());
296 dabs_coef_dx[ii] = 0.0;
297 }
298
299 abs_coef_per_species.resize(abs_xsec_per_species.nelem());
300
301 out3
302 << " Computing abs_coef and abs_coef_per_species from abs_xsec_per_species.\n";
303 // Loop through all tag groups
304 for (Index i = 0; i < abs_xsec_per_species.nelem(); ++i) {
305 out3 << " Tag group " << i << "\n";
306
307 // Make this element of abs_xsec_per_species the right size:
308 abs_coef_per_species[i].resize(abs_xsec_per_species[i].nrows(),
309 abs_xsec_per_species[i].ncols());
310 abs_coef_per_species[i] = 0; // Initialize all elements to 0.
311
312 // Loop through all altitudes
313 for (Index j = 0; j < abs_xsec_per_species[i].ncols(); j++) {
314 // Calculate total number density from pressure and temperature.
315 const Numeric n = number_density(abs_p[j], abs_t[j]);
316 const Numeric dn_dT = dnumber_density_dt(abs_p[j], abs_t[j]);
317 // Wasted calculations when Jacobians are not calculated...
318 // Though this is called seldom enough that it this fine? value is -1/t*n
319
320 // Loop through all frequencies
321 for (Index k = 0; k < abs_xsec_per_species[i].nrows(); k++) {
322 abs_coef_per_species[i](k, j) =
323 abs_xsec_per_species[i](k, j) * n * abs_vmrs(i, j);
324
325 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
326 const auto& deriv = jacobian_quantities[iq];
327
328 if (not deriv.propmattype()) continue;
329
330 if (deriv == Jacobian::Atm::Temperature) {
331 dabs_coef_dx[iq](k, j) +=
332 (dabs_xsec_per_species_dx[i][iq](k, j) * n +
333 abs_xsec_per_species[i](k, j) * dn_dT) *
334 abs_vmrs(i, j);
335 } else if (deriv == Jacobian::Line::VMR) {
336 bool seco = false, main = false;
337 for (const auto& s : abs_species[i]) {
338 if (species_match(
339 deriv, s.cia_2nd_species) or
340 s.type not_eq Species::TagType::Cia)
341 seco = true;
343 deriv,
344 s.Isotopologue()))
345 main = true;
346 }
347 if (main and seco) {
348 dabs_coef_dx[iq](k, j) +=
349 (dabs_xsec_per_species_dx[i][iq](k, j) * abs_vmrs(i, j) +
350 abs_xsec_per_species[i](k, j)) *
351 n;
352 } else if (main) {
353 dabs_coef_dx[iq](k, j) += abs_xsec_per_species[i](k, j) * n;
354 } else if (seco) {
355 dabs_coef_dx[iq](k, j) +=
356 dabs_xsec_per_species_dx[i][iq](k, j) * abs_vmrs(i, j) * n;
357 }
358 } else if (deriv == Jacobian::Special::ArrayOfSpeciesTagVMR and deriv == abs_species[i]) {
359 dabs_coef_dx[iq](k, j) += abs_xsec_per_species[i](k, j) * n;
360 } else {
361 dabs_coef_dx[iq](k, j) +=
362 dabs_xsec_per_species_dx[i][iq](k, j) * n * abs_vmrs(i, j);
363 }
364 }
365 }
366 }
367
368 // Add up to the total absorption:
369 abs_coef += abs_coef_per_species[i]; // In Matpack you can use the +=
370 // operator to do elementwise addition.
371 }
372}
373
374/* Workspace method: Doxygen documentation will be auto-generated */
375void abs_xsec_per_speciesInit( // WS Output:
376 ArrayOfMatrix& abs_xsec_per_species,
377 ArrayOfArrayOfMatrix& dabs_xsec_per_species_dx,
378 // WS Input:
379 const ArrayOfArrayOfSpeciesTag& tgs,
380 const ArrayOfRetrievalQuantity& jacobian_quantities,
381 const ArrayOfIndex& abs_species_active,
382 const Vector& f_grid,
383 const Vector& abs_p,
384 const Index& abs_xsec_agenda_checked,
385 const Verbosity& verbosity) {
387
388 ARTS_USER_ERROR_IF (!abs_xsec_agenda_checked,
389 "You must call *abs_xsec_agenda_checkedCalc* before calling this method.");
390
391 // Sizes
392 const Index nq = jacobian_quantities.nelem();
393 const Index nf = f_grid.nelem();
394 const Index np = abs_p.nelem();
395 const Index ns = tgs.nelem();
396
397 // We need to check that abs_species_active doesn't have more elements than
398 // abs_species (abs_xsec_agenda_checkedCalc doesn't know abs_species_active.
399 // Usually we come here through an agenda call, where abs_species_active has
400 // been properly created somewhere internally. But we might get here by
401 // direct call, and then need to be safe!).
402 ARTS_USER_ERROR_IF (ns < abs_species_active.nelem(),
403 "abs_species_active (n=", abs_species_active.nelem(),
404 ") not allowed to have more elements than abs_species (n=",
405 ns, ")!\n")
406
407 // Make elements the right size if they are not already the right size
408 if (abs_xsec_per_species.nelem() not_eq ns) abs_xsec_per_species.resize(ns);
409 if (dabs_xsec_per_species_dx.nelem() not_eq ns) dabs_xsec_per_species_dx.resize(ns);
410
411 // Loop abs_xsec_per_species and make each matrix the right size,
412 // initializing to zero.
413 // But skip inactive species, loop only over the active ones.
414 for (auto& i: abs_species_active) {
416 "*abs_species_active* contains an invalid species index.\n"
417 "Species index must be between 0 and ", ns - 1)
418
419 // Make elements the right size if they are not already the right size, then reset them
420 if (abs_xsec_per_species[i].nrows() == nf and abs_xsec_per_species[i].ncols() == np) {
421 abs_xsec_per_species[i] = 0.0;
422 } else {
423 abs_xsec_per_species[i] = Matrix(nf, np, 0.0);
424 }
425
426 // Make elements the right size if they are not already the right size, then reset them
427 if (dabs_xsec_per_species_dx[i].nelem() not_eq nq) {
428 dabs_xsec_per_species_dx[i] = ArrayOfMatrix(nq, Matrix(nf, np, 0.0));
429 } else {
430 for (Index j=0; j<nq; j++) {
431 if (dabs_xsec_per_species_dx[i][j].nrows() == nf and dabs_xsec_per_species_dx[i][j].ncols() == np) {
432 dabs_xsec_per_species_dx[i][j] = 0.0;
433 } else {
434 dabs_xsec_per_species_dx[i][j] = Matrix(nf, np, 0.0);
435 }
436 }
437 }
438 }
439
440 ostringstream os;
441 os << " Initialized abs_xsec_per_species.\n"
442 << " Number of frequencies : " << nf << "\n"
443 << " Number of pressure levels : " << np << "\n";
444 out3 << os.str();
445}
446
447
449 const ArrayOfVector& abs_cont_parameters,
450 const ArrayOfString& abs_cont_models) {
451 std::ostringstream os;
452
453 for (Index i = 0; i < abs_cont_names.nelem(); ++i)
454 os << "abs_xsec_per_speciesAddConts: " << i
455 << " name : " << abs_cont_names[i] << "\n";
456
457 for (Index i = 0; i < abs_cont_parameters.nelem(); ++i)
458 os << "abs_xsec_per_speciesAddConts: " << i
459 << " param: " << abs_cont_parameters[i] << "\n";
460
461 for (Index i = 0; i < abs_cont_models.nelem(); ++i)
462 os << "abs_xsec_per_speciesAddConts: " << i
463 << " option: " << abs_cont_models[i] << "\n";
464
465 os << "The following variables must have the same dimension:\n"
466 << "abs_cont_names: " << abs_cont_names.nelem() << "\n"
467 << "abs_cont_parameters: " << abs_cont_parameters.nelem();
468
469 return os.str();
470}
471
472
473/* Workspace method: Doxygen documentation will be auto-generated */
475 ArrayOfMatrix& abs_xsec_per_species,
476 ArrayOfArrayOfMatrix& dabs_xsec_per_species_dx,
477 // WS Input:
478 const ArrayOfArrayOfSpeciesTag& tgs,
479 const ArrayOfRetrievalQuantity& jacobian_quantities,
480 const ArrayOfIndex& abs_species_active,
481 const Vector& f_grid,
482 const Vector& abs_p,
483 const Vector& abs_t,
484 const Matrix& abs_vmrs,
485 const ArrayOfString& abs_cont_names,
486 const ArrayOfVector& abs_cont_parameters,
487 const ArrayOfString& abs_cont_models,
488 const Verbosity& verbosity) {
490
491 // Needed for some continua, and set here from abs_vmrs:
492 Vector abs_h2o, abs_n2, abs_o2;
493
494 // Check that all paramters that should have the number of tag
495 // groups as a dimension are consistent.
496 ARTS_USER_ERROR_IF (tgs.nelem() != abs_xsec_per_species.nelem() || tgs.nelem() != abs_vmrs.nrows(),
497 "The following variables must all have the same dimension:\n"
498 "tgs: ", tgs.nelem(), "\n"
499 "abs_xsec_per_species: ", abs_xsec_per_species.nelem(), "\n"
500 "abs_vmrs.nrows(): ", abs_vmrs.nrows())
501
502 // Jacobian overhead START
503 /* NOTE: if any of the functions inside continuum tags could
504 be made to give partial derivatives, then that would
505 speed things up. Also be aware that line specific
506 parameters cannot be retrieved while using these
507 models. */
508 const bool do_jac = supports_continuum(
509 jacobian_quantities); // Throws runtime error if line parameters are wanted since we cannot know if the line is in the Continuum...
510 const bool do_freq_jac = do_frequency_jacobian(jacobian_quantities);
511 const bool do_temp_jac = do_temperature_jacobian(jacobian_quantities);
512 Vector dfreq, dabs_t;
513 const Numeric df = frequency_perturbation(jacobian_quantities);
514 const Numeric dt = temperature_perturbation(jacobian_quantities);
515
516 if (do_freq_jac) {
517 dfreq.resize(f_grid.nelem());
518 for (Index iv = 0; iv < f_grid.nelem(); iv++) dfreq[iv] = f_grid[iv] + df;
519 }
520 if (do_temp_jac) {
521 dabs_t.resize(abs_t.nelem());
522 for (Index it = 0; it < abs_t.nelem(); it++) dabs_t[it] = abs_t[it] + dt;
523 }
524
525 Matrix jacs_df, jacs_dt, normal;
526 if (do_jac) {
527 if (do_freq_jac) jacs_df.resize(f_grid.nelem(), abs_p.nelem());
528 if (do_temp_jac) jacs_dt.resize(f_grid.nelem(), abs_p.nelem());
529 normal.resize(f_grid.nelem(), abs_p.nelem());
530 }
531 // Jacobian overhead END
532
533
534
535 // Check, that dimensions of abs_cont_names and
536 // abs_cont_parameters are consistent...
537 ARTS_USER_ERROR_IF (abs_cont_names.nelem() != abs_cont_parameters.nelem(),
538 continua_model_error_message(abs_cont_names, abs_cont_parameters, abs_cont_models))
539
540 // Check that abs_p, abs_t, and abs_vmrs have the same
541 // dimension. This could be a user error, so we throw a
542 // runtime_error.
543
544 ARTS_USER_ERROR_IF (abs_t.nelem() != abs_p.nelem(),
545 "Variable abs_t must have the same dimension as abs_p.\n"
546 "abs_t.nelem() = ", abs_t.nelem(), '\n',
547 "abs_p.nelem() = ", abs_p.nelem())
548
549 ARTS_USER_ERROR_IF (abs_vmrs.ncols() != abs_p.nelem(),
550 "Variable dimension abs_vmrs.ncols() must\n"
551 "be the same as abs_p.nelem().\n"
552 "abs_vmrs.ncols() = ", abs_vmrs.ncols(), '\n',
553 "abs_p.nelem() = ", abs_p.nelem())
554
555 // We set abs_h2o, abs_n2, and abs_o2 later, because we only want to
556 // do it if the parameters are really needed.
557
558 out3 << " Calculating continuum spectra.\n";
559
560 // Loop tag groups:
561 for (Index ii = 0; ii < abs_species_active.nelem(); ++ii) {
562 const Index i = abs_species_active[ii];
563
564 // Go through the tags in the current tag group to see if they
565 // are continuum tags:
566 for (Index s = 0; s < tgs[i].nelem(); ++s) {
567 // Continuum tags in the sense that we talk about here
568 // (including complete absorption models) are marked by a special type.
569 if (tgs[i][s].Type() == Species::TagType::PredefinedLegacy) {
570 // We have identified a continuum tag!
571
572 // Get only the continuum name. The full tag name is something like:
573 // H2O-HITRAN96Self-*-*. We want only the `H2O-HITRAN96Self' part:
574 const String name = tgs[i][s].Isotopologue().FullName();
575
576 // Check, if we have parameters for this model. For
577 // this, the model name must be listed in
578 // abs_cont_names.
579 const Index n =
580 find(abs_cont_names.begin(), abs_cont_names.end(), name) -
581 abs_cont_names.begin();
582
583 // n==abs_cont_names.nelem() indicates that
584 // the name was not found.
585 ARTS_USER_ERROR_IF (n == abs_cont_names.nelem(),
586 "Cannot find model ", name, " in abs_cont_names.")
587
588 // Ok, the tag specifies a valid continuum model and
589 // we have continuum parameters.
590
591 if (out3.sufficient_priority()) {
592 ostringstream os;
593 os << " Adding " << name << " to tag group " << i << ".\n";
594 out3 << os.str();
595 }
596
597 // find the options for this continuum tag from the input array
598 // of options. The actual field of the array is n:
599 const String ContOption = abs_cont_models[n];
600
601 // Set abs_h2o, abs_n2, and abs_o2 from the first matching species.
602 set_vmr_from_first_species(abs_h2o, "H2O", tgs, abs_vmrs);
603 set_vmr_from_first_species(abs_n2, "N2", tgs, abs_vmrs);
604 set_vmr_from_first_species(abs_o2, "O2", tgs, abs_vmrs);
605
606 // Add the continuum for this tag. The parameters in
607 // this call should be clear. The vmr is in
608 // abs_vmrs(i,Range(joker)). The other vmr variables,
609 // abs_h2o, abs_n2, and abs_o2 contains the real vmr of H2O,
610 // N2, nad O2, which are needed as additional information for
611 // certain continua:
612 // abs_h2o for
613 // O2-PWR88, O2-PWR93, O2-PWR98,
614 // O2-MPM85, O2-MPM87, O2-MPM89, O2-MPM92, O2-MPM93,
615 // O2-TRE05,
616 // O2-SelfContStandardType, O2-SelfContMPM93, O2-SelfContPWR93,
617 // N2-SelfContMPM93, N2-DryContATM01,
618 // N2-CIArotCKDMT252, N2-CIAfunCKDMT252
619 // abs_n2 for
620 // H2O-SelfContCKD24, H2O-ForeignContCKD24,
621 // O2-v0v0CKDMT100,
622 // CO2-ForeignContPWR93, CO2-ForeignContHo66
623 // abs_o2 for
624 // N2-CIArotCKDMT252, N2-CIAfunCKDMT252
625 if (!do_jac)
626 xsec_continuum_tag(abs_xsec_per_species[i],
627 name,
628 abs_cont_parameters[n],
629 abs_cont_models[n],
630 f_grid,
631 abs_p,
632 abs_t,
633 abs_n2,
634 abs_h2o,
635 abs_o2,
636 abs_vmrs(i, Range(joker)),
637 verbosity);
638 else // The Jacobian block
639 {
640 // Needs a reseted block here...
641 for (Index iv = 0; iv < f_grid.nelem(); iv++) {
642 for (Index ip = 0; ip < abs_p.nelem(); ip++) {
643 if (do_freq_jac) jacs_df(iv, ip) = 0.0;
644 if (do_temp_jac) jacs_dt(iv, ip) = 0.0;
645 normal(iv, ip) = 0.0;
646 }
647 }
648
649 // Normal calculations
650 xsec_continuum_tag(normal,
651 name,
652 abs_cont_parameters[n],
653 abs_cont_models[n],
654 f_grid,
655 abs_p,
656 abs_t,
657 abs_n2,
658 abs_h2o,
659 abs_o2,
660 abs_vmrs(i, Range(joker)),
661 verbosity);
662
663 // Frequency calculations
664 if (do_freq_jac)
665 xsec_continuum_tag(jacs_df,
666 name,
667 abs_cont_parameters[n],
668 abs_cont_models[n],
669 dfreq,
670 abs_p,
671 abs_t,
672 abs_n2,
673 abs_h2o,
674 abs_o2,
675 abs_vmrs(i, Range(joker)),
676 verbosity);
677
678 //Temperature calculations
679 if (do_temp_jac)
680 xsec_continuum_tag(jacs_dt,
681 name,
682 abs_cont_parameters[n],
683 abs_cont_models[n],
684 f_grid,
685 abs_p,
686 dabs_t,
687 abs_n2,
688 abs_h2o,
689 abs_o2,
690 abs_vmrs(i, Range(joker)),
691 verbosity);
692 for (Index iv = 0; iv < f_grid.nelem(); iv++) {
693 for (Index ip = 0; ip < abs_p.nelem(); ip++) {
694 abs_xsec_per_species[i](iv, ip) += normal(iv, ip);
695 for (Index iq = 0; iq < jacobian_quantities.nelem();
696 iq++) {
697 const auto& deriv = jacobian_quantities[iq];
698
699 if (not deriv.propmattype()) continue;
700
701 if (is_frequency_parameter(deriv))
702 dabs_xsec_per_species_dx[i][iq](iv, ip) +=
703 (jacs_df(iv, ip) - normal(iv, ip)) * (1. / df);
704 else if (deriv == Jacobian::Atm::Temperature)
705 dabs_xsec_per_species_dx[i][iq](iv, ip) +=
706 (jacs_dt(iv, ip) - normal(iv, ip)) * (1. / dt);
707 }
708 }
709 }
710 }
711 // Calling this function with a row of Matrix abs_vmrs
712 // is possible because it uses Views.
713 }
714 }
715 }
716}
717
718//======================================================================
719// Methods related to continua
720//======================================================================
721
722/* Workspace method: Doxygen documentation will be auto-generated */
723void abs_cont_descriptionInit( // WS Output:
724 ArrayOfString& abs_cont_names,
725 ArrayOfString& abs_cont_options,
726 ArrayOfVector& abs_cont_parameters,
727 const Verbosity& verbosity) {
729
730 abs_cont_names.resize(0);
731 abs_cont_options.resize(0);
732 abs_cont_parameters.resize(0);
733 out2 << " Initialized abs_cont_names \n"
734 " abs_cont_models\n"
735 " abs_cont_parameters.\n";
736}
737
738/* Workspace method: Doxygen documentation will be auto-generated */
739void abs_cont_descriptionAppend( // WS Output:
740 ArrayOfString& abs_cont_names,
741 ArrayOfString& abs_cont_models,
742 ArrayOfVector& abs_cont_parameters,
743 // Control Parameters:
744 const String& tagname,
745 const String& model,
746 const Vector& userparameters,
747 const Verbosity&) {
748 // First we have to check that name matches a continuum species tag.
749 check_continuum_model(tagname);
750
751 //cout << " + tagname: " << tagname << "\n";
752 //cout << " + model: " << model << "\n";
753 //cout << " + parameters: " << userparameters << "\n";
754
755 // Add name and parameters to the apropriate variables:
756 abs_cont_names.push_back(tagname);
757 abs_cont_models.push_back(model);
758 abs_cont_parameters.push_back(userparameters);
759}
760
761/* Workspace method: Doxygen documentation will be auto-generated */
763 StokesVector& nlte_source,
764 ArrayOfStokesVector& dnlte_source_dx,
765 // WS Input:
766 const ArrayOfMatrix& src_coef_per_species,
767 const ArrayOfMatrix& dsrc_coef_dx,
768 const ArrayOfRetrievalQuantity& jacobian_quantities,
769 const Vector& f_grid,
770 const Numeric& rtp_temperature,
771 const Verbosity&) {
772 // nlte_source has format
773 // [ abs_species, f_grid, stokes_dim ].
774 // src_coef_per_species has format ArrayOfMatrix (over species),
775 // where for each species the matrix has format [f_grid, abs_p].
776
777 Index n_species = src_coef_per_species.nelem(); // # species
778
779 ARTS_USER_ERROR_IF (not n_species,
780 "Must have at least one species.")
781
782 Index n_f = src_coef_per_species[0].nrows(); // # frequencies
783
784 // # pressures must be 1:
785 ARTS_USER_ERROR_IF (1 not_eq src_coef_per_species[0].ncols(),
786 "Must have exactly one pressure.")
787
788 // Check frequency dimension of propmat_clearsky
789 ARTS_USER_ERROR_IF (nlte_source.NumberOfFrequencies() not_eq n_f,
790 "Frequency dimension of nlte_source does not\n"
791 "match abs_coef_per_species.")
792
793 const Vector B = planck(f_grid, rtp_temperature);
794
795 StokesVector sv(n_f, nlte_source.StokesDimensions());
796 for (Index si = 0; si < n_species; ++si) {
797 sv.Kjj() = src_coef_per_species[si](joker, 0);
798 sv *= B;
799 nlte_source.Kjj() += sv.Kjj();
800 }
801
802 // Jacobian
803 for (Index ii = 0; ii < jacobian_quantities.nelem(); ii++) {
804 const auto& deriv = jacobian_quantities[ii];
805
806 if (not deriv.propmattype()) continue;
807
808 if (deriv == Jacobian::Atm::Temperature) {
809 const Vector dB = dplanck_dt(f_grid, rtp_temperature);
810
811 for (Index si = 0; si < n_species; ++si) {
812 sv.Kjj() = src_coef_per_species[si](joker, 0);
813 sv *= dB;
814 dnlte_source_dx[ii].Kjj() += sv.Kjj();
815 }
816
817 sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
818 sv *= B;
819 dnlte_source_dx[ii].Kjj() += sv.Kjj();
820 } else if (is_frequency_parameter(deriv)) {
821 const Vector dB = dplanck_df(f_grid, rtp_temperature);
822
823 for (Index si = 0; si < n_species; ++si) {
824 sv.Kjj() = src_coef_per_species[si](joker, 0);
825 sv *= dB;
826 dnlte_source_dx[ii].Kjj() += sv.Kjj();
827 }
828
829 sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
830 sv *= B;
831 dnlte_source_dx[ii].Kjj() += sv.Kjj();
832 } else {
833 sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
834 sv *= B;
835 dnlte_source_dx[ii].Kjj() += sv.Kjj();
836 }
837 }
838}
839
840/* Workspace method: Doxygen documentation will be auto-generated */
842 PropagationMatrix& propmat_clearsky,
843 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
844 // WS Input:
845 const ArrayOfMatrix& abs_coef_per_species,
846 const ArrayOfMatrix& dabs_coef_dx) {
847 // propmat_clearsky has format
848 // [ abs_species, f_grid, stokes_dim, stokes_dim ].
849 // abs_coef_per_species has format ArrayOfMatrix (over species),
850 // where for each species the matrix has format [f_grid, abs_p].
851
852 Index n_species = abs_coef_per_species.nelem(); // # species
853 Index n_jacs = dabs_coef_dx.nelem(); // # derivatives
854
855 ARTS_USER_ERROR_IF (0 == n_species,
856 "Must have at least one species.")
857
858 Index n_f = abs_coef_per_species[0].nrows(); // # frequencies
859
860 // # pressures must be 1:
861 ARTS_USER_ERROR_IF (1 not_eq abs_coef_per_species[0].ncols(),
862 "Must have exactly one pressure.")
863
864 // Check frequency dimension of propmat_clearsky
865 ARTS_USER_ERROR_IF (propmat_clearsky.NumberOfFrequencies() not_eq n_f,
866 "Frequency dimension of propmat_clearsky does not\n"
867 "match abs_coef_per_species.")
868
869 ARTS_USER_ERROR_IF (dpropmat_clearsky_dx.nelem() not_eq n_jacs,
870 "Must have the same dimension.")
871
872 // Loop species and stokes dimensions, and add to propmat_clearsky:
873 for (Index si = 0; si < n_species; ++si)
874 propmat_clearsky.Kjj() += abs_coef_per_species[si](joker, 0);
875
876 for (Index iqn = 0; iqn < n_jacs; iqn++) {
877 if (dabs_coef_dx[iqn].nrows() == n_f) {
878 ARTS_USER_ERROR_IF(dabs_coef_dx[iqn].ncols() not_eq 1, "Must have exactly one pressure.")
879 dpropmat_clearsky_dx[iqn].Kjj() += dabs_coef_dx[iqn](joker, 0);
880 }
881 }
882}
883
884/* Workspace method: Doxygen documentation will be auto-generated */
885void propmat_clearskyInit( //WS Output
886 PropagationMatrix& propmat_clearsky,
887 StokesVector& nlte_source,
888 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
889 ArrayOfStokesVector& dnlte_source_dx,
890 //WS Input
891 const ArrayOfRetrievalQuantity& jacobian_quantities,
892 const Vector& f_grid,
893 const Index& stokes_dim,
894 const Index& propmat_clearsky_agenda_checked,
895 const Verbosity&) {
896 const Index nf = f_grid.nelem();
897 const Index nq = jacobian_quantities.nelem();
898
899 ARTS_USER_ERROR_IF (!propmat_clearsky_agenda_checked,
900 "You must call *propmat_clearsky_agenda_checkedCalc* before calling this method.")
901
902 ARTS_USER_ERROR_IF (not nf, "No frequencies");
903
904 ARTS_USER_ERROR_IF (stokes_dim < 1 or stokes_dim > 4, "stokes_dim not in [1, 2, 3, 4]");
905
906 // Set size of propmat_clearsky or reset it's values
907 if (propmat_clearsky.StokesDimensions() == stokes_dim and propmat_clearsky.NumberOfFrequencies() == nf) {
908 propmat_clearsky.SetZero();
909 } else {
910 propmat_clearsky = PropagationMatrix(nf, stokes_dim);
911 }
912
913 // Set size of dpropmat_clearsky_dx or reset it's values
914 if (dpropmat_clearsky_dx.nelem() not_eq nq) {
915 dpropmat_clearsky_dx = ArrayOfPropagationMatrix(nq, PropagationMatrix(nf, stokes_dim));
916 } else {
917 for (auto& pm: dpropmat_clearsky_dx) {
918 if (pm.StokesDimensions() == stokes_dim and pm.NumberOfFrequencies() == nf) {
919 pm.SetZero();
920 } else {
921 pm = PropagationMatrix(nf, stokes_dim);
922 }
923 }
924 }
925
926 // Set size of nlte_source or reset it's values
927 if (nlte_source.StokesDimensions() == stokes_dim and nlte_source.NumberOfFrequencies() == nf) {
928 nlte_source.SetZero();
929 } else {
930 nlte_source = StokesVector(nf, stokes_dim);
931 }
932
933 // Set size of dnlte_source_dx or reset it's values
934 if (dnlte_source_dx.nelem() not_eq nq) {
935 dnlte_source_dx = ArrayOfStokesVector(nq, StokesVector(nf, stokes_dim));
936 } else {
937 for (auto& pm: dnlte_source_dx) {
938 if (pm.StokesDimensions() == stokes_dim and pm.NumberOfFrequencies() == nf) {
939 pm.SetZero();
940 } else {
941 pm = StokesVector(nf, stokes_dim);
942 }
943 }
944 }
945}
946
947/* Workspace method: Doxygen documentation will be auto-generated */
949 PropagationMatrix& propmat_clearsky,
950 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
951 const Index& stokes_dim,
952 const Index& atmosphere_dim,
953 const Vector& f_grid,
954 const ArrayOfArrayOfSpeciesTag& abs_species,
955 const ArrayOfRetrievalQuantity& jacobian_quantities,
956 const Vector& rtp_vmr,
957 const Vector& rtp_los,
958 const Vector& rtp_mag,
959 const Verbosity&) {
960 // All the physical constants joined into one static constant:
961 // (abs as e defined as negative)
962 static const Numeric FRconst =
966
967 ARTS_USER_ERROR_IF (stokes_dim < 3,
968 "To include Faraday rotation, stokes_dim >= 3 is required.")
969 ARTS_USER_ERROR_IF (atmosphere_dim == 1 && rtp_los.nelem() < 1,
970 "For applying propmat_clearskyAddFaraday, los needs to be specified\n"
971 "(at least zenith angle component for atmosphere_dim==1),\n"
972 "but it is not.\n")
973 ARTS_USER_ERROR_IF (atmosphere_dim > 1 && rtp_los.nelem() < 2,
974 "For applying propmat_clearskyAddFaraday, los needs to be specified\n"
975 "(both zenith and azimuth angle components for atmosphere_dim>1),\n"
976 "but it is not.\n")
977
978 const bool do_magn_jac = do_magnetic_jacobian(jacobian_quantities);
979 const Numeric dmag = magnetic_field_perturbation(jacobian_quantities);
980
981 Index ife = -1;
982 for (Index sp = 0; sp < abs_species.nelem() && ife < 0; sp++) {
983 if (abs_species[sp].FreeElectrons()) {
984 ife = sp;
985 }
986 }
987
988 ARTS_USER_ERROR_IF (ife < 0,
989 "Free electrons not found in *abs_species* and "
990 "Faraday rotation can not be calculated.");
991
992 const Numeric ne = rtp_vmr[ife];
993
994 if (ne != 0 && (rtp_mag[0] != 0 || rtp_mag[1] != 0 || rtp_mag[2] != 0)) {
995 // Include remaining terms, beside /f^2
996 const Numeric c1 =
997 2 * FRconst *
999 rtp_los, rtp_mag[0], rtp_mag[1], rtp_mag[2], atmosphere_dim);
1000
1001 Numeric dc1_u = 0.0, dc1_v = 0.0, dc1_w = 0.0;
1002 if (do_magn_jac) {
1003 dc1_u = (2 * FRconst *
1004 dotprod_with_los(rtp_los,
1005 rtp_mag[0] + dmag,
1006 rtp_mag[1],
1007 rtp_mag[2],
1008 atmosphere_dim) -
1009 c1) /
1010 dmag;
1011 dc1_v = (2 * FRconst *
1012 dotprod_with_los(rtp_los,
1013 rtp_mag[0],
1014 rtp_mag[1] + dmag,
1015 rtp_mag[2],
1016 atmosphere_dim) -
1017 c1) /
1018 dmag;
1019 dc1_w = (2 * FRconst *
1020 dotprod_with_los(rtp_los,
1021 rtp_mag[0],
1022 rtp_mag[1],
1023 rtp_mag[2] + dmag,
1024 atmosphere_dim) -
1025 c1) /
1026 dmag;
1027 }
1028
1029 for (Index iv = 0; iv < f_grid.nelem(); iv++) {
1030 const Numeric f2 = f_grid[iv] * f_grid[iv];
1031 const Numeric r = ne * c1 / f2;
1032 propmat_clearsky.AddFaraday(r, iv);
1033
1034 // The Jacobian loop
1035 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
1036 if (is_frequency_parameter(jacobian_quantities[iq]))
1037 dpropmat_clearsky_dx[iq].AddFaraday(-2.0 * ne * r / f_grid[iv], iv);
1038 else if (jacobian_quantities[iq] == Jacobian::Atm::MagneticU)
1039 dpropmat_clearsky_dx[iq].AddFaraday(ne * dc1_u / f2, iv);
1040 else if (jacobian_quantities[iq] == Jacobian::Atm::MagneticV)
1041 dpropmat_clearsky_dx[iq].AddFaraday(ne * dc1_v / f2, iv);
1042 else if (jacobian_quantities[iq] == Jacobian::Atm::MagneticW)
1043 dpropmat_clearsky_dx[iq].AddFaraday(ne * dc1_w / f2, iv);
1044 else if (jacobian_quantities[iq] == Jacobian::Atm::Electrons)
1045 dpropmat_clearsky_dx[iq].AddFaraday(r, iv);
1046 else if (jacobian_quantities[iq] == abs_species[ife])
1047 dpropmat_clearsky_dx[iq].AddFaraday(r, iv);
1048 }
1049 }
1050 }
1051}
1052
1053/* Workspace method: Doxygen documentation will be auto-generated */
1055 // WS Output:
1056 PropagationMatrix& propmat_clearsky,
1057 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
1058 // WS Input:
1059 const Index& stokes_dim,
1060 const Index& atmosphere_dim,
1061 const Vector& f_grid,
1062 const ArrayOfArrayOfSpeciesTag& abs_species,
1063 const ArrayOfRetrievalQuantity& jacobian_quantities,
1064 const Vector& rtp_vmr,
1065 const Vector& rtp_los,
1066 const Numeric& rtp_temperature,
1067 const ArrayOfArrayOfSingleScatteringData& scat_data,
1068 const Index& scat_data_checked,
1069 const Index& use_abs_as_ext,
1070 // Verbosity object:
1071 const Verbosity& verbosity) {
1073
1074 // (i)yCalc only checks scat_data_checked if cloudbox is on. It is off here,
1075 // though, i.e. we need to check it here explicitly. (Also, cloudboxOff sets
1076 // scat_data_checked=0 as it does not check it and as we ususally don't need
1077 // scat_data for clearsky cases, hence don't want to check them by
1078 // scat_data_checkedCalc in that case. This approach seems to be the more
1079 // handy compared to cloudboxOff setting scat_data_checked=1 without checking
1080 // it assuming we won't use it anyways.)
1081 ARTS_USER_ERROR_IF (scat_data_checked != 1,
1082 "The scat_data must be flagged to have "
1083 "passed a consistency check (scat_data_checked=1).")
1084
1085 const Index ns = TotalNumberOfElements(scat_data);
1086 Index np = 0;
1087 for (Index sp = 0; sp < abs_species.nelem(); sp++) {
1088 if (abs_species[sp].Particles()) {
1089 np++;
1090 }
1091 }
1092
1093 ARTS_USER_ERROR_IF (np == 0,
1094 "For applying propmat_clearskyAddParticles, *abs_species* needs to"
1095 "contain species 'particles', but it does not.\n")
1096
1097 ARTS_USER_ERROR_IF (ns != np,
1098 "Number of 'particles' entries in abs_species and of elements in\n"
1099 "*scat_data* needs to be identical. But you have " , np,
1100 " 'particles' entries\n"
1101 "and ", ns, " *scat_data* elements.\n")
1102
1103 ARTS_USER_ERROR_IF (atmosphere_dim == 1 && rtp_los.nelem() < 1,
1104 "For applying *propmat_clearskyAddParticles*, *rtp_los* needs to be specified\n"
1105 "(at least zenith angle component for atmosphere_dim==1),\n"
1106 "but it is not.\n")
1107 ARTS_USER_ERROR_IF (atmosphere_dim > 1 && rtp_los.nelem() < 2,
1108 "For applying *propmat_clearskyAddParticles*, *rtp_los* needs to be specified\n"
1109 "(both zenith and azimuth angle components for atmosphere_dim>1),\n"
1110 "but it is not.\n")
1111
1112 // Use for rescaling vmr of particulates
1113 Numeric rtp_vmr_sum = 0.0;
1114
1115 // Tests and setup partial derivatives
1116 const bool do_jac_temperature = do_temperature_jacobian(jacobian_quantities);
1117 const bool do_jac_frequencies = do_frequency_jacobian(jacobian_quantities);
1118 const Numeric dT = temperature_perturbation(jacobian_quantities);
1119
1120 const Index na = abs_species.nelem();
1121 Vector rtp_los_back;
1122 mirror_los(rtp_los_back, rtp_los, atmosphere_dim);
1123
1124 // 170918 JM: along with transition to use of new-type (aka
1125 // pre-f_grid-interpolated) scat_data, freq perturbation switched off. Typical
1126 // clear-sky freq perturbations yield insignificant effects in particle
1127 // properties. Hence, this feature is neglected here.
1128 if (do_jac_frequencies) {
1129 out1 << "WARNING:\n"
1130 << "Frequency perturbation not available for absorbing particles.\n";
1131 }
1132
1133 // creating temporary output containers
1134 ArrayOfArrayOfTensor5 ext_mat_Nse;
1135 ArrayOfArrayOfTensor4 abs_vec_Nse;
1136 ArrayOfArrayOfIndex ptypes_Nse;
1137 Matrix t_ok;
1138
1139 // preparing input in format needed
1140 Vector T_array;
1141 if (do_jac_temperature) {
1142 T_array.resize(2);
1143 T_array = rtp_temperature;
1144 T_array[1] += dT;
1145 } else {
1146 T_array.resize(1);
1147 T_array = rtp_temperature;
1148 }
1149 Matrix dir_array(1, 2);
1150 dir_array(0, joker) = rtp_los_back;
1151
1152 // ext/abs per scat element for all freqs at once
1153 opt_prop_NScatElems(ext_mat_Nse,
1154 abs_vec_Nse,
1155 ptypes_Nse,
1156 t_ok,
1157 scat_data,
1158 stokes_dim,
1159 T_array,
1160 dir_array,
1161 -1);
1162
1163 const Index nf = abs_vec_Nse[0][0].nbooks();
1164 Tensor3 tmp(nf, stokes_dim, stokes_dim);
1165
1166 // Internal computations necessary since it relies on zero start
1167 PropagationMatrix internal_propmat(propmat_clearsky.NumberOfFrequencies(), propmat_clearsky.StokesDimensions());
1168
1169 // loop over the scat_data and link them with correct vmr_field entry according
1170 // to the position of the particle type entries in abs_species.
1171 Index sp = 0;
1172 Index i_se_flat = 0;
1173 for (Index i_ss = 0; i_ss < scat_data.nelem(); i_ss++) {
1174 for (Index i_se = 0; i_se < scat_data[i_ss].nelem(); i_se++) {
1175 // forward to next particle entry in abs_species
1176 while (sp < na && not abs_species[sp].Particles())
1177 sp++;
1178 internal_propmat.SetZero();
1179
1180 // running beyond number of abs_species entries when looking for next
1181 // particle entry. shouldn't happen, though.
1182 ARTS_ASSERT(sp < na);
1183 ARTS_USER_ERROR_IF (rtp_vmr[sp] < 0.,
1184 "Negative absorbing particle 'vmr' (aka number density)"
1185 " encountered:\n"
1186 "scat species #", i_ss, ", scat elem #", i_se,
1187 " (vmr_field entry #", sp, ")\n")
1188
1189 if (rtp_vmr[sp] > 0.) {
1190 ARTS_USER_ERROR_IF (t_ok(i_se_flat, 0) < 0.,
1191 "Temperature interpolation error:\n"
1192 "scat species #", i_ss, ", scat elem #", i_se, "\n")
1193 if (use_abs_as_ext) {
1194 if (nf > 1)
1195 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1196 internal_propmat.AddAbsorptionVectorAtPosition(
1197 abs_vec_Nse[i_ss][i_se](iv, 0, 0, joker), iv);
1198 else
1199 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1200 internal_propmat.AddAbsorptionVectorAtPosition(
1201 abs_vec_Nse[i_ss][i_se](0, 0, 0, joker), iv);
1202 } else {
1203 if (nf > 1)
1204 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1205 internal_propmat.SetAtPosition(
1206 ext_mat_Nse[i_ss][i_se](iv, 0, 0, joker, joker), iv);
1207 else
1208 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1209 internal_propmat.SetAtPosition(
1210 ext_mat_Nse[i_ss][i_se](0, 0, 0, joker, joker), iv);
1211 }
1212 propmat_clearsky += rtp_vmr[sp] * internal_propmat;
1213 }
1214
1215 // For temperature derivatives (so we don't need to check it in jac loop)
1216 if (do_jac_temperature) {
1217 ARTS_USER_ERROR_IF (t_ok(i_se_flat, 1) < 0.,
1218 "Temperature interpolation error (in perturbation):\n"
1219 "scat species #", i_ss, ", scat elem #", i_se, "\n")
1220 }
1221
1222 // For number density derivatives
1223 if (jacobian_quantities.nelem()) rtp_vmr_sum += rtp_vmr[sp];
1224
1225 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
1226 const auto& deriv = jacobian_quantities[iq];
1227
1228 if (not deriv.propmattype()) continue;
1229
1230 if (deriv == Jacobian::Atm::Temperature) {
1231 if (use_abs_as_ext) {
1232 tmp(joker, joker, 0) =
1233 abs_vec_Nse[i_ss][i_se](joker, 1, 0, joker);
1234 tmp(joker, joker, 0) -=
1235 abs_vec_Nse[i_ss][i_se](joker, 0, 0, joker);
1236 } else {
1237 tmp = ext_mat_Nse[i_ss][i_se](joker, 1, 0, joker, joker);
1238 tmp -= ext_mat_Nse[i_ss][i_se](joker, 0, 0, joker, joker);
1239 }
1240
1241 tmp *= rtp_vmr[sp];
1242 tmp /= dT;
1243
1244 if (nf > 1)
1245 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1246 if (use_abs_as_ext)
1247 dpropmat_clearsky_dx[iq].AddAbsorptionVectorAtPosition(
1248 tmp(iv, joker, 0), iv);
1249 else
1250 dpropmat_clearsky_dx[iq].AddAtPosition(tmp(iv, joker, joker),
1251 iv);
1252 else
1253 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1254 if (use_abs_as_ext)
1255 dpropmat_clearsky_dx[iq].AddAbsorptionVectorAtPosition(
1256 tmp(0, joker, 0), iv);
1257 else
1258 dpropmat_clearsky_dx[iq].AddAtPosition(tmp(0, joker, joker),
1259 iv);
1260 }
1261
1262 else if (deriv == Jacobian::Atm::Particulates) {
1263 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1264 dpropmat_clearsky_dx[iq].AddAtPosition(internal_propmat, iv);
1265 }
1266
1267 else if (deriv == abs_species[sp]) {
1268 dpropmat_clearsky_dx[iq] += internal_propmat;
1269 }
1270 }
1271 sp++;
1272 i_se_flat++;
1273 }
1274 }
1275
1276 //checking that no further 'particle' entry left after all scat_data entries
1277 //are processes. this is basically not necessary. but checking it anyway to
1278 //really be safe. remove later, when more extensively tested.
1279 while (sp < na) {
1280 ARTS_ASSERT(abs_species[sp][0].Type() != Species::TagType::Particles);
1281 sp++;
1282 }
1283
1284 if (rtp_vmr_sum != 0.0) {
1285 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
1286 const auto& deriv = jacobian_quantities[iq];
1287
1288 if (not deriv.propmattype()) continue;
1289
1290 if (deriv == Jacobian::Atm::Particulates) {
1291 dpropmat_clearsky_dx[iq] /= rtp_vmr_sum;
1292 }
1293 }
1294 }
1295}
1296
1297
1299 const Vector& f_grid,
1300 const Numeric& sparse_df,
1301 const String& speedup_option,
1302 // Verbosity object:
1303 const Verbosity&)
1304{
1305 // Return empty for nothing
1306 if (not f_grid.nelem()) {
1307 sparse_f_grid.resize(0);
1308 return;
1309 };
1310
1311 switch (Options::toLblSpeedupOrThrow(speedup_option)) {
1312 case Options::LblSpeedup::LinearIndependent:
1313 sparse_f_grid = LineShape::linear_sparse_f_grid(f_grid, sparse_df);
1315 break;
1316 case Options::LblSpeedup::QuadraticIndependent:
1317 sparse_f_grid = LineShape::triple_sparse_f_grid(f_grid, sparse_df);
1318 break;
1319 case Options::LblSpeedup::None:
1320 sparse_f_grid.resize(0);
1321 break;
1322 case Options::LblSpeedup::FINAL: { /* Leave last */ }
1323 }
1324}
1325
1327 const Numeric& sparse_df,
1328 const String& speedup_option,
1329 // Verbosity object:
1330 const Verbosity& verbosity)
1331{
1332 Vector sparse_f_grid;
1333 sparse_f_gridFromFrequencyGrid(sparse_f_grid, f_grid, sparse_df, speedup_option, verbosity);
1334 return sparse_f_grid;
1335}
1336
1337/* Workspace method: Doxygen documentation will be auto-generated */
1338void propmat_clearskyAddLines( // Workspace reference:
1339 // WS Output:
1340 PropagationMatrix& propmat_clearsky,
1341 StokesVector& nlte_source,
1342 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
1343 ArrayOfStokesVector& dnlte_source_dx,
1344 // WS Input:
1345 const Vector& f_grid,
1346 const ArrayOfArrayOfSpeciesTag& abs_species,
1347 const ArrayOfRetrievalQuantity& jacobian_quantities,
1348 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1350 const Numeric& rtp_pressure,
1351 const Numeric& rtp_temperature,
1352 const EnergyLevelMap& rtp_nlte,
1353 const Vector& rtp_vmr,
1354 const Index& nlte_do,
1355 const Index& lbl_checked,
1356 // WS User Generic inputs
1357 const Numeric& sparse_df,
1358 const Numeric& sparse_lim,
1359 const String& speedup_option,
1360 const ArrayOfSpeciesTag& select_speciestags,
1361 const Index& robust,
1362 // Verbosity object:
1363 const Verbosity& verbosity) {
1364
1365 // Size of problem
1366 const Index nf = f_grid.nelem();
1367 const Index nq = jacobian_quantities.nelem();
1368 const Index ns = abs_species.nelem();
1369
1370 // Possible things that can go wrong in this code (excluding line parameters)
1371 ARTS_USER_ERROR_IF(not lbl_checked, "Must check LBL calculations")
1372 check_abs_species(abs_species);
1373 ARTS_USER_ERROR_IF(rtp_vmr.nelem() not_eq abs_species.nelem(),
1374 "*rtp_vmr* must match *abs_species*")
1375 ARTS_USER_ERROR_IF(propmat_clearsky.NumberOfFrequencies() not_eq nf,
1376 "*f_grid* must match *propmat_clearsky*")
1377 ARTS_USER_ERROR_IF(nlte_source.NumberOfFrequencies() not_eq nf,
1378 "*f_grid* must match *nlte_source*")
1379 ARTS_USER_ERROR_IF(not nq and (nq not_eq dpropmat_clearsky_dx.nelem()),
1380 "*dpropmat_clearsky_dx* must match derived form of *jacobian_quantities*")
1381 ARTS_USER_ERROR_IF(not nq and bad_propmat(dpropmat_clearsky_dx, f_grid),
1382 "*dpropmat_clearsky_dx* must have frequency dim same as *f_grid*")
1383 ARTS_USER_ERROR_IF(nlte_do and (nq not_eq dnlte_source_dx.nelem()),
1384 "*dnlte_source_dx* must match derived form of *jacobian_quantities* when non-LTE is on")
1385 ARTS_USER_ERROR_IF(nlte_do and bad_propmat(dnlte_source_dx, f_grid),
1386 "*dnlte_source_dx* must have frequency dim same as *f_grid* when non-LTE is on")
1387 ARTS_USER_ERROR_IF(any_negative(f_grid), "Negative frequency (at least one value).")
1388 ARTS_USER_ERROR_IF(not is_increasing(f_grid), "Must be sorted and increasing.")
1389 ARTS_USER_ERROR_IF(any_negative(rtp_vmr), "Negative VMR (at least one value).")
1390 ARTS_USER_ERROR_IF(any_negative(rtp_nlte.value), "Negative NLTE (at least one value).")
1391 ARTS_USER_ERROR_IF(rtp_temperature <= 0, "Non-positive temperature")
1392 ARTS_USER_ERROR_IF(rtp_pressure <= 0, "Non-positive pressure")
1393 ARTS_USER_ERROR_IF(sparse_lim > 0 and sparse_df > sparse_lim,
1394 "If sparse grids are to be used, the limit must be larger than the grid-spacing.\n"
1395 "The limit is ", sparse_lim, " Hz and the grid_spacing is ", sparse_df, " Hz")
1396
1397 if (not nf) return;
1398
1399 // Deal with sparse computational grid
1400 const Vector f_grid_sparse = create_sparse_f_grid_internal(f_grid, sparse_df, speedup_option, verbosity);
1401 const Options::LblSpeedup speedup_type = f_grid_sparse.nelem() ? Options::toLblSpeedupOrThrow(speedup_option) : Options::LblSpeedup::None;
1402 ARTS_USER_ERROR_IF(sparse_lim <= 0 and speedup_type not_eq Options::LblSpeedup::None,
1403 "Must have a sparse limit if you set speedup_option")
1404
1405 // Calculations data
1406 LineShape::ComputeData com(f_grid, jacobian_quantities, nlte_do);
1407 LineShape::ComputeData sparse_com(f_grid_sparse, jacobian_quantities, nlte_do);
1408
1409 for (Index ispecies = 0; ispecies < ns; ispecies++) {
1410 if (select_speciestags.nelem() and select_speciestags not_eq abs_species[ispecies]) continue;
1411
1412 // Skip it if there are no species or there is Zeeman requested
1413 if (not abs_species[ispecies].nelem() or abs_species[ispecies].Zeeman() or not abs_lines_per_species[ispecies].nelem())
1414 continue;
1415
1416 for (auto& band : abs_lines_per_species[ispecies]) {
1417 LineShape::compute(com, sparse_com, band, jacobian_quantities, rtp_nlte, band.BroadeningSpeciesVMR(rtp_vmr, abs_species), abs_species[ispecies], rtp_vmr[ispecies],
1418 isotopologue_ratios[band.Isotopologue()], rtp_pressure, rtp_temperature, 0, sparse_lim,
1419 Zeeman::Polarization::None, speedup_type, robust not_eq 0);
1420
1421 }
1422 }
1423
1424 switch (speedup_type) {
1425 case Options::LblSpeedup::LinearIndependent: com.interp_add_even(sparse_com); break;
1426 case Options::LblSpeedup::QuadraticIndependent: com.interp_add_triplequad(sparse_com); break;
1427 case Options::LblSpeedup::None: /* Do nothing */ break;
1428 case Options::LblSpeedup::FINAL: { /* Leave last */ }
1429 }
1430
1431 // Sum up the propagation matrix
1432 propmat_clearsky.Kjj() += com.F.real();
1433
1434 // Sum up the Jacobian
1435 for (Index j=0; j<nq; j++) {
1436 if (not jacobian_quantities[j].propmattype()) continue;
1437 dpropmat_clearsky_dx[j].Kjj() += com.dF.real()(joker, j);
1438 }
1439
1440 if (nlte_do) {
1441 // Sum up the source vector
1442 nlte_source.Kjj() += com.N.real();
1443
1444 // Sum up the Jacobian
1445 for (Index j=0; j<nq; j++) {
1446 if (not jacobian_quantities[j].propmattype()) continue;
1447 dnlte_source_dx[j].Kjj() += com.dN.real()(joker, j);
1448 }
1449 }
1450}
1451
1452/* Workspace method: Doxygen documentation will be auto-generated */
1453void propmat_clearskyAddXsecAgenda( // Workspace reference:
1454 Workspace& ws,
1455 // WS Output:
1456 PropagationMatrix& propmat_clearsky,
1457 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
1458 // WS Input:
1459 const Vector& f_grid,
1460 const ArrayOfArrayOfSpeciesTag& abs_species,
1461 const ArrayOfRetrievalQuantity& jacobian_quantities,
1462 const Numeric& rtp_pressure,
1463 const Numeric& rtp_temperature,
1464 const Vector& rtp_vmr,
1465 const Agenda& abs_xsec_agenda,
1466 // Verbosity object:
1467 const Verbosity& verbosity) {
1468
1469 // Output of AbsInputFromRteScalars:
1470 Vector abs_p;
1471 Vector abs_t;
1472 Matrix abs_vmrs;
1473 // Output of abs_h2oSet:
1474 Vector abs_h2o;
1475 // Output of abs_coefCalc:
1476 Matrix abs_coef;
1477 ArrayOfMatrix abs_coef_per_species, dabs_coef_dx;
1478
1480 abs_t,
1481 abs_vmrs,
1482 rtp_pressure,
1483 rtp_temperature,
1484 rtp_vmr,
1485 verbosity);
1486
1487 // Absorption cross sections per tag group.
1488 ArrayOfMatrix abs_xsec_per_species;
1489 ArrayOfArrayOfMatrix dabs_xsec_per_species_dx;
1490
1491 // Make all species active:
1492 ArrayOfIndex abs_species_active(abs_species.nelem());
1493 for (Index i = 0; i < abs_species.nelem(); ++i) abs_species_active[i] = i;
1494
1495 // Call agenda to calculate absorption:
1497 abs_xsec_per_species,
1498 dabs_xsec_per_species_dx,
1499 abs_species,
1500 jacobian_quantities,
1501 abs_species_active,
1502 f_grid,
1503 abs_p,
1504 abs_t,
1505 abs_vmrs,
1506 abs_xsec_agenda);
1507
1508 // Calculate absorption coefficients from cross sections:
1509 abs_coefCalcFromXsec(abs_coef,
1510 dabs_coef_dx,
1511 abs_coef_per_species,
1512 abs_xsec_per_species,
1513 dabs_xsec_per_species_dx,
1514 abs_species,
1515 jacobian_quantities,
1516 abs_vmrs,
1517 abs_p,
1518 abs_t,
1519 verbosity);
1520
1521 // Now add abs_coef_per_species to propmat_clearsky:
1523 dpropmat_clearsky_dx,
1524 abs_coef_per_species,
1525 dabs_coef_dx);
1526}
1527
1528/* Workspace method: Doxygen documentation will be auto-generated */
1530 const Vector& f_grid,
1531 const Index& stokes_dim,
1532 const Verbosity&) {
1533 propmat_clearsky = PropagationMatrix(f_grid.nelem(), stokes_dim);
1534}
1535
1536/* Workspace method: Doxygen documentation will be auto-generated */
1538 PropagationMatrix& propmat_clearsky, const Verbosity&) {
1539 for (Index i = 0; i < propmat_clearsky.NumberOfFrequencies(); i++)
1540 if (propmat_clearsky.Kjj()[i] < 0.0) propmat_clearsky.SetAtPosition(0.0, i);
1541}
1542
1543/* Workspace method: Doxygen documentation will be auto-generated */
1545 const Verbosity&) {
1547}
1548
1549/* Workspace method: Doxygen documentation will be auto-generated */
1551 const Verbosity&) {
1553}
1554
1555#ifdef ENABLE_NETCDF
1556/* Workspace method: Doxygen documentation will be auto-generated */
1557/* Included by Claudia Emde, 20100707 */
1558void WriteMolTau( //WS Input
1559 const Vector& f_grid,
1560 const Tensor3& z_field,
1561 const Tensor7& propmat_clearsky_field,
1562 const Index& atmosphere_dim,
1563 //Keyword
1564 const String& filename,
1565 const Verbosity&) {
1566 int retval, ncid;
1567 int nlev_dimid, nlyr_dimid, nwvl_dimid, stokes_dimid, none_dimid;
1568 int dimids[4];
1569 int wvlmin_varid, wvlmax_varid, z_varid, wvl_varid, tau_varid;
1570
1571 ARTS_USER_ERROR_IF (atmosphere_dim != 1,
1572 "WriteMolTau can only be used for atmosphere_dim=1")
1573
1574#pragma omp critical(netcdf__critical_region)
1575 {
1576 // Open file
1577 if ((retval = nc_create(filename.c_str(), NC_CLOBBER, &ncid)))
1578 nca_error(retval, "nc_create");
1579
1580 // Define dimensions
1581 if ((retval = nc_def_dim(ncid, "nlev", (int)z_field.npages(), &nlev_dimid)))
1582 nca_error(retval, "nc_def_dim");
1583
1584 if ((retval =
1585 nc_def_dim(ncid, "nlyr", (int)z_field.npages() - 1, &nlyr_dimid)))
1586 nca_error(retval, "nc_def_dim");
1587
1588 if ((retval = nc_def_dim(ncid, "nwvl", (int)f_grid.nelem(), &nwvl_dimid)))
1589 nca_error(retval, "nc_def_dim");
1590
1591 if ((retval = nc_def_dim(ncid, "none", 1, &none_dimid)))
1592 nca_error(retval, "nc_def_dim");
1593
1594 if ((retval = nc_def_dim(ncid,
1595 "nstk",
1596 (int)propmat_clearsky_field.nbooks(),
1597 &stokes_dimid)))
1598 nca_error(retval, "nc_def_dim");
1599
1600 // Define variables
1601 if ((retval = nc_def_var(
1602 ncid, "wvlmin", NC_DOUBLE, 1, &none_dimid, &wvlmin_varid)))
1603 nca_error(retval, "nc_def_var wvlmin");
1604
1605 if ((retval = nc_def_var(
1606 ncid, "wvlmax", NC_DOUBLE, 1, &none_dimid, &wvlmax_varid)))
1607 nca_error(retval, "nc_def_var wvlmax");
1608
1609 if ((retval = nc_def_var(ncid, "z", NC_DOUBLE, 1, &nlev_dimid, &z_varid)))
1610 nca_error(retval, "nc_def_var z");
1611
1612 if ((retval =
1613 nc_def_var(ncid, "wvl", NC_DOUBLE, 1, &nwvl_dimid, &wvl_varid)))
1614 nca_error(retval, "nc_def_var wvl");
1615
1616 dimids[0] = nlyr_dimid;
1617 dimids[1] = nwvl_dimid;
1618 dimids[2] = stokes_dimid;
1619 dimids[3] = stokes_dimid;
1620
1621 if ((retval =
1622 nc_def_var(ncid, "tau", NC_DOUBLE, 4, &dimids[0], &tau_varid)))
1623 nca_error(retval, "nc_def_var tau");
1624
1625 // Units
1626 if ((retval = nc_put_att_text(ncid, wvlmin_varid, "units", 2, "nm")))
1627 nca_error(retval, "nc_put_att_text");
1628
1629 if ((retval = nc_put_att_text(ncid, wvlmax_varid, "units", 2, "nm")))
1630 nca_error(retval, "nc_put_att_text");
1631
1632 if ((retval = nc_put_att_text(ncid, z_varid, "units", 2, "km")))
1633 nca_error(retval, "nc_put_att_text");
1634
1635 if ((retval = nc_put_att_text(ncid, wvl_varid, "units", 2, "nm")))
1636 nca_error(retval, "nc_put_att_text");
1637
1638 if ((retval = nc_put_att_text(ncid, tau_varid, "units", 1, "-")))
1639 nca_error(retval, "nc_put_att_text");
1640
1641 // End define mode. This tells netCDF we are done defining
1642 // metadata.
1643 if ((retval = nc_enddef(ncid))) nca_error(retval, "nc_enddef");
1644
1645 // Assign data
1646 double wvlmin[1];
1647 wvlmin[0] = SPEED_OF_LIGHT / f_grid[f_grid.nelem() - 1] * 1e9;
1648 if ((retval = nc_put_var_double(ncid, wvlmin_varid, &wvlmin[0])))
1649 nca_error(retval, "nc_put_var");
1650
1651 double wvlmax[1];
1652 wvlmax[0] = SPEED_OF_LIGHT / f_grid[0] * 1e9;
1653 if ((retval = nc_put_var_double(ncid, wvlmax_varid, &wvlmax[0])))
1654 nca_error(retval, "nc_put_var");
1655
1656 double z[z_field.npages()];
1657 for (int iz = 0; iz < z_field.npages(); iz++)
1658 z[iz] = z_field(z_field.npages() - 1 - iz, 0, 0) * 1e-3;
1659
1660 if ((retval = nc_put_var_double(ncid, z_varid, &z[0])))
1661 nca_error(retval, "nc_put_var");
1662
1663 double wvl[f_grid.nelem()];
1664 for (int iv = 0; iv < f_grid.nelem(); iv++)
1665 wvl[iv] = SPEED_OF_LIGHT / f_grid[f_grid.nelem() - 1 - iv] * 1e9;
1666
1667 if ((retval = nc_put_var_double(ncid, wvl_varid, &wvl[0])))
1668 nca_error(retval, "nc_put_var");
1669
1670 const Index zfnp = z_field.npages() - 1;
1671 const Index fgne = f_grid.nelem();
1672 const Index amfnb = propmat_clearsky_field.nbooks();
1673
1674 Tensor4 tau(zfnp, fgne, amfnb, amfnb, 0.);
1675
1676 // Calculate average tau for layers
1677 for (int is = 0; is < propmat_clearsky_field.nlibraries(); is++)
1678 for (int iz = 0; iz < zfnp; iz++)
1679 for (int iv = 0; iv < fgne; iv++)
1680 for (int is1 = 0; is1 < amfnb; is1++)
1681 for (int is2 = 0; is2 < amfnb; is2++)
1682 // sum up all species
1683 tau(iz, iv, is1, is2) +=
1684 0.5 *
1685 (propmat_clearsky_field(is,
1686 f_grid.nelem() - 1 - iv,
1687 is1,
1688 is2,
1689 z_field.npages() - 1 - iz,
1690 0,
1691 0) +
1692 propmat_clearsky_field(is,
1693 f_grid.nelem() - 1 - iv,
1694 is1,
1695 is2,
1696 z_field.npages() - 2 - iz,
1697 0,
1698 0)) *
1699 (z_field(z_field.npages() - 1 - iz, 0, 0) -
1700 z_field(z_field.npages() - 2 - iz, 0, 0));
1701
1702 if ((retval = nc_put_var_double(ncid, tau_varid, tau.get_c_array())))
1703 nca_error(retval, "nc_put_var");
1704
1705 // Close the file
1706 if ((retval = nc_close(ncid))) nca_error(retval, "nc_close");
1707 }
1708}
1709
1710#else
1711
1712void WriteMolTau( //WS Input
1713 const Vector& f_grid _U_,
1714 const Tensor3& z_field _U_,
1715 const Tensor7& propmat_clearsky_field _U_,
1716 const Index& atmosphere_dim _U_,
1717 //Keyword
1718 const String& filename _U_,
1719 const Verbosity&) {
1720 ARTS_USER_ERROR_IF(true,
1721 "The workspace method WriteMolTau is not available"
1722 "because ARTS was compiled without NetCDF support.");
1723}
1724
1725#endif /* ENABLE_NETCDF */
1726
1727/* Workspace method: Doxygen documentation will be auto-generated */
1729 // WS Output:
1730 ArrayOfMatrix& abs_xsec_per_species,
1732 // WS Input:
1733 const ArrayOfArrayOfSpeciesTag& abs_species,
1734 const ArrayOfRetrievalQuantity& jacobian_quantities,
1735 const ArrayOfIndex& abs_species_active,
1736 const Vector& f_grid,
1737 const Vector& abs_p,
1738 const Vector& abs_t,
1739 const Matrix& abs_vmrs,
1740 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1742 const Index& lbl_checked,
1743 const Verbosity&) {
1745 "abs_xsec_per_speciesAddLines",
1746 "2021-07-13",
1747 "This function is no longer up to date. It only exists to satisfy "
1748 "lookup table calculations before these are updated.\n"
1749 "Once the lookup table calculations are up-to-date, this function "
1750 "is fully replaced with propmat_clearskyAddLines, with better functionality\n")
1751
1752 ARTS_USER_ERROR_IF(not lbl_checked, "Must check LBL calculations")
1754 jacobian_quantities.nelem(),
1755 "There's a hard deprecation of derivatives using old style lbl-calculations "
1756 "with derivatives, switch to propmat_clearskyAddLines")
1758 abs_species.nelem() not_eq abs_xsec_per_species.nelem() or
1759 abs_species.nelem() not_eq abs_vmrs.nrows() or
1760 abs_species.nelem() not_eq abs_lines_per_species.nelem(),
1761 "The following variables must all have the same dimension:\n"
1762 "abs_species: ",
1763 abs_species.nelem(),
1764 '\n',
1765 "abs_xsec_per_species: ",
1766 abs_xsec_per_species.nelem(),
1767 '\n',
1768 "abs_vmrs: ",
1769 abs_vmrs.nrows(),
1770 '\n',
1771 "abs_lines_per_species: ",
1772 abs_lines_per_species.nelem(),
1773 '\n')
1775 min(abs_t) < 0,
1776 "Temperature must be at least 0 K. But you request an absorption\n"
1777 "calculation at ",
1778 min(abs_t),
1779 " K!")
1780
1781 // Size of problem
1782 const Index np = abs_p.nelem();
1783 const Index ns = abs_species_active.nelem();
1784
1785 // Calculations data
1786 constexpr Options::LblSpeedup speedup_type = Options::LblSpeedup::None;
1787 const EnergyLevelMap rtp_nlte;
1788
1789#pragma omp parallel for if (!arts_omp_in_parallel()) collapse(2) \
1790 schedule(dynamic)
1791 for (Index ip = 0; ip < np; ip++) {
1792 for (Index is = 0; is < ns; is++) {
1793 // Skip it if there are no species or there is Zeeman requested
1794 if (not abs_species[abs_species_active[is]].nelem() or
1795 abs_species[abs_species_active[is]].Zeeman() or
1796 not abs_lines_per_species[abs_species_active[is]].nelem())
1797 continue;
1798
1799 // Local compute values
1800 thread_local LineShape::ComputeData com(
1801 f_grid, jacobian_quantities, false);
1802 thread_local LineShape::ComputeData sparse_com(
1803 Vector(0), jacobian_quantities, false);
1804
1805 for (auto& band : abs_lines_per_species[abs_species_active[is]]) {
1807 com,
1808 sparse_com,
1809 band,
1810 jacobian_quantities,
1811 rtp_nlte,
1812 band.BroadeningSpeciesVMR(abs_vmrs(joker, ip), abs_species),
1813 abs_species[abs_species_active[is]],
1814 abs_vmrs(abs_species_active[is], ip),
1815 isotopologue_ratios[band.Isotopologue()],
1816 abs_p[ip],
1817 abs_t[ip],
1818 0,
1819 0,
1821 speedup_type,
1822 false);
1823 }
1824
1825 // Sum up the propagation matrix
1826 com.F /= number_density(abs_p[ip], abs_t[ip]) *
1827 abs_vmrs(abs_species_active[is], ip);
1828 abs_xsec_per_species[abs_species_active[is]](joker, ip) += com.F.real();
1829 }
1830 }
1831}
void set_vmr_from_first_species(Vector &vmr, const String &species_name, const ArrayOfArrayOfSpeciesTag &abs_species, const Matrix &abs_vmrs)
set_abs_from_first_species.
Definition: absorption.cc:110
Declarations required for the calculation of absorption coefficients.
This file contains the definition of Array.
Array< Matrix > ArrayOfMatrix
An array of matrices.
Definition: array.h:66
Index TotalNumberOfElements(const Array< Array< base > > &aa)
Determine total number of elements in an ArrayOfArray.
Definition: array.h:345
The global header file for ARTS.
int main(int argc, char **argv)
void abs_xsec_agendaExecute(Workspace &ws, ArrayOfMatrix &abs_xsec_per_species, ArrayOfArrayOfMatrix &dabs_xsec_per_species_dx, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfIndex &abs_species_active, const Vector &f_grid, const Vector &abs_p, const Vector &abs_t, const Matrix &abs_vmrs, const Agenda &input_agenda)
Definition: auto_md.cc:23646
void chk_size(const String &x_name, ConstVectorView x, const Index &c)
Runtime check for size of Vector.
Definition: check_input.cc:421
The Agenda class.
Definition: agenda_class.h:51
This can be used to make arrays out of anything.
Definition: array.h:108
Index nelem() const ARTS_NOEXCEPT
Number of elements.
Definition: array.h:197
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 nrows() const noexcept
Definition: matpackI.h:1061
Index ncols() const noexcept
Definition: matpackI.h:1062
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:541
The Matrix class.
Definition: matpackI.h:1270
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1063
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.
The range class.
Definition: matpackI.h:166
Stokes vector is as Propagation matrix but only has 4 possible values.
The Tensor3 class.
Definition: matpackIII.h:344
const Numeric * get_c_array() const ARTS_NOEXCEPT
Conversion to plain C-array.
Definition: matpackIV.cc:342
The Tensor4 class.
Definition: matpackIV.h:427
The Tensor7 class.
Definition: matpackVII.h:2397
The Vector class.
Definition: matpackI.h:908
void resize(Index n)
Resize function.
Definition: matpackI.cc:411
Workspace class.
Definition: workspace_ng.h:40
#define _U_
Definition: config.h:180
#define ARTS_ASSERT(condition,...)
Definition: debug.h:83
#define ARTS_USER_ERROR_IF(condition,...)
Definition: debug.h:134
#define DEPRECATED_FUNCTION(FUNCTION_NAME, DATE_STRING_OR_SILLY_REASON,...)
Definition: depr.h:9
void find_xml_file(String &filename, const Verbosity &verbosity)
Find an xml file.
Definition: file.cc:414
This file contains basic functions to handle ASCII files.
bool species_iso_match(const RetrievalQuantity &rq, const Species::IsotopeRecord &ir)
Returns if the Retrieval quantity is VMR derivative for all the species in the species tags.
Definition: jacobian.cc:1136
bool do_frequency_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants a frequency derivative.
Definition: jacobian.cc:1166
bool supports_continuum(const ArrayOfRetrievalQuantity &js)
Returns if the array supports continuum derivatives.
Definition: jacobian.cc:1092
bool species_match(const RetrievalQuantity &rq, const ArrayOfSpeciesTag &ast)
Returns if the Retrieval quantity is VMR derivative for all the species in the species tags.
Definition: jacobian.cc:1114
bool do_magnetic_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants a magnetic derivative.
Definition: jacobian.cc:1170
bool do_temperature_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants the temperature derivative.
Definition: jacobian.cc:1145
Numeric magnetic_field_perturbation(const ArrayOfRetrievalQuantity &js) noexcept
Returns the magnetic field perturbation if it exists.
Definition: jacobian.cc:1190
bool is_frequency_parameter(const RetrievalQuantity &t) noexcept
Returns if the Retrieval quantity is a frequency parameter in propagation matrix calculations.
Definition: jacobian.cc:1001
Numeric temperature_perturbation(const ArrayOfRetrievalQuantity &js) noexcept
Returns the temperature perturbation if it exists.
Definition: jacobian.cc:1174
Numeric frequency_perturbation(const ArrayOfRetrievalQuantity &js) noexcept
Returns the frequency perturbation if it exists.
Definition: jacobian.cc:1182
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
#define min(a, b)
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:215
const Numeric ELECTRON_CHARGE
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:1298
void isotopologue_ratiosInitFromHitran(SpeciesIsotopologueRatios &isotopologue_ratios, const Verbosity &)
WORKSPACE METHOD: isotopologue_ratiosInitFromHitran.
Definition: m_abs.cc:1550
void abs_xsec_per_speciesInit(ArrayOfMatrix &abs_xsec_per_species, ArrayOfArrayOfMatrix &dabs_xsec_per_species_dx, const ArrayOfArrayOfSpeciesTag &tgs, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfIndex &abs_species_active, const Vector &f_grid, const Vector &abs_p, const Index &abs_xsec_agenda_checked, const Verbosity &verbosity)
WORKSPACE METHOD: abs_xsec_per_speciesInit.
Definition: m_abs.cc:375
void abs_xsec_per_speciesAddConts(ArrayOfMatrix &abs_xsec_per_species, ArrayOfArrayOfMatrix &dabs_xsec_per_species_dx, const ArrayOfArrayOfSpeciesTag &tgs, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfIndex &abs_species_active, const Vector &f_grid, const Vector &abs_p, const Vector &abs_t, const Matrix &abs_vmrs, const ArrayOfString &abs_cont_names, const ArrayOfVector &abs_cont_parameters, const ArrayOfString &abs_cont_models, const Verbosity &verbosity)
WORKSPACE METHOD: abs_xsec_per_speciesAddConts.
Definition: m_abs.cc:474
void abs_speciesDefineAll(ArrayOfArrayOfSpeciesTag &abs_species, Index &propmat_clearsky_agenda_checked, Index &abs_xsec_agenda_checked, const Verbosity &verbosity)
WORKSPACE METHOD: abs_speciesDefineAll.
Definition: m_abs.cc:203
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:1326
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:885
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:762
void propmat_clearskyAddXsecAgenda(Workspace &ws, PropagationMatrix &propmat_clearsky, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, const Vector &f_grid, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfRetrievalQuantity &jacobian_quantities, const Numeric &rtp_pressure, const Numeric &rtp_temperature, const Vector &rtp_vmr, const Agenda &abs_xsec_agenda, const Verbosity &verbosity)
WORKSPACE METHOD: propmat_clearskyAddXsecAgenda.
Definition: m_abs.cc:1453
void abs_coefCalcFromXsec(Matrix &abs_coef, ArrayOfMatrix &dabs_coef_dx, ArrayOfMatrix &abs_coef_per_species, const ArrayOfMatrix &abs_xsec_per_species, const ArrayOfArrayOfMatrix &dabs_xsec_per_species_dx, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfRetrievalQuantity &jacobian_quantities, const Matrix &abs_vmrs, const Vector &abs_p, const Vector &abs_t, const Verbosity &verbosity)
Definition: m_abs.cc:245
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:739
String continua_model_error_message(const ArrayOfString &abs_cont_names, const ArrayOfVector &abs_cont_parameters, const ArrayOfString &abs_cont_models)
Definition: m_abs.cc:448
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:723
const Numeric VACUUM_PERMITTIVITY
void propmat_clearskyAddFromAbsCoefPerSpecies(PropagationMatrix &propmat_clearsky, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, const ArrayOfMatrix &abs_coef_per_species, const ArrayOfMatrix &dabs_coef_dx)
Definition: m_abs.cc:841
const Numeric ELECTRON_MASS
const Numeric PI
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:70
const Numeric SPEED_OF_LIGHT
void abs_speciesDefineAllInScenario(ArrayOfArrayOfSpeciesTag &tgs, Index &propmat_clearsky_agenda_checked, Index &abs_xsec_agenda_checked, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: abs_speciesDefineAllInScenario.
Definition: m_abs.cc:153
void abs_xsec_per_speciesAddLines(ArrayOfMatrix &abs_xsec_per_species, ArrayOfArrayOfMatrix &, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfIndex &abs_species_active, const Vector &f_grid, const Vector &abs_p, const Vector &abs_t, const Matrix &abs_vmrs, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const SpeciesIsotopologueRatios &isotopologue_ratios, const Index &lbl_checked, const Verbosity &)
WORKSPACE METHOD: abs_xsec_per_speciesAddLines.
Definition: m_abs.cc:1728
void propmat_clearskyForceNegativeToZero(PropagationMatrix &propmat_clearsky, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyForceNegativeToZero.
Definition: m_abs.cc:1537
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 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 ArrayOfSpeciesTag &select_speciestags, const Index &robust, const Verbosity &verbosity)
WORKSPACE METHOD: propmat_clearskyAddLines.
Definition: m_abs.cc:1338
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 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:1054
void propmat_clearskyZero(PropagationMatrix &propmat_clearsky, const Vector &f_grid, const Index &stokes_dim, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyZero.
Definition: m_abs.cc:1529
void abs_lines_per_speciesCreateFromLines(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfAbsorptionLines &abs_lines, const ArrayOfArrayOfSpeciesTag &tgs, const Verbosity &)
WORKSPACE METHOD: abs_lines_per_speciesCreateFromLines.
Definition: m_abs.cc:93
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:1558
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:224
void isotopologue_ratiosInitFromBuiltin(SpeciesIsotopologueRatios &isotopologue_ratios, const Verbosity &)
WORKSPACE METHOD: isotopologue_ratiosInitFromBuiltin.
Definition: m_abs.cc:1544
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 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:948
void abs_speciesSet(ArrayOfArrayOfSpeciesTag &abs_species, Index &abs_xsec_agenda_checked, 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:138
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:205
#define CREATE_OUT3
Definition: messages.h:207
#define CREATE_OUT2
Definition: messages.h:206
Index nelem(const Lines &l)
Number of lines.
SpeciesIsotopologueRatios isotopologue_ratios()
Particulates ENUMCLASS(Line, char, VMR, Strength, Center, ShapeG0X0, ShapeG0X1, ShapeG0X2, ShapeG0X3, ShapeD0X0, ShapeD0X1, ShapeD0X2, ShapeD0X3, ShapeG2X0, ShapeG2X1, ShapeG2X2, ShapeG2X3, ShapeD2X0, ShapeD2X1, ShapeD2X2, ShapeD2X3, ShapeFVCX0, ShapeFVCX1, ShapeFVCX2, ShapeFVCX3, ShapeETAX0, ShapeETAX1, ShapeETAX2, ShapeETAX3, ShapeYX0, ShapeYX1, ShapeYX2, ShapeYX3, ShapeGX0, ShapeGX1, ShapeGX2, ShapeGX3, ShapeDVX0, ShapeDVX1, ShapeDVX2, ShapeDVX3, ECS_SCALINGX0, ECS_SCALINGX1, ECS_SCALINGX2, ECS_SCALINGX3, ECS_BETAX0, ECS_BETAX1, ECS_BETAX2, ECS_BETAX3, ECS_LAMBDAX0, ECS_LAMBDAX1, ECS_LAMBDAX2, ECS_LAMBDAX3, ECS_DCX0, ECS_DCX1, ECS_DCX2, ECS_DCX3, NLTE) static_assert(Index(Line ArrayOfSpeciesTagVMR
Definition: jacobian.h:101
Temperature
Definition: jacobian.h:58
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:3090
bool good_linear_sparse_f_grid(const Vector &f_grid_dense, const Vector &f_grid_sparse) noexcept
Definition: lineshape.cc:3190
Vector triple_sparse_f_grid(const Vector &f_grid, const Numeric &sparse_df) noexcept
Definition: lineshape.cc:3203
Vector linear_sparse_f_grid(const Vector &f_grid, const Numeric &sparse_df) ARTS_NOEXCEPT
Definition: lineshape.cc:3168
constexpr bool same_or_joker(const IsotopeRecord &ir1, const IsotopeRecord &ir2) noexcept
constexpr IsotopologueRatios isotopologue_ratiosInitFromBuiltin()
Implements Zeeman modeling.
Definition: zeemandata.cc:331
invlib::Matrix< ArtsMatrix > Matrix
invlib wrapper type for ARTS matrices.
Definition: oem.h:33
invlib::Vector< ArtsVector > Vector
invlib wrapper type for ARTS vectors.
Definition: oem.h:31
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:668
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:1840
Declaration of functions in rte.cc.
void check_abs_species(const ArrayOfArrayOfSpeciesTag &abs_species)
Array< SingleLine > lines
A list of individual lines.
Index NumLines() const noexcept
Number of lines.
void ReverseLines() noexcept
Reverses the order of the internal lines.
void AppendSingleLine(SingleLine &&sl)
Appends a single line to the absorption lines.
Main computational data for the line shape and strength calculations.
Definition: lineshape.h:781
void interp_add_even(const ComputeData &sparse) ARTS_NOEXCEPT
Add a sparse grid to this grid via linear interpolation.
Definition: lineshape.cc:3227
ComplexVector N
Definition: lineshape.h:782
void interp_add_triplequad(const ComputeData &sparse) ARTS_NOEXCEPT
Add a sparse grid to this grid via square interpolation.
Definition: lineshape.cc:3268
ComplexVector F
Definition: lineshape.h:782
ComplexMatrix dF
Definition: lineshape.h:783
ComplexMatrix dN
Definition: lineshape.h:783
This file contains basic functions to handle XML data files.