ARTS 2.5.4 (git: 31ce4f0e)
absorptionlines.cc
Go to the documentation of this file.
1/* Copyright (C) 2019
2 Richard Larsson <ric.larsson@gmail.com>
3
4
5 This program is free software; you can redistribute it and/or modify it
6 under the terms of the GNU General Public License as published by the
7 Free Software Foundation; either version 2, or (at your option) any
8 later version.
9
10
11 This program is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16
17 You should have received a copy of the GNU General Public License
18 along with this program; if not, write to the Free Software
19 Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
20 USA. */
21
32#include "absorptionlines.h"
33
34#include <algorithm>
35#include <limits>
36#include <ostream>
37#include <string>
38
39#include "absorption.h"
40#include "arts_conversions.h"
41#include "debug.h"
42#include "enums.h"
43#include "file.h"
44#include "hitran_species.h"
45#include "jpl_species.h"
46#include "linescaling.h"
47#include "lineshapemodel.h"
48#include "quantum_numbers.h"
49#include "rational.h"
50#include "wigner_functions.h"
51
53 size_t k, Numeric T, Numeric P, const Vector& vmrs) const noexcept {
54 auto x = lines[k].lineshape.GetParams(T, T0, P, vmrs);
55
56 if (not DoLineMixing(P)) x.Y = x.G = x.DV = 0;
57
58 return x;
59}
60
61LineShape::Output Absorption::Lines::ShapeParameters(size_t k,
62 Numeric T,
63 Numeric P,
64 size_t m) const noexcept {
65 return lines[k].lineshape.GetParams(T, T0, P, m);
66}
67
69 size_t k, Numeric T, Numeric P, const Vector& vmrs) const noexcept {
70 auto x = lines[k].lineshape.GetTemperatureDerivs(T, T0, P, vmrs);
71
72 if (not DoLineMixing(P)) x.Y = x.G = x.DV = 0;
73
74 return x;
75}
76
78 const Species::Species spec) const noexcept {
79 // Is always first if this is self and self broadening exists
80 if (selfbroadening and spec == quantumidentity.Species()) {
81 return 0;
82 }
83
84 // First and last might be artificial so they should not be checked
85 const Index s = selfbroadening;
86 const Index e = broadeningspecies.nelem() - bathbroadening;
87 for (Index i = s; i < e; i++) {
88 if (spec == broadeningspecies[i]) {
89 return i;
90 }
91 }
92
93 // At this point, the ID is not explicitly among the broadeners, but bath broadening means its VMR still might matter
94 if (bathbroadening) return broadeningspecies.nelem() - 1;
95 return -1;
96}
97
99 size_t k,
100 Numeric T,
101 Numeric P,
102 const QuantumIdentifier& vmr_qid) const noexcept {
103 const Index pos = LineShapePos(vmr_qid);
104 if (pos >= 0) return lines[k].lineshape.GetVMRDerivs(T, T0, P, pos);
105 return LineShape::Output{};
106}
107
109 size_t k,
110 Numeric T,
111 Numeric P,
112 const Vector& vmrs,
113 const RetrievalQuantity& derivative) const noexcept {
114 const auto self = derivative.Mode() == LineShape::self_broadening;
115 const auto bath = derivative.Mode() == LineShape::bath_broadening;
116 const auto& ls = lines[k].lineshape;
117
118 if (derivative.QuantumIdentity().Species() not_eq Species() or
119 derivative.QuantumIdentity().Isotopologue() not_eq Isotopologue())
120 return 0;
121 if (self and selfbroadening)
122 return ls.GetInternalDeriv(T, T0, P, 0, vmrs, derivative.LineType());
123 if (self)
124 return ls.GetInternalDeriv(
125 T,
126 T0,
127 P,
128 LineShapePos(derivative.QuantumIdentity().Species()),
129 vmrs,
130 derivative.LineType());
131 if (bath and bathbroadening)
132 return ls.GetInternalDeriv(
133 T, T0, P, ls.nelem() - 1, vmrs, derivative.LineType());
134 if (bath) return 0;
135 return ls.GetInternalDeriv(
136 T,
137 T0,
138 P,
139 LineShapePos(derivative.QuantumIdentity().Species()),
140 vmrs,
141 derivative.LineType());
142}
143
145 // Default data and values for this type
147 data.selfbroadening = true;
148 data.bathbroadening = true;
149 data.lineshapetype = LineShape::Type::VP;
150 data.species.resize(2);
151
152 // This always contains the rest of the line to parse. At the
153 // beginning the entire line. Line gets shorter and shorter as we
154 // continue to extract stuff from the beginning.
155 String line;
156
157 // Look for more comments?
158 bool comment = true;
159
160 while (comment) {
161 // Return true if eof is reached:
162 if (is.eof()) return data;
163
164 // Throw runtime_error if stream is bad:
165 ARTS_USER_ERROR_IF(!is, "Stream bad.");
166
167 // Read line from file into linebuffer:
168 getline(is, line);
169
170 // It is possible that we were exactly at the end of the file before
171 // calling getline. In that case the previous eof() was still false
172 // because eof() evaluates only to true if one tries to read after the
173 // end of the file. The following check catches this.
174 if (line.nelem() == 0 && is.eof()) return data;
175
176 // @ as first character marks catalogue entry
177 char c;
178 extract(c, line, 1);
179
180 // check for empty line
181 if (c == '@') {
182 comment = false;
183 }
184 }
185
186 // read the arts identifier String
187 istringstream icecream(line);
188
189 String artsid;
190 icecream >> artsid;
191
192 if (artsid.length() != 0) {
193 // Set the species
194 const auto isotopologue = Species::Tag(Species::update_isot_name(artsid));
196 isotopologue.is_joker() or
197 isotopologue.type not_eq Species::TagType::Plain,
198 "A line catalog species can only be of the form \"Plain\", meaning it\nhas the form SPECIES-ISONUM.\n"
199 "Your input contains: ",
200 artsid,
201 ". which we cannot interpret as a plain species")
203 Species::find_species_index(isotopologue.Isotopologue())};
204
205 // Extract center frequency:
206 icecream >> double_imanip() >> data.line.F0;
207
208 Numeric psf;
209 // Extract pressure shift:
210 icecream >> double_imanip() >> psf;
211
212 // Extract intensity:
213 icecream >> double_imanip() >> data.line.I0;
214
215 // Extract reference temperature for Intensity in K:
216 icecream >> double_imanip() >> data.T0;
217
218 // Extract lower state energy:
219 icecream >> double_imanip() >> data.line.E0;
220
221 // Extract air broadening parameters:
222 Numeric agam, sgam;
223 icecream >> double_imanip() >> agam;
224 icecream >> double_imanip() >> sgam;
225
226 // Extract temperature coefficient of broadening parameters:
227 Numeric nair, nself;
228 icecream >> double_imanip() >> nair;
229 icecream >> double_imanip() >> nself;
230
231 // Extract reference temperature for broadening parameter in K:
232 Numeric tgam;
233 icecream >> double_imanip() >> tgam;
234
235 // Extract the aux parameters:
236 Index naux;
237 icecream >> naux;
238
239 // resize the aux array and read it
240 ArrayOfNumeric maux;
241 maux.resize(naux);
242
243 for (Index j = 0; j < naux; j++) {
244 icecream >> double_imanip() >> maux[j];
245 //cout << "maux" << j << " = " << maux[j] << "\n";
246 }
247
248 // Extract accuracies:
249 try {
250 Numeric unused_numeric;
251 icecream >> double_imanip() >> unused_numeric /*mdf*/;
252 icecream >> double_imanip() >> unused_numeric /*mdi0*/;
253 icecream >> double_imanip() >> unused_numeric /*dagam*/;
254 icecream >> double_imanip() >> unused_numeric /*dsgam*/;
255 icecream >> double_imanip() >> unused_numeric /*dnair*/;
256 icecream >> double_imanip() >> unused_numeric /*dnself*/;
257 icecream >> double_imanip() >> unused_numeric /*dpsf*/;
258 } catch (const std::runtime_error&) {
259 }
260
261 // Fix if tgam is different from ti0
262 if (tgam != data.T0) {
263 agam = agam * pow(tgam / data.T0, nair);
264 sgam = sgam * pow(tgam / data.T0, nself);
265 psf = psf * pow(tgam / data.T0, (Numeric).25 + (Numeric)1.5 * nair);
266 }
267
268 // Set line shape computer
269 data.line.lineshape = LineShape::Model(sgam, nself, agam, nair, psf);
270 }
271
272 // That's it!
273 data.bad = false;
274 return data;
275}
276
278 // Default data and values for this type
280 data.selfbroadening = true;
281 data.bathbroadening = false;
282 data.lineshapetype = LineShape::Type::VP;
283
284 // This always contains the rest of the line to parse. At the
285 // beginning the entire line. Line gets shorter and shorter as we
286 // continue to extract stuff from the beginning.
287 String line;
288
289 // Look for more comments?
290 bool comment = true;
291
292 while (comment) {
293 // Return true if eof is reached:
294 if (is.eof()) return data;
295
296 // Throw runtime_error if stream is bad:
297 ARTS_USER_ERROR_IF(!is, "Stream bad.");
298
299 // Read line from file into linebuffer:
300 getline(is, line);
301
302 // It is possible that we were exactly at the end of the file before
303 // calling getline. In that case the previous eof() was still false
304 // because eof() evaluates only to true if one tries to read after the
305 // end of the file. The following check catches this.
306 if (line.nelem() == 0 && is.eof()) return data;
307
308 // @ as first character marks catalogue entry
309 char c;
310 extract(c, line, 1);
311
312 // check for empty line
313 if (c == '@') {
314 comment = false;
315 }
316 }
317
318 // read the arts identifier String
319 istringstream icecream(line);
320
321 String artsid;
322 icecream >> artsid;
323
324 if (artsid.length() != 0) {
325 // Set line ID
326 const auto isotopologue = Species::Tag(Species::update_isot_name(artsid));
328 isotopologue.is_joker() or
329 isotopologue.type not_eq Species::TagType::Plain,
330 "A line catalog species can only be of the form \"Plain\", meaning it\nhas the form SPECIES-ISONUM.\n"
331 "Your input contains: ",
332 artsid,
333 ". which we cannot interpret as a plain species")
335 Species::find_species_index(isotopologue.Isotopologue())};
336
337 // Extract center frequency:
338 icecream >> double_imanip() >> data.line.F0;
339
340 // Extract intensity:
341 icecream >> double_imanip() >> data.line.I0;
342
343 // Extract reference temperature for Intensity in K:
344 icecream >> double_imanip() >> data.T0;
345
346 // Extract lower state energy:
347 icecream >> double_imanip() >> data.line.E0;
348
349 // Extract Einstein A-coefficient:
350 icecream >> double_imanip() >> data.line.A;
351
352 // Extract upper state stat. weight:
353 icecream >> double_imanip() >> data.line.gupp;
354
355 // Extract lower state stat. weight:
356 icecream >> double_imanip() >> data.line.glow;
357
359 data.lineshapetype,
360 data.selfbroadening,
361 data.bathbroadening,
362 data.line.lineshape,
363 data.species,
364 data.quantumidentity);
365 }
366
367 // That's it!
368 data.bad = false;
369 return data;
370}
371
373 // Default data and values for this type
375
376 LineShape::Model line_mixing_model;
377 bool lmd_found = false;
378
379 // This always contains the rest of the line to parse. At the
380 // beginning the entire line. Line gets shorter and shorter as we
381 // continue to extract stuff from the beginning.
382 String line;
383
384 // Look for more comments?
385 bool comment = true;
386
387 while (comment) {
388 // Return true if eof is reached:
389 if (is.eof()) return data;
390
391 // Throw runtime_error if stream is bad:
392 ARTS_USER_ERROR_IF(!is, "Stream bad.");
393
394 // Read line from file into linebuffer:
395 getline(is, line);
396
397 // It is possible that we were exactly at the end of the file before
398 // calling getline. In that case the previous eof() was still false
399 // because eof() evaluates only to true if one tries to read after the
400 // end of the file. The following check catches this.
401 if (line.nelem() == 0 && is.eof()) return data;
402
403 // @ as first character marks catalogue entry
404 char c;
405 extract(c, line, 1);
406
407 // check for empty line
408 if (c == '@') {
409 comment = false;
410 }
411 }
412
413 // read the arts identifier String
414 istringstream icecream(line);
415
416 try {
417 String artsid;
418 icecream >> artsid;
419
420 if (artsid.length() != 0) {
421 // Set line ID:
422 const auto isotopologue = Species::Tag(Species::update_isot_name(artsid));
424 isotopologue.is_joker() or
425 isotopologue.type not_eq Species::TagType::Plain,
426 "A line catalog species can only be of the form \"Plain\", meaning it\nhas the form SPECIES-ISONUM.\n"
427 "Your input contains: ",
428 artsid,
429 ". which we cannot interpret as a plain species")
431 Species::find_species_index(isotopologue.Isotopologue())};
432
433 // Extract center frequency:
434 icecream >> double_imanip() >> data.line.F0;
435
436 // Extract intensity:
437 icecream >> double_imanip() >> data.line.I0;
438
439 // Extract reference temperature for Intensity in K:
440 icecream >> double_imanip() >> data.T0;
441
442 // Extract lower state energy:
443 icecream >> double_imanip() >> data.line.E0;
444
445 // Extract Einstein A-coefficient:
446 icecream >> double_imanip() >> data.line.A;
447
448 // Extract upper state stat. weight:
449 icecream >> double_imanip() >> data.line.gupp;
450
451 // Extract lower state stat. weight:
452 icecream >> double_imanip() >> data.line.glow;
453
454 String token;
455 Index nelem;
456 icecream >> token;
457
458 while (icecream) {
459 // Read pressure broadening (LEGACY)
460 if (token == "PB") {
462 data.lineshapetype,
463 data.selfbroadening,
464 data.bathbroadening,
465 data.line.lineshape,
466 data.species,
467 data.quantumidentity);
468 icecream >> token;
469 } else if (token == "QN") {
470 // Quantum numbers
471
472 icecream >> token;
474 token != "UP", "Unknown quantum number tag: ", token)
475
476 icecream >> token;
477 String r;
478 while (icecream and token not_eq "LO") {
479 auto qn = Quantum::Number::toType(token);
480 icecream >> r;
481 icecream >> token;
482
483 if (good_enum(qn)) {
484 if (data.quantumidentity.val.has(qn)) {
486 val.set(r, true);
487 data.quantumidentity.val.set(val);
488 } else {
489 data.quantumidentity.val.add(qn).set(r, true);
490 }
491 }
492 }
493
495 !is || token != "LO",
496 "Error in catalog. Lower quantum number tag 'LO' not found.")
497
498 icecream >> token;
499 auto qn = Quantum::Number::toType(token);
500 while (icecream and
501 not(token == "LM" or token == "LF" or token == "ZM" or
502 token == "LSM" or token == "PB")) {
503 icecream >> r;
504 icecream >> token;
505
506 if (good_enum(qn)) {
507 if (data.quantumidentity.val.has(qn)) {
509 val.set(r, false);
510 data.quantumidentity.val.set(val);
511 } else {
512 data.quantumidentity.val.add(qn).set(r, false);
513 }
514 }
515
516 qn = Quantum::Number::toType(token);
517 }
518 } else if (token == "LM") {
519 LineShape::from_linemixingdata(icecream, line_mixing_model);
520 icecream >> token;
521 lmd_found = true;
522 } else if (token == "LF") {
524 data.lineshapetype,
525 data.selfbroadening,
526 data.bathbroadening,
527 data.line.lineshape,
528 data.species);
529 // if (data.selfbroadening) data.species[0] = data.quantumidentity.Species();
530 icecream >> token;
531 } else if (token == "ZM") {
532 // Zeeman effect
533 icecream >> data.line.zeeman;
534 icecream >> token;
535 } else if (token == "LSM") {
536 // Line shape modifications
537
538 // Starts with the number of modifications
539 icecream >> nelem;
540 for (Index lsm = 0; lsm < nelem; lsm++) {
541 icecream >> token;
542
543 // cutoff frequency
544 if (token == "CUT") {
545 icecream >> double_imanip() >> data.cutofffreq;
546 data.cutoff = CutoffType::ByLine;
547 }
548
549 // linemixing pressure limit
550 if (token == "LML") {
551 icecream >> double_imanip() >> data.linemixinglimit;
552 }
553
554 // mirroring
555 else if (token == "MTM") {
556 icecream >> data.mirroring;
557 }
558
559 // line normalization
560 else if (token == "LNT") {
561 icecream >> data.normalization;
562 }
563
564 else {
565 ARTS_USER_ERROR("Unknown line modifications given: ", token)
566 }
567 }
568 icecream >> token;
569 } else {
570 ARTS_USER_ERROR("Unknown line data tag in legacy reading routine: ",
571 token)
572 }
573 }
574 }
575 } catch (const std::runtime_error& e) {
576 ARTS_USER_ERROR("Parse error in catalog line: \n", line, '\n', e.what())
577 }
578
579 if (lmd_found)
580 data.line.lineshape.SetLineMixingModel(line_mixing_model.Data()[0]);
581
582 // That's it!
583 data.bad = false;
584 return data;
585}
586
588 istream& is) {
589 // Default data and values for this type
591 data.selfbroadening = true;
592 data.bathbroadening = true;
593 data.lineshapetype = LineShape::Type::VP;
594 data.species.resize(2);
595 data.species[1] = Species::Species::Bath;
596
597 // This contains the rest of the line to parse. At the beginning the
598 // entire line. Line gets shorter and shorter as we continue to
599 // extract stuff from the beginning.
600 String line;
601
602 // The first item is the molecule number:
603 Index mo;
604
605 // Look for more comments?
606 bool comment = true;
607
608 while (comment) {
609 // Return true if eof is reached:
610 if (is.eof()) return data;
611
612 // Throw runtime_error if stream is bad:
613 ARTS_USER_ERROR_IF(!is, "Stream bad.");
614
615 // Read line from file into linebuffer:
616 getline(is, line);
617
618 // It is possible that we were exactly at the end of the file before
619 // calling getline. In that case the previous eof() was still false
620 // because eof() evaluates only to true if one tries to read after the
621 // end of the file. The following check catches this.
622 if (line.nelem() == 0 && is.eof()) return data;
623
624 // If the catalogue is in dos encoding, throw away the
625 // additional carriage return
626 if (line[line.nelem() - 1] == 13) {
627 line.erase(line.nelem() - 1, 1);
628 }
629
630 mo = 0;
631 // Initialization of mo is important, because mo stays the same
632 // if line is empty.
633 extract(mo, line, 2);
634 comment = false;
635 }
636
637 // Extract isotopologue:
638 char iso;
639 extract(iso, line, 1);
640
641 // Set line data
643 data.species[0] = data.quantumidentity.Species();
644
645 // Position.
646 {
647 // HITRAN position in wavenumbers (cm^-1):
648 Numeric v;
649 // Conversion from wavenumber to Hz. If you multiply a line
650 // position in wavenumber (cm^-1) by this constant, you get the
651 // frequency in Hz.
652 constexpr Numeric w2Hz = Constant::c * 100.;
653
654 // Extract HITRAN postion:
655 extract(v, line, 12);
656
657 // ARTS position in Hz:
658 data.line.F0 = v * w2Hz;
659 }
660
661 // Intensity.
662 {
663 // HITRAN intensity is in cm-1/(molec * cm-2) at 296 Kelvin.
664 // It already includes the isotpic ratio.
665 // The first cm-1 is the frequency unit (it cancels with the
666 // 1/frequency unit of the line shape function).
667 //
668 // We need to do the following:
669 // 1. Convert frequency from wavenumber to Hz (factor 1e2 * c).
670 // 2. Convert [molec * cm-2] to [molec * m-2] (factor 1e-4).
671 // 3. Take out the isotopologue ratio.
672
673 constexpr Numeric hi2arts = 1e-2 * Constant::c;
674
675 Numeric s;
676
677 // Extract HITRAN intensity:
678 extract(s, line, 10);
679 // Convert to ARTS units (Hz / (molec * m-2) ), or shorter: Hz*m^2
680 data.line.I0 = s * hi2arts;
681 // Take out isotopologue ratio:
682 data.line.I0 /= Hitran::ratio_from_lookup(mo, iso);
683 }
684
685 // Einstein coefficient
686 { extract(data.line.A, line, 10); }
687
688 // Air broadening parameters.
689 Numeric agam, sgam;
690 {
691 // HITRAN parameter is in cm-1/atm at 296 Kelvin
692 // All parameters are HWHM (I hope this is true!)
693 Numeric gam;
694 // Conversion from wavenumber to Hz. If you multiply a value in
695 // wavenumber (cm^-1) by this constant, you get the value in Hz.
696 constexpr Numeric w2Hz = Constant::c * 1e2;
697 // Ok, put together the end-to-end conversion that we need:
698 constexpr Numeric hi2arts = w2Hz / Conversion::atm2pa(1);
699
700 // Extract HITRAN AGAM value:
701 extract(gam, line, 5);
702
703 // ARTS parameter in Hz/Pa:
704 agam = gam * hi2arts;
705
706 // Extract HITRAN SGAM value:
707 extract(gam, line, 5);
708
709 // ARTS parameter in Hz/Pa:
710 sgam = gam * hi2arts;
711
712 // If zero, set to agam:
713 if (0 == sgam) sgam = agam;
714
715 // cout << "agam, sgam = " << magam << ", " << msgam << endl;
716 }
717
718 // Lower state energy.
719 {
720 // HITRAN parameter is in wavenumbers (cm^-1).
721 // We have to convert this to the ARTS unit Joule.
722
723 // Extract from Catalogue line
724 extract(data.line.E0, line, 10);
725
726 // Convert to Joule:
727 data.line.E0 = wavenumber_to_joule(data.line.E0);
728 }
729
730 // Temperature coefficient of broadening parameters.
731 Numeric nair, nself;
732 {
733 // This is dimensionless, we can also extract directly.
734 extract(nair, line, 4);
735
736 // Set self broadening temperature coefficient to the same value:
737 nself = nair;
738 // cout << "mnair = " << mnair << endl;
739 }
740
741 // Pressure shift.
742 Numeric psf;
743 {
744 // HITRAN value in cm^-1 / atm. So the conversion goes exactly as
745 // for the broadening parameters.
746 Numeric d;
747 // Conversion from wavenumber to Hz. If you multiply a value in
748 // wavenumber (cm^-1) by this constant, you get the value in Hz.
749 constexpr Numeric w2Hz = Constant::c * 1e2;
750 // Ok, put together the end-to-end conversion that we need:
751 constexpr Numeric hi2arts = w2Hz / Conversion::atm2pa(1);
752
753 // Extract HITRAN value:
754 extract(d, line, 8);
755
756 // ARTS value in Hz/Pa
757 psf = d * hi2arts;
758 }
759 // Set the accuracies using the definition of HITRAN
760 // indices. If some are missing, they are set to -1.
761
762 // Upper state global quanta
763 {
764 Index eu;
765 extract(eu, line, 15);
766 }
767
768 // Lower state global quanta
769 {
770 Index el;
771 extract(el, line, 15);
772 }
773
774 // Upper state local quanta
775 {
776 Index eul;
777 extract(eul, line, 15);
778 }
779
780 // Lower state local quanta
781 {
782 Index ell;
783 extract(ell, line, 15);
784 }
785
786 // Accuracy index for frequency
787 {
788 Index df;
789 // Extract HITRAN value:
790 extract(df, line, 1);
791 }
792
793 // Accuracy index for intensity
794 {
795 Index di0;
796 // Extract HITRAN value:
797 extract(di0, line, 1);
798 }
799
800 // Accuracy index for air-broadened halfwidth
801 {
802 Index dgam;
803 // Extract HITRAN value:
804 extract(dgam, line, 1);
805 }
806
807 // Accuracy index for self-broadened half-width
808 {
809 Index dgam;
810 // Extract HITRAN value:
811 extract(dgam, line, 1);
812 }
813
814 // Accuracy index for temperature-dependence exponent for agam
815 {
816 Index dn;
817 // Extract HITRAN value:
818 extract(dn, line, 1);
819 }
820
821 // Accuracy index for pressure shift
822 {
823 Index dpsfi;
824 // Extract HITRAN value (given in cm-1):
825 extract(dpsfi, line, 1);
826 }
827
828 // These were all the parameters that we can extract from
829 // HITRAN 2004. However, we still have to set the reference temperatures
830 // to the appropriate value:
831
832 // Reference temperature for Intensity in K.
833 data.T0 = 296.0;
834
835 // Set line shape computer
836 data.line.lineshape = LineShape::hitran_model(sgam, nself, agam, nair, psf);
837 {
838 Index garbage;
839 extract(garbage, line, 13);
840
841 // The statistical weights
842 extract(data.line.gupp, line, 7);
843 extract(data.line.glow, line, 7);
844 }
845
846 // That's it!
847 data.bad = false;
848 return data;
849}
850
852 istream& is) {
853 // Default data and values for this type
855 data.selfbroadening = true;
856 data.bathbroadening = true;
857 data.lineshapetype = LineShape::Type::VP;
858 data.species.resize(2);
859 data.species[1] = Species::Species::Bath;
860
861 // This contains the rest of the line to parse. At the beginning the
862 // entire line. Line gets shorter and shorter as we continue to
863 // extract stuff from the beginning.
864 String line;
865
866 // The first item is the molecule number:
867 Index mo;
868
869 // Look for more comments?
870 bool comment = true;
871
872 while (comment) {
873 // Return true if eof is reached:
874 if (is.eof()) return data;
875
876 // Throw runtime_error if stream is bad:
877 ARTS_USER_ERROR_IF(!is, "Stream bad.");
878
879 // Read line from file into linebuffer:
880 getline(is, line);
881
882 // It is possible that we were exactly at the end of the file before
883 // calling getline. In that case the previous eof() was still false
884 // because eof() evaluates only to true if one tries to read after the
885 // end of the file. The following check catches this.
886 if (line.nelem() == 0 && is.eof()) return data;
887
888 // If the catalogue is in dos encoding, throw away the
889 // additional carriage return
890 if (line[line.nelem() - 1] == 13) {
891 line.erase(line.nelem() - 1, 1);
892 }
893
894 mo = 0;
895 // Initialization of mo is important, because mo stays the same
896 // if line is empty.
897 extract(mo, line, 2);
898 comment = false;
899 }
900
901 // Extract isotopologue:
902 char iso;
903 extract(iso, line, 1);
904
905 // Set line data
907 data.species[0] = data.quantumidentity.Species();
908
909 // Position.
910 {
911 // HITRAN position in wavenumbers (cm^-1):
912 Numeric v;
913 // Conversion from wavenumber to Hz. If you multiply a line
914 // position in wavenumber (cm^-1) by this constant, you get the
915 // frequency in Hz.
916 constexpr Numeric w2Hz = Constant::c * 100.;
917
918 // Extract HITRAN postion:
919 extract(v, line, 12);
920
921 // ARTS position in Hz:
922 data.line.F0 = v * w2Hz;
923 }
924
925 // Intensity.
926 {
927 // HITRAN intensity is in cm-1/(molec * cm-2) at 296 Kelvin.
928 // It already includes the isotpic ratio.
929 // The first cm-1 is the frequency unit (it cancels with the
930 // 1/frequency unit of the line shape function).
931 //
932 // We need to do the following:
933 // 1. Convert frequency from wavenumber to Hz (factor 1e2 * c).
934 // 2. Convert [molec * cm-2] to [molec * m-2] (factor 1e-4).
935 // 3. Take out the isotopologue ratio.
936
937 constexpr Numeric hi2arts = 1e-2 * Constant::c;
938
939 Numeric s;
940
941 // Extract HITRAN intensity:
942 extract(s, line, 10);
943 // Convert to ARTS units (Hz / (molec * m-2) ), or shorter: Hz*m^2
944 data.line.I0 = s * hi2arts;
945 // Take out isotopologue ratio:
946 data.line.I0 /= Hitran::ratio_from_lookup(mo, iso);
947 }
948
949 // Einstein coefficient
950 { extract(data.line.A, line, 10); }
951
952 // Air broadening parameters.
953 Numeric agam, sgam;
954 {
955 // HITRAN parameter is in cm-1/atm at 296 Kelvin
956 // All parameters are HWHM (I hope this is true!)
957 Numeric gam;
958 // Conversion from wavenumber to Hz. If you multiply a value in
959 // wavenumber (cm^-1) by this constant, you get the value in Hz.
960 constexpr Numeric w2Hz = Constant::c * 1e2;
961 // Ok, put together the end-to-end conversion that we need:
962 constexpr Numeric hi2arts = w2Hz / Conversion::atm2pa(1);
963
964 // Extract HITRAN AGAM value:
965 extract(gam, line, 5);
966
967 // ARTS parameter in Hz/Pa:
968 agam = gam * hi2arts;
969
970 // Extract HITRAN SGAM value:
971 extract(gam, line, 5);
972
973 // ARTS parameter in Hz/Pa:
974 sgam = gam * hi2arts;
975
976 // If zero, set to agam:
977 if (0 == sgam) sgam = agam;
978
979 // cout << "agam, sgam = " << magam << ", " << msgam << endl;
980 }
981
982 // Lower state energy.
983 {
984 // HITRAN parameter is in wavenumbers (cm^-1).
985 // We have to convert this to the ARTS unit Joule.
986
987 // Extract from Catalogue line
988 extract(data.line.E0, line, 10);
989
990 // Convert to Joule:
991 data.line.E0 = wavenumber_to_joule(data.line.E0);
992 }
993
994 // Temperature coefficient of broadening parameters.
995 Numeric nair, nself;
996 {
997 // This is dimensionless, we can also extract directly.
998 extract(nair, line, 4);
999
1000 // Set self broadening temperature coefficient to the same value:
1001 nself = nair;
1002 // cout << "mnair = " << mnair << endl;
1003 }
1004
1005 // Pressure shift.
1006 Numeric psf;
1007 {
1008 // HITRAN value in cm^-1 / atm. So the conversion goes exactly as
1009 // for the broadening parameters.
1010 Numeric d;
1011 // Conversion from wavenumber to Hz. If you multiply a value in
1012 // wavenumber (cm^-1) by this constant, you get the value in Hz.
1013 constexpr Numeric w2Hz = Constant::c * 1e2;
1014 // Ok, put together the end-to-end conversion that we need:
1015 constexpr Numeric hi2arts = w2Hz / Conversion::atm2pa(1);
1016
1017 // Extract HITRAN value:
1018 extract(d, line, 8);
1019
1020 // ARTS value in Hz/Pa
1021 psf = d * hi2arts;
1022 }
1023 // Set the accuracies using the definition of HITRAN
1024 // indices. If some are missing, they are set to -1.
1025
1026 // Upper state global quanta
1027 {
1028 Index eu;
1029 extract(eu, line, 15);
1030 }
1031
1032 // Lower state global quanta
1033 {
1034 Index el;
1035 extract(el, line, 15);
1036 }
1037
1038 // Upper state local quanta
1039 {
1040 Index eul;
1041 extract(eul, line, 15);
1042 }
1043
1044 // Lower state local quanta
1045 {
1046 Index ell;
1047 extract(ell, line, 15);
1048 }
1049
1050 // Accuracy index for frequency
1051 {
1052 Index df;
1053 // Extract HITRAN value:
1054 extract(df, line, 1);
1055 }
1056
1057 // Accuracy index for intensity
1058 {
1059 Index di0;
1060 // Extract HITRAN value:
1061 extract(di0, line, 1);
1062 }
1063
1064 // Accuracy index for air-broadened halfwidth
1065 {
1066 Index dgam;
1067 // Extract HITRAN value:
1068 extract(dgam, line, 1);
1069 }
1070
1071 // Accuracy index for self-broadened half-width
1072 {
1073 Index dgam;
1074 // Extract HITRAN value:
1075 extract(dgam, line, 1);
1076 }
1077
1078 // Accuracy index for temperature-dependence exponent for agam
1079 {
1080 Index dn;
1081 // Extract HITRAN value:
1082 extract(dn, line, 1);
1083 }
1084
1085 // Accuracy index for pressure shift
1086 {
1087 Index dpsfi;
1088 // Extract HITRAN value (given in cm-1):
1089 extract(dpsfi, line, 1);
1090 }
1091
1092 // These were all the parameters that we can extract from
1093 // HITRAN 2004. However, we still have to set the reference temperatures
1094 // to the appropriate value:
1095
1096 // Reference temperature for Intensity in K.
1097 data.T0 = 296.0;
1098
1099 // Set line shape computer
1100 data.line.lineshape = LineShape::hitran_model(sgam, nself, agam, nair, psf);
1101 {
1102 Index garbage;
1103 extract(garbage, line, 13);
1104
1105 // The statistical weights
1106 extract(data.line.gupp, line, 7);
1107 extract(data.line.glow, line, 7);
1108 }
1109
1110 // ADD QUANTUM NUMBER PARSING HERE!
1111 String upper, lower;
1112 std::stringstream ss;
1113 ss.str(line);
1114 ss >> upper >> lower;
1116
1117 // That's it!
1118 data.bad = false;
1119 return data;
1120}
1121
1123 istream& is) {
1124 // Default data and values for this type
1125 SingleLineExternal data;
1126 data.selfbroadening = true;
1127 data.bathbroadening = true;
1128 data.lineshapetype = LineShape::Type::VP;
1129 data.species.resize(2);
1130 data.species[1] = Species::Species::Bath;
1131
1132 // This contains the rest of the line to parse. At the beginning the
1133 // entire line. Line gets shorter and shorter as we continue to
1134 // extract stuff from the beginning.
1135 String line;
1136
1137 // The first item is the molecule number:
1138 Index mo;
1139
1140 // Look for more comments?
1141 bool comment = true;
1142
1143 while (comment) {
1144 // Return true if eof is reached:
1145 if (is.eof()) return data;
1146
1147 // Throw runtime_error if stream is bad:
1148 ARTS_USER_ERROR_IF(!is, "Stream bad.");
1149
1150 // Read line from file into linebuffer:
1151 getline(is, line);
1152
1153 // It is possible that we were exactly at the end of the file before
1154 // calling getline. In that case the previous eof() was still false
1155 // because eof() evaluates only to true if one tries to read after the
1156 // end of the file. The following check catches this.
1157 if (line.nelem() == 0 && is.eof()) return data;
1158
1159 // If the catalogue is in dos encoding, throw away the
1160 // additional carriage return
1161 if (line[line.nelem() - 1] == 13) {
1162 line.erase(line.nelem() - 1, 1);
1163 }
1164
1165 mo = 0;
1166 // Initialization of mo is important, because mo stays the same
1167 // if line is empty.
1168 extract(mo, line, 2);
1169 comment = false;
1170 }
1171
1172 // Extract isotopologue:
1173 char iso;
1174 extract(iso, line, 1);
1175
1176 // Set line data
1178 data.species[0] = data.quantumidentity.Species();
1179
1180 // Position.
1181 {
1182 // HITRAN position in wavenumbers (cm^-1):
1183 Numeric v;
1184 // Conversion from wavenumber to Hz. If you multiply a line
1185 // position in wavenumber (cm^-1) by this constant, you get the
1186 // frequency in Hz.
1187 constexpr Numeric w2Hz = Constant::c * 100.;
1188
1189 // Extract HITRAN postion:
1190 extract(v, line, 12);
1191
1192 // ARTS position in Hz:
1193 data.line.F0 = v * w2Hz;
1194 // cout << "mf = " << mf << endl;
1195 }
1196
1197 // Intensity.
1198 {
1199 // HITRAN intensity is in cm-1/(molec * cm-2) at 296 Kelvin.
1200 // It already includes the isotpic ratio.
1201 // The first cm-1 is the frequency unit (it cancels with the
1202 // 1/frequency unit of the line shape function).
1203 //
1204 // We need to do the following:
1205 // 1. Convert frequency from wavenumber to Hz (factor 1e2 * c).
1206 // 2. Convert [molec * cm-2] to [molec * m-2] (factor 1e-4).
1207 // 3. Take out the isotopologue ratio.
1208
1209 constexpr Numeric hi2arts = 1e-2 * Constant::c;
1210
1211 Numeric s;
1212
1213 // Extract HITRAN intensity:
1214 extract(s, line, 10);
1215 // Convert to ARTS units (Hz / (molec * m-2) ), or shorter: Hz*m^2
1216 data.line.I0 = s * hi2arts;
1217 // Take out isotopologue ratio:
1218 data.line.I0 /= Hitran::ratio_from_lookup(mo, iso);
1219 }
1220
1221 // Skip transition probability:
1222 {
1223 Numeric r;
1224 extract(r, line, 10);
1225 }
1226
1227 // Air broadening parameters.
1228 Numeric agam, sgam;
1229 {
1230 // HITRAN parameter is in cm-1/atm at 296 Kelvin
1231 // All parameters are HWHM (I hope this is true!)
1232 Numeric gam;
1233 // Conversion from wavenumber to Hz. If you multiply a value in
1234 // wavenumber (cm^-1) by this constant, you get the value in Hz.
1235 constexpr Numeric w2Hz = Constant::c * 1e2;
1236 // Ok, put together the end-to-end conversion that we need:
1237 constexpr Numeric hi2arts = w2Hz / Conversion::atm2pa(1);
1238
1239 // Extract HITRAN AGAM value:
1240 extract(gam, line, 5);
1241
1242 // ARTS parameter in Hz/Pa:
1243 agam = gam * hi2arts;
1244
1245 // Extract HITRAN SGAM value:
1246 extract(gam, line, 5);
1247
1248 // ARTS parameter in Hz/Pa:
1249 sgam = gam * hi2arts;
1250
1251 // If zero, set to agam:
1252 if (0 == sgam) sgam = agam;
1253
1254 // cout << "agam, sgam = " << magam << ", " << msgam << endl;
1255 }
1256
1257 // Lower state energy.
1258 {
1259 // HITRAN parameter is in wavenumbers (cm^-1).
1260 // We have to convert this to the ARTS unit Joule.
1261
1262 // Extract from Catalogue line
1263 extract(data.line.E0, line, 10);
1264
1265 // Convert to Joule:
1266 data.line.E0 = wavenumber_to_joule(data.line.E0);
1267 }
1268
1269 // Temperature coefficient of broadening parameters.
1270 Numeric nair, nself;
1271 {
1272 // This is dimensionless, we can also extract directly.
1273 extract(nair, line, 4);
1274
1275 // Set self broadening temperature coefficient to the same value:
1276 nself = nair;
1277 // cout << "mnair = " << mnair << endl;
1278 }
1279
1280 // Pressure shift.
1281 Numeric psf;
1282 {
1283 // HITRAN value in cm^-1 / atm. So the conversion goes exactly as
1284 // for the broadening parameters.
1285 Numeric d;
1286 // Conversion from wavenumber to Hz. If you multiply a value in
1287 // wavenumber (cm^-1) by this constant, you get the value in Hz.
1288 constexpr Numeric w2Hz = Constant::c * 1e2;
1289 // Ok, put together the end-to-end conversion that we need:
1290 constexpr Numeric hi2arts = w2Hz / Conversion::atm2pa(1);
1291
1292 // Extract HITRAN value:
1293 extract(d, line, 8);
1294
1295 // ARTS value in Hz/Pa
1296 psf = d * hi2arts;
1297 }
1298 // Set the accuracies using the definition of HITRAN
1299 // indices. If some are missing, they are set to -1.
1300
1301 //Skip upper state global quanta index
1302 {
1303 Index eu;
1304 extract(eu, line, 3);
1305 }
1306
1307 //Skip lower state global quanta index
1308 {
1309 Index el;
1310 extract(el, line, 3);
1311 }
1312
1313 //Skip upper state local quanta
1314 {
1315 Index eul;
1316 extract(eul, line, 9);
1317 }
1318
1319 //Skip lower state local quanta
1320 {
1321 Index ell;
1322 extract(ell, line, 9);
1323 }
1324
1325 // Accuracy index for frequency reference
1326 {
1327 Index df;
1328 // Extract HITRAN value:
1329 extract(df, line, 1);
1330 }
1331
1332 // Accuracy index for intensity reference
1333 {
1334 Index di0;
1335 // Extract HITRAN value:
1336 extract(di0, line, 1);
1337 }
1338
1339 // Accuracy index for halfwidth reference
1340 {
1341 Index dgam;
1342 // Extract HITRAN value:
1343 extract(dgam, line, 1);
1344 }
1345
1346 // These were all the parameters that we can extract from
1347 // HITRAN. However, we still have to set the reference temperatures
1348 // to the appropriate value:
1349
1350 // Reference temperature for Intensity in K.
1351 data.T0 = 296.0;
1352
1353 // Set line shape computer
1354 data.line.lineshape = LineShape::hitran_model(sgam, nself, agam, nair, psf);
1355
1356 // That's it!
1357 data.bad = false;
1358 return data;
1359}
1360
1362 // Default data and values for this type
1363 SingleLineExternal data;
1364 data.selfbroadening = true;
1365 data.bathbroadening = true;
1366 data.lineshapetype = LineShape::Type::VP;
1367 data.species.resize(2);
1368
1369 // This contains the rest of the line to parse. At the beginning the
1370 // entire line. Line gets shorter and shorter as we continue to
1371 // extract stuff from the beginning.
1372 String line;
1373
1374 // The first item is the molecule number:
1375 Index mo;
1376
1377 // Look for more comments?
1378 bool comment = true;
1379
1380 while (comment) {
1381 // Return true if eof is reached:
1382 if (is.eof()) return data;
1383
1384 // Throw runtime_error if stream is bad:
1385 ARTS_USER_ERROR_IF(!is, "Stream bad.");
1386
1387 // Read line from file into linebuffer:
1388 getline(is, line);
1389 if (line[0] == '>' or line[0] == '%') continue;
1390
1391 // It is possible that we were exactly at the end of the file before
1392 // calling getline. In that case the previous eof() was still false
1393 // because eof() evaluates only to true if one tries to read after the
1394 // end of the file. The following check catches this.
1395 if (line.nelem() == 0 && is.eof()) return data;
1396
1397 // If the catalogue is in dos encoding, throw away the
1398 // additional carriage return
1399 if (line[line.nelem() - 1] == 13) {
1400 line.erase(line.nelem() - 1, 1);
1401 }
1402
1403 mo = 0;
1404 // Initialization of mo is important, because mo stays the same
1405 // if line is empty.
1406 extract(mo, line, 2);
1407 comment = false;
1408 }
1409
1410 // Extract isotopologue:
1411 char iso;
1412 extract(iso, line, 1);
1413
1414 // Set line data
1416
1417 // Position.
1418 {
1419 // HITRAN position in wavenumbers (cm^-1):
1420 Numeric v;
1421 // Conversion from wavenumber to Hz. If you multiply a line
1422 // position in wavenumber (cm^-1) by this constant, you get the
1423 // frequency in Hz.
1424 constexpr Numeric w2Hz = Constant::c * 100.;
1425
1426 // Extract HITRAN postion:
1427 extract(v, line, 12);
1428
1429 // ARTS position in Hz:
1430 data.line.F0 = v * w2Hz;
1431 // cout << "mf = " << mf << endl;
1432 }
1433
1434 // Intensity.
1435 {
1436 // HITRAN intensity is in cm-1/(molec * cm-2) at 296 Kelvin.
1437 // It already includes the isotpic ratio.
1438 // The first cm-1 is the frequency unit (it cancels with the
1439 // 1/frequency unit of the line shape function).
1440 //
1441 // We need to do the following:
1442 // 1. Convert frequency from wavenumber to Hz (factor 1e2 * c).
1443 // 2. Convert [molec * cm-2] to [molec * m-2] (factor 1e-4).
1444 // 3. Take out the isotopologue ratio.
1445
1446 constexpr Numeric hi2arts = 1e-2 * Constant::c;
1447
1448 Numeric s;
1449 if (line[6] == 'D') line[6] = 'E';
1450 // Extract HITRAN intensity:
1451 extract(s,
1452 line,
1453 10); // NOTE: If error shooting, FORTRAN "D" is not read properly.
1454 // Convert to ARTS units (Hz / (molec * m-2) ), or shorter: Hz*m^2
1455 data.line.I0 = s * hi2arts;
1456 // Take out isotopologue ratio:
1457 data.line.I0 /= Hitran::ratio_from_lookup(mo, iso);
1458 }
1459
1460 // Skip transition probability:
1461 {
1462 Numeric r;
1463 extract(r, line, 10);
1464 }
1465
1466 // Air broadening parameters.
1467 Numeric sgam, agam;
1468 {
1469 // HITRAN parameter is in cm-1/atm at 296 Kelvin
1470 // All parameters are HWHM (I hope this is true!)
1471 Numeric gam;
1472 // Conversion from wavenumber to Hz. If you multiply a value in
1473 // wavenumber (cm^-1) by this constant, you get the value in Hz.
1474 constexpr Numeric w2Hz = Constant::c * 1e2;
1475 // Ok, put together the end-to-end conversion that we need:
1476 constexpr Numeric hi2arts = w2Hz / Conversion::atm2pa(1);
1477
1478 // Extract HITRAN AGAM value:
1479 extract(gam, line, 5);
1480
1481 // ARTS parameter in Hz/Pa:
1482 agam = gam * hi2arts;
1483
1484 // Extract HITRAN SGAM value:
1485 extract(gam, line, 5);
1486
1487 // ARTS parameter in Hz/Pa:
1488 sgam = gam * hi2arts;
1489
1490 // If zero, set to agam:
1491 if (0 == sgam) sgam = agam;
1492
1493 // cout << "agam, sgam = " << magam << ", " << msgam << endl;
1494 }
1495
1496 // Lower state energy.
1497 {
1498 // HITRAN parameter is in wavenumbers (cm^-1).
1499 // We have to convert this to the ARTS unit Joule.
1500
1501 // Extract from Catalogue line
1502 extract(data.line.E0, line, 10);
1503
1504 // Convert to Joule:
1505 data.line.E0 = wavenumber_to_joule(data.line.E0);
1506 }
1507
1508 // Temperature coefficient of broadening parameters.
1509 Numeric nair, nself;
1510 {
1511 // This is dimensionless, we can also extract directly.
1512 extract(nair, line, 4);
1513
1514 // Set self broadening temperature coefficient to the same value:
1515 nself = nair;
1516 // cout << "mnair = " << mnair << endl;
1517 }
1518
1519 // Pressure shift.
1520 Numeric psf;
1521 {
1522 // HITRAN value in cm^-1 / atm. So the conversion goes exactly as
1523 // for the broadening parameters.
1524 Numeric d;
1525 // Conversion from wavenumber to Hz. If you multiply a value in
1526 // wavenumber (cm^-1) by this constant, you get the value in Hz.
1527 constexpr Numeric w2Hz = Constant::c * 1e2;
1528 // Ok, put together the end-to-end conversion that we need:
1529 constexpr Numeric hi2arts = w2Hz / Conversion::atm2pa(1);
1530
1531 // Extract HITRAN value:
1532 extract(d, line, 8);
1533
1534 // ARTS value in Hz/Pa
1535 psf = d * hi2arts;
1536 }
1537 // Set the accuracies using the definition of HITRAN
1538 // indices. If some are missing, they are set to -1.
1539
1540 //Skip upper state global quanta index
1541 {
1542 Index eu;
1543 extract(eu, line, 3);
1544 }
1545
1546 //Skip lower state global quanta index
1547 {
1548 Index el;
1549 extract(el, line, 3);
1550 }
1551
1552 //Skip upper state local quanta
1553 {
1554 Index eul;
1555 extract(eul, line, 9);
1556 }
1557
1558 //Skip lower state local quanta
1559 {
1560 Index ell;
1561 extract(ell, line, 9);
1562 }
1563
1564 // Accuracy index for frequency reference
1565 {
1566 Index df;
1567 // Extract HITRAN value:
1568 extract(df, line, 1);
1569 }
1570
1571 // Accuracy index for intensity reference
1572 {
1573 Index di0;
1574 // Extract HITRAN value:
1575 extract(di0, line, 1);
1576 }
1577
1578 // Accuracy index for halfwidth reference
1579 {
1580 Index dgam;
1581 // Extract HITRAN value:
1582 extract(dgam, line, 1);
1583 }
1584
1585 // These were all the parameters that we can extract from
1586 // HITRAN. However, we still have to set the reference temperatures
1587 // to the appropriate value:
1588
1589 // Reference temperature for Intensity in K.
1590 // (This is fix for HITRAN)
1591 data.T0 = 296.0;
1592
1593 // Skip four
1594 {
1595 Index four;
1596 extract(four, line, 4);
1597 }
1598
1599 // This is the test for the last two characters of the line
1600 {
1601 /*
1602 * 0 is nothing,
1603 * -1 is linemixing on the next line,
1604 * -3 is the non-resonant line
1605 */
1606 Index test;
1607 extract(test, line, 2);
1608 //If the tag is as it should be, then a minus one means that more should be read
1609 if (test == -1 || test == -3)
1610 getline(is, line);
1611 else // the line is done and we are happy to leave
1612 {
1613 data.line.lineshape =
1614 LineShape::hitran_model(sgam, nself, agam, nair, psf);
1615
1616 data.bad = false;
1617 return data;
1618 }
1619 }
1620
1621 // In case we are unable to leave, the next line is a line mixing parameter line
1622
1623 // First is the molecular number. This should be the same as above.
1624 {
1625 Index mo2;
1626 extract(mo2, line, 2);
1627 // Skip one
1628
1629 ARTS_USER_ERROR_IF(mo != mo2, "There is an error in the line mixing\n");
1630 }
1631
1632 Vector Y(4), G(4), T(4);
1633
1634 // These are constants for AER but should be included because we need their grid.
1635 T[0] = 200;
1636 T[1] = 250;
1637 T[2] = 296;
1638 T[3] = 340;
1639
1640 // Next is the Y and G at various temperatures
1641 {
1642 Numeric Y_200K;
1643 extract(Y_200K, line, 13);
1644 Y[0] = Y_200K;
1645 }
1646 {
1647 Numeric G_200K;
1648 extract(G_200K, line, 11);
1649 G[0] = G_200K;
1650 }
1651 {
1652 Numeric Y_250K;
1653 extract(Y_250K, line, 13);
1654 Y[1] = Y_250K;
1655 }
1656 {
1657 Numeric G_250K;
1658 extract(G_250K, line, 11);
1659 G[1] = G_250K;
1660 }
1661 {
1662 Numeric Y_296K;
1663 extract(Y_296K, line, 13);
1664 Y[2] = Y_296K;
1665 }
1666 {
1667 Numeric G_296K;
1668 extract(G_296K, line, 11);
1669 G[2] = G_296K;
1670 }
1671 {
1672 Numeric Y_340K;
1673 extract(Y_340K, line, 13);
1674 Y[3] = Y_340K;
1675 }
1676 {
1677 Numeric G_340K;
1678 extract(G_340K, line, 11);
1679 G[3] = G_340K;
1680 }
1681
1682 // Cpnvert from per Atm and per Atm^2
1683 Y /= Conversion::atm2pa(1);
1685
1686 // ARTS uses (1-iY) as line-mixing factor, LBLRTM uses (1+iY), so we must change sign
1687 Y *= -1;
1688
1689 // Test that this is the end
1690 {
1691 Index test;
1692 extract(test, line, 2);
1693 if (test == -1) {
1695 nself,
1696 agam,
1697 nair,
1698 psf,
1699 {T[0],
1700 T[1],
1701 T[2],
1702 T[3],
1703 Y[0],
1704 Y[1],
1705 Y[2],
1706 Y[3],
1707 G[0],
1708 G[1],
1709 G[2],
1710 G[3]});
1711
1712 data.bad = false;
1713 return data;
1714 }
1715 if (test == -3) {
1717 nself,
1718 agam,
1719 nair,
1720 psf,
1721 {T[0],
1722 T[1],
1723 T[2],
1724 T[3],
1725 Y[0],
1726 Y[1],
1727 Y[2],
1728 Y[3],
1729 G[0],
1730 G[1],
1731 G[2],
1732 G[3]});
1733
1734 data.bad = false;
1735 return data;
1736 }
1737 return data;
1738 }
1739}
1740
1741std::vector<Absorption::Lines> Absorption::split_list_of_external_lines(
1742 std::vector<SingleLineExternal>& external_lines,
1743 const std::vector<QuantumNumberType>& localquantas,
1744 const std::vector<QuantumNumberType>& globalquantas) {
1745 std::vector<Lines> lines(0);
1746
1747 // Loop but make copies because we will need to modify some of the data
1748 while (external_lines.size()) {
1749 auto& sle = external_lines.back();
1750
1751 // Adapt broadening to fit with line catalog
1752 if (sle.selfbroadening) sle.species.front() = sle.quantumidentity.Species();
1753 if (sle.bathbroadening) sle.species.back() = Species::Species::Bath;
1754
1756 sle.quantumidentity.isotopologue_index};
1757 for (auto qn : globalquantas)
1758 if (sle.quantumidentity.val.has(qn))
1759 global_id.val.set(sle.quantumidentity.val[qn]);
1760
1761 Quantum::Number::LocalState local_id{};
1762 for (auto qn : localquantas)
1763 if (sle.quantumidentity.val.has(qn))
1764 local_id.val.set(sle.quantumidentity.val[qn]);
1765
1766 // Set the line
1767 auto line = sle.line;
1768 line.localquanta = local_id;
1769
1770 // Either find a line like this in the list of lines or start a new Lines
1771 auto band = std::find_if(lines.begin(), lines.end(), [&](const Lines& li) {
1772 return li.MatchWithExternal(sle, global_id);
1773 });
1774 if (band not_eq lines.end()) {
1775 band->AppendSingleLine(line);
1776 } else {
1777 lines.push_back(Lines(sle.selfbroadening,
1778 sle.bathbroadening,
1779 sle.cutoff,
1780 sle.mirroring,
1781 sle.population,
1782 sle.normalization,
1783 sle.lineshapetype,
1784 sle.T0,
1785 sle.cutofffreq,
1786 sle.linemixinglimit,
1787 global_id,
1788 sle.species,
1789 {line}));
1790 }
1791 external_lines.pop_back();
1792 }
1793
1794 return lines;
1795}
1796
1797namespace Absorption {
1798std::ostream& operator<<(std::ostream& os, const Absorption::Lines& lines) {
1799 for (auto& line : lines.lines) os << line << '\n';
1800 return os;
1801}
1802
1803std::istream& operator>>(std::istream& is, Lines& lines) {
1804 for (auto& line : lines.lines) is >> line;
1805 return is;
1806}
1807
1808std::ostream& operator<<(std::ostream& os, const Absorption::SingleLine& line) {
1809 os << line.F0 << ' ' << line.I0 << ' ' << line.E0 << ' ' << line.glow << ' '
1810 << line.gupp << ' ' << line.A << ' ' << line.zeeman << ' '
1811 << line.lineshape;
1812 if (line.localquanta.val.nelem()) os << ' ' << line.localquanta.values();
1813 return os;
1814}
1815
1816std::istream& operator>>(std::istream& is, Absorption::SingleLine& line) {
1817 is >> double_imanip() >> line.F0 >> line.I0 >> line.E0 >> line.glow >>
1818 line.gupp >> line.A;
1819 return is >> line.zeeman >> line.lineshape >> line.localquanta;
1820}
1821} // namespace Absorption
1822
1825}
1826
1828 std::ostringstream os;
1829
1830 os << "\nLines meta-data:\n";
1831 os << '\t' << "Species identity:\n";
1832 os << "\t\tSpecies: " << SpeciesName() << '\n';
1833 os << "\t\tIdentity: " << quantumidentity << '\n';
1834 os << '\t' << cutofftype2metadatastring(cutoff, cutofffreq);
1835 os << '\t' << populationtype2metadatastring(population);
1836 os << '\t' << normalizationtype2metadatastring(normalization);
1837 os << '\t' << LineShape::shapetype2metadatastring(lineshapetype);
1838 os << '\t' << mirroringtype2metadatastring(mirroring);
1839 os << '\t' << "The reference temperature for all line parameters is " << T0
1840 << " K.\n";
1841 if (linemixinglimit < 0)
1842 os << '\t' << "If applicable, there is no line mixing limit.\n";
1843 else
1844 os << '\t' << "If applicable, there is a line mixing limit at "
1845 << linemixinglimit << " Pa.\n";
1846
1847 if (not NumLines()) {
1848 os << "\tNo line data is available.\n";
1849 } else {
1850 os << "\tThere are " << NumLines() << " lines available.\n";
1851
1852 auto& line = lines.front();
1853 os << "\tThe front line has:\n";
1854 os << "\t\t"
1855 << "f0: " << line.F0 << " Hz\n";
1856 os << "\t\t"
1857 << "i0: " << line.I0 << " m^2/Hz\n";
1858 os << "\t\t"
1859 << "e0: " << line.E0 << " J\n";
1860 os << "\t\t"
1861 << "Lower stat. weight: " << line.glow << " [-]\n";
1862 os << "\t\t"
1863 << "Upper stat. weight: " << line.gupp << " [-]\n";
1864 os << "\t\t"
1865 << "A: " << line.A << " 1/s\n";
1866 os << "\t\t"
1867 << "Zeeman splitting of lower state: " << line.zeeman.gl() << " [-]\n";
1868 os << "\t\t"
1869 << "Zeeman splitting of upper state: " << line.zeeman.gu() << " [-]\n";
1870 os << "\t\t"
1871 << "Local quantum numbers: " << line.localquanta.val << "\n";
1872
1874 line.lineshape, selfbroadening, broadeningspecies, T0);
1875 os << "\t\t"
1876 << "Line shape parameters (are normalized by sum(VMR)):\n";
1877 for (auto& ls_form : ls_meta) os << "\t\t\t" << ls_form << "\n";
1878 }
1879
1880 return os.str();
1881}
1882
1884 lines.erase(lines.begin() + i);
1885}
1886
1888 auto line = lines[i];
1889 RemoveLine(i);
1890 return line;
1891}
1892
1894 std::reverse(lines.begin(), lines.end());
1895}
1896
1898 return quantumidentity.Isotopologue().mass;
1899}
1900
1902 const ConstVectorView& atm_vmrs,
1903 const ArrayOfArrayOfSpeciesTag& atm_spec) const {
1904 return LineShape::vmrs(atm_vmrs, atm_spec, broadeningspecies);
1905}
1906
1908 const ConstVectorView& atm_vmrs,
1909 const ArrayOfArrayOfSpeciesTag& atm_spec,
1910 const SpeciesIsotopologueRatios& ir,
1911 const Numeric& bath_mass) const {
1912 Vector mass = LineShape::mass(atm_vmrs, atm_spec, broadeningspecies, ir);
1913 if (bathbroadening and bath_mass > 0) mass[mass.nelem() - 1] = bath_mass;
1914 return mass;
1915}
1916
1918 const ConstVectorView& atm_vmrs,
1919 const ArrayOfArrayOfSpeciesTag& atm_spec) const {
1920 ARTS_USER_ERROR_IF(atm_vmrs.nelem() not_eq atm_spec.nelem(),
1921 "Bad species and vmr lists");
1922
1923 for (Index i = 0; i < atm_spec.nelem(); i++)
1924 if (atm_spec[i].nelem() and atm_spec[i].Species() == Species())
1925 return atm_vmrs[i];
1926 return 0;
1927}
1928
1930 Rational Jf, Rational Ji, Rational lf, Rational li, Rational k) {
1931 if (not iseven(Jf + lf + 1))
1932 return -sqrt(2 * Jf + 1) * wigner3j(Jf, k, Ji, li, lf - li, -lf);
1933 return +sqrt(2 * Jf + 1) * wigner3j(Jf, k, Ji, li, lf - li, -lf);
1934}
1935
1937 Rational Ji,
1938 Rational N) {
1939 if (not iseven(Jf + N))
1940 return -sqrt(6 * (2 * Jf + 1) * (2 * Ji + 1)) *
1941 wigner6j(1, 1, 1, Ji, Jf, N);
1942 return +sqrt(6 * (2 * Jf + 1) * (2 * Ji + 1)) * wigner6j(1, 1, 1, Ji, Jf, N);
1943}
1944
1946 // Default data and values for this type
1947 SingleLineExternal data;
1948 data.selfbroadening = true;
1949 data.bathbroadening = true;
1950 data.lineshapetype = LineShape::Type::VP;
1951 data.species.resize(2);
1952
1953 // This contains the rest of the line to parse. At the beginning the
1954 // entire line. Line gets shorter and shorter as we continue to
1955 // extract stuff from the beginning.
1956 String line;
1957
1958 // Look for more comments?
1959 bool comment = true;
1960
1961 while (comment) {
1962 // Return true if eof is reached:
1963 if (is.eof()) return data;
1964
1965 // Throw runtime_error if stream is bad:
1966 ARTS_USER_ERROR_IF(!is, "Stream bad.");
1967
1968 // Read line from file into linebuffer:
1969 getline(is, line);
1970
1971 // It is possible that we were exactly at the end of the file before
1972 // calling getline. In that case the previous eof() was still false
1973 // because eof() evaluates only to true if one tries to read after the
1974 // end of the file. The following check catches this.
1975 if (line.nelem() == 0 && is.eof()) return data;
1976
1977 // Because of the fixed FORTRAN format, we need to break up the line
1978 // explicitly in apropriate pieces. Not elegant, but works!
1979
1980 // Extract center frequency:
1981 // Initialization of v is important, because v stays the same
1982 // if line is empty.
1983 // JPL position in MHz:
1984 Numeric v = 0.0;
1985
1986 // Extract JPL position:
1987 extract(v, line, 13);
1988
1989 // check for empty line
1990 if (v != 0.0) {
1991 // ARTS position in Hz:
1992 data.line.F0 = v * 1E6;
1993
1994 comment = false;
1995 }
1996 }
1997
1998 // Accuracy for line position
1999 {
2000 Numeric df;
2001 extract(df, line, 8);
2002 }
2003
2004 // Intensity.
2005 {
2006 // JPL has log (10) of intensity in nm2 MHz at 300 Kelvin.
2007 //
2008 // We need to do the following:
2009 // 1. take 10^intensity
2010 // 2. convert to cm-1/(molecule * cm-2): devide by c * 1e10
2011 // 3. Convert frequency from wavenumber to Hz (factor 1e2 * c)
2012 // 4. Convert [molec * cm-2] to [molec * m-2] (factor 1e-4)
2013
2014 Numeric s;
2015
2016 // Extract JPL intensity:
2017 extract(s, line, 8);
2018
2019 // remove log
2020 s = pow((Numeric)10., s);
2021
2022 // Convert to ARTS units (Hz / (molec * m-2) ), or shorter: Hz*m^2
2023 data.line.I0 = s / 1E12;
2024 }
2025
2026 // Degrees of freedom
2027 {
2028 Index dr;
2029
2030 // Extract degrees of freedom
2031 extract(dr, line, 2);
2032 }
2033
2034 // Lower state energy.
2035 {
2036 // JPL parameter is in wavenumbers (cm^-1).
2037 // We have to convert this to the ARTS unit Joule.
2038
2039 // Extract from Catalogue line
2040 extract(data.line.E0, line, 10);
2041
2042 // Convert to Joule:
2043 data.line.E0 = wavenumber_to_joule(data.line.E0);
2044 }
2045
2046 // Upper state degeneracy
2047 {
2048 Index gup;
2049
2050 // Extract upper state degeneracy
2051 extract(gup, line, 3);
2052 }
2053
2054 // Tag number
2055 Index tag;
2056 {
2057 // Extract Tag number
2058 extract(tag, line, 7);
2059
2060 // make sure tag is not negative (damned jpl cat):
2061 tag = tag > 0 ? tag : -tag;
2062 }
2063
2064 // Set line ID
2066
2067 // Air broadening parameters: unknown to jpl, use old iup forward
2068 // model default values, which is mostly set to 0.0025 GHz/hPa, even
2069 // though for some lines the pressure broadening is given explicitly
2070 // in the program code. The explicitly given values are ignored and
2071 // only the default value is set. Self broadening was in general not
2072 // considered in the old forward model.
2073 Numeric agam, sgam;
2074 {
2075 // ARTS parameter in Hz/Pa:
2076 agam = 2.5E4;
2077
2078 // ARTS parameter in Hz/Pa:
2079 sgam = agam;
2080 }
2081
2082 // Temperature coefficient of broadening parameters. Was set to 0.75
2083 // in old forward model, even though for some lines the parameter is
2084 // given explicitly in the program code. The explicitly given values
2085 // are ignored and only the default value is set. Self broadening
2086 // not considered.
2087 Numeric nair, nself;
2088 {
2089 nair = 0.75;
2090 nself = 0.0;
2091 }
2092
2093 // Reference temperature for broadening parameter in K, was
2094 // generally set to 300 K in old forward model, with the exceptions
2095 // as already mentioned above: //DEPRECEATED but is same as for mti0 so moving on
2096 // {
2097 // mtgam = 300.0;
2098 // }
2099
2100 // Pressure shift: not given in JPL, set to 0
2101 Numeric psf;
2102 { psf = 0.0; }
2103
2104 // These were all the parameters that we can extract from
2105 // JPL. However, we still have to set the reference temperatures
2106 // to the appropriate value:
2107
2108 // Reference temperature for Intensity in K.
2109 data.T0 = 300.0;
2110
2111 // Set line shape computer
2112 data.line.lineshape = LineShape::Model(sgam, nself, agam, nair, psf);
2113
2114 // That's it!
2115 data.bad = false;
2116 return data;
2117}
2118
2119bool Absorption::Lines::OK() const noexcept {
2120 const Index nb = broadeningspecies.nelem();
2121
2122 // Check that the isotopologue is ok
2123 if (not Isotopologue().OK()) return false;
2124
2125 // Test that self and bath is covered by the range if set positive
2126 if (nb < (Index(selfbroadening) + Index(bathbroadening))) return false;
2127
2128 // Test that the temperature is physical
2129 if (T0 <= 0) return false;
2130
2131 // Test that all lines have the correct sized line shape model
2132 if (std::any_of(lines.cbegin(), lines.cend(), [nb](auto& line) {
2133 return line.LineShapeElems() != nb;
2134 }))
2135 return false;
2136
2137 // Test that all lines have the correct sized local quantum numbers
2138 if (std::any_of(lines.cbegin(), lines.cend(), [&](auto& line) {
2139 return line.LocalQuantumElems() != lines.front().LocalQuantumElems();
2140 }))
2141 return false;
2142
2143 // Otherwise everything is fine!
2144 return true;
2145}
2146
2149 SpeciesMass());
2150}
2151
2152namespace Absorption {
2154 std::ostringstream os;
2155 switch (in) {
2156 case CutoffType::None:
2157 os << "No cut-off will be applied.\n";
2158 break;
2159 case CutoffType::ByLine:
2160 os << "The lines will be cut-off " << cutoff
2161 << " Hz from the line center + D0.\n";
2162 break;
2163 case CutoffType::FINAL:
2164 break;
2165 }
2166 return os.str();
2167}
2168
2170 for (auto& val : localquanta.val) qid.val.set(val); // Copy to same struct
2171
2172 zeeman = Zeeman::Model(qid);
2173}
2174
2179 .Data()[0]);
2180}
2181
2184 LineShape::TemperatureModel::LM_AER, d[4], d[5], d[6], d[7]};
2186 LineShape::TemperatureModel::LM_AER, d[8], d[9], d[10], d[11]};
2187 for (auto& sm : lineshape.Data()) {
2188 sm.Y() = Y;
2189 sm.G() = G;
2190 }
2191}
2192
2195 bif >> F0 >> I0 >> E0 >> glow >> gupp >> A >> zeeman;
2196
2198 lineshape.read(bif);
2199
2201 for (auto& val : localquanta.val) val.read(bif);
2202
2203 return bif;
2204}
2205
2208 bof << F0 << I0 << E0 << glow << gupp << A << zeeman;
2209
2211 lineshape.write(bof);
2212
2214 for (auto& val : localquanta.val) val.write(bof);
2215
2216 return bof;
2217}
2218
2221 NumLines() not_eq 0 and NumLocalQuanta() not_eq sl.LocalQuantumElems(),
2222 "Error calling appending function, bad size of quantum numbers");
2223
2225 NumLines() not_eq 0 and
2226 sl.LineShapeElems() not_eq lines[0].LineShapeElems(),
2227 "Error calling appending function, bad size of broadening species");
2228
2229 lines.push_back(std::move(sl));
2230}
2231
2234 NumLines() not_eq 0 and NumLocalQuanta() not_eq sl.LocalQuantumElems(),
2235 "Error calling appending function, bad size of quantum numbers\n"
2236 "Type of quantum numbers in band: ",
2237 lines.front().localquanta.val,
2238 "\nType of quantum numbers in new line:",
2239 sl.localquanta.val);
2240
2242 NumLines() not_eq 0 and
2243 sl.LineShapeElems() not_eq lines.front().LineShapeElems(),
2244 "Error calling appending function, bad size of broadening species");
2245
2246 lines.push_back(sl);
2247}
2248
2250 const QuantumIdentifier& qid) const
2252 if (sle.bad) return false;
2253 if (sle.selfbroadening not_eq selfbroadening) return false;
2254 if (sle.bathbroadening not_eq bathbroadening) return false;
2255 if (sle.cutoff not_eq cutoff) return false;
2256 if (sle.mirroring not_eq mirroring) return false;
2257 if (sle.population not_eq population) return false;
2258 if (sle.normalization not_eq normalization) return false;
2259 if (sle.lineshapetype not_eq lineshapetype) return false;
2260 if (sle.T0 not_eq T0) return false;
2261 if (sle.cutofffreq not_eq cutofffreq) return false;
2262 if (sle.linemixinglimit not_eq linemixinglimit) return false;
2263 if (quantumidentity not_eq qid) return false;
2264 if (not std::equal(sle.species.cbegin(),
2265 sle.species.cend(),
2266 broadeningspecies.cbegin(),
2267 broadeningspecies.cend()))
2268 return false;
2269 if (NumLines() not_eq 0 and
2270 not lines.front().localquanta.same_types_as(sle.line.localquanta))
2271 return false;
2272 if (NumLines() not_eq 0 and
2273 not sle.line.lineshape.Match(lines.front().lineshape).first)
2274 return false;
2275 return true;
2276}
2277
2278std::pair<bool, bool> Lines::Match(const Lines& l) const noexcept {
2279 // Note: The pair here is first: matching and second: nullable
2280 if (l.selfbroadening not_eq selfbroadening) return {false, false};
2281 if (l.bathbroadening not_eq bathbroadening) return {false, false};
2282 if (l.cutoff not_eq cutoff) return {false, false};
2283 if (l.mirroring not_eq mirroring) return {false, false};
2284 if (l.population not_eq population) return {false, false};
2285 if (l.normalization not_eq normalization) return {false, false};
2286 if (l.lineshapetype not_eq lineshapetype) return {false, false};
2287 if (l.T0 not_eq T0) return {false, false};
2288 if (l.cutofffreq not_eq cutofffreq) return {false, false};
2289 if (l.linemixinglimit not_eq linemixinglimit) return {false, false};
2290 if (l.quantumidentity not_eq quantumidentity) return {false, false};
2291 if (not std::equal(l.broadeningspecies.cbegin(),
2292 l.broadeningspecies.cend(),
2293 broadeningspecies.cbegin(),
2294 broadeningspecies.cend()))
2295 return {false, false};
2296 if (NumLines() not_eq 0 and l.NumLines() not_eq 0 and
2297 not lines.front().localquanta.same_types_as(l.lines.front().localquanta))
2298 return {false, false};
2299 if (NumLines() not_eq 0 and l.NumLines() not_eq 0) {
2300 if (auto matchpair =
2301 l.lines.front().lineshape.Match(lines.front().lineshape);
2302 not matchpair.first)
2303 return matchpair;
2304 }
2305
2306 return {true, true};
2307}
2308
2310 std::sort(
2311 lines.begin(), lines.end(), [](const SingleLine& a, const SingleLine& b) {
2312 return a.F0 < b.F0;
2313 });
2314}
2315
2317 std::sort(lines.begin(),
2318 lines.end(),
2319 [](const SingleLine& a, const SingleLine& b) { return a.A < b.A; });
2320}
2321
2323 return NumLines() ? LineShape::ModelShape2MetaData(lines[0].lineshape) : "";
2324}
2325
2328 ARTS_ASSERT(qns.val.has(QuantumNumberType::F) or
2329 qns.val.has(QuantumNumberType::J))
2330 return qns.val.has(QuantumNumberType::F) ? qns.val[QuantumNumberType::F]
2331 : qns.val[QuantumNumberType::J];
2332}
2333
2334Index Lines::ZeemanCount(size_t k, Zeeman::Polarization type) const noexcept {
2335 if (type == Zeeman::Polarization::None) return 1;
2336
2337 // Select F before J but assume one of them exist
2338 auto& val = get(lines[k].localquanta);
2339 return Zeeman::nelem(val.upp(), val.low(), type);
2340}
2341
2344 Index i) const noexcept {
2345 if (type == Zeeman::Polarization::None) return 1.0;
2346
2347 // Select F before J but assume one of them exist
2348 auto& val = get(lines[k].localquanta);
2349 return lines[k].zeeman.Strength(val.upp(), val.low(), type, i);
2350}
2351
2354 Index i) const noexcept {
2355 if (type == Zeeman::Polarization::None) return 0.0;
2356
2357 // Select F before J but assume one of them exist
2358 auto& val = get(lines[k].localquanta);
2359 return lines[k].zeeman.Splitting(val.upp(), val.low(), type, i);
2360}
2361
2363 for (auto& line : lines) line.SetAutomaticZeeman(quantumidentity);
2364}
2365
2366Numeric Lines::F_mean(const ConstVectorView& wgts) const noexcept {
2367 const Numeric val =
2368 std::inner_product(lines.cbegin(),
2369 lines.cend(),
2370 wgts.begin(),
2371 0.0,
2372 std::plus<>(),
2373 [](const auto& a, const auto& b) { return a.F0 * b; });
2374 const Numeric div = wgts.sum();
2375 return val / div;
2376}
2377
2378Numeric Lines::F_mean(Numeric T) const noexcept {
2379 if (T <= 0) T = T0;
2380
2381 const Index n = NumLines();
2382 const Numeric QT = single_partition_function(T, Isotopologue());
2383 const Numeric QT0 = single_partition_function(T0, Isotopologue());
2384 const Numeric ratiopart = QT0 / QT;
2385
2386 Vector wgts(n);
2387 for (Index i = 0; i < n; i++) {
2388 const Numeric pop0 =
2389 (lines[i].gupp / QT0) * boltzman_factor(T0, lines[i].E0);
2390 const Numeric pop = pop0 * ratiopart * boltzman_ratio(T, T0, lines[i].E0);
2391 const Numeric dip_squared =
2392 -lines[i].I0 /
2393 (pop0 * lines[i].F0 *
2394 std::expm1(-(Constant::h * lines[i].F0) / (Constant::k * T0)));
2395 wgts[i] = pop * dip_squared;
2396 }
2397
2398 return F_mean(wgts);
2399}
2400
2401#pragma GCC diagnostic push
2402#pragma GCC diagnostic ignored "-Wreturn-type"
2403Numeric Lines::CutoffFreq(size_t k, Numeric shift) const noexcept {
2404 switch (cutoff) {
2405 case CutoffType::ByLine:
2406 return lines[k].F0 + cutofffreq +
2407 (mirroring == MirroringType::Manual ? -shift : shift);
2408 case CutoffType::None:
2410 case CutoffType::FINAL:
2411 break;
2412 }
2413}
2414#pragma GCC diagnostic pop
2415
2416#pragma GCC diagnostic push
2417#pragma GCC diagnostic ignored "-Wreturn-type"
2418Numeric Lines::CutoffFreqMinus(size_t k, Numeric shift) const noexcept {
2419 switch (cutoff) {
2420 case CutoffType::ByLine:
2421 return lines[k].F0 - cutofffreq +
2422 (mirroring == MirroringType::Manual ? -shift : shift);
2423 case CutoffType::None:
2424 return std::numeric_limits<Numeric>::lowest();
2425 case CutoffType::FINAL:
2426 break;
2427 }
2428}
2429#pragma GCC diagnostic pop
2430
2432 for (auto& line : lines) line.read(is);
2433 return is;
2434}
2435
2437 for (auto& line : lines) line.write(os);
2438 return os;
2439}
2440
2442 QuantumIdentifier qid = quantumidentity;
2443 for (auto& a : lines[k].localquanta.val) qid.val.set(a);
2444 return qid;
2445}
2446
2448 const Index n = NumLines();
2449 const Index m = NumBroadeners();
2450 if (not n) return;
2451
2452 for (Index j = 0; j < m; j++) {
2453 for (Index i = 0; i < LineShape::nVars; i++) {
2454 // Find a common type (same or none) or throw if there is none
2455 LineShape::TemperatureModel t = LineShape::TemperatureModel::None;
2456 for (auto& line : lines) {
2457 if (auto& data = line.lineshape[j].Data()[i];
2459 if (t == LineShape::TemperatureModel::None) t = data.type;
2461 t not_eq data.type,
2462 "Cannot make a common line shape model for the band as there are multiple non-empty types: ",
2463 data.type,
2464 " and ",
2465 t)
2466 }
2467 }
2468
2469 // Set the common type
2470 for (auto& line : lines) {
2471 if (auto& data = line.lineshape[j].Data()[i]; data.type not_eq t) {
2473 }
2474 }
2475 }
2476 }
2477}
2478} // namespace Absorption
2479
2481 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species) {
2482 for (auto& abs_lines : abs_lines_per_species) {
2483 for (auto& band : abs_lines) {
2484 switch (band.population) {
2485 case Absorption::PopulationType::LTE:
2486 LTE = true;
2487 break;
2488 case Absorption::PopulationType::NLTE:
2489 NLTE = true;
2490 break;
2491 case Absorption::PopulationType::VibTemps:
2492 VibTemps = true;
2493 break;
2494 case Absorption::PopulationType::ByHITRANFullRelmat:
2495 ByHITRANFullRelmat = true;
2496 break;
2497 case Absorption::PopulationType::ByHITRANRosenkranzRelmat:
2499 break;
2500 case Absorption::PopulationType::ByMakarovFullRelmat:
2501 ByMakarovFullRelmat = true;
2502 break;
2503 case Absorption::PopulationType::ByRovibLinearDipoleLineMixing:
2505 break;
2506 case Absorption::PopulationType::FINAL: { /* leave list */
2507 }
2508 }
2509 }
2510 }
2511}
2512
2514 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species) {
2515 for (auto& abs_lines : abs_lines_per_species) {
2516 for (auto& band : abs_lines) {
2517 switch (band.cutoff) {
2518 case Absorption::CutoffType::None:
2519 None = true;
2520 break;
2521 case Absorption::CutoffType::ByLine:
2522 ByLine = true;
2523 break;
2524 case Absorption::CutoffType::FINAL: { /* leave list */
2525 }
2526 }
2527 }
2528 }
2529}
2530
2532 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species) {
2533 for (auto& abs_lines : abs_lines_per_species) {
2534 for (auto& band : abs_lines) {
2535 switch (band.lineshapetype) {
2536 case LineShape::Type::VP:
2537 VP = true;
2538 break;
2539 case LineShape::Type::DP:
2540 DP = true;
2541 break;
2542 case LineShape::Type::LP:
2543 LP = true;
2544 break;
2545 case LineShape::Type::SDVP:
2546 SDVP = true;
2547 break;
2548 case LineShape::Type::HTP:
2549 HTP = true;
2550 break;
2551 case LineShape::Type::FINAL: { /* leave list */
2552 }
2553 }
2554 }
2555 }
2556}
2557
2559 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species) {
2560 for (auto& abs_lines : abs_lines_per_species) {
2561 for (auto& band : abs_lines) {
2562 switch (band.mirroring) {
2563 case Absorption::MirroringType::None:
2564 None = true;
2565 break;
2566 case Absorption::MirroringType::Lorentz:
2567 Lorentz = true;
2568 break;
2569 case Absorption::MirroringType::SameAsLineShape:
2570 SameAsLineShape = true;
2571 break;
2572 case Absorption::MirroringType::Manual:
2573 Manual = true;
2574 break;
2575 case Absorption::MirroringType::FINAL: { /* leave list */
2576 }
2577 }
2578 }
2579 }
2580}
2581
2583 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species) {
2584 for (auto& abs_lines : abs_lines_per_species) {
2585 for (auto& band : abs_lines) {
2586 switch (band.normalization) {
2587 case Absorption::NormalizationType::None:
2588 None = true;
2589 break;
2590 case Absorption::NormalizationType::VVH:
2591 VVH = true;
2592 break;
2593 case Absorption::NormalizationType::VVW:
2594 VVW = true;
2595 break;
2596 case Absorption::NormalizationType::RQ:
2597 RQ = true;
2598 break;
2599 case Absorption::NormalizationType::SFS:
2600 SFS = true;
2601 break;
2602 case Absorption::NormalizationType::FINAL: { /* leave list */
2603 }
2604 }
2605 }
2606 }
2607}
2608
2609namespace {
2610template <typename T>
2611constexpr std::size_t req_spaces(T my_enum) {
2612 constexpr std::size_t n = []() {
2613 std::size_t longest = 0;
2614 for (auto& x : Absorption::enumstrs::CutoffTypeNames) {
2615 longest = std::max(x.length(), longest);
2616 }
2617 for (auto& x : Absorption::enumstrs::MirroringTypeNames) {
2618 longest = std::max(x.length(), longest);
2619 }
2620 for (auto& x : Absorption::enumstrs::NormalizationTypeNames) {
2621 longest = std::max(x.length(), longest);
2622 }
2623 for (auto& x : Absorption::enumstrs::PopulationTypeNames) {
2624 longest = std::max(x.length(), longest);
2625 }
2626 for (auto& x : LineShape::enumstrs::TypeNames) {
2627 longest = std::max(x.length(), longest);
2628 }
2629 return longest + 1;
2630 }();
2631 return n - toString(my_enum).length();
2632}
2633
2634auto spaces(std::size_t n) { return std::basic_string(n, ' '); }
2635} // namespace
2636
2637std::ostream& operator<<(std::ostream& os, AbsorptionCutoffTagTypeStatus val) {
2638 // Trick to never forget to update
2639 Absorption::CutoffType x{Absorption::CutoffType::FINAL};
2640 switch (x) {
2641 case Absorption::CutoffType::FINAL:
2642 os << "Cutoff tag types:\n";
2643 [[fallthrough]];
2644 case AbsorptionCutoffType::ByLine:
2645 os << " ByLine:" << spaces(req_spaces(AbsorptionCutoffType::ByLine))
2646 << val.ByLine << '\n';
2647 [[fallthrough]];
2648 case AbsorptionCutoffType::None:
2649 os << " None:" << spaces(req_spaces(AbsorptionCutoffType::None))
2650 << val.None;
2651 }
2652 return os;
2653}
2654
2655std::ostream& operator<<(std::ostream& os,
2657 // Trick to never forget to update
2658 Absorption::MirroringType x{Absorption::MirroringType::FINAL};
2659 switch (x) {
2660 case Absorption::MirroringType::FINAL:
2661 os << "Mirroring tag types:\n";
2662 [[fallthrough]];
2663 case Absorption::MirroringType::None:
2664 os << " None:" << spaces(req_spaces(Absorption::MirroringType::None))
2665 << val.None << '\n';
2666 [[fallthrough]];
2667 case Absorption::MirroringType::Lorentz:
2668 os << " Lorentz:"
2669 << spaces(req_spaces(Absorption::MirroringType::Lorentz))
2670 << val.Lorentz << '\n';
2671 [[fallthrough]];
2672 case Absorption::MirroringType::SameAsLineShape:
2673 os << " SameAsLineShape:"
2674 << spaces(req_spaces(Absorption::MirroringType::SameAsLineShape))
2675 << val.SameAsLineShape << '\n';
2676 [[fallthrough]];
2677 case Absorption::MirroringType::Manual:
2678 os << " Manual:"
2679 << spaces(req_spaces(Absorption::MirroringType::Manual)) << val.Manual;
2680 }
2681 return os;
2682}
2683
2684std::ostream& operator<<(std::ostream& os,
2686 // Trick to never forget to update
2687 Absorption::NormalizationType x{Absorption::NormalizationType::FINAL};
2688 switch (x) {
2689 case Absorption::NormalizationType::FINAL:
2690 os << "Normalization tag types:\n";
2691 [[fallthrough]];
2692 case Absorption::NormalizationType::None:
2693 os << " None:"
2694 << spaces(req_spaces(Absorption::NormalizationType::None)) << val.None
2695 << '\n';
2696 [[fallthrough]];
2697 case Absorption::NormalizationType::VVH:
2698 os << " VVH:" << spaces(req_spaces(Absorption::NormalizationType::VVH))
2699 << val.VVH << '\n';
2700 [[fallthrough]];
2701 case Absorption::NormalizationType::VVW:
2702 os << " VVW:" << spaces(req_spaces(Absorption::NormalizationType::VVW))
2703 << val.VVW << '\n';
2704 [[fallthrough]];
2705 case Absorption::NormalizationType::RQ:
2706 os << " RQ:" << spaces(req_spaces(Absorption::NormalizationType::RQ))
2707 << val.RQ << '\n';
2708 [[fallthrough]];
2709 case Absorption::NormalizationType::SFS:
2710 os << " SFS:" << spaces(req_spaces(Absorption::NormalizationType::SFS))
2711 << val.SFS;
2712 }
2713 return os;
2714}
2715
2716std::ostream& operator<<(std::ostream& os,
2718 // Trick to never forget to update
2719 Absorption::PopulationType x{Absorption::PopulationType::FINAL};
2720 switch (x) {
2721 case Absorption::PopulationType::FINAL:
2722 os << "Population tag type:\n";
2723 [[fallthrough]];
2724 case Absorption::PopulationType::LTE:
2725 os << " LTE:" << spaces(req_spaces(Absorption::PopulationType::LTE))
2726 << val.LTE << '\n';
2727 [[fallthrough]];
2728 case Absorption::PopulationType::NLTE:
2729 os << " NLTE:" << spaces(req_spaces(Absorption::PopulationType::NLTE))
2730 << val.NLTE << '\n';
2731 [[fallthrough]];
2732 case Absorption::PopulationType::VibTemps:
2733 os << " VibTemps:"
2734 << spaces(req_spaces(Absorption::PopulationType::VibTemps))
2735 << val.VibTemps << '\n';
2736 [[fallthrough]];
2737 case Absorption::PopulationType::ByHITRANRosenkranzRelmat:
2738 os << " ByHITRANRosenkranzRelmat:"
2739 << spaces(req_spaces(
2740 Absorption::PopulationType::ByHITRANRosenkranzRelmat))
2741 << val.ByHITRANRosenkranzRelmat << '\n';
2742 [[fallthrough]];
2743 case Absorption::PopulationType::ByHITRANFullRelmat:
2744 os << " ByHITRANFullRelmat:"
2745 << spaces(req_spaces(Absorption::PopulationType::ByHITRANFullRelmat))
2746 << val.ByHITRANFullRelmat << '\n';
2747 [[fallthrough]];
2748 case Absorption::PopulationType::ByMakarovFullRelmat:
2749 os << " ByMakarovFullRelmat:"
2750 << spaces(req_spaces(Absorption::PopulationType::ByMakarovFullRelmat))
2751 << val.ByMakarovFullRelmat << '\n';
2752 [[fallthrough]];
2753 case Absorption::PopulationType::ByRovibLinearDipoleLineMixing:
2754 os << " ByRovibLinearDipoleLineMixing:"
2755 << spaces(req_spaces(
2756 Absorption::PopulationType::ByRovibLinearDipoleLineMixing))
2758 }
2759 return os;
2760}
2761
2762std::ostream& operator<<(std::ostream& os,
2764 // Trick to never forget to update
2765 LineShapeType x{LineShapeType::FINAL};
2766 switch (x) {
2767 case LineShapeType::FINAL:
2768 os << "Line shape tag type:\n";
2769 [[fallthrough]];
2770 case LineShapeType::DP:
2771 os << " DP:" << spaces(req_spaces(LineShapeType::DP)) << val.DP
2772 << '\n';
2773 [[fallthrough]];
2774 case LineShapeType::LP:
2775 os << " LP:" << spaces(req_spaces(LineShapeType::LP)) << val.LP
2776 << '\n';
2777 [[fallthrough]];
2778 case LineShapeType::VP:
2779 os << " VP:" << spaces(req_spaces(LineShapeType::VP)) << val.VP
2780 << '\n';
2781 [[fallthrough]];
2782 case LineShapeType::SDVP:
2783 os << " SDVP:" << spaces(req_spaces(LineShapeType::SDVP)) << val.SDVP
2784 << '\n';
2785 [[fallthrough]];
2786 case LineShapeType::HTP:
2787 os << " HTP:" << spaces(req_spaces(LineShapeType::HTP)) << val.HTP;
2788 }
2789 return os;
2790}
2791
2792std::ostream& operator<<(std::ostream& os, AbsorptionTagTypesStatus val) {
2793 return os << "Catalog tag summary:\n " << val.cutoff << "\n "
2794 << val.lineshapetype << "\n " << val.mirroring << "\n "
2795 << val.normalization << "\n " << val.population;
2796}
2797
2799 Index i,
2800 const ArrayOfArrayOfSpeciesTag& abs_species,
2801 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species) {
2802 const Index n = abs_species.nelem();
2803
2804 Index ispec{0};
2805 while (ispec < n and i >= abs_lines_per_species[ispec].nelem()) {
2806 i -= abs_lines_per_species[ispec].nelem();
2807 ++ispec;
2808 }
2809
2810 return {ispec, i};
2811}
2812
2813namespace Absorption {
2814bool any_cutoff(const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species) {
2815 return std::any_of(
2816 abs_lines_per_species.cbegin(),
2817 abs_lines_per_species.cend(),
2818 [](auto& abs_lines) {
2819 return std::any_of(
2820 abs_lines.cbegin(), abs_lines.cend(), [](auto& band) {
2821 return band.cutoff not_eq CutoffType::None;
2822 });
2823 });
2824}
2825
2827 ARTS_USER_ERROR_IF(Quantum::Number::common_value_type(x) == Quantum::Number::ValueType::S, "Cannot get a rational out from quantum number type ", x)
2828 if (quantumidentity.val.has(x)) {
2829 auto& val = quantumidentity.val[x];
2830 return std::max(val.low(), val.upp());
2831 }
2832
2833 Rational out{std::numeric_limits<Index>::lowest()};
2834 for (auto& line : lines) {
2835 ARTS_USER_ERROR_IF(not line.localquanta.val.has(x), "No ", x, " in some line(s)")
2836 auto& val = line.localquanta.val[x];
2837 out = std::max(std::max(val.low(), val.upp()), out);
2838 }
2839 return out;
2840}
2841} // namespace Absorption
Numeric wavenumber_to_joule(Numeric e)
A little helper function to convert energy from units of wavenumber (cm^-1) to Joule (J).
Definition: absorption.cc:87
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.
Common ARTS conversions.
Index nelem() const ARTS_NOEXCEPT
Definition: array.h:92
A constant view of a Vector.
Definition: matpackI.h:512
Index nelem() const noexcept
Returns the number of elements.
Definition: matpackI.h:536
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.
Index nelem() const
Return number of quantum numbers.
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.
Deals with internal derivatives, Jacobian definition, and OEM calculations.
Definition: jacobian.h:325
The Vector class.
Definition: matpackI.h:899
Main Zeeman Model.
Definition: zeemandata.h:296
Binary output file stream class.
Definition: bifstream.h:43
Binary output file stream class.
Definition: bofstream.h:42
Input manipulator class for doubles to enable nan and inf parsing.
Definition: double_imanip.h:42
Index nelem() const
Definition: mystring.h:189
Helper macros for debugging.
#define ARTS_NOEXCEPT
Definition: debug.h:80
#define ARTS_ASSERT(condition,...)
Definition: debug.h:83
#define ARTS_USER_ERROR(...)
Definition: debug.h:150
#define ARTS_USER_ERROR_IF(condition,...)
Definition: debug.h:134
constexpr std::string_view toString(EnergyLevelMapType x) noexcept
constexpr bool good_enum(EnumType x) noexcept
Checks if the enum number is good.
Definition: enums.h:22
This file contains basic functions to handle ASCII files.
#define max(a, b)
Numeric boltzman_factor(Numeric T, Numeric E0)
Computes exp(- E0/kT)
Definition: linescaling.cc:109
Numeric boltzman_ratio(const Numeric &T, const Numeric &T0, const Numeric &E0)
Computes exp(E0/c (T - T0) / (T * T0))
Definition: linescaling.cc:94
Numeric single_partition_function(const Numeric &T, const Species::IsotopeRecord &ir)
Computes the partition function at one temperature.
Definition: linescaling.cc:23
Constains various line scaling functions.
Contains the line shape namespace.
LineShape::Type LineShapeType
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
void extract(T &x, String &line, std::size_t n)
Extract something from the beginning of a string.
Definition: mystring.h:233
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)
Numeric cutoff
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 e
Elementary charge convenience name [C].
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.
constexpr Numeric l(const Index p0, const Index n, const Numeric x, const SortedVectorType &xi, const Index j, const std::pair< Numeric, Numeric > cycle={ -180, 180}) noexcept
QuantumIdentifier id_from_lookup(Index tag)
Finds the ID of the ARTS species from JPL.
Definition: jpl_species.cc:204
Model vector2modellm(Vector x, LegacyLineMixingData::TypeLM type)
LineShape::Model from legacy input vector.
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.
Definition: zeemandata.h:142
Polarization
Zeeman polarization selection.
Definition: zeemandata.h:44
Quantum::Number::Type QuantumNumberType
Contains the rational class definition.
Numeric pow(const Rational base, Numeric exp)
Power of.
Definition: rational.h:713
Numeric sqrt(const Rational r)
Square root.
Definition: rational.h:705
constexpr bool iseven(const Rational r) noexcept
Returns true if even integer.
Definition: rational.h:902
#define N
Definition: rng.cc:164
AbsorptionCutoffTagTypeStatus(const ArrayOfArrayOfAbsorptionLines &)
AbsorptionLineShapeTagTypeStatus(const ArrayOfArrayOfAbsorptionLines &)
AbsorptionMirroringTagTypeStatus(const ArrayOfArrayOfAbsorptionLines &)
AbsorptionNormalizationTagTypeStatus(const ArrayOfArrayOfAbsorptionLines &)
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.
Numeric ZeemanSplitting(size_t k, Zeeman::Polarization type, Index i) const noexcept
Returns the splitting of a Zeeman split line.
SingleLine PopLine(Index) noexcept
Pops a single line.
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.
LineShape::Output ShapeParameters_dT(size_t k, Numeric T, Numeric P, const Vector &vmrs) const noexcept
Line shape parameters temperature derivatives.
LineShape::Output ShapeParameters_dVMR(size_t k, Numeric T, Numeric P, const QuantumIdentifier &vmr_qid) const noexcept
Line shape parameters vmr derivative.
Array< SingleLine > lines
A list of individual lines.
bool OK() const noexcept
Numeric ShapeParameter_dInternal(size_t k, Numeric T, Numeric P, const Vector &vmrs, const RetrievalQuantity &derivative) const noexcept
Line shape parameter internal derivative.
Index NumLines() const noexcept
Number of lines.
bofstream & write(bofstream &os) const
Binary write for Lines.
Index NumBroadeners() const noexcept
Number of broadening species.
Numeric DopplerConstant(Numeric T) const noexcept
QuantumIdentifier QuantumIdentityOfLine(Index k) const noexcept
bifstream & read(bifstream &is)
Binary read for Lines.
LineShape::Output ShapeParameters(size_t k, Numeric T, Numeric P, const Vector &vmrs) const noexcept
Line shape parameters.
bool MatchWithExternal(const SingleLineExternal &sle, const QuantumIdentifier &quantumidentity) const ARTS_NOEXCEPT
Checks if an external line matches this structure.
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.
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.
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.
String MetaData() const
Returns a printable statement about the lines.
Index ZeemanCount(size_t k, Zeeman::Polarization type) const noexcept
Returns the number of Zeeman split lines.
Numeric SpeciesMass() const noexcept
Mass of the molecule.
Index LineShapePos(const Species::Species spec) const noexcept
Position among broadening species or -1.
QuantumIdentifier quantumidentity
Catalog ID.
Numeric ZeemanStrength(size_t k, Zeeman::Polarization type, Index i) const noexcept
Returns the strength of a Zeeman split line.
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.
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.
Index LocalQuantumElems() const noexcept
Number of lower quantum numbers.
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.
Index LineShapeElems() const noexcept
Number of lineshape elements.
void SetLineMixing2SecondOrderData(const Vector &d)
Set the line mixing model to 2nd order.
Numeric I0
Reference intensity.
Coefficients and temperature model for SingleSpeciesModel.
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.
Implements rational numbers to work with other ARTS types.
Definition: rational.h:43
String FullName() const noexcept
Definition: isotopologues.h:49
#define d
#define v
#define a
#define c
#define b
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.