28#include "matpack_math.h"
36 auto &lineshape = lines[k].lineshape;
39 return Output{lineshape.
G0(T, T0, P, vmrs), lineshape.D0(T, T0, P, vmrs),
40 lineshape.G2(T, T0, P, vmrs), lineshape.D2(T, T0, P, vmrs),
41 lineshape.FVC(T, T0, P, vmrs), lineshape.ETA(T, T0, P, vmrs),
42 lineshape.Y(T, T0, P, vmrs), lineshape.G(T, T0, P, vmrs),
43 lineshape.DV(T, T0, P, vmrs)}
44 .no_linemixing(not DoLineMixing(P));
51 auto &lineshape = lines[k].lineshape[pos];
53 return lineshape.at(T, T0, P).no_linemixing(not DoLineMixing(P));
57 size_t k, Numeric T, Numeric P,
const Vector& vmrs)
const ARTS_NOEXCEPT {
58 auto &lineshape = lines[k].lineshape;
61 return Output{lineshape.dG0dT(T, T0, P, vmrs), lineshape.dD0dT(T, T0, P, vmrs),
62 lineshape.dG2dT(T, T0, P, vmrs), lineshape.dD2dT(T, T0, P, vmrs),
63 lineshape.dFVCdT(T, T0, P, vmrs), lineshape.dETAdT(T, T0, P, vmrs),
64 lineshape.dYdT(T, T0, P, vmrs), lineshape.dGdT(T, T0, P, vmrs),
65 lineshape.dDVdT(T, T0, P, vmrs)}
70 size_t k, Numeric T, Numeric P,
size_t pos)
const ARTS_NOEXCEPT {
71 auto &lineshape = lines[k].lineshape[pos];
73 return lineshape.dT(T, T0, P).no_linemixing(not DoLineMixing(P));
79 if (selfbroadening and spec == quantumidentity.Species()) {
84 const Index s = selfbroadening;
85 const Index e = broadeningspecies.nelem() - bathbroadening;
86 for (Index i = s; i < e; i++) {
87 if (spec == broadeningspecies[i]) {
93 if (bathbroadening)
return broadeningspecies.nelem() - 1;
102 auto &lineshape = lines[k].lineshape;
104 const Index pos = LineShapePos(vmr_qid.Species());
108 out = lineshape[pos].at(T, T0, P);
109 const Index bath = lineshape.nelem() - 1;
110 if (bathbroadening and pos not_eq bath) {
111 out -= lineshape[bath].at(T, T0, P);
135 if (is.eof())
return data;
147 if (line.
nelem() == 0 && is.eof())
return data;
160 istringstream icecream(line);
165 if (artsid.length() != 0) {
169 isotopologue.is_joker() or
170 isotopologue.type not_eq Species::TagType::Plain,
171 "A line catalog species can only be of the form \"Plain\", meaning it\nhas the form SPECIES-ISONUM.\n"
172 "Your input contains: ",
174 ". which we cannot interpret as a plain species")
216 for (Index j = 0; j < naux; j++) {
223 Numeric unused_numeric;
231 }
catch (
const std::runtime_error&) {
235 if (tgam != data.
T0) {
236 agam = agam * pow(tgam / data.
T0, nair);
237 sgam = sgam * pow(tgam / data.
T0, nself);
238 psf = psf * pow(tgam / data.
T0, (Numeric).25 + (Numeric)1.5 * nair);
267 if (is.eof())
return data;
279 if (line.
nelem() == 0 && is.eof())
return data;
292 istringstream icecream(line);
297 if (artsid.length() != 0) {
301 isotopologue.is_joker() or
302 isotopologue.type not_eq Species::TagType::Plain,
303 "A line catalog species can only be of the form \"Plain\", meaning it\nhas the form SPECIES-ISONUM.\n"
304 "Your input contains: ",
306 ". which we cannot interpret as a plain species")
350 bool lmd_found =
false;
362 if (is.eof())
return data;
374 if (line.
nelem() == 0 && is.eof())
return data;
387 istringstream icecream(line);
393 if (artsid.length() != 0) {
397 isotopologue.is_joker() or
398 isotopologue.type not_eq Species::TagType::Plain,
399 "A line catalog species can only be of the form \"Plain\", meaning it\nhas the form SPECIES-ISONUM.\n"
400 "Your input contains: ",
402 ". which we cannot interpret as a plain species")
442 }
else if (token ==
"QN") {
447 token !=
"UP",
"Unknown quantum number tag: ", token)
451 while (icecream and token not_eq
"LO") {
452 auto qn = Quantum::Number::toType(token);
468 !is || token !=
"LO",
469 "Error in catalog. Lower quantum number tag 'LO' not found.")
472 auto qn = Quantum::Number::toType(token);
474 not(token ==
"LM" or token ==
"LF" or token ==
"ZM" or
475 token ==
"LSM" or token ==
"PB")) {
489 qn = Quantum::Number::toType(token);
491 }
else if (token ==
"LM") {
495 }
else if (token ==
"LF") {
504 }
else if (token ==
"ZM") {
508 }
else if (token ==
"LSM") {
513 for (Index lsm = 0; lsm <
nelem; lsm++) {
517 if (token ==
"CUT") {
519 data.
cutoff = CutoffType::ByLine;
523 if (token ==
"LML") {
528 else if (token ==
"MTM") {
533 else if (token ==
"LNT") {
548 }
catch (
const std::runtime_error& e) {
549 ARTS_USER_ERROR(
"Parse error in catalog line: \n", line,
'\n', e.what())
568 data.
species[1] = Species::Species::Bath;
583 if (is.eof())
return data;
595 if (line.
nelem() == 0 && is.eof())
return data;
599 if (line[line.
nelem() - 1] == 13) {
600 line.erase(line.
nelem() - 1, 1);
653 data.
line.
I0 = s * hi2arts;
677 agam = gam * hi2arts;
683 sgam = gam * hi2arts;
686 if (0 == sgam) sgam = agam;
832 data.
species[1] = Species::Species::Bath;
847 if (is.eof())
return data;
859 if (line.
nelem() == 0 && is.eof())
return data;
863 if (line[line.
nelem() - 1] == 13) {
864 line.erase(line.
nelem() - 1, 1);
917 data.
line.
I0 = s * hi2arts;
941 agam = gam * hi2arts;
947 sgam = gam * hi2arts;
950 if (0 == sgam) sgam = agam;
1085 std::stringstream ss;
1087 ss >> upper >> lower;
1103 data.
species[1] = Species::Species::Bath;
1114 bool comment =
true;
1118 if (is.eof())
return data;
1130 if (line.
nelem() == 0 && is.eof())
return data;
1134 if (line[line.
nelem() - 1] == 13) {
1135 line.erase(line.
nelem() - 1, 1);
1189 data.
line.
I0 = s * hi2arts;
1216 agam = gam * hi2arts;
1222 sgam = gam * hi2arts;
1225 if (0 == sgam) sgam = agam;
1243 Numeric nair, nself;
1351 bool comment =
true;
1355 if (is.eof())
return data;
1362 if (line[0] ==
'>' or line[0] ==
'%')
continue;
1368 if (line.
nelem() == 0 && is.eof())
return data;
1372 if (line[line.
nelem() - 1] == 13) {
1373 line.erase(line.
nelem() - 1, 1);
1422 if (line[6] ==
'D') line[6] =
'E';
1428 data.
line.
I0 = s * hi2arts;
1455 agam = gam * hi2arts;
1461 sgam = gam * hi2arts;
1464 if (0 == sgam) sgam = agam;
1482 Numeric nair, nself;
1582 if (test == -1 || test == -3)
1605 Vector Y(4), G(4), T(4);
1715 std::vector<SingleLineExternal>& external_lines,
1716 const std::vector<QuantumNumberType>& localquantas,
1717 const std::vector<QuantumNumberType>& globalquantas) {
1718 std::vector<Lines> lines(0);
1721 while (external_lines.size()) {
1722 auto& sle = external_lines.back();
1725 if (sle.selfbroadening) sle.species.front() = sle.quantumidentity.Species();
1726 if (sle.bathbroadening) sle.species.back() = Species::Species::Bath;
1730 for (
auto qn : globalquantas)
1731 if (sle.quantumidentity.val.has(qn))
1732 global_id.val.set(sle.quantumidentity.val[qn]);
1735 for (
auto qn : localquantas)
1736 if (sle.quantumidentity.val.has(qn))
1737 local_id.
val.
set(sle.quantumidentity.val[qn]);
1740 auto line = sle.line;
1741 line.localquanta = local_id;
1744 auto band = std::find_if(lines.begin(), lines.end(), [&](
const Lines& li) {
1745 return li.MatchWithExternal(sle, global_id);
1747 if (band not_eq lines.end()) {
1748 band->AppendSingleLine(line);
1750 lines.push_back(
Lines(sle.selfbroadening,
1759 sle.linemixinglimit,
1764 external_lines.pop_back();
1772 for (
auto& line : lines.
lines) os << line <<
'\n';
1777 for (
auto& line : lines.
lines) is >> line;
1782 os << line.
F0 <<
' ' << line.
I0 <<
' ' << line.
E0 <<
' ' << line.
glow <<
' '
1783 << line.
gupp <<
' ' << line.
A <<
' ' << line.
zeeman <<
' '
1791 line.
gupp >> line.
A;
1801 std::ostringstream os;
1803 os <<
"\nLines meta-data:\n";
1804 os <<
'\t' <<
"Species identity:\n";
1805 os <<
"\t\tSpecies: " << SpeciesName() <<
'\n';
1806 os <<
"\t\tIdentity: " << quantumidentity <<
'\n';
1808 os <<
'\t' << populationtype2metadatastring(population);
1809 os <<
'\t' << normalizationtype2metadatastring(normalization);
1810 os <<
'\t' << LineShape::shapetype2metadatastring(lineshapetype);
1811 os <<
'\t' << mirroringtype2metadatastring(mirroring);
1812 os <<
'\t' <<
"The reference temperature for all line parameters is " << T0
1814 if (linemixinglimit < 0)
1815 os <<
'\t' <<
"If applicable, there is no line mixing limit.\n";
1817 os <<
'\t' <<
"If applicable, there is a line mixing limit at "
1818 << linemixinglimit <<
" Pa.\n";
1820 if (not NumLines()) {
1821 os <<
"\tNo line data is available.\n";
1823 os <<
"\tThere are " << NumLines() <<
" lines available.\n";
1825 auto& line = lines.front();
1826 os <<
"\tThe front line has:\n";
1828 <<
"f0: " << line.F0 <<
" Hz\n";
1830 <<
"i0: " << line.I0 <<
" m^2/Hz\n";
1832 <<
"e0: " << line.E0 <<
" J\n";
1834 <<
"Lower stat. weight: " << line.glow <<
" [-]\n";
1836 <<
"Upper stat. weight: " << line.gupp <<
" [-]\n";
1838 <<
"A: " << line.A <<
" 1/s\n";
1840 <<
"Zeeman splitting of lower state: " << line.zeeman.gl() <<
" [-]\n";
1842 <<
"Zeeman splitting of upper state: " << line.zeeman.gu() <<
" [-]\n";
1844 <<
"Local quantum numbers: " << line.localquanta.val <<
"\n";
1847 line.lineshape, selfbroadening, broadeningspecies, T0);
1849 <<
"Line shape parameters (are normalized by sum(VMR)):\n";
1850 for (
auto& ls_form : ls_meta) os <<
"\t\t\t" << ls_form <<
"\n";
1857 lines.erase(lines.begin() + i);
1861 auto line = lines[i];
1867 std::reverse(lines.begin(), lines.end());
1871 return quantumidentity.Isotopologue().mass;
1875 const ConstVectorView& atm_vmrs,
1881 const ConstVectorView& atm_vmrs,
1884 const Numeric& bath_mass)
const {
1885 Vector mass =
LineShape::mass(atm_vmrs, atm_spec, broadeningspecies, ir);
1886 if (bathbroadening and bath_mass > 0) mass[mass.nelem() - 1] = bath_mass;
1891 const ConstVectorView& atm_vmrs,
1894 "Bad species and vmr lists");
1896 for (Index i = 0; i < atm_spec.
nelem(); i++)
1903 Rational Jf, Rational Ji, Rational lf, Rational li, Rational k) {
1904 if (not iseven(Jf + lf + 1))
1905 return -sqrt(2 * Jf + 1) *
wigner3j(Jf, k, Ji, li, lf - li, -lf);
1906 return +sqrt(2 * Jf + 1) *
wigner3j(Jf, k, Ji, li, lf - li, -lf);
1912 if (not iseven(Jf + N))
1913 return -sqrt(6 * (2 * Jf + 1) * (2 * Ji + 1)) *
1915 return +sqrt(6 * (2 * Jf + 1) * (2 * Ji + 1)) *
wigner6j(1, 1, 1, Ji, Jf, N);
1932 bool comment =
true;
1936 if (is.eof())
return data;
1948 if (line.
nelem() == 0 && is.eof())
return data;
1993 s = pow((Numeric)10., s);
2034 tag = tag > 0 ? tag : -tag;
2060 Numeric nair, nself;
2093 const Index nb = broadeningspecies.nelem();
2096 if (not Isotopologue().OK())
return false;
2099 if (nb < (Index(selfbroadening) + Index(bathbroadening)))
return false;
2102 if (T0 <= 0)
return false;
2105 if (std::any_of(lines.cbegin(), lines.cend(), [nb](
auto& line) {
2106 return line.LineShapeElems() != nb;
2111 if (std::any_of(lines.cbegin(), lines.cend(), [&](
auto& line) {
2112 return line.LocalQuantumElems() != lines.front().LocalQuantumElems();
2127 std::ostringstream os;
2129 case CutoffType::None:
2130 os <<
"No cut-off will be applied.\n";
2132 case CutoffType::ByLine:
2133 os <<
"The lines will be cut-off " <<
cutoff
2134 <<
" Hz from the line center + D0.\n";
2136 case CutoffType::FINAL:
2157 LineShape::TemperatureModel::LM_AER,
d[4],
d[5],
d[6],
d[7]};
2159 LineShape::TemperatureModel::LM_AER,
d[8],
d[9],
d[10],
d[11]};
2195 "Error calling appending function, bad size of quantum numbers\n"
2196 "Type of quantum numbers in band: ",
2197 lines.front().localquanta.val,
2198 "\nType of quantum numbers in new line:",
2199 sl.localquanta.val);
2203 sl.LineShapeElems() not_eq
lines[0].LineShapeElems(),
2204 "Error calling appending function, bad size of broadening species");
2206 lines.push_back(std::move(sl));
2216 if (sle.bad)
return false;
2217 if (sle.selfbroadening not_eq selfbroadening)
return false;
2218 if (sle.bathbroadening not_eq bathbroadening)
return false;
2219 if (sle.cutoff not_eq
cutoff)
return false;
2220 if (sle.mirroring not_eq mirroring)
return false;
2221 if (sle.population not_eq population)
return false;
2222 if (sle.normalization not_eq normalization)
return false;
2223 if (sle.lineshapetype not_eq lineshapetype)
return false;
2224 if (sle.T0 not_eq T0)
return false;
2225 if (sle.cutofffreq not_eq cutofffreq)
return false;
2226 if (sle.linemixinglimit not_eq linemixinglimit)
return false;
2227 if (quantumidentity not_eq qid)
return false;
2228 if (not std::equal(sle.species.cbegin(),
2230 broadeningspecies.cbegin(),
2231 broadeningspecies.cend()))
2233 if (NumLines() not_eq 0 and
2234 not lines.front().localquanta.same_types_as(sle.line.localquanta))
2236 if (NumLines() not_eq 0 and
2237 not sle.line.lineshape.Match(lines.front().lineshape).first)
2244 if (l.selfbroadening not_eq selfbroadening)
return {
false,
false};
2245 if (l.bathbroadening not_eq bathbroadening)
return {
false,
false};
2246 if (l.cutoff not_eq
cutoff)
return {
false,
false};
2247 if (l.mirroring not_eq mirroring)
return {
false,
false};
2248 if (l.population not_eq population)
return {
false,
false};
2249 if (l.normalization not_eq normalization)
return {
false,
false};
2250 if (l.lineshapetype not_eq lineshapetype)
return {
false,
false};
2251 if (l.T0 not_eq T0)
return {
false,
false};
2252 if (l.cutofffreq not_eq cutofffreq)
return {
false,
false};
2253 if (l.linemixinglimit not_eq linemixinglimit)
return {
false,
false};
2254 if (l.quantumidentity not_eq quantumidentity)
return {
false,
false};
2255 if (not std::equal(l.broadeningspecies.cbegin(),
2256 l.broadeningspecies.cend(),
2257 broadeningspecies.cbegin(),
2258 broadeningspecies.cend()))
2259 return {
false,
false};
2260 if (NumLines() not_eq 0 and l.NumLines() not_eq 0 and
2261 not lines.front().localquanta.same_types_as(l.lines.front().localquanta))
2262 return {
false,
false};
2263 if (NumLines() not_eq 0 and l.NumLines() not_eq 0) {
2264 if (
auto matchpair =
2265 l.lines.front().lineshape.Match(lines.front().lineshape);
2266 not matchpair.first)
2270 return {
true,
true};
2281 std::sort(
lines.begin(),
2299 return lines.size() ?
lines.front().localquanta.val.nelem() : 0;
2303 return population == PopulationType::ByMakarovFullRelmat or
2304 population == PopulationType::ByRovibLinearDipoleLineMixing;
2308 return linemixinglimit < 0 ? true : linemixinglimit > P;
2312 return qid.Isotopologue() == quantumidentity.Isotopologue() or
2313 (qid.Isotopologue().joker() and
2314 qid.Species() == quantumidentity.Species()) or
2315 std::any_of(broadeningspecies.begin(), broadeningspecies.end(),
2316 [s = qid.Species()](
auto &
a) { return a == s; });
2320 for (
auto &line :
lines) {
2321 for (
auto &shape : line.lineshape.Data()) {
2322 if (shape.Y().type not_eq LineShape::TemperatureModel::None or
2323 shape.G().type not_eq LineShape::TemperatureModel::None or
2324 shape.DV().type not_eq LineShape::TemperatureModel::None) {
2334 std::find(broadeningspecies.cbegin(), broadeningspecies.cend(), spec);
2335 ptr not_eq broadeningspecies.cend())
2336 return std::distance(broadeningspecies.cbegin(), ptr);
2343 qns.val.has(QuantumNumberType::J))
2344 return qns.val.has(QuantumNumberType::F) ? qns.val[QuantumNumberType::F]
2345 : qns.val[QuantumNumberType::J];
2352 auto& val =
get(lines[k].localquanta);
2362 auto& val =
get(lines[k].localquanta);
2363 return lines[k].zeeman.Strength(val.upp(), val.low(), type, i);
2372 auto& val =
get(lines[k].localquanta);
2373 return lines[k].zeeman.Splitting(val.upp(), val.low(), type, i);
2382 std::inner_product(lines.cbegin(),
2387 [](
const auto&
a,
const auto&
b) { return a.F0 * b; });
2388 const Numeric div = sum(wgts);
2395 const Index n = NumLines();
2398 const Numeric ratiopart = QT0 / QT;
2401 for (Index i = 0; i < n; i++) {
2402 const Numeric pop0 =
2404 const Numeric pop = pop0 * ratiopart *
boltzman_ratio(T, T0, lines[i].E0);
2405 const Numeric dip_squared =
2407 (pop0 * lines[i].F0 *
2409 wgts[i] = pop * dip_squared;
2412 return F_mean(wgts);
2417 case CutoffType::ByLine:
2418 return lines[k].F0 + cutofffreq +
2419 (mirroring == MirroringType::Manual ? -shift : shift);
2420 case CutoffType::None:
2421 return std::numeric_limits<Numeric>::max();
2422 case CutoffType::FINAL:
2425 return std::numeric_limits<Numeric>::max();
2430 case CutoffType::ByLine:
2431 return lines[k].F0 - cutofffreq +
2432 (mirroring == MirroringType::Manual ? -shift : shift);
2433 case CutoffType::None:
2434 return std::numeric_limits<Numeric>::lowest();
2435 case CutoffType::FINAL:
2439 return std::numeric_limits<Numeric>::lowest();
2443 for (
auto& line :
lines) line.read(is);
2448 for (
auto& line :
lines) line.write(os);
2454 for (
auto&
a : lines[k].localquanta.val) qid.
val.
set(
a);
2463 for (Index j = 0; j < m; j++) {
2466 LineShape::TemperatureModel t = LineShape::TemperatureModel::None;
2467 for (
auto& line :
lines) {
2468 if (
auto& data = line.lineshape[j].Data()[i];
2470 if (t == LineShape::TemperatureModel::None) t = data.type;
2473 "Cannot make a common line shape model for the band as there are multiple non-empty types: ",
2481 for (
auto& line :
lines) {
2482 if (
auto& data = line.lineshape[j].Data()[i]; data.type not_eq t) {
2493 for (
auto& abs_lines : abs_lines_per_species) {
2494 for (
auto& band : abs_lines) {
2495 switch (band.population) {
2496 case Absorption::PopulationType::LTE:
2499 case Absorption::PopulationType::NLTE:
2502 case Absorption::PopulationType::VibTemps:
2505 case Absorption::PopulationType::ByHITRANFullRelmat:
2508 case Absorption::PopulationType::ByHITRANRosenkranzRelmat:
2511 case Absorption::PopulationType::ByMakarovFullRelmat:
2514 case Absorption::PopulationType::ByRovibLinearDipoleLineMixing:
2517 case Absorption::PopulationType::FINAL: {
2526 for (
auto& abs_lines : abs_lines_per_species) {
2527 for (
auto& band : abs_lines) {
2528 switch (band.cutoff) {
2529 case Absorption::CutoffType::None:
2532 case Absorption::CutoffType::ByLine:
2535 case Absorption::CutoffType::FINAL: {
2544 for (
auto& abs_lines : abs_lines_per_species) {
2545 for (
auto& band : abs_lines) {
2546 switch (band.lineshapetype) {
2547 case LineShape::Type::DP:
2550 case LineShape::Type::LP:
2553 case LineShape::Type::VP:
2556 case LineShape::Type::SDVP:
2559 case LineShape::Type::HTP:
2562 case LineShape::Type::SplitLP:
2565 case LineShape::Type::SplitVP:
2568 case LineShape::Type::SplitSDVP:
2571 case LineShape::Type::SplitHTP:
2574 case LineShape::Type::FINAL: {
2583 for (
auto& abs_lines : abs_lines_per_species) {
2584 for (
auto& band : abs_lines) {
2585 switch (band.mirroring) {
2586 case Absorption::MirroringType::None:
2589 case Absorption::MirroringType::Lorentz:
2592 case Absorption::MirroringType::SameAsLineShape:
2595 case Absorption::MirroringType::Manual:
2598 case Absorption::MirroringType::FINAL: {
2607 for (
auto& abs_lines : abs_lines_per_species) {
2608 for (
auto& band : abs_lines) {
2609 switch (band.normalization) {
2610 case Absorption::NormalizationType::None:
2613 case Absorption::NormalizationType::VVH:
2616 case Absorption::NormalizationType::VVW:
2619 case Absorption::NormalizationType::RQ:
2622 case Absorption::NormalizationType::SFS:
2625 case Absorption::NormalizationType::FINAL: {
2633template <
typename T>
2634constexpr std::size_t req_spaces(T my_enum) {
2635 constexpr std::size_t n = []() {
2636 std::size_t longest = 0;
2637 for (
auto& x : Absorption::enumstrs::CutoffTypeNames) {
2638 longest = std::max(x.length(), longest);
2640 for (
auto& x : Absorption::enumstrs::MirroringTypeNames) {
2641 longest = std::max(x.length(), longest);
2643 for (
auto& x : Absorption::enumstrs::NormalizationTypeNames) {
2644 longest = std::max(x.length(), longest);
2646 for (
auto& x : Absorption::enumstrs::PopulationTypeNames) {
2647 longest = std::max(x.length(), longest);
2649 for (
auto& x : LineShape::enumstrs::TypeNames) {
2650 longest = std::max(x.length(), longest);
2654 return n -
toString(my_enum).length();
2657auto spaces(std::size_t n) {
return std::basic_string(n,
' '); }
2662 Absorption::CutoffType x{Absorption::CutoffType::FINAL};
2664 case Absorption::CutoffType::FINAL:
2665 os <<
"Cutoff tag types:\n";
2667 case AbsorptionCutoffType::ByLine:
2668 os <<
" ByLine:" << spaces(req_spaces(AbsorptionCutoffType::ByLine))
2671 case AbsorptionCutoffType::None:
2672 os <<
" None:" << spaces(req_spaces(AbsorptionCutoffType::None))
2681 Absorption::MirroringType x{Absorption::MirroringType::FINAL};
2683 case Absorption::MirroringType::FINAL:
2684 os <<
"Mirroring tag types:\n";
2686 case Absorption::MirroringType::None:
2687 os <<
" None:" << spaces(req_spaces(Absorption::MirroringType::None))
2688 << val.
None <<
'\n';
2690 case Absorption::MirroringType::Lorentz:
2692 << spaces(req_spaces(Absorption::MirroringType::Lorentz))
2695 case Absorption::MirroringType::SameAsLineShape:
2696 os <<
" SameAsLineShape:"
2697 << spaces(req_spaces(Absorption::MirroringType::SameAsLineShape))
2700 case Absorption::MirroringType::Manual:
2702 << spaces(req_spaces(Absorption::MirroringType::Manual)) << val.
Manual;
2710 Absorption::NormalizationType x{Absorption::NormalizationType::FINAL};
2712 case Absorption::NormalizationType::FINAL:
2713 os <<
"Normalization tag types:\n";
2715 case Absorption::NormalizationType::None:
2717 << spaces(req_spaces(Absorption::NormalizationType::None)) << val.
None
2720 case Absorption::NormalizationType::VVH:
2721 os <<
" VVH:" << spaces(req_spaces(Absorption::NormalizationType::VVH))
2724 case Absorption::NormalizationType::VVW:
2725 os <<
" VVW:" << spaces(req_spaces(Absorption::NormalizationType::VVW))
2728 case Absorption::NormalizationType::RQ:
2729 os <<
" RQ:" << spaces(req_spaces(Absorption::NormalizationType::RQ))
2732 case Absorption::NormalizationType::SFS:
2733 os <<
" SFS:" << spaces(req_spaces(Absorption::NormalizationType::SFS))
2742 Absorption::PopulationType x{Absorption::PopulationType::FINAL};
2744 case Absorption::PopulationType::FINAL:
2745 os <<
"Population tag type:\n";
2747 case Absorption::PopulationType::LTE:
2748 os <<
" LTE:" << spaces(req_spaces(Absorption::PopulationType::LTE))
2751 case Absorption::PopulationType::NLTE:
2752 os <<
" NLTE:" << spaces(req_spaces(Absorption::PopulationType::NLTE))
2753 << val.
NLTE <<
'\n';
2755 case Absorption::PopulationType::VibTemps:
2757 << spaces(req_spaces(Absorption::PopulationType::VibTemps))
2760 case Absorption::PopulationType::ByHITRANRosenkranzRelmat:
2761 os <<
" ByHITRANRosenkranzRelmat:"
2762 << spaces(req_spaces(
2763 Absorption::PopulationType::ByHITRANRosenkranzRelmat))
2766 case Absorption::PopulationType::ByHITRANFullRelmat:
2767 os <<
" ByHITRANFullRelmat:"
2768 << spaces(req_spaces(Absorption::PopulationType::ByHITRANFullRelmat))
2771 case Absorption::PopulationType::ByMakarovFullRelmat:
2772 os <<
" ByMakarovFullRelmat:"
2773 << spaces(req_spaces(Absorption::PopulationType::ByMakarovFullRelmat))
2776 case Absorption::PopulationType::ByRovibLinearDipoleLineMixing:
2777 os <<
" ByRovibLinearDipoleLineMixing:"
2778 << spaces(req_spaces(
2779 Absorption::PopulationType::ByRovibLinearDipoleLineMixing))
2790 case LineShapeType::FINAL:
2791 os <<
"Line shape tag type:\n";
2793 case LineShapeType::DP:
2794 os <<
" DP:" << spaces(req_spaces(LineShapeType::DP)) << val.
DP
2797 case LineShapeType::LP:
2798 os <<
" LP:" << spaces(req_spaces(LineShapeType::LP)) << val.
LP
2801 case LineShapeType::VP:
2802 os <<
" VP:" << spaces(req_spaces(LineShapeType::VP)) << val.
VP
2805 case LineShapeType::SDVP:
2806 os <<
" SDVP:" << spaces(req_spaces(LineShapeType::SDVP)) << val.
SDVP
2809 case LineShapeType::HTP:
2810 os <<
" HTP:" << spaces(req_spaces(LineShapeType::HTP)) << val.
HTP
2813 case LineShapeType::SplitLP:
2814 os <<
" SplitLP:" << spaces(req_spaces(LineShapeType::SplitLP)) << val.
SplitLP
2817 case LineShapeType::SplitVP:
2818 os <<
" SplitVP:" << spaces(req_spaces(LineShapeType::SplitVP)) << val.
SplitVP
2821 case LineShapeType::SplitSDVP:
2822 os <<
" SplitSDVP:" << spaces(req_spaces(LineShapeType::SplitSDVP)) << val.
SplitSDVP
2825 case LineShapeType::SplitHTP:
2826 os <<
" SplitHTP:" << spaces(req_spaces(LineShapeType::SplitHTP)) << val.
SplitHTP;
2832 return os <<
"Catalog tag summary:\n " << val.
cutoff <<
"\n "
2841 const Index n = abs_species.
nelem();
2844 while (ispec < n and i >= abs_lines_per_species[ispec].nelem()) {
2845 i -= abs_lines_per_species[ispec].
nelem();
2855 abs_lines_per_species.cbegin(),
2856 abs_lines_per_species.cend(),
2857 [](
auto& abs_lines) {
2859 abs_lines.cbegin(), abs_lines.cend(), [](auto& band) {
2860 return band.cutoff not_eq CutoffType::None;
2867 if (quantumidentity.val.has(x)) {
2868 auto& val = quantumidentity.val[x];
2869 return std::max(val.low(), val.upp());
2872 Rational out{std::numeric_limits<Index>::lowest()};
2873 for (
auto& line : lines) {
2875 auto& val = line.localquanta.val[x];
2876 out = std::max(std::max(val.low(), val.upp()), out);
Numeric wavenumber_to_joule(Numeric e)
A little helper function to convert energy from units of wavenumber (cm^-1) to Joule (J).
Declarations required for the calculation of absorption coefficients.
std::ostream & operator<<(std::ostream &os, AbsorptionCutoffTagTypeStatus val)
AbsorptionSpeciesBandIndex flat_index(Index i, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species)
Get a flat index pair for species and band.
Contains the absorption namespace.
This can be used to make arrays out of anything.
Index nelem() const ARTS_NOEXCEPT
Main line shape model class.
void SetLineMixingModel(SingleSpeciesModel x)
Sets the same line mixing model to all species.
bifstream & read(bifstream &bif)
Binary read for Model.
bofstream & write(bofstream &bof) const
Binary write for Model.
const std::vector< SingleSpeciesModel > & Data() const noexcept
The line shape model data.
void set(Value v)
Sets the value if it exists or adds it otherwise.
bool has(Types... ts) const ARTS_NOEXCEPT
Returns whether all the Types are part of the list, the types must be sorted.
Value & add(Type t)
Add for manipulation.
Index nelem() const ARTS_NOEXCEPT
Return number of quantum numbers.
Binary output file stream class.
Binary output file stream class.
Input manipulator class for doubles to enable nan and inf parsing.
Helper macros for debugging.
#define ARTS_ASSERT(condition,...)
#define ARTS_USER_ERROR(...)
#define ARTS_USER_ERROR_IF(condition,...)
constexpr std::string_view toString(EnergyLevelMapType x) noexcept
constexpr bool good_enum(EnumType x) noexcept
Checks if the enum number is good.
This file contains basic functions to handle ASCII files.
Numeric boltzman_factor(Numeric T, Numeric E0)
Computes exp(- E0/kT)
Numeric boltzman_ratio(const Numeric &T, const Numeric &T0, const Numeric &E0)
Computes exp(E0/c (T - T0) / (T * T0))
Numeric single_partition_function(const Numeric &T, const Species::IsotopeRecord &ir)
Computes the partition function at one temperature.
Constains various line scaling functions.
Contains the line shape namespace.
LineShape::Type LineShapeType
void extract(T &x, String &line, std::size_t n)
Extract something from the beginning of a string.
Namespace to contain things required for absorption calculations.
SingleLineExternal ReadFromHitranOnlineStream(istream &is)
Read from HITRAN online.
String cutofftype2metadatastring(CutoffType in, Numeric cutoff)
SingleLineExternal ReadFromJplStream(istream &is)
Read from JPL.
const Quantum::Number::Value & get(const Quantum::Number::LocalState &qns) ARTS_NOEXCEPT
Numeric reduced_rovibrational_dipole(Rational Jf, Rational Ji, Rational lf, Rational li, Rational k=Rational(1))
Compute the reduced rovibrational dipole moment.
SingleLineExternal ReadFromArtscat3Stream(istream &is)
Read from ARTSCAT-3.
Index nelem(const Lines &l)
Number of lines.
SingleLineExternal ReadFromLBLRTMStream(istream &is)
Read from LBLRTM.
std::ostream & operator<<(std::ostream &os, const Absorption::Lines &lines)
std::istream & operator>>(std::istream &is, Lines &lines)
SingleLineExternal ReadFromArtscat5Stream(istream &is)
Read from ARTSCAT-5.
bool any_cutoff(const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species)
Numeric reduced_magnetic_quadrapole(Rational Jf, Rational Ji, Rational N)
Compute the reduced magnetic quadrapole moment.
SingleLineExternal ReadFromHitran2004Stream(istream &is)
Read from newer HITRAN.
SingleLineExternal ReadFromHitran2001Stream(istream &is)
Read from HITRAN before 2004.
SingleLineExternal ReadFromArtscat4Stream(istream &is)
Read from ARTSCAT-4.
std::vector< Lines > split_list_of_external_lines(std::vector< SingleLineExternal > &external_lines, const std::vector< QuantumNumberType > &localquantas={}, const std::vector< QuantumNumberType > &globalquantas={})
Splits a list of lines into proper Lines.
constexpr Numeric doppler_broadening_const_squared
Doppler broadening constant squared [kg/T]^2.
constexpr Numeric k
Boltzmann constant convenience name [J/K].
constexpr Numeric c
Speed of light convenience name [m/s].
constexpr Numeric h
Planck constant convenience name [J s].
constexpr auto atm2pa(auto x) noexcept
Conversion from Atm to Pa.
QuantumIdentifier id_from_lookup(Index mol, char isochar)
Finds the ID of the ARTS species from HITRAN.
Numeric ratio_from_lookup(Index mol, char isochar)
Finds the isotopologue ratio of the species from HITRAN.
QuantumIdentifier id_from_lookup(Index tag)
Finds the ID of the ARTS species from JPL.
Model vector2modellm(Vector x, LegacyLineMixingData::TypeLM type)
LineShape::Model from legacy input vector.
Computations of line shape derived parameters.
Model lblrtm_model(Numeric sgam, Numeric nself, Numeric agam, Numeric nair, Numeric psf, std::array< Numeric, 12 > aer_interp)
constexpr ModelParameters modelparameterGetEmpty(const TemperatureModel t) noexcept
std::istream & from_linefunctiondata(std::istream &data, Type &type, bool &self, bool &bath, Model &m, ArrayOfSpecies &species)
constexpr bool modelparameterEmpty(const ModelParameters mp) noexcept
Model hitran_model(Numeric sgam, Numeric nself, Numeric agam, Numeric nair, Numeric psf)
Vector mass(const ConstVectorView &atmospheric_vmrs, const ArrayOfArrayOfSpeciesTag &atmospheric_species, const ArrayOfSpecies &lineshape_species, const SpeciesIsotopologueRatios &ir) ARTS_NOEXCEPT
Returns a mass vector for this model's main calculations.
std::istream & from_pressurebroadeningdata(std::istream &data, LineShape::Type &type, bool &self, bool &bath, Model &m, ArrayOfSpecies &species, const QuantumIdentifier &qid)
Legacy reading of old deprecated PressureBroadeningData class.
ArrayOfString ModelMetaDataArray(const LineShape::Model &m, const bool self, const ArrayOfSpecies &sts, const Numeric T0)
String ModelShape2MetaData(const Model &m)
std::istream & from_artscat4(std::istream &is, Type &type, bool &self, bool &bath, Model &m, ArrayOfSpecies &species, const QuantumIdentifier &qid)
std::istream & from_linemixingdata(std::istream &data, Model &lsc)
Legacy reading of old deprecated LineMixingData class.
constexpr Index nVars
Current max number of line shape variables.
Vector vmrs(const ConstVectorView &atmospheric_vmrs, const ArrayOfArrayOfSpeciesTag &atmospheric_species, const ArrayOfSpecies &lineshape_species) ARTS_NOEXCEPT
Returns a VMR vector for this model's main calculations.
constexpr auto pow2(auto x) noexcept
power of two
constexpr ValueType common_value_type(ValueType a, ValueType b) noexcept
Return a common type between a and b.
ValueList from_hitran(std::string_view upp, std::string_view low)
String update_isot_name(const String &old_name)
Updates the name of the isotopologue based on updates of the isotopologues.
constexpr Index find_species_index(const Species spec, const std::string_view isot) noexcept
constexpr Index nelem(Rational Ju, Rational Jl, Polarization type) noexcept
Gives the number of elements of the polarization type of this transition.
Polarization
Zeeman polarization selection.
Quantum::Number::Type QuantumNumberType
AbsorptionCutoffTagTypeStatus(const ArrayOfArrayOfAbsorptionLines &)
AbsorptionLineShapeTagTypeStatus(const ArrayOfArrayOfAbsorptionLines &)
AbsorptionMirroringTagTypeStatus(const ArrayOfArrayOfAbsorptionLines &)
AbsorptionNormalizationTagTypeStatus(const ArrayOfArrayOfAbsorptionLines &)
bool ByHITRANRosenkranzRelmat
bool ByRovibLinearDipoleLineMixing
AbsorptionPopulationTagTypeStatus(const ArrayOfArrayOfAbsorptionLines &)
Helper struct for flat_index.
AbsorptionNormalizationTagTypeStatus normalization
AbsorptionPopulationTagTypeStatus population
AbsorptionLineShapeTagTypeStatus lineshapetype
AbsorptionCutoffTagTypeStatus cutoff
AbsorptionMirroringTagTypeStatus mirroring
void sort_by_einstein()
Sort inner line list by Einstein coefficient.
void sort_by_frequency()
Sort inner line list by frequency.
Index NumBroadeners() const ARTS_NOEXCEPT
Number of broadening species.
SingleLine PopLine(Index) noexcept
Pops a single line.
Index BroadeningSpeciesPosition(Species::Species spec) const noexcept
Position of species if available or -1 else.
PopulationType population
Line population distribution.
Index NumLocalQuanta() const noexcept
Number of broadening species.
void RemoveLine(Index) noexcept
Removes a single line.
Numeric F_mean(Numeric T=0) const noexcept
Mean frequency by weight of line strength.
void SetAutomaticZeeman() noexcept
Set Zeeman effect for all lines that have the correct quantum numbers.
Vector BroadeningSpeciesVMR(const ConstVectorView &, const ArrayOfArrayOfSpeciesTag &) const
Returns the VMRs of the broadening species.
Numeric ZeemanSplitting(size_t k, Zeeman::Polarization type, Index i) const ARTS_NOEXCEPT
Returns the splitting of a Zeeman split line.
Array< SingleLine > lines
A list of individual lines.
Index NumLines() const noexcept
Number of lines.
bofstream & write(bofstream &os) const
Binary write for Lines.
Species::Species Species() const noexcept
Species Enum.
Numeric DopplerConstant(Numeric T) const noexcept
bool OnTheFlyLineMixing() const noexcept
On-the-fly line mixing.
Numeric ZeemanStrength(size_t k, Zeeman::Polarization type, Index i) const ARTS_NOEXCEPT
Returns the strength of a Zeeman split line.
QuantumIdentifier QuantumIdentityOfLine(Index k) const noexcept
bool OK() const ARTS_NOEXCEPT
bifstream & read(bifstream &is)
Binary read for Lines.
bool MatchWithExternal(const SingleLineExternal &sle, const QuantumIdentifier &quantumidentity) const ARTS_NOEXCEPT
Checks if an external line matches this structure.
Species::IsotopeRecord Isotopologue() const noexcept
Isotopologue Index.
String LineShapeMetaData() const noexcept
Meta data for the line shape if it exists.
void MakeLineShapeModelCommon()
Make a common line shape if possible.
String SpeciesName() const noexcept
Species Name.
Index LineShapePos(const Species::Species spec) const ARTS_NOEXCEPT
Position among broadening species or -1.
LineShape::Output ShapeParameters(size_t k, Numeric T, Numeric P, const Vector &vmrs) const ARTS_NOEXCEPT
Line shape parameters.
Index ZeemanCount(size_t k, Zeeman::Polarization type) const ARTS_NOEXCEPT
Returns the number of Zeeman split lines.
Numeric CutoffFreq(size_t k, Numeric shift=0) const noexcept
Returns cutoff frequency or maximum value.
void ReverseLines() noexcept
Reverses the order of the internal lines.
std::pair< bool, bool > Match(const Lines &l) const noexcept
Checks if another line list matches this structure.
void AppendSingleLine(SingleLine &&sl)
Appends a single line to the absorption lines.
ArrayOfSpecies broadeningspecies
A list of broadening species.
Numeric SelfVMR(const ConstVectorView &, const ArrayOfArrayOfSpeciesTag &) const
Returns the VMR of the species.
Numeric CutoffFreqMinus(size_t k, Numeric shift=0) const noexcept
Returns negative cutoff frequency or lowest value.
bool AnyLinemixing() const noexcept
String MetaData() const
Returns a printable statement about the lines.
bool DoLineMixing(Numeric P) const noexcept
Returns if the pressure should do line mixing.
Numeric SpeciesMass() const noexcept
Mass of the molecule.
QuantumIdentifier quantumidentity
Catalog ID.
bool DoVmrDerivative(const QuantumIdentifier &qid) const noexcept
LineShape::Output ShapeParameters_dVMR(size_t k, Numeric T, Numeric P, const QuantumIdentifier &vmr_qid) const ARTS_NOEXCEPT
Line shape parameters vmr derivative.
LineShape::Output ShapeParameters_dT(size_t k, Numeric T, Numeric P, const Vector &vmrs) const ARTS_NOEXCEPT
Line shape parameters temperature derivatives.
Vector BroadeningSpeciesMass(const ConstVectorView &, const ArrayOfArrayOfSpeciesTag &, const SpeciesIsotopologueRatios &, const Numeric &bath_mass=0) const
Returns the mass of the broadening species.
Single line reading output.
QuantumIdentifier quantumidentity
NormalizationType normalization
LineShape::Type lineshapetype
Computations and data for a single absorption line.
Numeric E0
Lower state energy level.
bofstream & write(bofstream &bof) const
Binary write for AbsorptionLines.
bifstream & read(bifstream &bif)
Binary read for AbsorptionLines.
Numeric F0
Central frequency.
Quantum::Number::LocalState localquanta
Local quantum numbers.
LineShape::Model lineshape
Line shape model.
void SetAutomaticZeeman(QuantumIdentifier qid)
Set Zeeman effect by automatic detection.
Zeeman::Model zeeman
Zeeman model.
Numeric A
Einstein spontaneous emission coefficient.
Numeric glow
Lower level statistical weight.
void SetLineMixing2AER(const Vector &d)
Set the line mixing model to AER kind.
Numeric gupp
Upper level statistical weight.
void SetLineMixing2SecondOrderData(const Vector &d)
Set the line mixing model to 2nd order.
Numeric I0
Reference intensity.
Coefficients and temperature model for SingleSpeciesModel.
constexpr Output & no_linemixing(bool do_no_linemixing)
Turns of line mixing if true. Return *this.
A logical struct for global quantum numbers with species identifiers.
Species::Species Species() const noexcept
Species::IsotopeRecord Isotopologue() const noexcept
A logical struct for local quantum numbers.
A complete quantum number value with type information.
constexpr void set(std::string_view s, bool upp)
Set level value.
Struct containing all information needed about one isotope.
String FullName() const noexcept
Numeric wigner6j(const Rational j1, const Rational j2, const Rational j3, const Rational l1, const Rational l2, const Rational l3)
Wigner 6J symbol.
Numeric wigner3j(const Rational j1, const Rational j2, const Rational j3, const Rational m1, const Rational m2, const Rational m3)
Wigner 3J symbol.
Wigner symbol interactions.