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