ARTS 2.5.0 (git: 9ee3ac6c)
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 = createEmptyCopy(lines);
127 for (Index k=lines.NumLines()-1; k>=0; k--)
128 if (low <= lines.F0(k) and upp >= lines.F0(k))
129 these_lines.AppendSingleLine(lines.PopLine(k));
130
131 // Append these lines after sorting them if there are any of them
132 if (these_lines.NumLines()) {
133 these_lines.ReverseLines();
134 abs_lines_per_species[i].push_back(these_lines);
135 }
136
137 // If this means we have deleted all lines, then we leave
138 if (lines.NumLines() == 0)
139 goto leave_inner_loop;
140 }
141 else {
142 abs_lines_per_species[i].push_back(lines);
143 goto leave_inner_loop;
144 }
145 }
146 }
147 leave_inner_loop: {}
148 }
149}
150
151/* Workspace method: Doxygen documentation will be auto-generated */
154 Index& propmat_clearsky_agenda_checked,
155 Index& abs_xsec_agenda_checked,
156 // Control Parameters:
157 const String& basename,
158 const Verbosity& verbosity) {
160
161 // Invalidate agenda check flags
162 propmat_clearsky_agenda_checked = false;
163 abs_xsec_agenda_checked = false;
164
165 // We want to make lists of included and excluded species:
166 ArrayOfString included(0), excluded(0);
167
168 tgs.resize(0);
169
170 for (Index i = 0; i < Index(Species::Species::FINAL); ++i) {
171 const String specname = Species::toShortName(Species::Species(i));
172
173 String filename = basename;
174 if (basename.length() && basename[basename.length() - 1] != '/')
175 filename += ".";
176 filename += specname;
177
178 try {
179 find_xml_file(filename, verbosity);
180 // Add to included list:
181 included.push_back(specname);
182
183 // Add this tag group to tgs:
184 tgs.emplace_back(ArrayOfSpeciesTag(specname));
185 } catch (const std::runtime_error& e) {
186 // The file for the species could not be found.
187 excluded.push_back(specname);
188 }
189 }
190
191 // Some nice output:
192 out2 << " Included Species (" << included.nelem() << "):\n";
193 for (Index i = 0; i < included.nelem(); ++i)
194 out2 << " " << included[i] << "\n";
195
196 out2 << " Excluded Species (" << excluded.nelem() << "):\n";
197 for (Index i = 0; i < excluded.nelem(); ++i)
198 out2 << " " << excluded[i] << "\n";
199}
200
201/* Workspace method: Doxygen documentation will be auto-generated */
202void abs_speciesDefineAll( // WS Output:
203 ArrayOfArrayOfSpeciesTag& abs_species,
204 Index& propmat_clearsky_agenda_checked,
205 Index& abs_xsec_agenda_checked,
206 // Control Parameters:
207 const Verbosity& verbosity) {
208 // Species lookup data:
209
210 // We want to make lists of all species
211 ArrayOfString specs(0);
212 for (Index i = 0; i < Index(Species::Species::FINAL); ++i) {
213 if (Species::Species(i) not_eq Species::Species::Bath) {
214 specs.emplace_back(Species::toShortName(Species::Species(i)));
215 }
216 }
217
218 // Set the values
219 abs_speciesSet(abs_species, abs_xsec_agenda_checked, propmat_clearsky_agenda_checked, specs, verbosity);
220}
221
222/* Workspace method: Doxygen documentation will be auto-generated */
223void AbsInputFromAtmFields( // WS Output:
224 Vector& abs_p,
225 Vector& abs_t,
226 Matrix& abs_vmrs,
227 // WS Input:
228 const Index& atmosphere_dim,
229 const Vector& p_grid,
230 const Tensor3& t_field,
231 const Tensor4& vmr_field,
232 const Verbosity&) {
233 // First, make sure that we really have a 1D atmosphere:
234 ARTS_USER_ERROR_IF (1 != atmosphere_dim,
235 "Atmospheric dimension must be 1D, but atmosphere_dim is ",
236 atmosphere_dim, ".")
237
238 abs_p = p_grid;
239 abs_t = t_field(joker, 0, 0);
240 abs_vmrs = vmr_field(joker, joker, 0, 0);
241}
242
243/* Workspace method: Doxygen documentation will be auto-generated */
244void abs_coefCalcFromXsec( // WS Output:
245 Matrix& abs_coef,
246 ArrayOfMatrix& dabs_coef_dx,
247 ArrayOfMatrix& abs_coef_per_species,
248 // WS Input:
249 const ArrayOfMatrix& abs_xsec_per_species,
250 const ArrayOfArrayOfMatrix& dabs_xsec_per_species_dx,
251 const ArrayOfArrayOfSpeciesTag& abs_species,
252 const ArrayOfRetrievalQuantity& jacobian_quantities,
253 const Matrix& abs_vmrs,
254 const Vector& abs_p,
255 const Vector& abs_t,
256 const Verbosity& verbosity) {
258
259 // Check that abs_vmrs and abs_xsec_per_species really have compatible
260 // dimensions. In abs_vmrs there should be one row for each tg:
261 ARTS_USER_ERROR_IF(abs_vmrs.nrows() != abs_xsec_per_species.nelem(),
262 "Variable abs_vmrs must have compatible dimension to abs_xsec_per_species.\n"
263 "abs_vmrs.nrows() = ", abs_vmrs.nrows(), "\n"
264 "abs_xsec_per_species.nelem() = ", abs_xsec_per_species.nelem())
265
266 // Check that number of altitudes are compatible. We only check the
267 // first element, this is possilble because within arts all elements
268 // are on the same altitude grid.
269 ARTS_USER_ERROR_IF(abs_vmrs.ncols() != abs_xsec_per_species[0].ncols(),
270 "Variable abs_vmrs must have same numbers of altitudes as abs_xsec_per_species.\n"
271 "abs_vmrs.ncols() = ", abs_vmrs.ncols(), "\n"
272 "abs_xsec_per_species[0].ncols() = ", abs_xsec_per_species[0].ncols());
273
274 // Check dimensions of abs_p and abs_t:
275 chk_size("abs_p", abs_p, abs_vmrs.ncols());
276 chk_size("abs_t", abs_t, abs_vmrs.ncols());
277
278 // Initialize abs_coef and abs_coef_per_species. The array dimension of abs_coef_per_species
279 // is the same as that of abs_xsec_per_species. The dimension of abs_coef should
280 // be equal to one of the abs_xsec_per_species enries.
281
282 abs_coef.resize(abs_xsec_per_species[0].nrows(),
283 abs_xsec_per_species[0].ncols());
284 abs_coef = 0;
285
286 dabs_coef_dx.resize(jacobian_quantities.nelem());
287
288 for (Index ii = 0; ii < jacobian_quantities.nelem(); ii++) {
289 const auto& deriv = jacobian_quantities[ii];
290
291 if (not deriv.propmattype()) continue;
292
293 dabs_coef_dx[ii].resize(abs_xsec_per_species[0].nrows(),
294 abs_xsec_per_species[0].ncols());
295 dabs_coef_dx[ii] = 0.0;
296 }
297
298 abs_coef_per_species.resize(abs_xsec_per_species.nelem());
299
300 out3
301 << " Computing abs_coef and abs_coef_per_species from abs_xsec_per_species.\n";
302 // Loop through all tag groups
303 for (Index i = 0; i < abs_xsec_per_species.nelem(); ++i) {
304 out3 << " Tag group " << i << "\n";
305
306 // Make this element of abs_xsec_per_species the right size:
307 abs_coef_per_species[i].resize(abs_xsec_per_species[i].nrows(),
308 abs_xsec_per_species[i].ncols());
309 abs_coef_per_species[i] = 0; // Initialize all elements to 0.
310
311 // Loop through all altitudes
312 for (Index j = 0; j < abs_xsec_per_species[i].ncols(); j++) {
313 // Calculate total number density from pressure and temperature.
314 const Numeric n = number_density(abs_p[j], abs_t[j]);
315 const Numeric dn_dT = dnumber_density_dt(abs_p[j], abs_t[j]);
316 // Wasted calculations when Jacobians are not calculated...
317 // Though this is called seldom enough that it this fine? value is -1/t*n
318
319 // Loop through all frequencies
320 for (Index k = 0; k < abs_xsec_per_species[i].nrows(); k++) {
321 abs_coef_per_species[i](k, j) =
322 abs_xsec_per_species[i](k, j) * n * abs_vmrs(i, j);
323
324 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
325 const auto& deriv = jacobian_quantities[iq];
326
327 if (not deriv.propmattype()) continue;
328
329 if (deriv == Jacobian::Atm::Temperature) {
330 dabs_coef_dx[iq](k, j) +=
331 (dabs_xsec_per_species_dx[i][iq](k, j) * n +
332 abs_xsec_per_species[i](k, j) * dn_dT) *
333 abs_vmrs(i, j);
334 } else if (deriv == Jacobian::Line::VMR) {
335 bool seco = false, main = false;
336 for (const auto& s : abs_species[i]) {
337 if (species_match(
338 deriv, s.cia_2nd_species) or
339 s.type not_eq Species::TagType::Cia)
340 seco = true;
342 deriv,
343 s.Isotopologue()))
344 main = true;
345 }
346 if (main and seco) {
347 dabs_coef_dx[iq](k, j) +=
348 (dabs_xsec_per_species_dx[i][iq](k, j) * abs_vmrs(i, j) +
349 abs_xsec_per_species[i](k, j)) *
350 n;
351 } else if (main) {
352 dabs_coef_dx[iq](k, j) += abs_xsec_per_species[i](k, j) * n;
353 } else if (seco) {
354 dabs_coef_dx[iq](k, j) +=
355 dabs_xsec_per_species_dx[i][iq](k, j) * abs_vmrs(i, j) * n;
356 }
357 } else if (deriv == Jacobian::Special::ArrayOfSpeciesTagVMR) {
358 dabs_coef_dx[iq](k, j) +=
359 dabs_xsec_per_species_dx[i][iq](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::Predefined) {
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 if (name == "O2-MPM2020") continue;
577
578 // Check, if we have parameters for this model. For
579 // this, the model name must be listed in
580 // abs_cont_names.
581 const Index n =
582 find(abs_cont_names.begin(), abs_cont_names.end(), name) -
583 abs_cont_names.begin();
584
585 // n==abs_cont_names.nelem() indicates that
586 // the name was not found.
587 ARTS_USER_ERROR_IF (n == abs_cont_names.nelem(),
588 "Cannot find model ", name, " in abs_cont_names.")
589
590 // Ok, the tag specifies a valid continuum model and
591 // we have continuum parameters.
592
593 if (out3.sufficient_priority()) {
594 ostringstream os;
595 os << " Adding " << name << " to tag group " << i << ".\n";
596 out3 << os.str();
597 }
598
599 // find the options for this continuum tag from the input array
600 // of options. The actual field of the array is n:
601 const String ContOption = abs_cont_models[n];
602
603 // Set abs_h2o, abs_n2, and abs_o2 from the first matching species.
604 set_vmr_from_first_species(abs_h2o, "H2O", tgs, abs_vmrs);
605 set_vmr_from_first_species(abs_n2, "N2", tgs, abs_vmrs);
606 set_vmr_from_first_species(abs_o2, "O2", tgs, abs_vmrs);
607
608 // Add the continuum for this tag. The parameters in
609 // this call should be clear. The vmr is in
610 // abs_vmrs(i,Range(joker)). The other vmr variables,
611 // abs_h2o, abs_n2, and abs_o2 contains the real vmr of H2O,
612 // N2, nad O2, which are needed as additional information for
613 // certain continua:
614 // abs_h2o for
615 // O2-PWR88, O2-PWR93, O2-PWR98,
616 // O2-MPM85, O2-MPM87, O2-MPM89, O2-MPM92, O2-MPM93,
617 // O2-TRE05,
618 // O2-SelfContStandardType, O2-SelfContMPM93, O2-SelfContPWR93,
619 // N2-SelfContMPM93, N2-DryContATM01,
620 // N2-CIArotCKDMT252, N2-CIAfunCKDMT252
621 // abs_n2 for
622 // H2O-SelfContCKD24, H2O-ForeignContCKD24,
623 // O2-v0v0CKDMT100,
624 // CO2-ForeignContPWR93, CO2-ForeignContHo66
625 // abs_o2 for
626 // N2-CIArotCKDMT252, N2-CIAfunCKDMT252
627 if (!do_jac)
628 xsec_continuum_tag(abs_xsec_per_species[i],
629 name,
630 abs_cont_parameters[n],
631 abs_cont_models[n],
632 f_grid,
633 abs_p,
634 abs_t,
635 abs_n2,
636 abs_h2o,
637 abs_o2,
638 abs_vmrs(i, Range(joker)),
639 verbosity);
640 else // The Jacobian block
641 {
642 // Needs a reseted block here...
643 for (Index iv = 0; iv < f_grid.nelem(); iv++) {
644 for (Index ip = 0; ip < abs_p.nelem(); ip++) {
645 if (do_freq_jac) jacs_df(iv, ip) = 0.0;
646 if (do_temp_jac) jacs_dt(iv, ip) = 0.0;
647 normal(iv, ip) = 0.0;
648 }
649 }
650
651 // Normal calculations
652 xsec_continuum_tag(normal,
653 name,
654 abs_cont_parameters[n],
655 abs_cont_models[n],
656 f_grid,
657 abs_p,
658 abs_t,
659 abs_n2,
660 abs_h2o,
661 abs_o2,
662 abs_vmrs(i, Range(joker)),
663 verbosity);
664
665 // Frequency calculations
666 if (do_freq_jac)
667 xsec_continuum_tag(jacs_df,
668 name,
669 abs_cont_parameters[n],
670 abs_cont_models[n],
671 dfreq,
672 abs_p,
673 abs_t,
674 abs_n2,
675 abs_h2o,
676 abs_o2,
677 abs_vmrs(i, Range(joker)),
678 verbosity);
679
680 //Temperature calculations
681 if (do_temp_jac)
682 xsec_continuum_tag(jacs_dt,
683 name,
684 abs_cont_parameters[n],
685 abs_cont_models[n],
686 f_grid,
687 abs_p,
688 dabs_t,
689 abs_n2,
690 abs_h2o,
691 abs_o2,
692 abs_vmrs(i, Range(joker)),
693 verbosity);
694 for (Index iv = 0; iv < f_grid.nelem(); iv++) {
695 for (Index ip = 0; ip < abs_p.nelem(); ip++) {
696 abs_xsec_per_species[i](iv, ip) += normal(iv, ip);
697 for (Index iq = 0; iq < jacobian_quantities.nelem();
698 iq++) {
699 const auto& deriv = jacobian_quantities[iq];
700
701 if (not deriv.propmattype()) continue;
702
703 if (is_frequency_parameter(deriv))
704 dabs_xsec_per_species_dx[i][iq](iv, ip) +=
705 (jacs_df(iv, ip) - normal(iv, ip)) * (1. / df);
706 else if (deriv == Jacobian::Atm::Temperature)
707 dabs_xsec_per_species_dx[i][iq](iv, ip) +=
708 (jacs_dt(iv, ip) - normal(iv, ip)) * (1. / dt);
709 }
710 }
711 }
712 }
713 // Calling this function with a row of Matrix abs_vmrs
714 // is possible because it uses Views.
715 }
716 }
717 }
718}
719
720//======================================================================
721// Methods related to continua
722//======================================================================
723
724/* Workspace method: Doxygen documentation will be auto-generated */
725void abs_cont_descriptionInit( // WS Output:
726 ArrayOfString& abs_cont_names,
727 ArrayOfString& abs_cont_options,
728 ArrayOfVector& abs_cont_parameters,
729 const Verbosity& verbosity) {
731
732 abs_cont_names.resize(0);
733 abs_cont_options.resize(0);
734 abs_cont_parameters.resize(0);
735 out2 << " Initialized abs_cont_names \n"
736 " abs_cont_models\n"
737 " abs_cont_parameters.\n";
738}
739
740/* Workspace method: Doxygen documentation will be auto-generated */
741void abs_cont_descriptionAppend( // WS Output:
742 ArrayOfString& abs_cont_names,
743 ArrayOfString& abs_cont_models,
744 ArrayOfVector& abs_cont_parameters,
745 // Control Parameters:
746 const String& tagname,
747 const String& model,
748 const Vector& userparameters,
749 const Verbosity&) {
750 // First we have to check that name matches a continuum species tag.
751 check_continuum_model(tagname);
752
753 //cout << " + tagname: " << tagname << "\n";
754 //cout << " + model: " << model << "\n";
755 //cout << " + parameters: " << userparameters << "\n";
756
757 // Add name and parameters to the apropriate variables:
758 abs_cont_names.push_back(tagname);
759 abs_cont_models.push_back(model);
760 abs_cont_parameters.push_back(userparameters);
761}
762
763/* Workspace method: Doxygen documentation will be auto-generated */
765 StokesVector& nlte_source,
766 ArrayOfStokesVector& dnlte_source_dx,
767 // WS Input:
768 const ArrayOfMatrix& src_coef_per_species,
769 const ArrayOfMatrix& dsrc_coef_dx,
770 const ArrayOfRetrievalQuantity& jacobian_quantities,
771 const Vector& f_grid,
772 const Numeric& rtp_temperature,
773 const Verbosity&) {
774 // nlte_source has format
775 // [ abs_species, f_grid, stokes_dim ].
776 // src_coef_per_species has format ArrayOfMatrix (over species),
777 // where for each species the matrix has format [f_grid, abs_p].
778
779 Index n_species = src_coef_per_species.nelem(); // # species
780
781 ARTS_USER_ERROR_IF (not n_species,
782 "Must have at least one species.")
783
784 Index n_f = src_coef_per_species[0].nrows(); // # frequencies
785
786 // # pressures must be 1:
787 ARTS_USER_ERROR_IF (1 not_eq src_coef_per_species[0].ncols(),
788 "Must have exactly one pressure.")
789
790 // Check frequency dimension of propmat_clearsky
791 ARTS_USER_ERROR_IF (nlte_source.NumberOfFrequencies() not_eq n_f,
792 "Frequency dimension of nlte_source does not\n"
793 "match abs_coef_per_species.")
794
795 const Vector B = planck(f_grid, rtp_temperature);
796
797 StokesVector sv(n_f, nlte_source.StokesDimensions());
798 for (Index si = 0; si < n_species; ++si) {
799 sv.Kjj() = src_coef_per_species[si](joker, 0);
800 sv *= B;
801 nlte_source.Kjj() += sv.Kjj();
802 }
803
804 // Jacobian
805 for (Index ii = 0; ii < jacobian_quantities.nelem(); ii++) {
806 const auto& deriv = jacobian_quantities[ii];
807
808 if (not deriv.propmattype()) continue;
809
810 if (deriv == Jacobian::Atm::Temperature) {
811 const Vector dB = dplanck_dt(f_grid, rtp_temperature);
812
813 for (Index si = 0; si < n_species; ++si) {
814 sv.Kjj() = src_coef_per_species[si](joker, 0);
815 sv *= dB;
816 dnlte_source_dx[ii].Kjj() += sv.Kjj();
817 }
818
819 sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
820 sv *= B;
821 dnlte_source_dx[ii].Kjj() += sv.Kjj();
822 } else if (is_frequency_parameter(deriv)) {
823 const Vector dB = dplanck_df(f_grid, rtp_temperature);
824
825 for (Index si = 0; si < n_species; ++si) {
826 sv.Kjj() = src_coef_per_species[si](joker, 0);
827 sv *= dB;
828 dnlte_source_dx[ii].Kjj() += sv.Kjj();
829 }
830
831 sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
832 sv *= B;
833 dnlte_source_dx[ii].Kjj() += sv.Kjj();
834 } else {
835 sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
836 sv *= B;
837 dnlte_source_dx[ii].Kjj() += sv.Kjj();
838 }
839 }
840}
841
842/* Workspace method: Doxygen documentation will be auto-generated */
844 PropagationMatrix& propmat_clearsky,
845 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
846 // WS Input:
847 const ArrayOfMatrix& abs_coef_per_species,
848 const ArrayOfMatrix& dabs_coef_dx,
849 const ArrayOfRetrievalQuantity& jacobian_quantities,
850 const ArrayOfArrayOfSpeciesTag& abs_species) {
851 // propmat_clearsky has format
852 // [ abs_species, f_grid, stokes_dim, stokes_dim ].
853 // abs_coef_per_species has format ArrayOfMatrix (over species),
854 // where for each species the matrix has format [f_grid, abs_p].
855
856 Index n_species = abs_coef_per_species.nelem(); // # species
857
858 ARTS_USER_ERROR_IF (0 == n_species,
859 "Must have at least one species.")
860
861 Index n_f = abs_coef_per_species[0].nrows(); // # frequencies
862
863 // # pressures must be 1:
864 ARTS_USER_ERROR_IF (1 not_eq abs_coef_per_species[0].ncols(),
865 "Must have exactly one pressure.")
866
867 // Check frequency dimension of propmat_clearsky
868 ARTS_USER_ERROR_IF (propmat_clearsky.NumberOfFrequencies() not_eq n_f,
869 "Frequency dimension of propmat_clearsky does not\n"
870 "match abs_coef_per_species.")
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 < dabs_coef_dx.nelem(); iqn++) {
877 bool found_special = false;
878 for (Index isp=0; isp<abs_species.nelem(); isp++) {
879 if (jacobian_quantities[iqn] == abs_species[isp]) {
880 dpropmat_clearsky_dx[iqn].Kjj() += abs_coef_per_species[isp](joker, 0);
881 found_special = true;
882 break;
883 }
884 }
885
886 if (not found_special) {
887 if (dabs_coef_dx[iqn].nrows() == n_f) {
888 if (dabs_coef_dx[iqn].ncols() == 1) {
889 dpropmat_clearsky_dx[iqn].Kjj() += dabs_coef_dx[iqn](joker, 0);
890 } else {
891 ARTS_USER_ERROR("Must have exactly one pressure.");
892 }
893 }
894 }
895 }
896}
897
898/* Workspace method: Doxygen documentation will be auto-generated */
899void propmat_clearskyInit( //WS Output
900 PropagationMatrix& propmat_clearsky,
901 StokesVector& nlte_source,
902 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
903 ArrayOfStokesVector& dnlte_source_dx,
904 //WS Input
905 const ArrayOfRetrievalQuantity& jacobian_quantities,
906 const Vector& f_grid,
907 const Index& stokes_dim,
908 const Index& propmat_clearsky_agenda_checked,
909 const Verbosity&) {
910 const Index nf = f_grid.nelem();
911 const Index nq = jacobian_quantities.nelem();
912
913 ARTS_USER_ERROR_IF (!propmat_clearsky_agenda_checked,
914 "You must call *propmat_clearsky_agenda_checkedCalc* before calling this method.")
915
916 ARTS_USER_ERROR_IF (not nf, "No frequencies");
917
918 ARTS_USER_ERROR_IF (stokes_dim < 1 or stokes_dim > 4, "stokes_dim not in [1, 2, 3, 4]");
919
920 // Set size of propmat_clearsky or reset it's values
921 if (propmat_clearsky.StokesDimensions() == stokes_dim and propmat_clearsky.NumberOfFrequencies() == nf) {
922 propmat_clearsky.SetZero();
923 } else {
924 propmat_clearsky = PropagationMatrix(nf, stokes_dim);
925 }
926
927 // Set size of dpropmat_clearsky_dx or reset it's values
928 if (dpropmat_clearsky_dx.nelem() not_eq nq) {
929 dpropmat_clearsky_dx = ArrayOfPropagationMatrix(nq, PropagationMatrix(nf, stokes_dim));
930 } else {
931 for (auto& pm: dpropmat_clearsky_dx) {
932 if (pm.StokesDimensions() == stokes_dim and pm.NumberOfFrequencies() == nf) {
933 pm.SetZero();
934 } else {
935 pm = PropagationMatrix(nf, stokes_dim);
936 }
937 }
938 }
939
940 // Set size of nlte_source or reset it's values
941 if (nlte_source.StokesDimensions() == stokes_dim and nlte_source.NumberOfFrequencies() == nf) {
942 nlte_source.SetZero();
943 } else {
944 nlte_source = StokesVector(nf, stokes_dim);
945 }
946
947 // Set size of dnlte_source_dx or reset it's values
948 if (dnlte_source_dx.nelem() not_eq nq) {
949 dnlte_source_dx = ArrayOfStokesVector(nq, StokesVector(nf, stokes_dim));
950 } else {
951 for (auto& pm: dnlte_source_dx) {
952 if (pm.StokesDimensions() == stokes_dim and pm.NumberOfFrequencies() == nf) {
953 pm.SetZero();
954 } else {
955 pm = StokesVector(nf, stokes_dim);
956 }
957 }
958 }
959}
960
961/* Workspace method: Doxygen documentation will be auto-generated */
963 PropagationMatrix& propmat_clearsky,
964 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
965 const Index& stokes_dim,
966 const Index& atmosphere_dim,
967 const Vector& f_grid,
968 const ArrayOfArrayOfSpeciesTag& abs_species,
969 const ArrayOfRetrievalQuantity& jacobian_quantities,
970 const Vector& rtp_vmr,
971 const Vector& rtp_los,
972 const Vector& rtp_mag,
973 const Verbosity&) {
974 // All the physical constants joined into one static constant:
975 // (abs as e defined as negative)
976 static const Numeric FRconst =
980
981 ARTS_USER_ERROR_IF (stokes_dim < 3,
982 "To include Faraday rotation, stokes_dim >= 3 is required.")
983 ARTS_USER_ERROR_IF (atmosphere_dim == 1 && rtp_los.nelem() < 1,
984 "For applying propmat_clearskyAddFaraday, los needs to be specified\n"
985 "(at least zenith angle component for atmosphere_dim==1),\n"
986 "but it is not.\n")
987 ARTS_USER_ERROR_IF (atmosphere_dim > 1 && rtp_los.nelem() < 2,
988 "For applying propmat_clearskyAddFaraday, los needs to be specified\n"
989 "(both zenith and azimuth angle components for atmosphere_dim>1),\n"
990 "but it is not.\n")
991
992 const bool do_magn_jac = do_magnetic_jacobian(jacobian_quantities);
993 const Numeric dmag = magnetic_field_perturbation(jacobian_quantities);
994
995 Index ife = -1;
996 for (Index sp = 0; sp < abs_species.nelem() && ife < 0; sp++) {
997 if (abs_species[sp].FreeElectrons()) {
998 ife = sp;
999 }
1000 }
1001
1002 ARTS_USER_ERROR_IF (ife < 0,
1003 "Free electrons not found in *abs_species* and "
1004 "Faraday rotation can not be calculated.");
1005
1006 const Numeric ne = rtp_vmr[ife];
1007
1008 if (ne != 0 && (rtp_mag[0] != 0 || rtp_mag[1] != 0 || rtp_mag[2] != 0)) {
1009 // Include remaining terms, beside /f^2
1010 const Numeric c1 =
1011 2 * FRconst *
1013 rtp_los, rtp_mag[0], rtp_mag[1], rtp_mag[2], atmosphere_dim);
1014
1015 Numeric dc1_u = 0.0, dc1_v = 0.0, dc1_w = 0.0;
1016 if (do_magn_jac) {
1017 dc1_u = (2 * FRconst *
1018 dotprod_with_los(rtp_los,
1019 rtp_mag[0] + dmag,
1020 rtp_mag[1],
1021 rtp_mag[2],
1022 atmosphere_dim) -
1023 c1) /
1024 dmag;
1025 dc1_v = (2 * FRconst *
1026 dotprod_with_los(rtp_los,
1027 rtp_mag[0],
1028 rtp_mag[1] + dmag,
1029 rtp_mag[2],
1030 atmosphere_dim) -
1031 c1) /
1032 dmag;
1033 dc1_w = (2 * FRconst *
1034 dotprod_with_los(rtp_los,
1035 rtp_mag[0],
1036 rtp_mag[1],
1037 rtp_mag[2] + dmag,
1038 atmosphere_dim) -
1039 c1) /
1040 dmag;
1041 }
1042
1043 for (Index iv = 0; iv < f_grid.nelem(); iv++) {
1044 const Numeric f2 = f_grid[iv] * f_grid[iv];
1045 const Numeric r = ne * c1 / f2;
1046 propmat_clearsky.AddFaraday(r, iv);
1047
1048 // The Jacobian loop
1049 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
1050 if (is_frequency_parameter(jacobian_quantities[iq]))
1051 dpropmat_clearsky_dx[iq].AddFaraday(-2.0 * ne * r / f_grid[iv], iv);
1052 else if (jacobian_quantities[iq] == Jacobian::Atm::MagneticU)
1053 dpropmat_clearsky_dx[iq].AddFaraday(ne * dc1_u / f2, iv);
1054 else if (jacobian_quantities[iq] == Jacobian::Atm::MagneticV)
1055 dpropmat_clearsky_dx[iq].AddFaraday(ne * dc1_v / f2, iv);
1056 else if (jacobian_quantities[iq] == Jacobian::Atm::MagneticW)
1057 dpropmat_clearsky_dx[iq].AddFaraday(ne * dc1_w / f2, iv);
1058 else if (jacobian_quantities[iq] == Jacobian::Atm::Electrons)
1059 dpropmat_clearsky_dx[iq].AddFaraday(r, iv);
1060 else if (jacobian_quantities[iq] == abs_species[ife])
1061 dpropmat_clearsky_dx[iq].AddFaraday(r, iv);
1062 }
1063 }
1064 }
1065}
1066
1067/* Workspace method: Doxygen documentation will be auto-generated */
1069 // WS Output:
1070 PropagationMatrix& propmat_clearsky,
1071 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
1072 // WS Input:
1073 const Index& stokes_dim,
1074 const Index& atmosphere_dim,
1075 const Vector& f_grid,
1076 const ArrayOfArrayOfSpeciesTag& abs_species,
1077 const ArrayOfRetrievalQuantity& jacobian_quantities,
1078 const Vector& rtp_vmr,
1079 const Vector& rtp_los,
1080 const Numeric& rtp_temperature,
1081 const ArrayOfArrayOfSingleScatteringData& scat_data,
1082 const Index& scat_data_checked,
1083 const Index& use_abs_as_ext,
1084 // Verbosity object:
1085 const Verbosity& verbosity) {
1087
1088 // (i)yCalc only checks scat_data_checked if cloudbox is on. It is off here,
1089 // though, i.e. we need to check it here explicitly. (Also, cloudboxOff sets
1090 // scat_data_checked=0 as it does not check it and as we ususally don't need
1091 // scat_data for clearsky cases, hence don't want to check them by
1092 // scat_data_checkedCalc in that case. This approach seems to be the more
1093 // handy compared to cloudboxOff setting scat_data_checked=1 without checking
1094 // it assuming we won't use it anyways.)
1095 ARTS_USER_ERROR_IF (scat_data_checked != 1,
1096 "The scat_data must be flagged to have "
1097 "passed a consistency check (scat_data_checked=1).")
1098
1099 const Index ns = TotalNumberOfElements(scat_data);
1100 Index np = 0;
1101 for (Index sp = 0; sp < abs_species.nelem(); sp++) {
1102 if (abs_species[sp].Particles()) {
1103 np++;
1104 }
1105 }
1106
1107 ARTS_USER_ERROR_IF (np == 0,
1108 "For applying propmat_clearskyAddParticles, *abs_species* needs to"
1109 "contain species 'particles', but it does not.\n")
1110
1111 ARTS_USER_ERROR_IF (ns != np,
1112 "Number of 'particles' entries in abs_species and of elements in\n"
1113 "*scat_data* needs to be identical. But you have " , np,
1114 " 'particles' entries\n"
1115 "and ", ns, " *scat_data* elements.\n")
1116
1117 ARTS_USER_ERROR_IF (atmosphere_dim == 1 && rtp_los.nelem() < 1,
1118 "For applying *propmat_clearskyAddParticles*, *rtp_los* needs to be specified\n"
1119 "(at least zenith angle component for atmosphere_dim==1),\n"
1120 "but it is not.\n")
1121 ARTS_USER_ERROR_IF (atmosphere_dim > 1 && rtp_los.nelem() < 2,
1122 "For applying *propmat_clearskyAddParticles*, *rtp_los* needs to be specified\n"
1123 "(both zenith and azimuth angle components for atmosphere_dim>1),\n"
1124 "but it is not.\n")
1125
1126 // Use for rescaling vmr of particulates
1127 Numeric rtp_vmr_sum = 0.0;
1128
1129 // Tests and setup partial derivatives
1130 const bool do_jac_temperature = do_temperature_jacobian(jacobian_quantities);
1131 const bool do_jac_frequencies = do_frequency_jacobian(jacobian_quantities);
1132 const Numeric dT = temperature_perturbation(jacobian_quantities);
1133
1134 const Index na = abs_species.nelem();
1135 Vector rtp_los_back;
1136 mirror_los(rtp_los_back, rtp_los, atmosphere_dim);
1137
1138 // 170918 JM: along with transition to use of new-type (aka
1139 // pre-f_grid-interpolated) scat_data, freq perturbation switched off. Typical
1140 // clear-sky freq perturbations yield insignificant effects in particle
1141 // properties. Hence, this feature is neglected here.
1142 if (do_jac_frequencies) {
1143 out1 << "WARNING:\n"
1144 << "Frequency perturbation not available for absorbing particles.\n";
1145 }
1146
1147 // creating temporary output containers
1148 ArrayOfArrayOfTensor5 ext_mat_Nse;
1149 ArrayOfArrayOfTensor4 abs_vec_Nse;
1150 ArrayOfArrayOfIndex ptypes_Nse;
1151 Matrix t_ok;
1152
1153 // preparing input in format needed
1154 Vector T_array;
1155 if (do_jac_temperature) {
1156 T_array.resize(2);
1157 T_array = rtp_temperature;
1158 T_array[1] += dT;
1159 } else {
1160 T_array.resize(1);
1161 T_array = rtp_temperature;
1162 }
1163 Matrix dir_array(1, 2);
1164 dir_array(0, joker) = rtp_los_back;
1165
1166 // ext/abs per scat element for all freqs at once
1167 opt_prop_NScatElems(ext_mat_Nse,
1168 abs_vec_Nse,
1169 ptypes_Nse,
1170 t_ok,
1171 scat_data,
1172 stokes_dim,
1173 T_array,
1174 dir_array,
1175 -1);
1176
1177 const Index nf = abs_vec_Nse[0][0].nbooks();
1178 Tensor3 tmp(nf, stokes_dim, stokes_dim);
1179
1180 // Internal computations necessary since it relies on zero start
1181 PropagationMatrix internal_propmat(propmat_clearsky.NumberOfFrequencies(), propmat_clearsky.StokesDimensions());
1182
1183 // loop over the scat_data and link them with correct vmr_field entry according
1184 // to the position of the particle type entries in abs_species.
1185 Index sp = 0;
1186 Index i_se_flat = 0;
1187 for (Index i_ss = 0; i_ss < scat_data.nelem(); i_ss++) {
1188 for (Index i_se = 0; i_se < scat_data[i_ss].nelem(); i_se++) {
1189 // forward to next particle entry in abs_species
1190 while (sp < na && not abs_species[sp].Particles())
1191 sp++;
1192 internal_propmat.SetZero();
1193
1194 // running beyond number of abs_species entries when looking for next
1195 // particle entry. shouldn't happen, though.
1196 ARTS_ASSERT(sp < na);
1197 ARTS_USER_ERROR_IF (rtp_vmr[sp] < 0.,
1198 "Negative absorbing particle 'vmr' (aka number density)"
1199 " encountered:\n"
1200 "scat species #", i_ss, ", scat elem #", i_se,
1201 " (vmr_field entry #", sp, ")\n")
1202
1203 if (rtp_vmr[sp] > 0.) {
1204 ARTS_USER_ERROR_IF (t_ok(i_se_flat, 0) < 0.,
1205 "Temperature interpolation error:\n"
1206 "scat species #", i_ss, ", scat elem #", i_se, "\n")
1207 if (use_abs_as_ext) {
1208 if (nf > 1)
1209 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1210 internal_propmat.AddAbsorptionVectorAtPosition(
1211 abs_vec_Nse[i_ss][i_se](iv, 0, 0, joker), iv);
1212 else
1213 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1214 internal_propmat.AddAbsorptionVectorAtPosition(
1215 abs_vec_Nse[i_ss][i_se](0, 0, 0, joker), iv);
1216 } else {
1217 if (nf > 1)
1218 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1219 internal_propmat.SetAtPosition(
1220 ext_mat_Nse[i_ss][i_se](iv, 0, 0, joker, joker), iv);
1221 else
1222 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1223 internal_propmat.SetAtPosition(
1224 ext_mat_Nse[i_ss][i_se](0, 0, 0, joker, joker), iv);
1225 }
1226 propmat_clearsky += rtp_vmr[sp] * internal_propmat;
1227 }
1228
1229 // For temperature derivatives (so we don't need to check it in jac loop)
1230 if (do_jac_temperature) {
1231 ARTS_USER_ERROR_IF (t_ok(i_se_flat, 1) < 0.,
1232 "Temperature interpolation error (in perturbation):\n"
1233 "scat species #", i_ss, ", scat elem #", i_se, "\n")
1234 }
1235
1236 // For number density derivatives
1237 if (jacobian_quantities.nelem()) rtp_vmr_sum += rtp_vmr[sp];
1238
1239 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
1240 const auto& deriv = jacobian_quantities[iq];
1241
1242 if (not deriv.propmattype()) continue;
1243
1244 if (deriv == Jacobian::Atm::Temperature) {
1245 if (use_abs_as_ext) {
1246 tmp(joker, joker, 0) =
1247 abs_vec_Nse[i_ss][i_se](joker, 1, 0, joker);
1248 tmp(joker, joker, 0) -=
1249 abs_vec_Nse[i_ss][i_se](joker, 0, 0, joker);
1250 } else {
1251 tmp = ext_mat_Nse[i_ss][i_se](joker, 1, 0, joker, joker);
1252 tmp -= ext_mat_Nse[i_ss][i_se](joker, 0, 0, joker, joker);
1253 }
1254
1255 tmp *= rtp_vmr[sp];
1256 tmp /= dT;
1257
1258 if (nf > 1)
1259 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1260 if (use_abs_as_ext)
1261 dpropmat_clearsky_dx[iq].AddAbsorptionVectorAtPosition(
1262 tmp(iv, joker, 0), iv);
1263 else
1264 dpropmat_clearsky_dx[iq].AddAtPosition(tmp(iv, joker, joker),
1265 iv);
1266 else
1267 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1268 if (use_abs_as_ext)
1269 dpropmat_clearsky_dx[iq].AddAbsorptionVectorAtPosition(
1270 tmp(0, joker, 0), iv);
1271 else
1272 dpropmat_clearsky_dx[iq].AddAtPosition(tmp(0, joker, joker),
1273 iv);
1274 }
1275
1276 else if (deriv == Jacobian::Atm::Particulates) {
1277 for (Index iv = 0; iv < f_grid.nelem(); iv++)
1278 dpropmat_clearsky_dx[iq].AddAtPosition(internal_propmat, iv);
1279 }
1280
1281 else if (deriv == abs_species[sp]) {
1282 dpropmat_clearsky_dx[iq] += internal_propmat;
1283 }
1284 }
1285 sp++;
1286 i_se_flat++;
1287 }
1288 }
1289
1290 //checking that no further 'particle' entry left after all scat_data entries
1291 //are processes. this is basically not necessary. but checking it anyway to
1292 //really be safe. remove later, when more extensively tested.
1293 while (sp < na) {
1294 ARTS_ASSERT(abs_species[sp][0].Type() != Species::TagType::Particles);
1295 sp++;
1296 }
1297
1298 if (rtp_vmr_sum != 0.0) {
1299 for (Index iq = 0; iq < jacobian_quantities.nelem(); iq++) {
1300 const auto& deriv = jacobian_quantities[iq];
1301
1302 if (not deriv.propmattype()) continue;
1303
1304 if (deriv == Jacobian::Atm::Particulates) {
1305 dpropmat_clearsky_dx[iq] /= rtp_vmr_sum;
1306 }
1307 }
1308 }
1309}
1310
1311
1313 const Vector& f_grid,
1314 const Numeric& sparse_df,
1315 const String& speedup_option,
1316 // Verbosity object:
1317 const Verbosity&)
1318{
1319 // Return empty for nothing
1320 if (not f_grid.nelem()) {
1321 sparse_f_grid.resize(0);
1322 return;
1323 };
1324
1325 switch (Options::toLblSpeedupOrThrow(speedup_option)) {
1326 case Options::LblSpeedup::LinearIndependent:
1327 sparse_f_grid = LineShape::linear_sparse_f_grid(f_grid, sparse_df);
1329 break;
1331 sparse_f_grid = LineShape::triple_sparse_f_grid(f_grid, sparse_df);
1332 break;
1334 sparse_f_grid.resize(0);
1335 break;
1336 case Options::LblSpeedup::FINAL: { /* Leave last */ }
1337 }
1338}
1339
1341 const Numeric& sparse_df,
1342 const String& speedup_option,
1343 // Verbosity object:
1344 const Verbosity& verbosity)
1345{
1346 Vector sparse_f_grid;
1347 sparse_f_gridFromFrequencyGrid(sparse_f_grid, f_grid, sparse_df, speedup_option, verbosity);
1348 return sparse_f_grid;
1349}
1350
1351/* Workspace method: Doxygen documentation will be auto-generated */
1352void propmat_clearskyAddLines( // Workspace reference:
1353 // WS Output:
1354 PropagationMatrix& propmat_clearsky,
1355 StokesVector& nlte_source,
1356 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
1357 ArrayOfStokesVector& dnlte_source_dx,
1358 // WS Input:
1359 const Vector& f_grid,
1360 const ArrayOfArrayOfSpeciesTag& abs_species,
1361 const ArrayOfRetrievalQuantity& jacobian_quantities,
1362 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1364 const Numeric& rtp_pressure,
1365 const Numeric& rtp_temperature,
1366 const EnergyLevelMap& rtp_nlte,
1367 const Vector& rtp_vmr,
1368 const Index& nlte_do,
1369 const Index& lbl_checked,
1370 // WS User Generic inputs
1371 const Numeric& sparse_df,
1372 const Numeric& sparse_lim,
1373 const String& speedup_option,
1374 const ArrayOfSpeciesTag& select_speciestags,
1375 // Verbosity object:
1376 const Verbosity& verbosity) {
1377
1378 // Size of problem
1379 const Index nf = f_grid.nelem();
1380 const Index nq = jacobian_quantities.nelem();
1381 const Index ns = abs_species.nelem();
1382
1383 // Possible things that can go wrong in this code (excluding line parameters)
1384 ARTS_USER_ERROR_IF(not lbl_checked, "Must check LBL calculations")
1385 check_abs_species(abs_species);
1386 ARTS_USER_ERROR_IF(rtp_vmr.nelem() not_eq abs_species.nelem(),
1387 "*rtp_vmr* must match *abs_species*")
1388 ARTS_USER_ERROR_IF(propmat_clearsky.NumberOfFrequencies() not_eq nf,
1389 "*f_grid* must match *propmat_clearsky*")
1390 ARTS_USER_ERROR_IF(nlte_source.NumberOfFrequencies() not_eq nf,
1391 "*f_grid* must match *nlte_source*")
1392 ARTS_USER_ERROR_IF(not nq and (nq not_eq dpropmat_clearsky_dx.nelem()),
1393 "*dpropmat_clearsky_dx* must match derived form of *jacobian_quantities*")
1394 ARTS_USER_ERROR_IF(not nq and bad_propmat(dpropmat_clearsky_dx, f_grid),
1395 "*dpropmat_clearsky_dx* must have frequency dim same as *f_grid*")
1396 ARTS_USER_ERROR_IF(nlte_do and (nq not_eq dnlte_source_dx.nelem()),
1397 "*dnlte_source_dx* must match derived form of *jacobian_quantities* when non-LTE is on")
1398 ARTS_USER_ERROR_IF(nlte_do and bad_propmat(dnlte_source_dx, f_grid),
1399 "*dnlte_source_dx* must have frequency dim same as *f_grid* when non-LTE is on")
1400 ARTS_USER_ERROR_IF(any_negative(f_grid), "Negative frequency (at least one value).")
1401 ARTS_USER_ERROR_IF(not is_increasing(f_grid), "Must be sorted and increasing.")
1402 ARTS_USER_ERROR_IF(any_negative(rtp_vmr), "Negative VMR (at least one value).")
1403 ARTS_USER_ERROR_IF(any_negative(rtp_nlte.Data()), "Negative NLTE (at least one value).")
1404 ARTS_USER_ERROR_IF(rtp_temperature <= 0, "Non-positive temperature")
1405 ARTS_USER_ERROR_IF(rtp_pressure <= 0, "Non-positive pressure")
1406 ARTS_USER_ERROR_IF(sparse_lim > 0 and sparse_df > sparse_lim,
1407 "If sparse grids are to be used, the limit must be larger than the grid-spacing.\n"
1408 "The limit is ", sparse_lim, " Hz and the grid_spacing is ", sparse_df, " Hz")
1409
1410 if (not nf) return;
1411
1412 // Deal with sparse computational grid
1413 const Vector f_grid_sparse = create_sparse_f_grid_internal(f_grid, sparse_df, speedup_option, verbosity);
1414 const Options::LblSpeedup speedup_type = f_grid_sparse.nelem() ? Options::toLblSpeedupOrThrow(speedup_option) : Options::LblSpeedup::None;
1415 ARTS_USER_ERROR_IF(sparse_lim <= 0 and speedup_type not_eq Options::LblSpeedup::None,
1416 "Must have a sparse limit if you set speedup_option")
1417
1418 // Calculations data
1419 LineShape::ComputeData com(f_grid, jacobian_quantities, nlte_do);
1420 LineShape::ComputeData sparse_com(f_grid_sparse, jacobian_quantities, nlte_do);
1421
1422 for (Index ispecies = 0; ispecies < ns; ispecies++) {
1423 if (select_speciestags.nelem() and select_speciestags not_eq abs_species[ispecies]) continue;
1424
1425 // Skip it if there are no species or there is Zeeman requested
1426 if (not abs_species[ispecies].nelem() or abs_species[ispecies].Zeeman() or not abs_lines_per_species[ispecies].nelem())
1427 continue;
1428
1429 for (auto& band : abs_lines_per_species[ispecies]) {
1430 LineShape::compute(com, sparse_com, band, jacobian_quantities, rtp_nlte, band.BroadeningSpeciesVMR(rtp_vmr, abs_species), abs_species[ispecies], rtp_vmr[ispecies],
1431 isotopologue_ratios[band.Isotopologue()], rtp_pressure, rtp_temperature, 0, sparse_lim,
1432 false, Zeeman::Polarization::Pi, speedup_type);
1433
1434 }
1435 }
1436
1437 switch (speedup_type) {
1438 case Options::LblSpeedup::LinearIndependent: com.interp_add_even(sparse_com); break;
1440 case Options::LblSpeedup::None: /* Do nothing */ break;
1441 case Options::LblSpeedup::FINAL: { /* Leave last */ }
1442 }
1443
1444 // Sum up the propagation matrix
1445 propmat_clearsky.Kjj() += com.F.real();
1446
1447 // Sum up the Jacobian
1448 for (Index j=0; j<nq; j++) {
1449 if (not jacobian_quantities[j].propmattype()) continue;
1450 dpropmat_clearsky_dx[j].Kjj() += com.dF.real()(joker, j);
1451 }
1452
1453 if (nlte_do) {
1454 // Sum up the source vector
1455 nlte_source.Kjj() += com.N.real();
1456
1457 // Sum up the Jacobian
1458 for (Index j=0; j<nq; j++) {
1459 if (not jacobian_quantities[j].propmattype()) continue;
1460 dnlte_source_dx[j].Kjj() += com.dN.real()(joker, j);
1461 }
1462 }
1463}
1464
1465/* Workspace method: Doxygen documentation will be auto-generated */
1466void propmat_clearskyAddXsecAgenda( // Workspace reference:
1467 Workspace& ws,
1468 // WS Output:
1469 PropagationMatrix& propmat_clearsky,
1470 ArrayOfPropagationMatrix& dpropmat_clearsky_dx,
1471 // WS Input:
1472 const Vector& f_grid,
1473 const ArrayOfArrayOfSpeciesTag& abs_species,
1474 const ArrayOfRetrievalQuantity& jacobian_quantities,
1475 const Numeric& rtp_pressure,
1476 const Numeric& rtp_temperature,
1477 const Vector& rtp_vmr,
1478 const Agenda& abs_xsec_agenda,
1479 // Verbosity object:
1480 const Verbosity& verbosity) {
1481
1482 // Output of AbsInputFromRteScalars:
1483 Vector abs_p;
1484 Vector abs_t;
1485 Matrix abs_vmrs;
1486 // Output of abs_h2oSet:
1487 Vector abs_h2o;
1488 // Output of abs_coefCalc:
1489 Matrix abs_coef;
1490 ArrayOfMatrix abs_coef_per_species, dabs_coef_dx;
1491
1493 abs_t,
1494 abs_vmrs,
1495 rtp_pressure,
1496 rtp_temperature,
1497 rtp_vmr,
1498 verbosity);
1499
1500 // Absorption cross sections per tag group.
1501 ArrayOfMatrix abs_xsec_per_species;
1502 ArrayOfArrayOfMatrix dabs_xsec_per_species_dx;
1503
1504 // Make all species active:
1505 ArrayOfIndex abs_species_active(abs_species.nelem());
1506 for (Index i = 0; i < abs_species.nelem(); ++i) abs_species_active[i] = i;
1507
1508 // Call agenda to calculate absorption:
1510 abs_xsec_per_species,
1511 dabs_xsec_per_species_dx,
1512 abs_species,
1513 jacobian_quantities,
1514 abs_species_active,
1515 f_grid,
1516 abs_p,
1517 abs_t,
1518 abs_vmrs,
1519 abs_xsec_agenda);
1520
1521 // Calculate absorption coefficients from cross sections:
1522 abs_coefCalcFromXsec(abs_coef,
1523 dabs_coef_dx,
1524 abs_coef_per_species,
1525 abs_xsec_per_species,
1526 dabs_xsec_per_species_dx,
1527 abs_species,
1528 jacobian_quantities,
1529 abs_vmrs,
1530 abs_p,
1531 abs_t,
1532 verbosity);
1533
1534 // Now add abs_coef_per_species to propmat_clearsky:
1536 dpropmat_clearsky_dx,
1537 abs_coef_per_species,
1538 dabs_coef_dx,
1539 jacobian_quantities, abs_species);
1540}
1541
1542/* Workspace method: Doxygen documentation will be auto-generated */
1544 const Vector& f_grid,
1545 const Index& stokes_dim,
1546 const Verbosity&) {
1547 propmat_clearsky = PropagationMatrix(f_grid.nelem(), stokes_dim);
1548}
1549
1550/* Workspace method: Doxygen documentation will be auto-generated */
1552 PropagationMatrix& propmat_clearsky, const Verbosity&) {
1553 for (Index i = 0; i < propmat_clearsky.NumberOfFrequencies(); i++)
1554 if (propmat_clearsky.Kjj()[i] < 0.0) propmat_clearsky.SetAtPosition(0.0, i);
1555}
1556
1557/* Workspace method: Doxygen documentation will be auto-generated */
1559 const Verbosity&) {
1561}
1562
1563/* Workspace method: Doxygen documentation will be auto-generated */
1565 const String& option,
1566 const Verbosity&) {
1567 isotopologue_ratios = Hitran::isotopologue_ratios(Hitran::toTypeOrThrow(option));
1568}
1569
1570#ifdef ENABLE_NETCDF
1571/* Workspace method: Doxygen documentation will be auto-generated */
1572/* Included by Claudia Emde, 20100707 */
1573void WriteMolTau( //WS Input
1574 const Vector& f_grid,
1575 const Tensor3& z_field,
1576 const Tensor7& propmat_clearsky_field,
1577 const Index& atmosphere_dim,
1578 //Keyword
1579 const String& filename,
1580 const Verbosity&) {
1581 int retval, ncid;
1582 int nlev_dimid, nlyr_dimid, nwvl_dimid, stokes_dimid, none_dimid;
1583 int dimids[4];
1584 int wvlmin_varid, wvlmax_varid, z_varid, wvl_varid, tau_varid;
1585
1586 ARTS_USER_ERROR_IF (atmosphere_dim != 1,
1587 "WriteMolTau can only be used for atmosphere_dim=1")
1588
1589#pragma omp critical(netcdf__critical_region)
1590 {
1591 // Open file
1592 if ((retval = nc_create(filename.c_str(), NC_CLOBBER, &ncid)))
1593 nca_error(retval, "nc_create");
1594
1595 // Define dimensions
1596 if ((retval = nc_def_dim(ncid, "nlev", (int)z_field.npages(), &nlev_dimid)))
1597 nca_error(retval, "nc_def_dim");
1598
1599 if ((retval =
1600 nc_def_dim(ncid, "nlyr", (int)z_field.npages() - 1, &nlyr_dimid)))
1601 nca_error(retval, "nc_def_dim");
1602
1603 if ((retval = nc_def_dim(ncid, "nwvl", (int)f_grid.nelem(), &nwvl_dimid)))
1604 nca_error(retval, "nc_def_dim");
1605
1606 if ((retval = nc_def_dim(ncid, "none", 1, &none_dimid)))
1607 nca_error(retval, "nc_def_dim");
1608
1609 if ((retval = nc_def_dim(ncid,
1610 "nstk",
1611 (int)propmat_clearsky_field.nbooks(),
1612 &stokes_dimid)))
1613 nca_error(retval, "nc_def_dim");
1614
1615 // Define variables
1616 if ((retval = nc_def_var(
1617 ncid, "wvlmin", NC_DOUBLE, 1, &none_dimid, &wvlmin_varid)))
1618 nca_error(retval, "nc_def_var wvlmin");
1619
1620 if ((retval = nc_def_var(
1621 ncid, "wvlmax", NC_DOUBLE, 1, &none_dimid, &wvlmax_varid)))
1622 nca_error(retval, "nc_def_var wvlmax");
1623
1624 if ((retval = nc_def_var(ncid, "z", NC_DOUBLE, 1, &nlev_dimid, &z_varid)))
1625 nca_error(retval, "nc_def_var z");
1626
1627 if ((retval =
1628 nc_def_var(ncid, "wvl", NC_DOUBLE, 1, &nwvl_dimid, &wvl_varid)))
1629 nca_error(retval, "nc_def_var wvl");
1630
1631 dimids[0] = nlyr_dimid;
1632 dimids[1] = nwvl_dimid;
1633 dimids[2] = stokes_dimid;
1634 dimids[3] = stokes_dimid;
1635
1636 if ((retval =
1637 nc_def_var(ncid, "tau", NC_DOUBLE, 4, &dimids[0], &tau_varid)))
1638 nca_error(retval, "nc_def_var tau");
1639
1640 // Units
1641 if ((retval = nc_put_att_text(ncid, wvlmin_varid, "units", 2, "nm")))
1642 nca_error(retval, "nc_put_att_text");
1643
1644 if ((retval = nc_put_att_text(ncid, wvlmax_varid, "units", 2, "nm")))
1645 nca_error(retval, "nc_put_att_text");
1646
1647 if ((retval = nc_put_att_text(ncid, z_varid, "units", 2, "km")))
1648 nca_error(retval, "nc_put_att_text");
1649
1650 if ((retval = nc_put_att_text(ncid, wvl_varid, "units", 2, "nm")))
1651 nca_error(retval, "nc_put_att_text");
1652
1653 if ((retval = nc_put_att_text(ncid, tau_varid, "units", 1, "-")))
1654 nca_error(retval, "nc_put_att_text");
1655
1656 // End define mode. This tells netCDF we are done defining
1657 // metadata.
1658 if ((retval = nc_enddef(ncid))) nca_error(retval, "nc_enddef");
1659
1660 // Assign data
1661 double wvlmin[1];
1662 wvlmin[0] = SPEED_OF_LIGHT / f_grid[f_grid.nelem() - 1] * 1e9;
1663 if ((retval = nc_put_var_double(ncid, wvlmin_varid, &wvlmin[0])))
1664 nca_error(retval, "nc_put_var");
1665
1666 double wvlmax[1];
1667 wvlmax[0] = SPEED_OF_LIGHT / f_grid[0] * 1e9;
1668 if ((retval = nc_put_var_double(ncid, wvlmax_varid, &wvlmax[0])))
1669 nca_error(retval, "nc_put_var");
1670
1671 double z[z_field.npages()];
1672 for (int iz = 0; iz < z_field.npages(); iz++)
1673 z[iz] = z_field(z_field.npages() - 1 - iz, 0, 0) * 1e-3;
1674
1675 if ((retval = nc_put_var_double(ncid, z_varid, &z[0])))
1676 nca_error(retval, "nc_put_var");
1677
1678 double wvl[f_grid.nelem()];
1679 for (int iv = 0; iv < f_grid.nelem(); iv++)
1680 wvl[iv] = SPEED_OF_LIGHT / f_grid[f_grid.nelem() - 1 - iv] * 1e9;
1681
1682 if ((retval = nc_put_var_double(ncid, wvl_varid, &wvl[0])))
1683 nca_error(retval, "nc_put_var");
1684
1685 const Index zfnp = z_field.npages() - 1;
1686 const Index fgne = f_grid.nelem();
1687 const Index amfnb = propmat_clearsky_field.nbooks();
1688
1689 Tensor4 tau(zfnp, fgne, amfnb, amfnb, 0.);
1690
1691 // Calculate average tau for layers
1692 for (int is = 0; is < propmat_clearsky_field.nlibraries(); is++)
1693 for (int iz = 0; iz < zfnp; iz++)
1694 for (int iv = 0; iv < fgne; iv++)
1695 for (int is1 = 0; is1 < amfnb; is1++)
1696 for (int is2 = 0; is2 < amfnb; is2++)
1697 // sum up all species
1698 tau(iz, iv, is1, is2) +=
1699 0.5 *
1700 (propmat_clearsky_field(is,
1701 f_grid.nelem() - 1 - iv,
1702 is1,
1703 is2,
1704 z_field.npages() - 1 - iz,
1705 0,
1706 0) +
1707 propmat_clearsky_field(is,
1708 f_grid.nelem() - 1 - iv,
1709 is1,
1710 is2,
1711 z_field.npages() - 2 - iz,
1712 0,
1713 0)) *
1714 (z_field(z_field.npages() - 1 - iz, 0, 0) -
1715 z_field(z_field.npages() - 2 - iz, 0, 0));
1716
1717 if ((retval = nc_put_var_double(ncid, tau_varid, tau.get_c_array())))
1718 nca_error(retval, "nc_put_var");
1719
1720 // Close the file
1721 if ((retval = nc_close(ncid))) nca_error(retval, "nc_close");
1722 }
1723}
1724
1725#else
1726
1727void WriteMolTau( //WS Input
1728 const Vector& f_grid _U_,
1729 const Tensor3& z_field _U_,
1730 const Tensor7& propmat_clearsky_field _U_,
1731 const Index& atmosphere_dim _U_,
1732 //Keyword
1733 const String& filename _U_,
1734 const Verbosity&) {
1735 ARTS_USER_ERROR_IF(true,
1736 "The workspace method WriteMolTau is not available"
1737 "because ARTS was compiled without NetCDF support.");
1738}
1739
1740#endif /* ENABLE_NETCDF */
1741
1742/* Workspace method: Doxygen documentation will be auto-generated */
1744 // WS Output:
1745 ArrayOfMatrix& abs_xsec_per_species,
1747 // WS Input:
1748 const ArrayOfArrayOfSpeciesTag& abs_species,
1749 const ArrayOfRetrievalQuantity& jacobian_quantities,
1750 const ArrayOfIndex& abs_species_active,
1751 const Vector& f_grid,
1752 const Vector& abs_p,
1753 const Vector& abs_t,
1754 const Matrix& abs_vmrs,
1755 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1757 const Index& lbl_checked,
1758 const Verbosity&) {
1759 DEPRECATED_FUNCTION("abs_xsec_per_speciesAddLines", "2021-07-13",
1760 "This function is no longer up to date. It only exists to satisfy "
1761 "lookup table calculations before these are updated.\n"
1762 "Once the lookup table calculations are up-to-date, this function "
1763 "is fully replaced with propmat_clearskyAddLines, with better functionality\n" )
1764
1765 ARTS_USER_ERROR_IF(not lbl_checked, "Must check LBL calculations")
1766 ARTS_USER_ERROR_IF(jacobian_quantities.nelem(), "There's a hard deprecation of derivatives using old style lbl-calculations with derivatives, switch to propmat_clearskyAddLines")
1767 ARTS_USER_ERROR_IF (abs_species.nelem() not_eq abs_xsec_per_species.nelem() or
1768 abs_species.nelem() not_eq abs_vmrs.nrows() or
1769 abs_species.nelem() not_eq abs_lines_per_species.nelem(),
1770 "The following variables must all have the same dimension:\n"
1771 "abs_species: ", abs_species.nelem(), '\n',
1772 "abs_xsec_per_species: ", abs_xsec_per_species.nelem(), '\n',
1773 "abs_vmrs: ", abs_vmrs.nrows(), '\n',
1774 "abs_lines_per_species: ", abs_lines_per_species.nelem(), '\n')
1775 ARTS_USER_ERROR_IF (min(abs_t) < 0,
1776 "Temperature must be at least 0 K. But you request an absorption\n"
1777 "calculation at ", min(abs_t), " K!")
1778
1779 // Size of problem
1780 const Index np = abs_p.nelem();
1781
1782 // Calculations data
1783 LineShape::ComputeData com(f_grid, jacobian_quantities, false);
1784 LineShape::ComputeData sparse_com(Vector(0), jacobian_quantities, false);
1785 constexpr Options::LblSpeedup speedup_type = Options::LblSpeedup::None;
1786 const EnergyLevelMap rtp_nlte;
1787
1788 for (Index ip=0; ip<np; ip++) {
1789 for (Index ispecies: abs_species_active) {
1790 // Skip it if there are no species or there is Zeeman requested
1791 if (not abs_species[ispecies].nelem() or abs_species[ispecies].Zeeman() or not abs_lines_per_species[ispecies].nelem())
1792 continue;
1793
1794 // Reset for legacy VMR jacobian
1795 com.reset();
1796 sparse_com.reset();
1797
1798 for (auto& band : abs_lines_per_species[ispecies]) {
1799 LineShape::compute(com, sparse_com, band, jacobian_quantities, rtp_nlte, band.BroadeningSpeciesVMR(abs_vmrs(joker, ip), abs_species), abs_species[ispecies], abs_vmrs(ispecies, ip),
1800 isotopologue_ratios[band.Isotopologue()], abs_p[ip], abs_t[ip], 0, 0, false, Zeeman::Polarization::Pi, speedup_type);
1801
1802 }
1803
1804 // Sum up the propagation matrix
1805 com.F /= number_density(abs_p[ip], abs_t[ip]) * abs_vmrs(ispecies, ip);
1806 abs_xsec_per_species[ispecies](joker, ip) += com.F.real();
1807 }
1808 }
1809}
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:65
Index TotalNumberOfElements(const Array< Array< base > > &aa)
Determine total number of elements in an ArrayOfArray.
Definition: array.h:343
The global header file for ARTS.
Index nrows
Index ncols
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:22890
void abs_speciesSet(ArrayOfArrayOfSpeciesTag &abs_species, Index &abs_xsec_agenda_checked, Index &propmat_clearsky_agenda_checked, const ArrayOfString &species, const Verbosity &verbosity)
WORKSPACE METHOD: abs_speciesSet.
void chk_size(const String &x_name, ConstVectorView x, const Index &c)
Runtime check for size of Vector.
Definition: check_input.cc:421
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.
The Agenda class.
Definition: agenda_class.h:44
This can be used to make arrays out of anything.
Definition: array.h:107
Index nelem() const ARTS_NOEXCEPT
Number of elements.
Definition: array.h:195
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 ARTS_NOEXCEPT
Returns the number of rows.
Definition: matpackI.cc:433
Index ncols() const ARTS_NOEXCEPT
Returns the number of columns.
Definition: matpackI.cc:436
Index npages() const
Returns the number of pages.
Definition: matpackIII.h:144
Index nlibraries() const
Returns the number of libraries.
Definition: matpackVII.cc:42
Index nbooks() const
Returns the number of books.
Definition: matpackVII.cc:51
Index nelem() const ARTS_NOEXCEPT
Returns the number of elements.
Definition: matpackI.cc:48
const Tensor4 & Data() const noexcept
Energy level type.
The Matrix class.
Definition: matpackI.h:1225
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1056
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:165
Stokes vector is as Propagation matrix but only has 4 possible values.
The Tensor3 class.
Definition: matpackIII.h:339
const Numeric * get_c_array() const ARTS_NOEXCEPT
Conversion to plain C-array.
Definition: matpackIV.cc:354
The Tensor4 class.
Definition: matpackIV.h:421
The Tensor7 class.
Definition: matpackVII.h:2382
The Vector class.
Definition: matpackI.h:876
void resize(Index n)
Resize function.
Definition: matpackI.cc:408
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(...)
Definition: debug.h:150
#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.
int main(int argc, char **argv)
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:1138
bool do_frequency_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants a frequency derivative.
Definition: jacobian.cc:1168
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:1172
bool do_temperature_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants the temperature derivative.
Definition: jacobian.cc:1147
Numeric magnetic_field_perturbation(const ArrayOfRetrievalQuantity &js) noexcept
Returns the magnetic field perturbation if it exists.
Definition: jacobian.cc:1192
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:1176
Numeric frequency_perturbation(const ArrayOfRetrievalQuantity &js) noexcept
Returns the frequency perturbation if it exists.
Definition: jacobian.cc:1184
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:1312
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:202
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:1340
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:899
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:764
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:1466
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:244
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:741
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:725
const Numeric VACUUM_PERMITTIVITY
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 Verbosity &verbosity)
WORKSPACE METHOD: propmat_clearskyAddLines.
Definition: m_abs.cc:1352
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:152
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:1743
void propmat_clearskyForceNegativeToZero(PropagationMatrix &propmat_clearsky, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyForceNegativeToZero.
Definition: m_abs.cc:1551
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:1068
void propmat_clearskyZero(PropagationMatrix &propmat_clearsky, const Vector &f_grid, const Index &stokes_dim, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyZero.
Definition: m_abs.cc:1543
void propmat_clearskyAddFromAbsCoefPerSpecies(PropagationMatrix &propmat_clearsky, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, const ArrayOfMatrix &abs_coef_per_species, const ArrayOfMatrix &dabs_coef_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfArrayOfSpeciesTag &abs_species)
Definition: m_abs.cc:843
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:1573
void isotopologue_ratiosInitFromHitran(SpeciesIsotopologueRatios &isotopologue_ratios, const String &option, const Verbosity &)
WORKSPACE METHOD: isotopologue_ratiosInitFromHitran.
Definition: m_abs.cc:1564
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:223
void isotopologue_ratiosInitFromBuiltin(SpeciesIsotopologueRatios &isotopologue_ratios, const Verbosity &)
WORKSPACE METHOD: isotopologue_ratiosInitFromBuiltin.
Definition: m_abs.cc:1558
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:962
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:135
Implementation of Matrix, Vector, and such stuff.
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
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.
Lines createEmptyCopy(const Lines &al) noexcept
Creates a copy of the input lines structure.
SpeciesIsotopologueRatios isotopologue_ratios(Type type)
Temperature
Definition: jacobian.h:57
bool good_linear_sparse_f_grid(const Vector &f_grid_dense, const Vector &f_grid_sparse) noexcept
Definition: lineshape.cc:3196
Vector triple_sparse_f_grid(const Vector &f_grid, const Numeric &sparse_df) noexcept
Definition: lineshape.cc:3209
Vector linear_sparse_f_grid(const Vector &f_grid, const Numeric &sparse_df) ARTS_NOEXCEPT
Definition: lineshape.cc:3174
void compute(ComputeData &com, ComputeData &sparse_com, const AbsorptionLines &band, const ArrayOfRetrievalQuantity &jacobian_quantities, const EnergyLevelMap &nlte, const Vector &vmrs, const ArrayOfSpeciesTag &self_tag, const Numeric &self_vmr, const Numeric &isot_ratio, const Numeric &P, const Numeric &T, const Numeric &H, const Numeric &sparse_lim, const bool do_zeeman, const Zeeman::Polarization zeeman_polarization, const Options::LblSpeedup speedup_type) ARTS_NOEXCEPT
Compute the line shape in its entirety.
Definition: lineshape.cc:3116
X3 LineCenter QuadraticIndependent
Definition: constants.h:577
X3 LineCenter None
Definition: constants.h:576
constexpr bool same_or_joker(const IsotopeRecord &ir1, const IsotopeRecord &ir2) noexcept
constexpr IsotopologueRatios isotopologue_ratiosInitFromBuiltin()
Implements Zeeman modeling.
Definition: zeemandata.cc:281
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:2010
Declaration of functions in rte.cc.
void check_abs_species(const ArrayOfArrayOfSpeciesTag &abs_species)
void interp_add_even(const ComputeData &sparse) ARTS_NOEXCEPT
Definition: lineshape.cc:3233
ComplexVector N
Definition: lineshape.h:687
void interp_add_triplequad(const ComputeData &sparse) ARTS_NOEXCEPT
Definition: lineshape.cc:3274
ComplexVector F
Definition: lineshape.h:687
void reset() noexcept
Definition: lineshape.h:700
ComplexMatrix dF
Definition: lineshape.h:688
ComplexMatrix dN
Definition: lineshape.h:688
This file contains basic functions to handle XML data files.