ARTS  2.4.0(git:4fb77825)
m_abs.cc
Go to the documentation of this file.
1 
2 /* Copyright (C) 2000-2012
3  Stefan Buehler <sbuehler@ltu.se>
4  Patrick Eriksson <patrick.eriksson@chalmers.se>
5  Axel von Engeln <engeln@uni-bremen.de>
6  Thomas Kuhn <tkuhn@uni-bremen.de>
7 
8  This program is free software; you can redistribute it and/or modify it
9  under the terms of the GNU General Public License as published by the
10  Free Software Foundation; either version 2, or (at your option) any
11  later version.
12 
13  This program is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with this program; if not, write to the Free Software
20  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
21  USA. */
22 
23 //
24 
33 #include <algorithm>
34 #include <cmath>
35 #include "absorption.h"
36 #include "array.h"
37 #include "arts.h"
38 #include "auto_md.h"
39 #include "check_input.h"
40 #include "legacy_continua.h"
41 #include "file.h"
42 #include "global_data.h"
43 #include "jacobian.h"
44 #include "m_xml.h"
45 #include "math_funcs.h"
46 #include "matpackI.h"
47 #include "messages.h"
48 #include "montecarlo.h"
49 #include "optproperties.h"
50 #include "parameters.h"
51 #include "physics_funcs.h"
52 #include "rte.h"
53 #include "xml_io.h"
54 
55 #ifdef ENABLE_NETCDF
56 #include <netcdf.h>
57 #include "nc_io.h"
58 #endif
59 
60 extern const Numeric ELECTRON_CHARGE;
61 extern const Numeric ELECTRON_MASS;
62 extern const Numeric PI;
63 extern const Numeric SPEED_OF_LIGHT;
64 extern const Numeric VACUUM_PERMITTIVITY;
65 
66 /* Workspace method: Doxygen documentation will be auto-generated */
67 void AbsInputFromRteScalars( // WS Output:
68  Vector& abs_p,
69  Vector& abs_t,
71  // WS Input:
72  const Numeric& rtp_pressure,
73  const Numeric& rtp_temperature,
74  const Vector& rtp_vmr,
75  const Verbosity&) {
76  // Prepare abs_p:
77  abs_p.resize(1);
79 
80  // Prepare abs_t:
81  abs_t.resize(1);
83 
84  // Prepare abs_vmrs:
85  abs_vmrs.resize(rtp_vmr.nelem(), 1);
86  abs_vmrs = rtp_vmr;
87 }
88 
89 /* Workspace method: Doxygen documentation will be auto-generated */
92  // WS Input:
94  const ArrayOfArrayOfSpeciesTag& tgs,
95  const Verbosity&) {
96  // Size is set but inner size will now change from the original definition of species tags...
97  abs_lines_per_species.resize(tgs.nelem());
98 
99  // The inner arrays need to be emptied, because they may contain lines
100  // from a previous calculation
101  for (auto &tg : abs_lines_per_species)
102  tg.resize(0);
103 
104  // Take copies because we have to support frequency ranges, so might have to delete
105  for (AbsorptionLines lines: abs_lines) {
106 
107  // Skip empty lines
108  if (lines.NumLines() == 0) continue;
109 
110  // Loop all the tags
111  for (Index i=0; i<tgs.nelem() and lines.NumLines(); i++) {
112  for (auto& this_tag: tgs[i]) {
113  // Test species first, we might leave as we leave
114  if (this_tag.Species() not_eq lines.Species()) break;
115 
116  // Test isotopologue, we have to hit the end of the list for no isotopologue or the exact value
117  if (this_tag.Isotopologue() not_eq SpeciesDataOfBand(lines).Isotopologue().nelem() and
118  this_tag.Isotopologue() not_eq lines.Isotopologue())
119  continue;
120 
121  // If there is a frequency range, we have to check so that only selected lines are included
122  if (this_tag.Lf() >= 0 or this_tag.Uf() >= 0) {
123  const Numeric low = (this_tag.Lf() >= 0) ? this_tag.Lf() : std::numeric_limits<Numeric>::lowest();
124  const Numeric upp = (this_tag.Uf() >= 0) ? this_tag.Uf() : std::numeric_limits<Numeric>::max();
125 
126  // Fill up a copy of the line record to match with the wished frequency criteria
127  AbsorptionLines these_lines = createEmptyCopy(lines);
128  for (Index k=lines.NumLines()-1; k>=0; k--)
129  if (low <= lines.F0(k) and upp >= lines.F0(k))
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 */
157  // Control Parameters:
158  const String& basename,
159  const Verbosity& verbosity) {
160  CREATE_OUT2;
161 
162  // Invalidate agenda check flags
164  abs_xsec_agenda_checked = false;
165 
166  // Species lookup data:
168 
169  // We want to make lists of included and excluded species:
170  ArrayOfString included(0), excluded(0);
171 
172  tgs.resize(0);
173 
174  for (Index i = 0; i < species_data.nelem(); ++i) {
175  const String specname = species_data[i].Name();
176 
177  String filename = basename;
178  if (basename.length() && basename[basename.length() - 1] != '/')
179  filename += ".";
180  filename += specname;
181 
182  try {
183  find_xml_file(filename, verbosity);
184  // Add to included list:
185  included.push_back(specname);
186 
187  // Convert name of species to a SpeciesTag object:
188  SpeciesTag this_tag(specname);
189 
190  // Create Array of SpeciesTags with length 1
191  // (our tag group has only one tag):
192  ArrayOfSpeciesTag this_group(1);
193  this_group[0] = this_tag;
194 
195  // Add this tag group to tgs:
196  tgs.push_back(this_group);
197  } catch (const std::runtime_error& e) {
198  // The file for the species could not be found.
199  excluded.push_back(specname);
200  }
201  }
202 
203  // Some nice output:
204  out2 << " Included Species (" << included.nelem() << "):\n";
205  for (Index i = 0; i < included.nelem(); ++i)
206  out2 << " " << included[i] << "\n";
207 
208  out2 << " Excluded Species (" << excluded.nelem() << "):\n";
209  for (Index i = 0; i < excluded.nelem(); ++i)
210  out2 << " " << excluded[i] << "\n";
211 }
212 
213 /* Workspace method: Doxygen documentation will be auto-generated */
214 void abs_speciesDefineAll( // WS Output:
218  // Control Parameters:
219  const Verbosity& verbosity) {
220  // Species lookup data:
222 
223  // We want to make lists of all species
224  ArrayOfString specs(0);
225  for (auto& spec: species_data) {
226  specs.push_back(spec.Name());
227  }
228 
229  // Set the values
231 }
232 
233 /* Workspace method: Doxygen documentation will be auto-generated */
234 void AbsInputFromAtmFields( // WS Output:
235  Vector& abs_p,
236  Vector& abs_t,
237  Matrix& abs_vmrs,
238  // WS Input:
239  const Index& atmosphere_dim,
240  const Vector& p_grid,
241  const Tensor3& t_field,
242  const Tensor4& vmr_field,
243  const Verbosity&) {
244  // First, make sure that we really have a 1D atmosphere:
245  if (1 != atmosphere_dim) {
246  ostringstream os;
247  os << "Atmospheric dimension must be 1D, but atmosphere_dim is "
248  << atmosphere_dim << ".";
249  throw runtime_error(os.str());
250  }
251 
252  abs_p = p_grid;
253  abs_t = t_field(joker, 0, 0);
254  abs_vmrs = vmr_field(joker, joker, 0, 0);
255 }
256 
257 /* Workspace method: Doxygen documentation will be auto-generated */
258 void abs_coefCalcFromXsec( // WS Output:
259  Matrix& abs_coef,
260  Matrix& src_coef,
265  // WS Input:
272  const Matrix& abs_vmrs,
273  const Vector& abs_p,
274  const Vector& abs_t,
275  const Verbosity& verbosity) {
276  CREATE_OUT3;
277 
278  // Check that abs_vmrs and abs_xsec_per_species really have compatible
279  // dimensions. In abs_vmrs there should be one row for each tg:
280  if (abs_vmrs.nrows() != abs_xsec_per_species.nelem()) {
281  ostringstream os;
282  os << "Variable abs_vmrs must have compatible dimension to abs_xsec_per_species.\n"
283  << "abs_vmrs.nrows() = " << abs_vmrs.nrows() << "\n"
284  << "abs_xsec_per_species.nelem() = " << abs_xsec_per_species.nelem();
285  throw runtime_error(os.str());
286  }
287 
288  // Check that number of altitudes are compatible. We only check the
289  // first element, this is possilble because within arts all elements
290  // are on the same altitude grid.
291  if (abs_vmrs.ncols() != abs_xsec_per_species[0].ncols()) {
292  ostringstream os;
293  os << "Variable abs_vmrs must have same numbers of altitudes as abs_xsec_per_species.\n"
294  << "abs_vmrs.ncols() = " << abs_vmrs.ncols() << "\n"
295  << "abs_xsec_per_species[0].ncols() = "
296  << abs_xsec_per_species[0].ncols();
297  throw runtime_error(os.str());
298  }
299 
300  // Check dimensions of abs_p and abs_t:
301  chk_size("abs_p", abs_p, abs_vmrs.ncols());
302  chk_size("abs_t", abs_t, abs_vmrs.ncols());
303 
304  // Will these calculations deal with nlte?
305  const bool do_src =
306  (src_xsec_per_species.nelem() == abs_xsec_per_species.nelem())
307  ? not src_xsec_per_species[0].empty() ? true : false
308  : false;
309  const Index src_rows = do_src ? src_xsec_per_species[0].nrows() : 0;
310 
311  // Initialize abs_coef and abs_coef_per_species. The array dimension of abs_coef_per_species
312  // is the same as that of abs_xsec_per_species. The dimension of abs_coef should
313  // be equal to one of the abs_xsec_per_species enries.
314 
315  abs_coef.resize(abs_xsec_per_species[0].nrows(),
316  abs_xsec_per_species[0].ncols());
317  abs_coef = 0;
318  src_coef.resize(src_rows, src_xsec_per_species[0].ncols());
319  src_coef = 0;
320 
321  const ArrayOfIndex jacobian_quantities_position =
323  dabs_coef_dx.resize(jacobian_quantities_position.nelem());
324  dsrc_coef_dx.resize(do_src ? jacobian_quantities_position.nelem() : 0);
325 
326  for (Index ii = 0; ii < jacobian_quantities_position.nelem(); ii++) {
327  if (jacobian_quantities[jacobian_quantities_position[ii]] not_eq
329  dabs_coef_dx[ii].resize(abs_xsec_per_species[0].nrows(),
331  dabs_coef_dx[ii] = 0.0;
332  if (do_src) {
333  dsrc_coef_dx[ii].resize(src_rows, src_xsec_per_species[0].ncols());
334  dsrc_coef_dx[ii] = 0.0;
335  }
336  }
337  }
338 
341 
342  out3
343  << " Computing abs_coef and abs_coef_per_species from abs_xsec_per_species.\n";
344  // Loop through all tag groups
345  for (Index i = 0; i < abs_xsec_per_species.nelem(); ++i) {
346  out3 << " Tag group " << i << "\n";
347 
348  // Make this element of abs_xsec_per_species the right size:
351  abs_coef_per_species[i] = 0; // Initialize all elements to 0.
352 
353  if (do_src) {
354  src_coef_per_species[i].resize(src_rows, src_xsec_per_species[i].ncols());
355  src_coef_per_species[i] = 0; // Initialize all elements to 0.
356  }
357 
358  // Loop through all altitudes
359  for (Index j = 0; j < abs_xsec_per_species[i].ncols(); j++) {
360  // Calculate total number density from pressure and temperature.
361  const Numeric n = number_density(abs_p[j], abs_t[j]);
362  const Numeric dn_dT = dnumber_density_dt(abs_p[j], abs_t[j]);
363  // Wasted calculations when Jacobians are not calculated...
364  // Though this is called seldom enough that it this fine? value is -1/t*n
365 
366  // Loop through all frequencies
367  for (Index k = 0; k < abs_xsec_per_species[i].nrows(); k++) {
368  abs_coef_per_species[i](k, j) =
369  abs_xsec_per_species[i](k, j) * n * abs_vmrs(i, j);
370  if (do_src)
371  src_coef_per_species[i](k, j) =
372  src_xsec_per_species[i](k, j) * n * abs_vmrs(i, j);
373 
374  for (Index iq = 0; iq < jacobian_quantities_position.nelem(); iq++) {
375  if (jacobian_quantities[jacobian_quantities_position[iq]] ==
376  JacPropMatType::Temperature) {
377  dabs_coef_dx[iq](k, j) +=
378  (dabs_xsec_per_species_dx[i][iq](k, j) * n +
379  abs_xsec_per_species[i](k, j) * dn_dT) *
380  abs_vmrs(i, j);
381  if (do_src)
382  dsrc_coef_dx[iq](k, j) +=
383  (dsrc_xsec_per_species_dx[i][iq](k, j) * n +
384  src_xsec_per_species[i](k, j) * dn_dT) *
385  abs_vmrs(i, j);
386  } else if (jacobian_quantities[jacobian_quantities_position[iq]] ==
387  JacPropMatType::VMR) {
388  bool seco = false, main = false;
389  for (const auto& s : abs_species[i]) {
390  if (species_match(
391  jacobian_quantities[jacobian_quantities_position[iq]],
392  s.BathSpecies()) or
393  s.Type() not_eq SpeciesTag::TYPE_CIA)
394  seco = true;
395  if (species_iso_match(
396  jacobian_quantities[jacobian_quantities_position[iq]],
397  s.Species(),
398  s.Isotopologue()))
399  main = true;
400  }
401  if (main and seco) {
402  dabs_coef_dx[iq](k, j) +=
403  (dabs_xsec_per_species_dx[i][iq](k, j) * abs_vmrs(i, j) +
404  abs_xsec_per_species[i](k, j)) *
405  n;
406  if (do_src)
407  dsrc_coef_dx[iq](k, j) +=
408  (dsrc_xsec_per_species_dx[i][iq](k, j) * abs_vmrs(i, j) +
409  src_xsec_per_species[i](k, j)) *
410  n;
411  } else if (main) {
412  dabs_coef_dx[iq](k, j) += abs_xsec_per_species[i](k, j) * n;
413  if (do_src)
414  dsrc_coef_dx[iq](k, j) += src_xsec_per_species[i](k, j) * n;
415  } else if (seco) {
416  dabs_coef_dx[iq](k, j) +=
417  dabs_xsec_per_species_dx[i][iq](k, j) * abs_vmrs(i, j) * n;
418  if (do_src)
419  dsrc_coef_dx[iq](k, j) +=
420  dsrc_xsec_per_species_dx[i][iq](k, j) * abs_vmrs(i, j) * n;
421  }
422  } else if (jacobian_quantities
423  [jacobian_quantities_position[iq]] not_eq
425  dabs_coef_dx[iq](k, j) +=
426  dabs_xsec_per_species_dx[i][iq](k, j) * n * abs_vmrs(i, j);
427  if (do_src)
428  dsrc_coef_dx[iq](k, j) +=
429  dsrc_xsec_per_species_dx[i][iq](k, j) * n * abs_vmrs(i, j);
430  }
431  }
432  }
433  }
434 
435  // Add up to the total absorption:
436  abs_coef += abs_coef_per_species[i]; // In Matpack you can use the +=
437  // operator to do elementwise addition.
438  if (do_src) src_coef += src_coef_per_species[i];
439  }
440 }
441 
442 /* Workspace method: Doxygen documentation will be auto-generated */
443 void abs_xsec_per_speciesInit( // WS Output:
448  // WS Input:
449  const ArrayOfArrayOfSpeciesTag& tgs,
452  const Vector& f_grid,
453  const Vector& abs_p,
455  const Index& nlte_do,
456  const Verbosity& verbosity) {
457  CREATE_OUT3;
458 
460  throw runtime_error(
461  "You must call *abs_xsec_agenda_checkedCalc* before calling this method.");
462 
463  // We need to check that abs_species_active doesn't have more elements than
464  // abs_species (abs_xsec_agenda_checkedCalc doesn't know abs_species_active.
465  // Usually we come here through an agenda call, where abs_species_active has
466  // been properly created somewhere internally. But we might get here by
467  // direct call, and then need to be safe!).
468  if (tgs.nelem() < abs_species_active.nelem()) {
469  ostringstream os;
470  os << "abs_species_active (n=" << abs_species_active.nelem()
471  << ") not allowed to have more elements than abs_species (n="
472  << tgs.nelem() << ")!\n";
473  throw runtime_error(os.str());
474  }
475 
476  // Initialize abs_xsec_per_species. The array dimension of abs_xsec_per_species
477  // is the same as that of abs_lines_per_species.
478  abs_xsec_per_species.resize(tgs.nelem());
479  src_xsec_per_species.resize(tgs.nelem());
480 
481  const bool do_jac = supports_propmat_clearsky(
482  jacobian_quantities); //minus one is a flag for a non-species...
483  const ArrayOfIndex jacobian_quantities_position =
485 
486  dabs_xsec_per_species_dx.resize(do_jac ? tgs.nelem() : 0);
487  dsrc_xsec_per_species_dx.resize(do_jac ? tgs.nelem() : 0);
488 
489  // Loop abs_xsec_per_species and make each matrix the right size,
490  // initializing to zero.
491  // But skip inactive species, loop only over the active ones.
492  for (Index ii = 0; ii < abs_species_active.nelem(); ++ii) {
493  const Index i = abs_species_active[ii];
494  // Check that abs_species_active index is not higher than the number
495  // of species
496  if (i >= tgs.nelem()) {
497  ostringstream os;
498  os << "*abs_species_active* contains an invalid species index.\n"
499  << "Species index must be between 0 and " << tgs.nelem() - 1;
500  throw std::runtime_error(os.str());
501  }
502  // Make this element of abs_xsec_per_species the right size:
503  abs_xsec_per_species[i].resize(f_grid.nelem(), abs_p.nelem());
504  abs_xsec_per_species[i] = 0; // Matpack can set all elements like this.
505  if (nlte_do) {
506  src_xsec_per_species[i].resize(f_grid.nelem(), abs_p.nelem());
507  src_xsec_per_species[i] = 0;
508  } else {
509  src_xsec_per_species[i].resize(0, 0);
510  }
511 
512  if (do_jac) {
514  ArrayOfMatrix(jacobian_quantities_position.nelem(),
515  Matrix(f_grid.nelem(), abs_p.nelem(), 0.0));
516  if (nlte_do)
518  ArrayOfMatrix(jacobian_quantities_position.nelem(),
519  Matrix(f_grid.nelem(), abs_p.nelem(), 0.0));
520  }
521  }
522 
523  ostringstream os;
524  os << " Initialized abs_xsec_per_species.\n"
525  << " Number of frequencies : " << f_grid.nelem() << "\n"
526  << " Number of pressure levels : " << abs_p.nelem() << "\n";
527  out3 << os.str();
528 }
529 
530 /* Workspace method: Doxygen documentation will be auto-generated */
531 void abs_xsec_per_speciesAddConts( // WS Output:
534  // WS Input:
535  const ArrayOfArrayOfSpeciesTag& tgs,
538  const Vector& f_grid,
539  const Vector& abs_p,
540  const Vector& abs_t,
541  const Matrix& abs_vmrs,
545  const Verbosity& verbosity) {
546  CREATE_OUT3;
547 
548  // Needed for some continua, and set here from abs_vmrs:
549  Vector abs_h2o, abs_n2, abs_o2;
550 
551  // Check that all paramters that should have the number of tag
552  // groups as a dimension are consistent.
553  {
554  const Index n_tgs = tgs.nelem();
555  const Index n_xsec = abs_xsec_per_species.nelem();
556  const Index n_vmrs = abs_vmrs.nrows();
557 
558  if (n_tgs != n_xsec || n_tgs != n_vmrs) {
559  ostringstream os;
560  os << "The following variables must all have the same dimension:\n"
561  << "tgs: " << tgs.nelem() << "\n"
562  << "abs_xsec_per_species: " << abs_xsec_per_species.nelem() << "\n"
563  << "abs_vmrs.nrows(): " << abs_vmrs.nrows();
564  throw runtime_error(os.str());
565  }
566  }
567 
568  // Jacobian overhead START
569  /* NOTE: if any of the functions inside continuum tags could
570  be made to give partial derivatives, then that would
571  speed things up. Also be aware that line specific
572  parameters cannot be retrieved while using these
573  models. */
574  const bool do_jac = supports_continuum(
575  jacobian_quantities); // Throws runtime error if line parameters are wanted since we cannot know if the line is in the Continuum...
576  const bool do_freq_jac = do_frequency_jacobian(jacobian_quantities);
577  const bool do_temp_jac = do_temperature_jacobian(jacobian_quantities);
578  Vector dfreq, dabs_t;
581  const ArrayOfIndex jacobian_quantities_position =
583 
584  if (do_freq_jac) {
585  dfreq.resize(f_grid.nelem());
586  for (Index iv = 0; iv < f_grid.nelem(); iv++) dfreq[iv] = f_grid[iv] + df;
587  }
588  if (do_temp_jac) {
589  dabs_t.resize(abs_t.nelem());
590  for (Index it = 0; it < abs_t.nelem(); it++) dabs_t[it] = abs_t[it] + dt;
591  }
592 
593  Matrix jacs_df, jacs_dt, normal;
594  if (do_jac) {
595  if (do_freq_jac) jacs_df.resize(f_grid.nelem(), abs_p.nelem());
596  if (do_temp_jac) jacs_dt.resize(f_grid.nelem(), abs_p.nelem());
597  normal.resize(f_grid.nelem(), abs_p.nelem());
598  }
599  // Jacobian overhead END
600 
601  // Check, that dimensions of abs_cont_names and
602  // abs_cont_parameters are consistent...
603  if (abs_cont_names.nelem() != abs_cont_parameters.nelem()) {
604  ostringstream os;
605 
606  for (Index i = 0; i < abs_cont_names.nelem(); ++i)
607  os << "abs_xsec_per_speciesAddConts: " << i
608  << " name : " << abs_cont_names[i] << "\n";
609 
610  for (Index i = 0; i < abs_cont_parameters.nelem(); ++i)
611  os << "abs_xsec_per_speciesAddConts: " << i
612  << " param: " << abs_cont_parameters[i] << "\n";
613 
614  for (Index i = 0; i < abs_cont_models.nelem(); ++i)
615  os << "abs_xsec_per_speciesAddConts: " << i
616  << " option: " << abs_cont_models[i] << "\n";
617 
618  os << "The following variables must have the same dimension:\n"
619  << "abs_cont_names: " << abs_cont_names.nelem() << "\n"
620  << "abs_cont_parameters: " << abs_cont_parameters.nelem();
621 
622  throw runtime_error(os.str());
623  }
624 
625  // Check that abs_p, abs_t, and abs_vmrs have the same
626  // dimension. This could be a user error, so we throw a
627  // runtime_error.
628 
629  if (abs_t.nelem() != abs_p.nelem()) {
630  ostringstream os;
631  os << "Variable abs_t must have the same dimension as abs_p.\n"
632  << "abs_t.nelem() = " << abs_t.nelem() << '\n'
633  << "abs_p.nelem() = " << abs_p.nelem();
634  throw runtime_error(os.str());
635  }
636 
637  if (abs_vmrs.ncols() != abs_p.nelem()) {
638  ostringstream os;
639  os << "Variable dimension abs_vmrs.ncols() must\n"
640  << "be the same as abs_p.nelem().\n"
641  << "abs_vmrs.ncols() = " << abs_vmrs.ncols() << '\n'
642  << "abs_p.nelem() = " << abs_p.nelem();
643  throw runtime_error(os.str());
644  }
645 
646  // We set abs_h2o, abs_n2, and abs_o2 later, because we only want to
647  // do it if the parameters are really needed.
648 
649  out3 << " Calculating continuum spectra.\n";
650 
651  // Loop tag groups:
652  for (Index ii = 0; ii < abs_species_active.nelem(); ++ii) {
653  const Index i = abs_species_active[ii];
654 
656 
657  // Go through the tags in the current tag group to see if they
658  // are continuum tags:
659  for (Index s = 0; s < tgs[i].nelem(); ++s) {
660  // Continuum tags in the sense that we talk about here
661  // (including complete absorption models) are marked by a special type.
662  if (tgs[i][s].Type() == SpeciesTag::TYPE_PREDEF) {
663  // We have identified a continuum tag!
664 
665  // Get only the continuum name. The full tag name is something like:
666  // H2O-HITRAN96Self-*-*. We want only the `H2O-HITRAN96Self' part:
667  const String name = species_data[tgs[i][s].Species()].Name() + "-" +
668  species_data[tgs[i][s].Species()]
669  .Isotopologue()[tgs[i][s].Isotopologue()]
670  .Name();
671 
672  if (name == "O2-MPM2020") continue;
673 
674  // Check, if we have parameters for this model. For
675  // this, the model name must be listed in
676  // abs_cont_names.
677  const Index n =
678  find(abs_cont_names.begin(), abs_cont_names.end(), name) -
679  abs_cont_names.begin();
680 
681  // n==abs_cont_names.nelem() indicates that
682  // the name was not found.
683  if (n == abs_cont_names.nelem()) {
684  ostringstream os;
685  os << "Cannot find model " << name << " in abs_cont_names.";
686  throw runtime_error(os.str());
687  }
688 
689  // Ok, the tag specifies a valid continuum model and
690  // we have continuum parameters.
691 
692  if (out3.sufficient_priority()) {
693  ostringstream os;
694  os << " Adding " << name << " to tag group " << i << ".\n";
695  out3 << os.str();
696  }
697 
698  // find the options for this continuum tag from the input array
699  // of options. The actual field of the array is n:
700  const String ContOption = abs_cont_models[n];
701 
702  // Set abs_h2o, abs_n2, and abs_o2 from the first matching species.
703  set_vmr_from_first_species(abs_h2o, "H2O", tgs, abs_vmrs);
704  set_vmr_from_first_species(abs_n2, "N2", tgs, abs_vmrs);
705  set_vmr_from_first_species(abs_o2, "O2", tgs, abs_vmrs);
706 
707  // Add the continuum for this tag. The parameters in
708  // this call should be clear. The vmr is in
709  // abs_vmrs(i,Range(joker)). The other vmr variables,
710  // abs_h2o, abs_n2, and abs_o2 contains the real vmr of H2O,
711  // N2, nad O2, which are needed as additional information for
712  // certain continua:
713  // abs_h2o for
714  // O2-PWR88, O2-PWR93, O2-PWR98,
715  // O2-MPM85, O2-MPM87, O2-MPM89, O2-MPM92, O2-MPM93,
716  // O2-TRE05,
717  // O2-SelfContStandardType, O2-SelfContMPM93, O2-SelfContPWR93,
718  // N2-SelfContMPM93, N2-DryContATM01,
719  // N2-CIArotCKDMT252, N2-CIAfunCKDMT252
720  // abs_n2 for
721  // H2O-SelfContCKD24, H2O-ForeignContCKD24,
722  // O2-v0v0CKDMT100,
723  // CO2-ForeignContPWR93, CO2-ForeignContHo66
724  // abs_o2 for
725  // N2-CIArotCKDMT252, N2-CIAfunCKDMT252
726  if (!do_jac)
728  name,
730  abs_cont_models[n],
731  f_grid,
732  abs_p,
733  abs_t,
734  abs_n2,
735  abs_h2o,
736  abs_o2,
737  abs_vmrs(i, Range(joker)),
738  verbosity);
739  else // The Jacobian block
740  {
741  // Needs a reseted block here...
742  for (Index iv = 0; iv < f_grid.nelem(); iv++) {
743  for (Index ip = 0; ip < abs_p.nelem(); ip++) {
744  if (do_freq_jac) jacs_df(iv, ip) = 0.0;
745  if (do_temp_jac) jacs_dt(iv, ip) = 0.0;
746  normal(iv, ip) = 0.0;
747  }
748  }
749 
750  // Normal calculations
751  xsec_continuum_tag(normal,
752  name,
754  abs_cont_models[n],
755  f_grid,
756  abs_p,
757  abs_t,
758  abs_n2,
759  abs_h2o,
760  abs_o2,
761  abs_vmrs(i, Range(joker)),
762  verbosity);
763 
764  // Frequency calculations
765  if (do_freq_jac)
766  xsec_continuum_tag(jacs_df,
767  name,
769  abs_cont_models[n],
770  dfreq,
771  abs_p,
772  abs_t,
773  abs_n2,
774  abs_h2o,
775  abs_o2,
776  abs_vmrs(i, Range(joker)),
777  verbosity);
778 
779  //Temperature calculations
780  if (do_temp_jac)
781  xsec_continuum_tag(jacs_dt,
782  name,
784  abs_cont_models[n],
785  f_grid,
786  abs_p,
787  dabs_t,
788  abs_n2,
789  abs_h2o,
790  abs_o2,
791  abs_vmrs(i, Range(joker)),
792  verbosity);
793  for (Index iv = 0; iv < f_grid.nelem(); iv++) {
794  for (Index ip = 0; ip < abs_p.nelem(); ip++) {
795  abs_xsec_per_species[i](iv, ip) += normal(iv, ip);
796  for (Index iq = 0; iq < jacobian_quantities_position.nelem();
797  iq++) {
799  jacobian_quantities[jacobian_quantities_position[iq]]))
800  dabs_xsec_per_species_dx[i][iq](iv, ip) +=
801  (jacs_df(iv, ip) - normal(iv, ip)) * (1. / df);
802  else if (jacobian_quantities
803  [jacobian_quantities_position[iq]] ==
804  JacPropMatType::Temperature)
805  dabs_xsec_per_species_dx[i][iq](iv, ip) +=
806  (jacs_dt(iv, ip) - normal(iv, ip)) * (1. / dt);
807  }
808  }
809  }
810  }
811  // Calling this function with a row of Matrix abs_vmrs
812  // is possible because it uses Views.
813  }
814  }
815  }
816 }
817 
818 //======================================================================
819 // Methods related to continua
820 //======================================================================
821 
822 /* Workspace method: Doxygen documentation will be auto-generated */
823 void abs_cont_descriptionInit( // WS Output:
825  ArrayOfString& abs_cont_options,
827  const Verbosity& verbosity) {
828  CREATE_OUT2;
829 
830  abs_cont_names.resize(0);
831  abs_cont_options.resize(0);
832  abs_cont_parameters.resize(0);
833  out2 << " Initialized abs_cont_names \n"
834  " abs_cont_models\n"
835  " abs_cont_parameters.\n";
836 }
837 
838 /* Workspace method: Doxygen documentation will be auto-generated */
839 void abs_cont_descriptionAppend( // WS Output:
843  // Control Parameters:
844  const String& tagname,
845  const String& model,
846  const Vector& userparameters,
847  const Verbosity&) {
848  // First we have to check that name matches a continuum species tag.
849  check_continuum_model(tagname);
850 
851  //cout << " + tagname: " << tagname << "\n";
852  //cout << " + model: " << model << "\n";
853  //cout << " + parameters: " << userparameters << "\n";
854 
855  // Add name and parameters to the apropriate variables:
856  abs_cont_names.push_back(tagname);
857  abs_cont_models.push_back(model);
858  abs_cont_parameters.push_back(userparameters);
859 }
860 
861 /* Workspace method: Doxygen documentation will be auto-generated */
866  // WS Input:
870  const Vector& f_grid,
871  const Numeric& rtp_temperature,
872  const Verbosity&) {
873  // nlte_source has format
874  // [ abs_species, f_grid, stokes_dim ].
875  // src_coef_per_species has format ArrayOfMatrix (over species),
876  // where for each species the matrix has format [f_grid, abs_p].
877 
878  Index n_species = src_coef_per_species.nelem(); // # species
879 
880  if (not n_species) {
881  ostringstream os;
882  os << "Must have at least one species.";
883  throw runtime_error(os.str());
884  }
885 
886  Index n_f = src_coef_per_species[0].nrows(); // # frequencies
887 
888  // # pressures must be 1:
889  if (1 not_eq src_coef_per_species[0].ncols()) {
890  ostringstream os;
891  os << "Must have exactly one pressure.";
892  throw runtime_error(os.str());
893  }
894 
895  // Check species dimension of propmat_clearsky
896  if (nlte_source.nelem() not_eq n_species) {
897  ostringstream os;
898  os << "Species dimension of propmat_clearsky does not\n"
899  << "match src_coef_per_species.";
900  throw std::runtime_error(os.str());
901  }
902 
903  // Check frequency dimension of propmat_clearsky
904  if (nlte_source[0].NumberOfFrequencies() not_eq n_f) {
905  ostringstream os;
906  os << "Frequency dimension of propmat_clearsky does not\n"
907  << "match abs_coef_per_species.";
908  throw runtime_error(os.str());
909  }
910 
911  const ArrayOfIndex jacobian_quantities_position =
913  Vector B(n_f);
914 
915  for (Index iv = 0; iv < n_f; iv++)
916  B[iv] = planck(f_grid[iv], rtp_temperature);
917 
918  StokesVector sv(n_f, nlte_source[0].StokesDimensions());
919  for (Index si = 0; si < n_species; ++si) {
920  sv.Kjj() = src_coef_per_species[si](joker, 0);
921  sv *= B;
922  nlte_source[si].Kjj() += sv.Kjj();
923  }
924 
925  // Jacobian
926  for (Index ii = 0; ii < jacobian_quantities_position.nelem(); ii++) {
927  if (jacobian_quantities[jacobian_quantities_position[ii]] ==
928  JacPropMatType::Temperature) {
929  Vector dB(n_f);
930  for (Index iv = 0; iv < n_f; iv++)
931  dB[iv] = dplanck_dt(f_grid[iv], rtp_temperature);
932 
933  for (Index si = 0; si < n_species; ++si) {
934  sv.Kjj() = src_coef_per_species[si](joker, 0);
935  sv *= dB;
936  nlte_dsource_dx[ii].Kjj() += sv.Kjj();
937  }
938 
939  sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
940  sv *= B;
941  dnlte_dx_source[ii].Kjj() += sv.Kjj();
942  } else if (is_frequency_parameter(
943  jacobian_quantities[jacobian_quantities_position[ii]])) {
944  Vector dB(n_f);
945  for (Index iv = 0; iv < n_f; iv++)
946  dB[iv] = dplanck_df(f_grid[iv], rtp_temperature);
947 
948  for (Index si = 0; si < n_species; ++si) {
949  sv.Kjj() = src_coef_per_species[si](joker, 0);
950  sv *= dB;
951  nlte_dsource_dx[ii].Kjj() += sv.Kjj();
952  }
953 
954  sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
955  sv *= B;
956  dnlte_dx_source[ii].Kjj() += sv.Kjj();
957  } else if (jacobian_quantities[jacobian_quantities_position[ii]] not_eq
959  sv.Kjj() = dsrc_coef_dx[ii](joker, 0);
960  sv *= B;
961  dnlte_dx_source[ii].Kjj() += sv.Kjj();
962  } else { /* All is fine! */
963  }
964  }
965 }
966 
967 /* Workspace method: Doxygen documentation will be auto-generated */
971  // WS Input:
974  const Verbosity&) {
975  // propmat_clearsky has format
976  // [ abs_species, f_grid, stokes_dim, stokes_dim ].
977  // abs_coef_per_species has format ArrayOfMatrix (over species),
978  // where for each species the matrix has format [f_grid, abs_p].
979 
980  Index n_species = abs_coef_per_species.nelem(); // # species
981 
982  if (0 == n_species) {
983  ostringstream os;
984  os << "Must have at least one species.";
985  throw runtime_error(os.str());
986  }
987 
988  Index n_f = abs_coef_per_species[0].nrows(); // # frequencies
989 
990  // # pressures must be 1:
991  if (1 not_eq abs_coef_per_species[0].ncols()) {
992  ostringstream os;
993  os << "Must have exactly one pressure.";
994  throw runtime_error(os.str());
995  }
996 
997  // Check species dimension of propmat_clearsky
998  if (propmat_clearsky.nelem() not_eq n_species) {
999  ostringstream os;
1000  os << "Species dimension of propmat_clearsky does not\n"
1001  << "match abs_coef_per_species.";
1002  throw runtime_error(os.str());
1003  }
1004 
1005  // Check frequency dimension of propmat_clearsky
1006  if (propmat_clearsky[0].NumberOfFrequencies() not_eq n_f) {
1007  ostringstream os;
1008  os << "Frequency dimension of propmat_clearsky does not\n"
1009  << "match abs_coef_per_species.";
1010  throw runtime_error(os.str());
1011  }
1012 
1013  // Loop species and stokes dimensions, and add to propmat_clearsky:
1014  for (Index si = 0; si < n_species; ++si)
1015  propmat_clearsky[si].Kjj() += abs_coef_per_species[si](joker, 0);
1016 
1017  for (Index iqn = 0; iqn < dabs_coef_dx.nelem(); iqn++) {
1018  if (dabs_coef_dx[iqn].nrows() == n_f) {
1019  if (dabs_coef_dx[iqn].ncols() == 1) {
1020  dpropmat_clearsky_dx[iqn].Kjj() += dabs_coef_dx[iqn](joker, 0);
1021  } else
1022  throw std::runtime_error("Must have exactly one pressure.");
1023  }
1024  }
1025 }
1026 
1027 /* Workspace method: Doxygen documentation will be auto-generated */
1028 void propmat_clearskyInit( //WS Output
1034  //WS Input
1037  const Vector& f_grid,
1038  const Index& stokes_dim,
1040  const Index& nlte_do,
1041  const Verbosity&) {
1043  throw runtime_error(
1044  "You must call *propmat_clearsky_agenda_checkedCalc* before calling this method.");
1045 
1046  Index nf = f_grid.nelem();
1047 
1048  if (not abs_species.nelem())
1049  throw std::runtime_error("abs_species.nelem() = 0");
1050 
1051  if (not nf) throw runtime_error("nf = 0");
1052 
1053  if (not stokes_dim) throw runtime_error("stokes_dim = 0");
1054 
1056 
1058  abs_species.nelem(), PropagationMatrix(nf, stokes_dim));
1061 
1063  StokesVector(nf, stokes_dim))
1064  : ArrayOfStokesVector(0);
1067  : ArrayOfStokesVector(0);
1070  : ArrayOfStokesVector(0);
1071 }
1072 
1073 /* Workspace method: Doxygen documentation will be auto-generated */
1077  const Index& stokes_dim,
1078  const Index& atmosphere_dim,
1079  const Vector& f_grid,
1082  const Vector& rtp_vmr,
1083  const Vector& rtp_los,
1084  const Vector& rtp_mag,
1085  const Verbosity&) {
1086  // All the physical constants joined into one static constant:
1087  // (abs as e defined as negative)
1088  static const Numeric FRconst =
1091  ELECTRON_MASS));
1092 
1093  if (stokes_dim < 3)
1094  throw runtime_error(
1095  "To include Faraday rotation, stokes_dim >= 3 is required.");
1096 
1097  if (atmosphere_dim == 1 && rtp_los.nelem() < 1) {
1098  ostringstream os;
1099  os << "For applying propmat_clearskyAddFaraday, los needs to be specified\n"
1100  << "(at least zenith angle component for atmosphere_dim==1),\n"
1101  << "but it is not.\n";
1102  throw runtime_error(os.str());
1103  } else if (atmosphere_dim > 1 && rtp_los.nelem() < 2) {
1104  ostringstream os;
1105  os << "For applying propmat_clearskyAddFaraday, los needs to be specified\n"
1106  << "(both zenith and azimuth angle components for atmosphere_dim>1),\n"
1107  << "but it is not.\n";
1108  throw runtime_error(os.str());
1109  }
1110 
1111  const bool do_jac = supports_faraday(jacobian_quantities);
1112  const bool do_magn_jac = do_magnetic_jacobian(jacobian_quantities);
1114  const ArrayOfIndex jacobian_quantities_position =
1116 
1117  Index ife = -1;
1118  for (Index sp = 0; sp < abs_species.nelem() && ife < 0; sp++) {
1120  ife = sp;
1121  }
1122  }
1123 
1124  if (ife < 0) {
1125  throw runtime_error(
1126  "Free electrons not found in *abs_species* and "
1127  "Faraday rotation can not be calculated.");
1128  } else {
1129  const Numeric ne = rtp_vmr[ife];
1130 
1131  if (ne != 0 && (rtp_mag[0] != 0 || rtp_mag[1] != 0 || rtp_mag[2] != 0)) {
1132  // Include remaining terms, beside /f^2
1133  const Numeric c1 =
1134  2 * FRconst * ne *
1137 
1138  Numeric dc1_u = 0.0, dc1_v = 0.0, dc1_w = 0.0;
1139  if (do_magn_jac) {
1140  dc1_u = (2 * FRconst * ne *
1142  rtp_mag[0] + dmag,
1143  rtp_mag[1],
1144  rtp_mag[2],
1145  atmosphere_dim) -
1146  c1) /
1147  dmag;
1148  dc1_v = (2 * FRconst * ne *
1150  rtp_mag[0],
1151  rtp_mag[1] + dmag,
1152  rtp_mag[2],
1153  atmosphere_dim) -
1154  c1) /
1155  dmag;
1156  dc1_w = (2 * FRconst * ne *
1158  rtp_mag[0],
1159  rtp_mag[1],
1160  rtp_mag[2] + dmag,
1161  atmosphere_dim) -
1162  c1) /
1163  dmag;
1164  }
1165 
1166  if (not do_jac) {
1167  for (Index iv = 0; iv < f_grid.nelem(); iv++) {
1168  const Numeric r = c1 / (f_grid[iv] * f_grid[iv]);
1169  propmat_clearsky[ife].SetFaraday(r, iv);
1170  }
1171  } else {
1172  for (Index iv = 0; iv < f_grid.nelem(); iv++) {
1173  const Numeric f2 = f_grid[iv] * f_grid[iv];
1174  const Numeric r = c1 / f2;
1175  propmat_clearsky[ife].SetFaraday(r, iv);
1176 
1177  // The Jacobian loop
1178  for (Index iq = 0; iq < jacobian_quantities_position.nelem(); iq++) {
1180  jacobian_quantities[jacobian_quantities_position[iq]]))
1181  dpropmat_clearsky_dx[iq].AddFaraday(-2.0 * r / f_grid[iv], iv);
1182  else if (jacobian_quantities[jacobian_quantities_position[iq]] ==
1183  JacPropMatType::MagneticU)
1184  dpropmat_clearsky_dx[iq].AddFaraday(dc1_u / f2, iv);
1185  else if (jacobian_quantities[jacobian_quantities_position[iq]] ==
1186  JacPropMatType::MagneticV)
1187  dpropmat_clearsky_dx[iq].AddFaraday(dc1_v / f2, iv);
1188  else if (jacobian_quantities[jacobian_quantities_position[iq]] ==
1189  JacPropMatType::MagneticW)
1190  dpropmat_clearsky_dx[iq].AddFaraday(dc1_w / f2, iv);
1191  else if (jacobian_quantities[jacobian_quantities_position[iq]] ==
1192  JacPropMatType::Electrons)
1193  dpropmat_clearsky_dx[iq].AddFaraday(r / ne, iv);
1194  }
1195  }
1196  }
1197  }
1198  }
1199 }
1200 
1201 /* Workspace method: Doxygen documentation will be auto-generated */
1203  // WS Output:
1206  // WS Input:
1207  const Index& stokes_dim,
1208  const Index& atmosphere_dim,
1209  const Vector& f_grid,
1212  const Vector& rtp_vmr,
1213  const Vector& rtp_los,
1214  const Numeric& rtp_temperature,
1216  const Index& scat_data_checked,
1217  const Index& use_abs_as_ext,
1218  // Verbosity object:
1219  const Verbosity& verbosity) {
1220  CREATE_OUT1;
1221 
1222  // (i)yCalc only checks scat_data_checked if cloudbox is on. It is off here,
1223  // though, i.e. we need to check it here explicitly. (Also, cloudboxOff sets
1224  // scat_data_checked=0 as it does not check it and as we ususally don't need
1225  // scat_data for clearsky cases, hence don't want to check them by
1226  // scat_data_checkedCalc in that case. This approach seems to be the more
1227  // handy compared to cloudboxOff setting scat_data_checked=1 without checking
1228  // it assuming we won't use it anyways.)
1229  if (scat_data_checked != 1)
1230  throw runtime_error(
1231  "The scat_data must be flagged to have "
1232  "passed a consistency check (scat_data_checked=1).");
1233 
1235  Index np = 0;
1236  for (Index sp = 0; sp < abs_species.nelem(); sp++) {
1237  if (abs_species[sp][0].Type() == SpeciesTag::TYPE_PARTICLES) {
1238  np++;
1239  }
1240  }
1241 
1242  if (np == 0) {
1243  ostringstream os;
1244  os << "For applying propmat_clearskyAddParticles, *abs_species* needs to"
1245  << "contain species 'particles', but it does not.\n";
1246  throw runtime_error(os.str());
1247  }
1248 
1249  if (ns != np) {
1250  ostringstream os;
1251  os << "Number of 'particles' entries in abs_species and of elements in\n"
1252  << "*scat_data* needs to be identical. But you have " << np
1253  << " 'particles' entries\n"
1254  << "and " << ns << " *scat_data* elements.\n";
1255  throw runtime_error(os.str());
1256  }
1257 
1258  if (atmosphere_dim == 1 && rtp_los.nelem() < 1) {
1259  ostringstream os;
1260  os << "For applying *propmat_clearskyAddParticles*, *rtp_los* needs to be specified\n"
1261  << "(at least zenith angle component for atmosphere_dim==1),\n"
1262  << "but it is not.\n";
1263  throw runtime_error(os.str());
1264  } else if (atmosphere_dim > 1 && rtp_los.nelem() < 2) {
1265  ostringstream os;
1266  os << "For applying *propmat_clearskyAddParticles*, *rtp_los* needs to be specified\n"
1267  << "(both zenith and azimuth angle components for atmosphere_dim>1),\n"
1268  << "but it is not.\n";
1269  throw runtime_error(os.str());
1270  }
1271 
1272  // Use for rescaling vmr of particulates
1273  Numeric rtp_vmr_sum = 0.0;
1274 
1275  // Tests and setup partial derivatives
1276  const bool do_jac = supports_particles(jacobian_quantities);
1277  const bool do_jac_temperature = do_temperature_jacobian(jacobian_quantities);
1278  const bool do_jac_frequencies = do_frequency_jacobian(jacobian_quantities);
1279  const ArrayOfIndex jacobian_quantities_position =
1282 
1283  const Index na = abs_species.nelem();
1284  Vector rtp_los_back;
1285  mirror_los(rtp_los_back, rtp_los, atmosphere_dim);
1286 
1287  // 170918 JM: along with transition to use of new-type (aka
1288  // pre-f_grid-interpolated) scat_data, freq perturbation switched off. Typical
1289  // clear-sky freq perturbations yield insignificant effects in particle
1290  // properties. Hence, this feature is neglected here.
1291  if (do_jac_frequencies) {
1292  out1 << "WARNING:\n"
1293  << "Frequency perturbation not available for absorbing particles.\n";
1294  }
1295 
1296  // creating temporary output containers
1297  ArrayOfArrayOfTensor5 ext_mat_Nse;
1298  ArrayOfArrayOfTensor4 abs_vec_Nse;
1299  ArrayOfArrayOfIndex ptypes_Nse;
1300  Matrix t_ok;
1301 
1302  // preparing input in format needed
1303  Vector T_array;
1304  if (do_jac_temperature) {
1305  T_array.resize(2);
1306  T_array = rtp_temperature;
1307  T_array[1] += dT;
1308  } else {
1309  T_array.resize(1);
1310  T_array = rtp_temperature;
1311  }
1312  Matrix dir_array(1, 2);
1313  dir_array(0, joker) = rtp_los_back;
1314 
1315  // ext/abs per scat element for all freqs at once
1316  opt_prop_NScatElems(ext_mat_Nse,
1317  abs_vec_Nse,
1318  ptypes_Nse,
1319  t_ok,
1320  scat_data,
1321  stokes_dim,
1322  T_array,
1323  dir_array,
1324  -1);
1325 
1326  const Index nf = abs_vec_Nse[0][0].nbooks();
1327  Tensor3 tmp(nf, stokes_dim, stokes_dim);
1328 
1329  // loop over the scat_data and link them with correct vmr_field entry according
1330  // to the position of the particle type entries in abs_species.
1331  Index sp = 0;
1332  Index i_se_flat = 0;
1333  for (Index i_ss = 0; i_ss < scat_data.nelem(); i_ss++) {
1334  for (Index i_se = 0; i_se < scat_data[i_ss].nelem(); i_se++) {
1335  // forward to next particle entry in abs_species
1336  while (sp < na && abs_species[sp][0].Type() != SpeciesTag::TYPE_PARTICLES)
1337  sp++;
1338 
1339  // running beyond number of abs_species entries when looking for next
1340  // particle entry. shouldn't happen, though.
1341  assert(sp < na);
1342  if (rtp_vmr[sp] < 0.) {
1343  ostringstream os;
1344  os << "Negative absorbing particle 'vmr' (aka number density)"
1345  << " encountered:\n"
1346  << "scat species #" << i_ss << ", scat elem #" << i_se
1347  << " (vmr_field entry #" << sp << ")\n";
1348  throw runtime_error(os.str());
1349  } else if (rtp_vmr[sp] > 0.) {
1350  if (t_ok(i_se_flat, 0) < 0.) {
1351  ostringstream os;
1352  os << "Temperature interpolation error:\n"
1353  << "scat species #" << i_ss << ", scat elem #" << i_se << "\n";
1354  throw runtime_error(os.str());
1355  } else {
1356  if (use_abs_as_ext) {
1357  if (nf > 1)
1358  for (Index iv = 0; iv < f_grid.nelem(); iv++)
1359  propmat_clearsky[sp].AddAbsorptionVectorAtPosition(
1360  abs_vec_Nse[i_ss][i_se](iv, 0, 0, joker), iv);
1361  else
1362  for (Index iv = 0; iv < f_grid.nelem(); iv++)
1363  propmat_clearsky[sp].AddAbsorptionVectorAtPosition(
1364  abs_vec_Nse[i_ss][i_se](0, 0, 0, joker), iv);
1365  } else {
1366  if (nf > 1)
1367  for (Index iv = 0; iv < f_grid.nelem(); iv++)
1368  propmat_clearsky[sp].SetAtPosition(
1369  ext_mat_Nse[i_ss][i_se](iv, 0, 0, joker, joker), iv);
1370  else
1371  for (Index iv = 0; iv < f_grid.nelem(); iv++)
1372  propmat_clearsky[sp].SetAtPosition(
1373  ext_mat_Nse[i_ss][i_se](0, 0, 0, joker, joker), iv);
1374  }
1375  propmat_clearsky[sp] *= rtp_vmr[sp];
1376  }
1377 
1378  // For temperature derivatives (so we don't need to check it in jac loop)
1379  if (do_jac_temperature) {
1380  if (t_ok(i_se_flat, 1) < 0.) {
1381  ostringstream os;
1382  os << "Temperature interpolation error (in perturbation):\n"
1383  << "scat species #" << i_ss << ", scat elem #" << i_se << "\n";
1384  throw runtime_error(os.str());
1385  }
1386  }
1387 
1388  // For number density derivatives
1389  if (do_jac) rtp_vmr_sum += rtp_vmr[sp];
1390 
1391  for (Index iq = 0; iq < jacobian_quantities_position.nelem(); iq++) {
1392  // we don't do freq perturbations here, i.e. nothing to do.
1393  /*
1394  if( ppd(iq)==JQT_frequency || ppd(iq)==JQT_wind_magnitude ||
1395  ppd(iq)==JQT_wind_u || ppd(iq)==JQT_wind_v || ppd(iq)==JQT_wind_w )
1396  {
1397  dpropmat_clearsky_dx[iq] = 0.;
1398  }
1399 
1400  else if( ppd(iq) == JQT_temperature )
1401  */
1402  if (jacobian_quantities[jacobian_quantities_position[iq]] ==
1403  JacPropMatType::Temperature) {
1404  if (use_abs_as_ext) {
1405  tmp(joker, joker, 0) =
1406  abs_vec_Nse[i_ss][i_se](joker, 1, 0, joker);
1407  tmp(joker, joker, 0) -=
1408  abs_vec_Nse[i_ss][i_se](joker, 0, 0, joker);
1409  } else {
1410  tmp = ext_mat_Nse[i_ss][i_se](joker, 1, 0, joker, joker);
1411  tmp -= ext_mat_Nse[i_ss][i_se](joker, 0, 0, joker, joker);
1412  }
1413 
1414  tmp *= rtp_vmr[sp];
1415  tmp /= dT;
1416 
1417  if (nf > 1)
1418  for (Index iv = 0; iv < f_grid.nelem(); iv++)
1419  if (use_abs_as_ext)
1420  dpropmat_clearsky_dx[iq].AddAbsorptionVectorAtPosition(
1421  tmp(iv, joker, 0), iv);
1422  else
1423  dpropmat_clearsky_dx[iq].AddAtPosition(tmp(iv, joker, joker),
1424  iv);
1425  else
1426  for (Index iv = 0; iv < f_grid.nelem(); iv++)
1427  if (use_abs_as_ext)
1428  dpropmat_clearsky_dx[iq].AddAbsorptionVectorAtPosition(
1429  tmp(0, joker, 0), iv);
1430  else
1431  dpropmat_clearsky_dx[iq].AddAtPosition(tmp(0, joker, joker),
1432  iv);
1433  }
1434  /* // alternative version
1435  else if( ppd(iq) == JQT_temperature )
1436  {
1437  if( nf > 1 )
1438  for( Index iv=0; iv<f_grid.nelem(); iv++ )
1439  if( use_abs_as_ext )
1440  {
1441  atmp = abs_vec_Nse[i_ss][i_se](iv,1,0,joker);
1442  atmp -= abs_vec_Nse[i_ss][i_se](iv,0,0,joker);
1443  atmp *= (rtp_vmr[sp] / ppd.Temperature_Perturbation());
1444  dpropmat_clearsky_dx[iq].AddAbsorptionVectorAtPosition(atmp, iv);
1445  }
1446  else
1447  {
1448  etmp = ext_mat_Nse[i_ss][i_se](iv,1,0,joker,joker) -
1449  ext_mat_Nse[i_ss][i_se](iv,0,0,joker,joker);
1450  etmp *= (rtp_vmr[sp] / ppd.Temperature_Perturbation());
1451  dpropmat_clearsky_dx[iq].AddAtPosition(etmp, iv);
1452  }
1453  else
1454  if( use_abs_as_ext )
1455  {
1456  atmp = (abs_vec_Nse[i_ss][i_se](0,1,0,joker) -
1457  abs_vec_Nse[i_ss][i_se](0,0,0,joker));
1458  atmp *= (rtp_vmr[sp] / ppd.Temperature_Perturbation());
1459  for( Index iv=0; iv<f_grid.nelem(); iv++ )
1460  dpropmat_clearsky_dx[iq].AddAbsorptionVectorAtPosition(atmp, iv);
1461  }
1462  else
1463  {
1464  etmp = ext_mat_Nse[i_ss][i_se](0,1,0,joker,joker) -
1465  ext_mat_Nse[i_ss][i_se](0,0,0,joker,joker);
1466  etmp *= (rtp_vmr[sp] / ppd.Temperature_Perturbation());
1467  for( Index iv=0; iv<f_grid.nelem(); iv++ )
1468  dpropmat_clearsky_dx[iq].AddAtPosition(etmp, iv);
1469  }
1470  }
1471  */
1472 
1473  else if (jacobian_quantities[jacobian_quantities_position[iq]] ==
1474  JacPropMatType::Particulates) {
1475  for (Index iv = 0; iv < f_grid.nelem(); iv++)
1476  dpropmat_clearsky_dx[iq].AddAtPosition(propmat_clearsky[sp], iv);
1477  }
1478  }
1479  }
1480  sp++;
1481  i_se_flat++;
1482  }
1483  }
1484 
1485  //checking that no further 'particle' entry left after all scat_data entries
1486  //are processes. this is basically not necessary. but checking it anyway to
1487  //really be safe. remove later, when more extensively tested.
1488  while (sp < na) {
1489  assert(abs_species[sp][0].Type() != SpeciesTag::TYPE_PARTICLES);
1490  sp++;
1491  }
1492 
1493  if (rtp_vmr_sum != 0.0) {
1494  for (Index iq = 0; iq < jacobian_quantities_position.nelem(); iq++) {
1495  if (jacobian_quantities[jacobian_quantities_position[iq]] ==
1496  JacPropMatType::Particulates) {
1497  dpropmat_clearsky_dx[iq] /= rtp_vmr_sum;
1498  }
1499  }
1500  }
1501 }
1502 
1503 /* Workspace method: Doxygen documentation will be auto-generated */
1504 void propmat_clearskyAddOnTheFly( // Workspace reference:
1505  Workspace& ws,
1506  // WS Output:
1512  // WS Input:
1513  const Vector& f_grid,
1516  const Numeric& rtp_pressure,
1517  const Numeric& rtp_temperature,
1518  const EnergyLevelMap& rtp_nlte,
1519  const Vector& rtp_vmr,
1520  const Agenda& abs_xsec_agenda,
1521  // Verbosity object:
1522  const Verbosity& verbosity) {
1523 
1524  // Output of AbsInputFromRteScalars:
1525  Vector abs_p;
1526  Vector abs_t;
1527  Matrix abs_vmrs;
1528  // Output of abs_h2oSet:
1529  Vector abs_h2o;
1530  // Output of abs_coefCalc:
1533  dsrc_coef_dx;
1534 
1536  abs_t,
1537  abs_vmrs,
1538  rtp_pressure,
1540  rtp_vmr,
1541  verbosity);
1542 
1543  // Absorption cross sections per tag group.
1547 
1548  // Make all species active:
1550  for (Index i = 0; i < abs_species.nelem(); ++i) abs_species_active[i] = i;
1551 
1552  // Call agenda to calculate absorption:
1558  abs_species,
1561  f_grid,
1562  abs_p,
1563  abs_t,
1564  rtp_nlte,
1565  abs_vmrs,
1566  abs_xsec_agenda);
1567 
1568  // Calculate absorption coefficients from cross sections:
1570  src_coef,
1571  dabs_coef_dx,
1572  dsrc_coef_dx,
1579  abs_species,
1581  abs_vmrs,
1582  abs_p,
1583  abs_t,
1584  verbosity);
1585 
1586  // Now add abs_coef_per_species to propmat_clearsky:
1590  dabs_coef_dx,
1591  verbosity);
1592 
1593  // Now turn nlte_source from absorption into a proper source function
1594  if (not nlte_source.empty())
1599  dsrc_coef_dx,
1601  f_grid,
1603  verbosity);
1604 }
1605 
1606 /* Workspace method: Doxygen documentation will be auto-generated */
1608  const Vector& f_grid,
1609  const Index& stokes_dim,
1610  const Verbosity&) {
1611  propmat_clearsky.resize(1);
1613  propmat_clearsky[0].SetZero();
1614 }
1615 
1616 /* Workspace method: Doxygen documentation will be auto-generated */
1619  for (auto& pm : propmat_clearsky)
1620  for (Index i = 0; i < pm.NumberOfFrequencies(); i++)
1621  if (pm.Kjj()[i] < 0.0) pm.SetAtPosition(0.0, i);
1622 }
1623 
1624 /* Workspace method: Doxygen documentation will be auto-generated */
1626  const Verbosity&) {
1628 }
1629 
1630 /* Workspace method: Doxygen documentation will be auto-generated */
1632  const Verbosity&) {
1634 }
1635 
1636 #ifdef ENABLE_NETCDF
1637 /* Workspace method: Doxygen documentation will be auto-generated */
1638 /* Included by Claudia Emde, 20100707 */
1639 void WriteMolTau( //WS Input
1640  const Vector& f_grid,
1641  const Tensor3& z_field,
1643  const Index& atmosphere_dim,
1644  //Keyword
1645  const String& filename,
1646  const Verbosity&) {
1647  int retval, ncid;
1648  int nlev_dimid, nlyr_dimid, nwvl_dimid, stokes_dimid, none_dimid;
1649  int dimids[4];
1650  int wvlmin_varid, wvlmax_varid, z_varid, wvl_varid, tau_varid;
1651 
1652  if (atmosphere_dim != 1)
1653  throw runtime_error("WriteMolTau can only be used for atmosphere_dim=1");
1654 
1655 #pragma omp critical(netcdf__critical_region)
1656  {
1657  // Open file
1658  if ((retval = nc_create(filename.c_str(), NC_CLOBBER, &ncid)))
1659  nca_error(retval, "nc_create");
1660 
1661  // Define dimensions
1662  if ((retval = nc_def_dim(ncid, "nlev", (int)z_field.npages(), &nlev_dimid)))
1663  nca_error(retval, "nc_def_dim");
1664 
1665  if ((retval =
1666  nc_def_dim(ncid, "nlyr", (int)z_field.npages() - 1, &nlyr_dimid)))
1667  nca_error(retval, "nc_def_dim");
1668 
1669  if ((retval = nc_def_dim(ncid, "nwvl", (int)f_grid.nelem(), &nwvl_dimid)))
1670  nca_error(retval, "nc_def_dim");
1671 
1672  if ((retval = nc_def_dim(ncid, "none", 1, &none_dimid)))
1673  nca_error(retval, "nc_def_dim");
1674 
1675  if ((retval = nc_def_dim(ncid,
1676  "nstk",
1677  (int)propmat_clearsky_field.nbooks(),
1678  &stokes_dimid)))
1679  nca_error(retval, "nc_def_dim");
1680 
1681  // Define variables
1682  if ((retval = nc_def_var(
1683  ncid, "wvlmin", NC_DOUBLE, 1, &none_dimid, &wvlmin_varid)))
1684  nca_error(retval, "nc_def_var wvlmin");
1685 
1686  if ((retval = nc_def_var(
1687  ncid, "wvlmax", NC_DOUBLE, 1, &none_dimid, &wvlmax_varid)))
1688  nca_error(retval, "nc_def_var wvlmax");
1689 
1690  if ((retval = nc_def_var(ncid, "z", NC_DOUBLE, 1, &nlev_dimid, &z_varid)))
1691  nca_error(retval, "nc_def_var z");
1692 
1693  if ((retval =
1694  nc_def_var(ncid, "wvl", NC_DOUBLE, 1, &nwvl_dimid, &wvl_varid)))
1695  nca_error(retval, "nc_def_var wvl");
1696 
1697  dimids[0] = nlyr_dimid;
1698  dimids[1] = nwvl_dimid;
1699  dimids[2] = stokes_dimid;
1700  dimids[3] = stokes_dimid;
1701 
1702  if ((retval =
1703  nc_def_var(ncid, "tau", NC_DOUBLE, 4, &dimids[0], &tau_varid)))
1704  nca_error(retval, "nc_def_var tau");
1705 
1706  // Units
1707  if ((retval = nc_put_att_text(ncid, wvlmin_varid, "units", 2, "nm")))
1708  nca_error(retval, "nc_put_att_text");
1709 
1710  if ((retval = nc_put_att_text(ncid, wvlmax_varid, "units", 2, "nm")))
1711  nca_error(retval, "nc_put_att_text");
1712 
1713  if ((retval = nc_put_att_text(ncid, z_varid, "units", 2, "km")))
1714  nca_error(retval, "nc_put_att_text");
1715 
1716  if ((retval = nc_put_att_text(ncid, wvl_varid, "units", 2, "nm")))
1717  nca_error(retval, "nc_put_att_text");
1718 
1719  if ((retval = nc_put_att_text(ncid, tau_varid, "units", 1, "-")))
1720  nca_error(retval, "nc_put_att_text");
1721 
1722  // End define mode. This tells netCDF we are done defining
1723  // metadata.
1724  if ((retval = nc_enddef(ncid))) nca_error(retval, "nc_enddef");
1725 
1726  // Assign data
1727  double wvlmin[1];
1728  wvlmin[0] = SPEED_OF_LIGHT / f_grid[f_grid.nelem() - 1] * 1e9;
1729  if ((retval = nc_put_var_double(ncid, wvlmin_varid, &wvlmin[0])))
1730  nca_error(retval, "nc_put_var");
1731 
1732  double wvlmax[1];
1733  wvlmax[0] = SPEED_OF_LIGHT / f_grid[0] * 1e9;
1734  if ((retval = nc_put_var_double(ncid, wvlmax_varid, &wvlmax[0])))
1735  nca_error(retval, "nc_put_var");
1736 
1737  double z[z_field.npages()];
1738  for (int iz = 0; iz < z_field.npages(); iz++)
1739  z[iz] = z_field(z_field.npages() - 1 - iz, 0, 0) * 1e-3;
1740 
1741  if ((retval = nc_put_var_double(ncid, z_varid, &z[0])))
1742  nca_error(retval, "nc_put_var");
1743 
1744  double wvl[f_grid.nelem()];
1745  for (int iv = 0; iv < f_grid.nelem(); iv++)
1746  wvl[iv] = SPEED_OF_LIGHT / f_grid[f_grid.nelem() - 1 - iv] * 1e9;
1747 
1748  if ((retval = nc_put_var_double(ncid, wvl_varid, &wvl[0])))
1749  nca_error(retval, "nc_put_var");
1750 
1751  const Index zfnp = z_field.npages() - 1;
1752  const Index fgne = f_grid.nelem();
1753  const Index amfnb = propmat_clearsky_field.nbooks();
1754 
1755  Tensor4 tau(zfnp, fgne, amfnb, amfnb, 0.);
1756 
1757  // Calculate average tau for layers
1758  for (int is = 0; is < propmat_clearsky_field.nlibraries(); is++)
1759  for (int iz = 0; iz < zfnp; iz++)
1760  for (int iv = 0; iv < fgne; iv++)
1761  for (int is1 = 0; is1 < amfnb; is1++)
1762  for (int is2 = 0; is2 < amfnb; is2++)
1763  // sum up all species
1764  tau(iz, iv, is1, is2) +=
1765  0.5 *
1767  f_grid.nelem() - 1 - iv,
1768  is1,
1769  is2,
1770  z_field.npages() - 1 - iz,
1771  0,
1772  0) +
1774  f_grid.nelem() - 1 - iv,
1775  is1,
1776  is2,
1777  z_field.npages() - 2 - iz,
1778  0,
1779  0)) *
1780  (z_field(z_field.npages() - 1 - iz, 0, 0) -
1781  z_field(z_field.npages() - 2 - iz, 0, 0));
1782 
1783  if ((retval = nc_put_var_double(ncid, tau_varid, tau.get_c_array())))
1784  nca_error(retval, "nc_put_var");
1785 
1786  // Close the file
1787  if ((retval = nc_close(ncid))) nca_error(retval, "nc_close");
1788  }
1789 }
1790 
1791 #else
1792 
1793 void WriteMolTau( //WS Input
1794  const Vector& f_grid _U_,
1795  const Tensor3& z_field _U_,
1797  const Index& atmosphere_dim _U_,
1798  //Keyword
1799  const String& filename _U_,
1800  const Verbosity&) {
1801  throw runtime_error(
1802  "The workspace method WriteMolTau is not available"
1803  "because ARTS was compiled without NetCDF support.");
1804 }
1805 
1806 #endif /* ENABLE_NETCDF */
1807 
1808 /* Workspace method: Doxygen documentation will be auto-generated */
1810  // WS Output:
1815  // WS Input:
1819  const Vector& f_grid,
1820  const Vector& abs_p,
1821  const Vector& abs_t,
1822  const EnergyLevelMap& abs_nlte,
1823  const Matrix& abs_vmrs,
1827  const Index& lbl_checked,
1828  const Verbosity&) {
1829  if (not abs_lines_per_species.nelem()) return;
1830 
1831  if (not lbl_checked)
1832  throw std::runtime_error("Please set lbl_checked true to use this function");
1833 
1834  // Check that all temperatures are above 0 K
1835  if (min(abs_t) < 0) {
1836  std::ostringstream os;
1837  os << "Temperature must be at least 0 K. But you request an absorption\n"
1838  << "calculation at " << min(abs_t) << " K!";
1839  throw std::runtime_error(os.str());
1840  }
1841 
1842  // Check that all parameters that should have the number of tag
1843  // groups as a dimension are consistent.
1844  {
1845  const Index n_tgs = abs_species.nelem();
1846  const Index n_xsec = abs_xsec_per_species.nelem();
1847  const Index n_vmrs = abs_vmrs.nrows();
1848  const Index n_lines = abs_lines_per_species.nelem();
1849 
1850  if (n_tgs not_eq n_xsec or n_tgs not_eq n_vmrs or n_tgs not_eq n_lines) {
1851  std::ostringstream os;
1852  os << "The following variables must all have the same dimension:\n"
1853  << "abs_species: " << abs_species.nelem() << '\n'
1854  << "abs_xsec_per_species: " << abs_xsec_per_species.nelem() << '\n'
1855  << "abs_vmrs: " << abs_vmrs.nrows() << '\n'
1856  << "abs_lines_per_species: " << abs_lines_per_species.nelem() << '\n';
1857  throw std::runtime_error(os.str());
1858  }
1859  }
1860 
1861  // Meta variables that explain the calculations required
1862  const bool do_jac = supports_propmat_clearsky(jacobian_quantities);
1863  const bool do_lte = abs_nlte.Data().empty();
1865 
1866  // Skipping uninteresting data
1867  static Matrix dummy1(0, 0);
1868  static ArrayOfMatrix dummy2(0);
1869 
1870  // Call xsec_species for each tag group.
1871  for (Index ii = 0; ii < abs_species_active.nelem(); ++ii) {
1872  const Index i = abs_species_active[ii];
1873 
1874  if (not abs_species[i].nelem() or is_zeeman(abs_species[i]))
1875  continue;
1876 
1877  for (auto& lines: abs_lines_per_species[i]) {
1878  xsec_species(
1881  dummy1,
1882  do_jac ? dabs_xsec_per_species_dx[i] : dummy2,
1883  (do_jac and not do_lte) ? dsrc_xsec_per_species_dx[i] : dummy2,
1884  dummy2,
1886  jac_pos,
1887  f_grid,
1888  abs_p,
1889  abs_t,
1890  abs_nlte,
1891  abs_vmrs,
1892  abs_species,
1893  lines,
1894  isotopologue_ratios.getIsotopologueRatio(lines.QuantumIdentity()),
1895  partition_functions.getParamType(lines.QuantumIdentity()),
1896  partition_functions.getParam(lines.QuantumIdentity()));
1897  }
1898  } // End of species for loop.
1899 }
Matrix
The Matrix class.
Definition: matpackI.h:1193
legacy_continua.h
supports_continuum
bool supports_continuum(const ArrayOfRetrievalQuantity &js)
Returns if the array supports continuum derivatives.
Definition: jacobian.cc:1208
ARTS::Var::atmosphere_dim
Index atmosphere_dim(Workspace &ws) noexcept
Definition: autoarts.h:2510
ARTS::Var::z_field
Tensor3 z_field(Workspace &ws) noexcept
Definition: autoarts.h:7690
propmat_clearskyAddFaraday
void propmat_clearskyAddFaraday(ArrayOfPropagationMatrix &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:1074
magnetic_field_perturbation
Numeric magnetic_field_perturbation(const ArrayOfRetrievalQuantity &js) noexcept
Returns the magnetic field perturbation if it exists.
Definition: jacobian.cc:1320
ARTS::Var::rtp_mag
Vector rtp_mag(Workspace &ws) noexcept
Definition: autoarts.h:5736
ARTS::Group::PropagationMatrix
PropagationMatrix PropagationMatrix
Definition: autoarts.h:93
ARTS::Var::rtp_los
Vector rtp_los(Workspace &ws) noexcept
Definition: autoarts.h:5715
propmat_clearskyAddFromAbsCoefPerSpecies
void propmat_clearskyAddFromAbsCoefPerSpecies(ArrayOfPropagationMatrix &propmat_clearsky, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, const ArrayOfMatrix &abs_coef_per_species, const ArrayOfMatrix &dabs_coef_dx, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyAddFromAbsCoefPerSpecies.
Definition: m_abs.cc:968
StokesVector
Stokes vector is as Propagation matrix but only has 4 possible values.
Definition: propagationmatrix.h:1075
auto_md.h
ARTS::Var::jacobian_quantities
ArrayOfRetrievalQuantity jacobian_quantities(Workspace &ws) noexcept
Definition: autoarts.h:3924
absorption.h
Declarations required for the calculation of absorption coefficients.
Absorption::Lines::AppendSingleLine
void AppendSingleLine(SingleLine &&sl)
Appends a single line to the absorption lines.
Definition: absorptionlines.h:700
Tensor3
The Tensor3 class.
Definition: matpackIII.h:339
species_iso_match
bool species_iso_match(const RetrievalQuantity &rq, const Index species, const Index iso)
Returns if the Retrieval quantity is VMR derivative for all the species in the species tags.
Definition: jacobian.cc:1268
VACUUM_PERMITTIVITY
const Numeric VACUUM_PERMITTIVITY
fillSpeciesAuxDataWithPartitionFunctionsFromSpeciesData
void fillSpeciesAuxDataWithPartitionFunctionsFromSpeciesData(SpeciesAuxData &sad)
Fill SpeciesAuxData with default partition functions from species data.
Definition: absorption.cc:417
abs_speciesDefineAll
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:214
SpeciesDataOfBand
const SpeciesRecord & SpeciesDataOfBand(const AbsorptionLines &band)
Returns the species data.
Definition: absorption.cc:611
Absorption::Lines::NumLines
Index NumLines() const noexcept
Number of lines.
Definition: absorptionlines.h:852
joker
const Joker joker
nc_io.h
This file contains basic functions to handle NetCDF data files.
dplanck_dt
Numeric dplanck_dt(const Numeric &f, const Numeric &t)
dplanck_dt
Definition: physics_funcs.cc:337
abs
#define abs(x)
Definition: legacy_continua.cc:20626
ARTS::Var::verbosity
Verbosity verbosity(Workspace &ws) noexcept
Definition: autoarts.h:7112
ARTS::Var::dnlte_dx_source
ArrayOfStokesVector dnlte_dx_source(Workspace &ws) noexcept
Definition: autoarts.h:3021
ARTS::Var::dsrc_xsec_per_species_dx
ArrayOfArrayOfMatrix dsrc_xsec_per_species_dx(Workspace &ws) noexcept
Definition: autoarts.h:3321
temperature_perturbation
Numeric temperature_perturbation(const ArrayOfRetrievalQuantity &js) noexcept
Returns the temperature perturbation if it exists.
Definition: jacobian.cc:1304
check_continuum_model
void check_continuum_model(const String &name)
An auxiliary functions that checks if a given continuum model is listed in species_data....
Definition: legacy_continua.cc:20406
PropagationMatrix
Definition: propagationmatrix.h:87
Vector::resize
void resize(Index n)
Resize function.
Definition: matpackI.cc:404
ARTS::Var::abs_coef_per_species
ArrayOfMatrix abs_coef_per_species(Workspace &ws) noexcept
Definition: autoarts.h:1800
ARTS::Var::abs_lines_per_species
ArrayOfArrayOfAbsorptionLines abs_lines_per_species(Workspace &ws) noexcept
Definition: autoarts.h:2022
propmat_clearskyForceNegativeToZero
void propmat_clearskyForceNegativeToZero(ArrayOfPropagationMatrix &propmat_clearsky, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyForceNegativeToZero.
Definition: m_abs.cc:1617
ARTS::Var::rtp_vmr
Vector rtp_vmr(Workspace &ws) noexcept
Definition: autoarts.h:5826
propmat_clearskyAddParticles
void propmat_clearskyAddParticles(ArrayOfPropagationMatrix &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:1202
ARTS::Var::stokes_dim
Index stokes_dim(Workspace &ws) noexcept
Definition: autoarts.h:6650
supports_particles
bool supports_particles(const ArrayOfRetrievalQuantity &js)
Returns if the array supports particulate derivatives.
Definition: jacobian.cc:1236
ELECTRON_CHARGE
const Numeric ELECTRON_CHARGE
do_temperature_jacobian
bool do_temperature_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants the temperature derivative.
Definition: jacobian.cc:1279
global_data::species_data
const Array< SpeciesRecord > species_data
Species Data.
Definition: partition_function_data.cc:42
ARTS::Var::abs_xsec_per_species
ArrayOfMatrix abs_xsec_per_species(Workspace &ws) noexcept
Definition: autoarts.h:2299
CREATE_OUT2
#define CREATE_OUT2
Definition: messages.h:206
Absorption::Lines::ReverseLines
void ReverseLines() noexcept
Reverses the order of the internal lines.
Definition: absorptionlines.cc:2720
montecarlo.h
TotalNumberOfElements
Index TotalNumberOfElements(const Array< Array< base > > &aa)
Determine total number of elements in an ArrayOfArray.
Definition: array.h:343
Tensor4
The Tensor4 class.
Definition: matpackIV.h:421
array.h
This file contains the definition of Array.
abs_xsec_per_speciesAddLines
void abs_xsec_per_speciesAddLines(ArrayOfMatrix &abs_xsec_per_species, ArrayOfMatrix &src_xsec_per_species, ArrayOfArrayOfMatrix &dabs_xsec_per_species_dx, ArrayOfArrayOfMatrix &dsrc_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 EnergyLevelMap &abs_nlte, const Matrix &abs_vmrs, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const SpeciesAuxData &isotopologue_ratios, const SpeciesAuxData &partition_functions, const Index &lbl_checked, const Verbosity &)
WORKSPACE METHOD: abs_xsec_per_speciesAddLines.
Definition: m_abs.cc:1809
SpeciesAuxData
Auxiliary data for isotopologues.
Definition: absorption.h:217
frequency_perturbation
Numeric frequency_perturbation(const ArrayOfRetrievalQuantity &js) noexcept
Returns the frequency perturbation if it exists.
Definition: jacobian.cc:1312
CREATE_OUT1
#define CREATE_OUT1
Definition: messages.h:205
Agenda
The Agenda class.
Definition: agenda_class.h:44
ARTS::Var::abs_vmrs
Matrix abs_vmrs(Workspace &ws) noexcept
Definition: autoarts.h:2267
ARTS::Var::abs_coef
Matrix abs_coef(Workspace &ws) noexcept
Definition: autoarts.h:1790
ARTS::Var::propmat_clearsky_field
Tensor7 propmat_clearsky_field(Workspace &ws) noexcept
Definition: autoarts.h:5433
matpackI.h
Implementation of Matrix, Vector, and such stuff.
xsec_continuum_tag
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.
Definition: legacy_continua.cc:16340
opt_prop_NScatElems
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.
Definition: optproperties.cc:248
Array
This can be used to make arrays out of anything.
Definition: array.h:108
SpeciesTag
A tag group can consist of the sum of several of these.
Definition: abs_species_tags.h:44
number_density
Numeric number_density(const Numeric &p, const Numeric &t)
number_density
Definition: physics_funcs.cc:219
Absorption::nelem
Index nelem(const Lines &l)
Number of lines.
Definition: absorptionlines.h:1820
is_zeeman
bool is_zeeman(const ArrayOfSpeciesTag &tg)
Is this a Zeeman tag group?
Definition: abs_species_tags.cc:702
min
#define min(a, b)
Definition: legacy_continua.cc:20628
set_vmr_from_first_species
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:597
SPEED_OF_LIGHT
const Numeric SPEED_OF_LIGHT
CREATE_OUT3
#define CREATE_OUT3
Definition: messages.h:207
abs_speciesDefineAllInScenario
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
messages.h
Declarations having to do with the four output streams.
propmat_clearskyAddOnTheFly
void propmat_clearskyAddOnTheFly(Workspace &ws, ArrayOfPropagationMatrix &propmat_clearsky, ArrayOfStokesVector &nlte_source, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, ArrayOfStokesVector &dnlte_dx_source, ArrayOfStokesVector &nlte_dsource_dx, const Vector &f_grid, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfRetrievalQuantity &jacobian_quantities, const Numeric &rtp_pressure, const Numeric &rtp_temperature, const EnergyLevelMap &rtp_nlte, const Vector &rtp_vmr, const Agenda &abs_xsec_agenda, const Verbosity &verbosity)
WORKSPACE METHOD: propmat_clearskyAddOnTheFly.
Definition: m_abs.cc:1504
ARTS::Var::rtp_nlte
EnergyLevelMap rtp_nlte(Workspace &ws) noexcept
Definition: autoarts.h:5754
_U_
#define _U_
Definition: config.h:183
my_basic_string< char >
abs_xsec_agendaExecute
void abs_xsec_agendaExecute(Workspace &ws, ArrayOfMatrix &abs_xsec_per_species, ArrayOfMatrix &src_xsec_per_species, ArrayOfArrayOfMatrix &dabs_xsec_per_species_dx, ArrayOfArrayOfMatrix &dsrc_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 EnergyLevelMap &abs_nlte, const Matrix &abs_vmrs, const Agenda &input_agenda)
Definition: auto_md.cc:23582
equivalent_propmattype_indexes
ArrayOfIndex equivalent_propmattype_indexes(const ArrayOfRetrievalQuantity &js)
Returns a list of positions for the derivatives in Propagation Matrix calculations.
Definition: jacobian.cc:1099
ARTS::Var::abs_cont_parameters
ArrayOfVector abs_cont_parameters(Workspace &ws) noexcept
Definition: autoarts.h:1978
ARTS::Var::nlte_source
ArrayOfStokesVector nlte_source(Workspace &ws) noexcept
Definition: autoarts.h:4648
WriteMolTau
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:1639
ARTS::Var::vmr_field
Tensor4 vmr_field(Workspace &ws) noexcept
Definition: autoarts.h:7130
physics_funcs.h
This file contains declerations of functions of physical character.
abs_xsec_per_speciesAddConts
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:531
ARTS::Var::abs_species
ArrayOfArrayOfSpeciesTag abs_species(Workspace &ws) noexcept
Definition: autoarts.h:2157
ArrayOfPropagationMatrix
Array< PropagationMatrix > ArrayOfPropagationMatrix
Definition: propagationmatrix.h:978
ARTS::Var::p_grid
Vector p_grid(Workspace &ws) noexcept
Definition: autoarts.h:4763
optproperties.h
Scattering database structure and functions.
SpeciesTag::TYPE_FREE_ELECTRONS
@ TYPE_FREE_ELECTRONS
Definition: abs_species_tags.h:152
find_xml_file
void find_xml_file(String &filename, const Verbosity &verbosity)
Find an xml file.
Definition: file.cc:414
jacobian.h
Routines for setting up the jacobian.
AbsInputFromRteScalars
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 &)
WORKSPACE METHOD: AbsInputFromRteScalars.
Definition: m_abs.cc:67
dplanck_df
Numeric dplanck_df(const Numeric &f, const Numeric &t)
dplanck_df
Definition: physics_funcs.cc:407
Numeric
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
Verbosity
Definition: messages.h:49
abs_cont_descriptionInit
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:823
math_funcs.h
main
int main(int argc, char **argv)
This is the main function of ARTS.
Definition: main.cc:612
global_data.h
do_frequency_jacobian
bool do_frequency_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants a frequency derivative.
Definition: jacobian.cc:1296
ARTS::Var::f_grid
Vector f_grid(Workspace &ws) noexcept
Definition: autoarts.h:3449
Absorption::Lines
Definition: absorptionlines.h:547
partition_functionsInitFromBuiltin
void partition_functionsInitFromBuiltin(SpeciesAuxData &partition_functions, const Verbosity &)
WORKSPACE METHOD: partition_functionsInitFromBuiltin.
Definition: m_abs.cc:1631
ARTS::Var::src_coef_per_species
ArrayOfMatrix src_coef_per_species(Workspace &ws) noexcept
Definition: autoarts.h:6625
EnergyLevelMap
Definition: energylevelmap.h:60
ARTS::Var::abs_xsec_agenda
Agenda abs_xsec_agenda(Workspace &ws) noexcept
Definition: autoarts.h:2274
ArrayOfMatrix
Array< Matrix > ArrayOfMatrix
An array of matrices.
Definition: array.h:66
abs_speciesSet
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.
Definition: m_abs_lookup.cc:1974
ARTS::Var::dabs_coef_dx
ArrayOfMatrix dabs_coef_dx(Workspace &ws) noexcept
Definition: autoarts.h:2968
max
#define max(a, b)
Definition: legacy_continua.cc:20629
ARTS::Var::scat_data
ArrayOfArrayOfSingleScatteringData scat_data(Workspace &ws) noexcept
Definition: autoarts.h:5840
Matrix::resize
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1056
AbsInputFromAtmFields
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:234
LineShape::Type
Type
Type of line shape to compute.
Definition: lineshapemodel.h:788
abs_coefCalcFromXsec
void abs_coefCalcFromXsec(Matrix &abs_coef, Matrix &src_coef, ArrayOfMatrix &dabs_coef_dx, ArrayOfMatrix &dsrc_coef_dx, ArrayOfMatrix &abs_coef_per_species, ArrayOfMatrix &src_coef_per_species, const ArrayOfMatrix &abs_xsec_per_species, const ArrayOfMatrix &src_xsec_per_species, const ArrayOfArrayOfMatrix &dabs_xsec_per_species_dx, const ArrayOfArrayOfMatrix &dsrc_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)
WORKSPACE METHOD: abs_coefCalcFromXsec.
Definition: m_abs.cc:258
PI
const Numeric PI
parameters.h
This file contains header information for the dealing with command line parameters.
m_xml.h
Workspace methods and template functions for supergeneric XML IO.
ARTS::Var::dsrc_coef_dx
ArrayOfMatrix dsrc_coef_dx(Workspace &ws) noexcept
Definition: autoarts.h:3303
ARTS::Var::rtp_pressure
Numeric rtp_pressure(Workspace &ws) noexcept
Definition: autoarts.h:5790
supports_faraday
bool supports_faraday(const ArrayOfRetrievalQuantity &js)
Returns if the array supports Faraday derivatives.
Definition: jacobian.cc:1230
nlte_sourceFromTemperatureAndSrcCoefPerSpecies
void nlte_sourceFromTemperatureAndSrcCoefPerSpecies(ArrayOfStokesVector &nlte_source, ArrayOfStokesVector &dnlte_dx_source, ArrayOfStokesVector &nlte_dsource_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 &)
WORKSPACE METHOD: nlte_sourceFromTemperatureAndSrcCoefPerSpecies.
Definition: m_abs.cc:862
ARTS::Var::nrows
Index nrows(Workspace &ws) noexcept
Definition: autoarts.h:4682
ARTS::Var::abs_t
Vector abs_t(Workspace &ws) noexcept
Definition: autoarts.h:2186
ARTS::Var::nlte_do
Index nlte_do(Workspace &ws) noexcept
Definition: autoarts.h:4567
ARTS::Var::scat_data_checked
Index scat_data_checked(Workspace &ws) noexcept
Definition: autoarts.h:5850
ArrayOfStokesVector
Array< StokesVector > ArrayOfStokesVector
Definition: propagationmatrix.h:1354
abs_cont_descriptionAppend
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:839
Isotopologue
QuantumIdentifier::QType Isotopologue
Definition: arts_api_classes.cc:242
planck
Numeric planck(const Numeric &f, const Numeric &t)
planck
Definition: physics_funcs.cc:269
Range
The range class.
Definition: matpackI.h:160
ARTS::Var::abs_species_active
ArrayOfIndex abs_species_active(Workspace &ws) noexcept
Definition: autoarts.h:2170
ELECTRON_MASS
const Numeric ELECTRON_MASS
ARTS::Var::abs_lines
ArrayOfAbsorptionLines abs_lines(Workspace &ws) noexcept
Definition: autoarts.h:2013
ARTS::Var::lbl_checked
Index lbl_checked(Workspace &ws) noexcept
Definition: autoarts.h:3996
dnumber_density_dt
Numeric dnumber_density_dt(const Numeric &p, const Numeric &t)
dnumber_density_dT
Definition: physics_funcs.cc:243
fillSpeciesAuxDataWithIsotopologueRatiosFromSpeciesData
void fillSpeciesAuxDataWithIsotopologueRatiosFromSpeciesData(SpeciesAuxData &sad)
Fill SpeciesAuxData with default isotopologue ratios from species data.
Definition: absorption.cc:394
ARTS::Var::propmat_clearsky_agenda_checked
Index propmat_clearsky_agenda_checked(Workspace &ws) noexcept
Definition: autoarts.h:5414
ARTS::Var::dpropmat_clearsky_dx
ArrayOfPropagationMatrix dpropmat_clearsky_dx(Workspace &ws) noexcept
Definition: autoarts.h:3269
ARTS::Var::rtp_temperature
Numeric rtp_temperature(Workspace &ws) noexcept
Definition: autoarts.h:5807
abs_xsec_per_speciesInit
void abs_xsec_per_speciesInit(ArrayOfMatrix &abs_xsec_per_species, ArrayOfMatrix &src_xsec_per_species, ArrayOfArrayOfMatrix &dabs_xsec_per_species_dx, ArrayOfArrayOfMatrix &dsrc_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 Index &nlte_do, const Verbosity &verbosity)
WORKSPACE METHOD: abs_xsec_per_speciesInit.
Definition: m_abs.cc:443
do_magnetic_jacobian
bool do_magnetic_jacobian(const ArrayOfRetrievalQuantity &js) noexcept
Returns if the array wants a magnetic derivative.
Definition: jacobian.cc:1300
ARTS::Var::src_xsec_per_species
ArrayOfMatrix src_xsec_per_species(Workspace &ws) noexcept
Definition: autoarts.h:6641
mirror_los
void mirror_los(Vector &los_mirrored, ConstVectorView los, const Index &atmosphere_dim)
Determines the backward direction for a given line-of-sight.
Definition: rte.cc:2290
dotprod_with_los
Numeric dotprod_with_los(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:730
propmat_clearskyInit
void propmat_clearskyInit(ArrayOfPropagationMatrix &propmat_clearsky, ArrayOfStokesVector &nlte_source, ArrayOfPropagationMatrix &dpropmat_clearsky_dx, ArrayOfStokesVector &dnlte_dx_source, ArrayOfStokesVector &nlte_dsource_dx, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfRetrievalQuantity &jacobian_quantities, const Vector &f_grid, const Index &stokes_dim, const Index &propmat_clearsky_agenda_checked, const Index &nlte_do, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyInit.
Definition: m_abs.cc:1028
rte.h
Declaration of functions in rte.cc.
ARTS::Var::t_field
Tensor3 t_field(Workspace &ws) noexcept
Definition: autoarts.h:6947
Workspace
Workspace class.
Definition: workspace_ng.h:40
PropagationMatrix::Kjj
VectorView Kjj(const Index iz=0, const Index ia=0)
Vector view to diagonal elements.
Definition: propagationmatrix.h:799
ARTS::Var::src_coef
Matrix src_coef(Workspace &ws) noexcept
Definition: autoarts.h:6615
ARTS::Var::abs_xsec_agenda_checked
Index abs_xsec_agenda_checked(Workspace &ws) noexcept
Definition: autoarts.h:2283
ARTS::Var::partition_functions
SpeciesAuxData partition_functions(Workspace &ws) noexcept
Definition: autoarts.h:4860
ARTS::Var::isotopologue_ratios
SpeciesAuxData isotopologue_ratios(Workspace &ws) noexcept
Definition: autoarts.h:3669
SpeciesTag::TYPE_PREDEF
@ TYPE_PREDEF
Definition: abs_species_tags.h:150
file.h
This file contains basic functions to handle ASCII files.
ARTS::Var::abs_cont_models
ArrayOfString abs_cont_models(Workspace &ws) noexcept
Definition: autoarts.h:1812
supports_propmat_clearsky
bool supports_propmat_clearsky(const ArrayOfRetrievalQuantity &js)
Returns if the array supports propagation matrix derivatives.
Definition: jacobian.cc:1240
Absorption::createEmptyCopy
Lines createEmptyCopy(const Lines &al) noexcept
Creates a copy of the input lines structure.
Definition: absorptionlines.cc:2701
chk_size
void chk_size(const String &x_name, ConstVectorView x, const Index &c)
Runtime check for size of Vector.
Definition: check_input.cc:466
Index
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
ARTS::Var::abs_p
Vector abs_p(Workspace &ws) noexcept
Definition: autoarts.h:2129
ARTS::Var::dabs_xsec_per_species_dx
ArrayOfArrayOfMatrix dabs_xsec_per_species_dx(Workspace &ws) noexcept
Definition: autoarts.h:2986
check_input.h
ARTS::Var::nlte_dsource_dx
ArrayOfStokesVector nlte_dsource_dx(Workspace &ws) noexcept
Definition: autoarts.h:4581
SpeciesTag::TYPE_CIA
@ TYPE_CIA
Definition: abs_species_tags.h:151
ARTS::Var::ncols
Index ncols(Workspace &ws) noexcept
Definition: autoarts.h:4546
xsec_species
void xsec_species(Matrix &xsec, Matrix &source, Matrix &phase, ArrayOfMatrix &dxsec_dx, ArrayOfMatrix &dsource_dx, ArrayOfMatrix &dphase_dx, const ArrayOfRetrievalQuantity &jacobian_quantities, const ArrayOfIndex &jacobian_propmat_positions, const Vector &f_grid, const Vector &abs_p, const Vector &abs_t, const EnergyLevelMap &abs_nlte, const Matrix &abs_vmrs, const ArrayOfArrayOfSpeciesTag &abs_species, const AbsorptionLines &band, const Numeric &isot_ratio, const SpeciesAuxData::AuxType &partfun_type, const ArrayOfGriddedField1 &partfun_data)
Cross-section algorithm.
Definition: absorption.cc:616
is_frequency_parameter
bool is_frequency_parameter(const RetrievalQuantity &t) noexcept
Returns if the Retrieval quantity is a frequency parameter in propagation matrix calculations.
Definition: jacobian.cc:1120
Vector
The Vector class.
Definition: matpackI.h:860
ARTS::Var::propmat_clearsky
ArrayOfPropagationMatrix propmat_clearsky(Workspace &ws) noexcept
Definition: autoarts.h:5398
ARTS::Var::abs_cont_names
ArrayOfString abs_cont_names(Workspace &ws) noexcept
Definition: autoarts.h:1966
abs_lines_per_speciesCreateFromLines
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:90
nca_error
void nca_error(const int e, const String s)
Throws a runtime error for the given NetCDF error code.
Definition: nc_io.cc:622
isotopologue_ratiosInitFromBuiltin
void isotopologue_ratiosInitFromBuiltin(SpeciesAuxData &isotopologue_ratios, const Verbosity &)
WORKSPACE METHOD: isotopologue_ratiosInitFromBuiltin.
Definition: m_abs.cc:1625
SpeciesTag::TYPE_PARTICLES
@ TYPE_PARTICLES
Definition: abs_species_tags.h:153
Tensor4View::get_c_array
const Numeric * get_c_array() const
Conversion to plain C-array.
Definition: matpackIV.cc:358
species_match
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:1244
Array::nelem
Index nelem() const
Number of elements.
Definition: array.h:195
ARTS::Group::StokesVector
StokesVector StokesVector
Definition: autoarts.h:101
JacPropMatType::VMR
@ VMR
propmat_clearskyZero
void propmat_clearskyZero(ArrayOfPropagationMatrix &propmat_clearsky, const Vector &f_grid, const Index &stokes_dim, const Verbosity &)
WORKSPACE METHOD: propmat_clearskyZero.
Definition: m_abs.cc:1607
Tensor7
The Tensor7 class.
Definition: matpackVII.h:2382
ns
#define ns
Definition: legacy_continua.cc:22183
arts.h
The global header file for ARTS.
xml_io.h
This file contains basic functions to handle XML data files.
ARTS::Var::abs_nlte
EnergyLevelMap abs_nlte(Workspace &ws) noexcept
Definition: autoarts.h:2111