ARTS 2.5.4 (git: bcd8c674)
m_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
30#include "absorptionlines.h"
31#include "array.h"
32#include "arts_omp.h"
33#include "artstime.h"
34#include "auto_md.h"
35#include "debug.h"
36#include "enums.h"
37#include "file.h"
38#include "global_data.h"
39#include "m_xml.h"
40#include "quantum_numbers.h"
41#include <algorithm>
42#include <atomic>
43#include <exception>
44#include <iterator>
45
49
51 const Verbosity&)
52{
53 abs_lines.erase(std::remove_if(abs_lines.begin(), abs_lines.end(),
54 [](auto& x){return x.NumLines() == 0;}),
55 abs_lines.end());
56}
57
59 const Verbosity& verbosity)
60{
61 const Index n = abs_lines.nelem();
62
63 for (Index i=0; i<n; i++) {
64 AbsorptionLines& band = abs_lines[i];
65 if (band.NumLines()) {
66 for (Index j=i+1; j<n; j++) {
67 if (band.Match(abs_lines[j]).first) {
68 for (auto& line: abs_lines[j].lines) {
69 band.AppendSingleLine(line);
70 }
71 abs_lines[j].lines.clear();
72 }
73 }
74 }
75 }
76
77 abs_linesRemoveEmptyBands(abs_lines, verbosity);
78}
79
81 const Verbosity& verbosity)
82{
83 for (auto& abs_lines: abs_lines_per_species) abs_linesFlatten(abs_lines, verbosity);
84}
85
89
90
91
100 const QuantumIdentifier& global_qid) {
101 auto band =
102 std::find_if(abs_lines.begin(), abs_lines.end(), [&](const auto& li) {
103 return li.MatchWithExternal(sline, global_qid);
104 });
105 if (band not_eq abs_lines.end()) {
106 band->AppendSingleLine(sline.line);
107 } else {
108 abs_lines.emplace_back(sline.selfbroadening,
109 sline.bathbroadening,
110 sline.cutoff,
111 sline.mirroring,
112 sline.population,
113 sline.normalization,
114 sline.lineshapetype,
115 sline.T0,
116 sline.cutofffreq,
117 sline.linemixinglimit,
118 global_qid,
119 sline.species,
120 Array<AbsorptionSingleLine>{sline.line});
121 }
122}
123
125
132 const ArrayOfAbsorptionLines& local_lines) {
133 for (auto& band : local_lines) {
134 if (auto ptr = std::find_first_of(
135 abs_lines.begin(),
136 abs_lines.end(),
137 &band,
138 &band,
139 [](const auto& abs_band, const auto& local_band) {
140 return abs_band.Match(local_band).first;
141 });
142 ptr == abs_lines.end())
143 abs_lines.push_back(band);
144 else
145 for (auto& line : band.lines) ptr->lines.push_back(line);
146 }
147}
148
157 for(auto qn: qns) {
158 if (qid.val.has(qn)) {
159 out.val.set(qid.val[qn]);
160 }
161 }
162 return out;
163}
164
171Array<QuantumNumberType> string2vecqn(std::string_view qnstr) {
172 using namespace Quantum::Number;
173
174 if (qnstr == "DEFAULT_GLOBAL") {
176 }
177 if (qnstr == "DEFAULT_LOCAL") {
178 return Array<Type>(local_types);
179 }
180
181 const Index N = count_items(qnstr);
182 Array<Type> nums(N);
183 for (Index i = 0; i < N; i++) {
184 nums[i] = toTypeOrThrow(items(qnstr, i));
185 }
186
187 return nums;
188}
189
190bool check_local(const Array<QuantumNumberType>& local_state) {
191 using namespace Quantum::Number;
192 for (auto qn: local_state) if (common_value_type(common_value_type(qn), ValueType::H) not_eq ValueType::H) return false;
193 return true;
194}
195
196/* Workspace method: Doxygen documentation will be auto-generated */
198 const String& artscat_file,
199 const Numeric& fmin,
200 const Numeric& fmax,
201 const String& globalquantumnumbers,
202 const String& localquantumnumbers,
203 const String& normalization_option,
204 const String& mirroring_option,
205 const String& population_option,
206 const String& lineshapetype_option,
207 const String& cutoff_option,
208 const Numeric& cutoff_value,
209 const Numeric& linemixinglimit_value,
210 const Verbosity& verbosity)
211{
212 abs_lines.resize(0);
213
214 // Global numbers
215 const Array<QuantumNumberType> global_nums = string2vecqn(globalquantumnumbers);
216
217 // Local numbers
218 const Array<QuantumNumberType> local_nums = string2vecqn(localquantumnumbers);
219 ARTS_USER_ERROR_IF(not check_local(local_nums), "Can only have non-string values in the local state")
220
222
223 ArtsXMLTag tag(verbosity);
224 Index nelem;
225
226 // ARTSCAT data
227 shared_ptr<istream> ifs = nullptr;
228 xml_find_and_open_input_file(ifs, artscat_file, verbosity);
229 istream& is_xml = *ifs;
230 auto a = FILE_TYPE_ASCII;
231 auto b = NUMERIC_TYPE_DOUBLE;
232 auto c = ENDIAN_TYPE_LITTLE;
234 verbosity);
235
236 tag.read_from_stream(is_xml);
237 tag.check_name("Array");
238
239 Index num_arrays;
240 tag.get_attribute_value("nelem", num_arrays);
241
242 ArrayOfAbsorptionLines local_bands(0);
243 for (Index i=0; i<num_arrays; i++) {
244 tag.read_from_stream(is_xml);
245 tag.check_name("ArrayOfLineRecord");
246
247 tag.get_attribute_value("nelem", nelem);
248
249 String version;
250 tag.get_attribute_value("version", version);
251
252 Index artscat_version;
253
254 if (version == "3") {
255 artscat_version = 3;
256 } else if (version.substr(0, 8) != "ARTSCAT-") {
258 "The ARTS line file you are trying to read does not contain a valid version tag.\n"
259 "Probably it was created with an older version of ARTS that used different units.")
260 } else {
261 istringstream is(version.substr(8));
262 is >> artscat_version;
263 }
264
265 ARTS_USER_ERROR_IF (artscat_version < 3 or artscat_version > 5,
266 "Unknown ARTS line file version: ", version)
267
268 bool go_on = true;
269 Index n = 0;
270 while (n<nelem) {
271 n++;
272
273 if (go_on) {
275 switch(artscat_version) {
276 case 3:
278 break;
279 case 4:
281 break;
282 case 5:
284 break;
285 default:
286 ARTS_ASSERT (false, "Bad version!");
287 }
288
289 ARTS_USER_ERROR_IF(sline.bad, "Cannot read line ", n)
290 if (sline.line.F0 < fmin) continue;
291 if (sline.line.F0 > fmax) {go_on = false; continue;}
292
294
295 // Get the global quantum number identifier
296 const QuantumIdentifier global_qid = global_quantumidentifier(global_nums, sline.quantumidentity);
297
298 // Get local quantum numbers into the line
299 for(auto qn: local_nums) {
300 if (sline.quantumidentity.val.has(qn)) {
301 sline.line.localquanta.val.set(sline.quantumidentity.val[qn]);
302 }
303 }
304
305 merge_external_line(local_bands, sline, global_qid);
306 if (local_bands.nelem() > merge_local_lines_size) {
307 merge_local_lines(abs_lines, local_bands);
308 local_bands.resize(0);
309 }
310 } else {
311 String line;
312 getline(is_xml, line);
313 }
314 }
315
316 tag.read_from_stream(is_xml);
317 tag.check_name("/ArrayOfLineRecord");
318 }
319
320 merge_local_lines(abs_lines, local_bands);
321
322 abs_linesNormalization(abs_lines, normalization_option, verbosity);
323 abs_linesMirroring(abs_lines, mirroring_option, verbosity);
324 abs_linesPopulation(abs_lines, population_option, verbosity);
325 abs_linesLineShapeType(abs_lines, lineshapetype_option, verbosity);
326 abs_linesCutoff(abs_lines, cutoff_option, cutoff_value, verbosity);
327 abs_linesLinemixingLimit(abs_lines, linemixinglimit_value, verbosity);
328
329 tag.read_from_stream(is_xml);
330 tag.check_name("/Array");
331}
332
333/* Workspace method: Doxygen documentation will be auto-generated */
335 const String& artscat_file,
336 const Numeric& fmin,
337 const Numeric& fmax,
338 const String& globalquantumnumbers,
339 const String& localquantumnumbers,
340 const String& normalization_option,
341 const String& mirroring_option,
342 const String& population_option,
343 const String& lineshapetype_option,
344 const String& cutoff_option,
345 const Numeric& cutoff_value,
346 const Numeric& linemixinglimit_value,
347 const Verbosity& verbosity)
348{
349 abs_lines.resize(0);
350
351 // Global numbers
352 const Array<QuantumNumberType> global_nums = string2vecqn(globalquantumnumbers);
353
354 // Local numbers
355 const Array<QuantumNumberType> local_nums = string2vecqn(localquantumnumbers);
356 ARTS_USER_ERROR_IF(not check_local(local_nums), "Can only have non-string values in the local state")
357
359
360 ArtsXMLTag tag(verbosity);
361 Index nelem;
362
363 // ARTSCAT data
364 shared_ptr<istream> ifs = nullptr;
365 xml_find_and_open_input_file(ifs, artscat_file, verbosity);
366 istream& is_xml = *ifs;
367 auto a = FILE_TYPE_ASCII;
368 auto b = NUMERIC_TYPE_DOUBLE;
369 auto c = ENDIAN_TYPE_LITTLE;
371 verbosity);
372
373 tag.read_from_stream(is_xml);
374 tag.check_name("ArrayOfLineRecord");
375
376 tag.get_attribute_value("nelem", nelem);
377
378 String version;
379 tag.get_attribute_value("version", version);
380
381 Index artscat_version;
382
383 if (version == "3") {
384 artscat_version = 3;
385 } else if (version.substr(0, 8) != "ARTSCAT-") {
387 "The ARTS line file you are trying to read does not contain a valid version tag.\n"
388 "Probably it was created with an older version of ARTS that used different units.")
389 } else {
390 istringstream is(version.substr(8));
391 is >> artscat_version;
392 }
393
394 ARTS_USER_ERROR_IF (artscat_version < 3 or artscat_version > 5,
395 "Unknown ARTS line file version: ", version)
396
397 bool go_on = true;
398 Index n = 0;
399 ArrayOfAbsorptionLines local_bands(0);
400 while (n<nelem) {
401 n++;
402
403 if (go_on) {
405 switch(artscat_version) {
406 case 3:
408 break;
409 case 4:
411 break;
412 case 5:
414 break;
415 default:
416 ARTS_ASSERT (false, "Bad version!");
417 }
418
419 ARTS_USER_ERROR_IF(sline.bad, "Cannot read line ", n)
420 if (sline.line.F0 < fmin) continue;
421 if (sline.line.F0 > fmax) {go_on = false; continue;}
422
424
425 // Get the global quantum number identifier
426 const QuantumIdentifier global_qid = global_quantumidentifier(global_nums, sline.quantumidentity);
427
428 // Get local quantum numbers into the line
429 for(auto qn: local_nums) {
430 if (sline.quantumidentity.val.has(qn)) {
431 sline.line.localquanta.val.set(sline.quantumidentity.val[qn]);
432 }
433 }
434
435 merge_external_line(local_bands, sline, global_qid);
436 if (local_bands.nelem() > merge_local_lines_size) {
437 merge_local_lines(abs_lines, local_bands);
438 local_bands.resize(0);
439 }
440 } else {
441 String line;
442 getline(is_xml, line);
443 }
444 }
445
446 merge_local_lines(abs_lines, local_bands);
447
448 abs_linesNormalization(abs_lines, normalization_option, verbosity);
449 abs_linesMirroring(abs_lines, mirroring_option, verbosity);
450 abs_linesPopulation(abs_lines, population_option, verbosity);
451 abs_linesLineShapeType(abs_lines, lineshapetype_option, verbosity);
452 abs_linesCutoff(abs_lines, cutoff_option, cutoff_value, verbosity);
453 abs_linesLinemixingLimit(abs_lines, linemixinglimit_value, verbosity);
454
455 tag.read_from_stream(is_xml);
456 tag.check_name("/ArrayOfLineRecord");
457}
458
459/* Workspace method: Doxygen documentation will be auto-generated */
461 const ArrayOfArrayOfSpeciesTag& abs_species,
462 const String& basename,
463 const Numeric& fmin,
464 const Numeric& fmax,
465 const String& globalquantumnumbers,
466 const String& localquantumnumbers,
467 const Index& ignore_missing,
468 const String& normalization_option,
469 const String& mirroring_option,
470 const String& population_option,
471 const String& lineshapetype_option,
472 const String& cutoff_option,
473 const Numeric& cutoff_value,
474 const Numeric& linemixinglimit_value,
475 const Verbosity& verbosity)
476{
477 abs_lines.resize(0);
478
479 // Build a set of species indices. Duplicates are ignored.
480 const std::set<Species::Species> unique_species = lbl_species(abs_species);
481
482 String tmpbasename = basename;
483 if (basename.length() && basename[basename.length() - 1] != '/') {
484 tmpbasename += '.';
485 }
486
487 // Read catalogs for each identified species and put them all into
488 // abs_lines.
489 abs_lines.resize(0);
490 for (auto& spec: unique_species) {
491 ArrayOfAbsorptionLines more_abs_lines;
492
493 try {
494 ReadARTSCAT(more_abs_lines,
495 tmpbasename + String(Species::toShortName(spec)) + ".xml",
496 fmin,
497 fmax,
498 globalquantumnumbers,
499 localquantumnumbers,
500 normalization_option,
501 mirroring_option,
502 population_option,
503 lineshapetype_option,
504 cutoff_option,
505 cutoff_value,
506 linemixinglimit_value,
507 verbosity);
508
509 // Either find a line like this in the list of lines or start a new Lines
510 for (auto& newband: more_abs_lines) {
511 bool found = false;
512 for (auto& band: abs_lines) {
513 if (band.Match(newband).first) {
514 for (Index k=0; k<newband.NumLines(); k++) {
515 band.AppendSingleLine(newband.lines[k]);
516 found = true;
517 }
518 }
519 }
520 if (not found) {
521 abs_lines.push_back(newband);
522 }
523 }
524 } catch (const std::exception& e) {
525 ARTS_USER_ERROR_IF (not ignore_missing,
526 "Errors in calls by *propmat_clearskyAddZeeman*:\n",
527 e.what())
528 continue;
529 }
530 }
531
532 for (auto& band: abs_lines)
533 band.sort_by_frequency();
534}
535
536/* Workspace method: Doxygen documentation will be auto-generated */
538 const String& hitran_file,
539 const Numeric& fmin,
540 const Numeric& fmax,
541 const String& globalquantumnumbers,
542 const String& localquantumnumbers,
543 const String& hitran_type,
544 const String& normalization_option,
545 const String& mirroring_option,
546 const String& population_option,
547 const String& lineshapetype_option,
548 const String& cutoff_option,
549 const Numeric& cutoff_value,
550 const Numeric& linemixinglimit_value,
551 const Verbosity& verbosity)
552{
553 abs_lines.resize(0);
554
555 // Global numbers
556 const Array<QuantumNumberType> global_nums = string2vecqn(globalquantumnumbers);
557
558 // Local numbers
559 const Array<QuantumNumberType> local_nums = string2vecqn(localquantumnumbers);
560 ARTS_USER_ERROR_IF(not check_local(local_nums), "Can only have non-string values in the local state")
561
562 // HITRAN type
563 const Options::HitranType hitran_version = Options::toHitranTypeOrThrow(hitran_type);
564
565 // Hitran data
566 ifstream is;
567 open_input_file(is, hitran_file);
568
569 ArrayOfAbsorptionLines local_bands(0);
570 bool go_on = true;
571 while (go_on) {
572
574 switch (hitran_version) {
577 break;
580 break;
581 case Options::HitranType::Online:
583 break;
584 default:
585 ARTS_ASSERT(false, "The HitranType enum class has to be fully updated!\n");
586 }
587
588 if (sline.bad) {
589 if (is.eof())
590 break;
591 ARTS_USER_ERROR("Cannot read line ", nelem(abs_lines) + nelem(local_bands) + 1);
592 }
593 if (sline.line.F0 < fmin)
594 continue; // Skip this line
595 if (sline.line.F0 > fmax)
596 break; // We assume sorted so quit here
597
598 // Set Zeeman if implemented
600
601 // Get the global quantum number identifier
602 const QuantumIdentifier global_qid = global_quantumidentifier(global_nums, sline.quantumidentity);
603
604 // Get local quantum numbers into the line
605 for(auto qn: local_nums) {
606 if (sline.quantumidentity.val.has(qn)) {
607 sline.line.localquanta.val.set(sline.quantumidentity.val[qn]);
608 }
609 }
610
611 // Either find a line like this in the list of lines or start a new Lines
612 merge_external_line(local_bands, sline, global_qid);
613 if (local_bands.nelem() > merge_local_lines_size) {
614 merge_local_lines(abs_lines, local_bands);
615 local_bands.resize(0);
616 }
617 }
618
619 merge_local_lines(abs_lines, local_bands);
620
621 abs_linesNormalization(abs_lines, normalization_option, verbosity);
622 abs_linesMirroring(abs_lines, mirroring_option, verbosity);
623 abs_linesPopulation(abs_lines, population_option, verbosity);
624 abs_linesLineShapeType(abs_lines, lineshapetype_option, verbosity);
625 abs_linesCutoff(abs_lines, cutoff_option, cutoff_value, verbosity);
626 abs_linesLinemixingLimit(abs_lines, linemixinglimit_value, verbosity);
627}
628
629/* Workspace method: Doxygen documentation will be auto-generated */
631 const String& lblrtm_file,
632 const Numeric& fmin,
633 const Numeric& fmax,
634 const String& globalquantumnumbers,
635 const String& localquantumnumbers,
636 const String& normalization_option,
637 const String& mirroring_option,
638 const String& population_option,
639 const String& lineshapetype_option,
640 const String& cutoff_option,
641 const Numeric& cutoff_value,
642 const Numeric& linemixinglimit_value,
643 const Verbosity& verbosity)
644{
645 abs_lines.resize(0);
646
647 // Global numbers
648 const Array<QuantumNumberType> global_nums = string2vecqn(globalquantumnumbers);
649
650 // Local numbers
651 const Array<QuantumNumberType> local_nums = string2vecqn(localquantumnumbers);
652 ARTS_USER_ERROR_IF(not check_local(local_nums), "Can only have non-string values in the local state")
653
654 // LBLRTM data
655 ifstream is;
656 open_input_file(is, lblrtm_file);
657
658 std::vector<Absorption::SingleLineExternal> v(0);
659
660 bool go_on = true;
661 while (go_on) {
662 v.push_back(Absorption::ReadFromLBLRTMStream(is));
663
664 if (v.back().bad) {
665 v.pop_back();
666 go_on = false;
667 } else if (v.back().line.F0 < fmin) {
668 v.pop_back();
669 } else if (v.back().line.F0 > fmax) {
670 v.pop_back();
671 go_on = false;
672 }
673 }
674
675 for (auto& x: v)
676 x.line.zeeman = Zeeman::GetAdvancedModel(x.quantumidentity);
677
678 auto x = Absorption::split_list_of_external_lines(v, local_nums, global_nums);
679 abs_lines.resize(0);
680 abs_lines.reserve(x.size());
681 while (x.size()) {
682 abs_lines.push_back(x.back());
683 abs_lines.back().sort_by_frequency();
684 x.pop_back();
685 }
686
687 abs_linesNormalization(abs_lines, normalization_option, verbosity);
688 abs_linesMirroring(abs_lines, mirroring_option, verbosity);
689 abs_linesPopulation(abs_lines, population_option, verbosity);
690 abs_linesLineShapeType(abs_lines, lineshapetype_option, verbosity);
691 abs_linesCutoff(abs_lines, cutoff_option, cutoff_value, verbosity);
692 abs_linesLinemixingLimit(abs_lines, linemixinglimit_value, verbosity);
693}
694
695/* Workspace method: Doxygen documentation will be auto-generated */
697 const String& jpl_file,
698 const Numeric& fmin,
699 const Numeric& fmax,
700 const String& globalquantumnumbers,
701 const String& localquantumnumbers,
702 const String& normalization_option,
703 const String& mirroring_option,
704 const String& population_option,
705 const String& lineshapetype_option,
706 const String& cutoff_option,
707 const Numeric& cutoff_value,
708 const Numeric& linemixinglimit_value,
709 const Verbosity& verbosity)
710{
711 abs_lines.resize(0);
712
713 // Global numbers
714 const Array<QuantumNumberType> global_nums = string2vecqn(globalquantumnumbers);
715
716 // Local numbers
717 const Array<QuantumNumberType> local_nums = string2vecqn(localquantumnumbers);
718 ARTS_USER_ERROR_IF(not check_local(local_nums), "Can only have non-string values in the local state")
719
720 // LBLRTM data
721 ifstream is;
722 open_input_file(is, jpl_file);
723
724 std::vector<Absorption::SingleLineExternal> v(0);
725
726 bool go_on = true;
727 while (go_on) {
728 v.push_back(Absorption::ReadFromJplStream(is));
729
730 if (v.back().bad) {
731 v.pop_back();
732 go_on = false;
733 } else if (v.back().line.F0 < fmin) {
734 v.pop_back();
735 } else if (v.back().line.F0 > fmax) {
736 v.pop_back();
737 go_on = false;
738 }
739 }
740
741 for (auto& x: v)
742 x.line.zeeman = Zeeman::GetAdvancedModel(x.quantumidentity);
743
744 auto x = Absorption::split_list_of_external_lines(v, local_nums, global_nums);
745 abs_lines.resize(0);
746 abs_lines.reserve(x.size());
747 while (x.size()) {
748 abs_lines.push_back(x.back());
749 abs_lines.back().sort_by_frequency();
750 x.pop_back();
751 }
752
753 abs_linesNormalization(abs_lines, normalization_option, verbosity);
754 abs_linesMirroring(abs_lines, mirroring_option, verbosity);
755 abs_linesPopulation(abs_lines, population_option, verbosity);
756 abs_linesLineShapeType(abs_lines, lineshapetype_option, verbosity);
757 abs_linesCutoff(abs_lines, cutoff_option, cutoff_value, verbosity);
758 abs_linesLinemixingLimit(abs_lines, linemixinglimit_value, verbosity);
759}
760
764
765/* Workspace method: Doxygen documentation will be auto-generated */
767 const ArrayOfAbsorptionLines& abs_lines,
768 const String& basename,
769 const Verbosity& verbosity) {
770 // Set the true name of the saving
771 String true_basename = basename;
772 if (not(true_basename.back() == '.' or true_basename.back() == '/'))
773 true_basename += '.';
774
775 // Find all species
776 ArrayOfString specs(0);
777 for (auto& band: abs_lines) {
778 auto specname = band.SpeciesName();
779
780 bool any = false;
781 for (auto& thisname: specs) {
782 if (any) break;
783 if (thisname == specname) any = true;
784 }
785
786 if (not any)
787 specs.push_back(specname);
788 }
789
790 // Make all species into a species tag array
791 Index throwaway;
793 abs_speciesSet(as, throwaway, specs, verbosity);
794
795 // Split lines by species
797 abs_lines_per_speciesCreateFromLines(alps, abs_lines, as, verbosity);
798
799 // Save the arrays
800 for (Index i=0; i<specs.nelem(); i++) {
801 auto& name = specs[i];
802 auto& lines = alps[i];
803
804 WriteXML(output_format, lines,
805 true_basename + name + ".xml",
806 0, "", "", "", verbosity);
807 }
808}
809
810/* Workspace method: Doxygen documentation will be auto-generated */
812 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
813 const String& basename,
814 const Verbosity& verbosity) {
815 // Compact to abs_lines
816 ArrayOfAbsorptionLines abs_lines(0);
817 for (auto& lines: abs_lines_per_species) {
818 for (auto& band: lines) {
819 abs_lines.push_back(band);
820 }
821 }
822
823 // Save using the other function
824 abs_linesWriteSpeciesSplitCatalog(output_format, abs_lines, basename, verbosity);
825}
826
827/* Workspace method: Doxygen documentation will be auto-generated */
829 const String& basename,
830 const Index& robust,
831 const Verbosity& verbosity) {
832 abs_lines.resize(0);
833
835 std::size_t bands_found{0};
836
837 String tmpbasename = basename;
838 if (basename.length() && basename[basename.length() - 1] != '/') {
839 tmpbasename += '.';
840 }
841
842 // Read catalogs for each identified species and put them all into
843 // abs_lines
844 for (auto& ir: Species::Isotopologues) {
845 String filename = tmpbasename + ir.FullName() + ".xml";
846 if (find_xml_file_existence(filename)) {
847 ArrayOfAbsorptionLines speclines;
848 xml_read_from_file(filename, speclines, verbosity);
849 for (auto& band: speclines) {
850 abs_lines.push_back(band);
851 bands_found++;
852 }
853 }
854 }
855
856 ARTS_USER_ERROR_IF (not bands_found and not robust,
857 "Cannot find any bands in the directory you are reading");
858 out3 << "Found " << bands_found << " bands\n";
859}
860
861/* Workspace method: Doxygen documentation will be auto-generated */
863 const ArrayOfArrayOfSpeciesTag& abs_species,
864 const String& basename,
865 const Index& robust,
866 const Verbosity& verbosity) {
867 abs_lines_per_species.resize(0);
868
870
871 // Build a set of species indices. Duplicates are ignored.
872 const std::set<Species::Species> unique_species = lbl_species(abs_species);
873
874 String tmpbasename = basename;
875 if (basename.length() && basename[basename.length() - 1] != '/') {
876 tmpbasename += '.';
877 }
878
879 // Read catalogs for each identified species and put them all into
880 // abs_lines in aa parallel manner
882 ArrayOfArrayOfAbsorptionLines par_abs_lines(n);
883 ArrayOfString errors;
884 std::atomic_bool error{false};
885
886 #pragma omp parallel for schedule(dynamic) if (!arts_omp_in_parallel())
887 for (std::size_t ispec=0; ispec<unique_species.size(); ispec++) {
888
889 const auto i = arts_omp_get_thread_num();
890 const auto& spec = *std::next(unique_species.cbegin(), ispec);
891 const auto isots = Species::isotopologues(spec);
892
893 if (not error.load()) {
894 try {
895 for (const auto& isot: isots) {
896 String filename = tmpbasename + isot.FullName() + ".xml";
897 if (find_xml_file_existence(filename)) {
898 ArrayOfAbsorptionLines speclines;
899 xml_read_from_file(filename, speclines, verbosity);
900 for (auto& band: speclines) {
901 par_abs_lines[i].push_back(std::move(band));
902 }
903 }
904 }
905 } catch (std::exception& e) {
906 error.store(true);
907
908 #pragma omp critical
909 errors.push_back(var_string("\n\nError in thread ", i, " for species ", spec, ":\n", e.what(), "\n\n"));
910 }
911 }
912 }
913
914 ARTS_USER_ERROR_IF(error.load(), "Encountered ", errors.nelem(), ':', errors)
915
916 const Index bands_found = nelem(par_abs_lines);
917 ARTS_USER_ERROR_IF (not bands_found and not robust,
918 "Cannot find any bands in the directory you are reading");
919 out3 << "Found " << bands_found << " bands\n";
920
921 ArrayOfAbsorptionLines abs_lines;
922 abs_lines.reserve(bands_found);
923 for (auto& lines: par_abs_lines) {
924 for (auto& band: lines) {
925 abs_lines.push_back(std::move(band));
926 }
927 }
928
929 abs_lines_per_speciesCreateFromLines(abs_lines_per_species, abs_lines, abs_species, verbosity);
930}
931
935
936/* Workspace method: Doxygen documentation will be auto-generated */
938{
939 for (auto& rlines: replacing_lines) {
940 Index number_of_matching_bands = 0;
941 for (auto& tlines: abs_lines) {
942 if (tlines.Match(rlines).first) {
943 number_of_matching_bands++;
944 for (auto& rline: rlines.lines) {
945 Index number_of_matching_single_lines = 0;
946 for (auto& tline: tlines.lines) {
947 if (tline.localquanta.val == rline.localquanta.val) {
948 number_of_matching_single_lines++;
949 tline = rline;
950 }
951 }
952
953 ARTS_USER_ERROR_IF (number_of_matching_single_lines not_eq 1,
954 "Error: Did not match to a single single line. "
955 "This means the input data has not been understood. "
956 "This function needs exactly one match.");
957 }
958 tlines.sort_by_frequency();
959 }
960 }
961
962 ARTS_USER_ERROR_IF (number_of_matching_bands not_eq 1,
963 "Error: Did not match to a single set of absorption lines. "
964 "This means the input data has not been understood. "
965 "This function needs exactly one match.");
966 }
967}
968
969/* Workspace method: Doxygen documentation will be auto-generated */
970void abs_linesAppendWithLines(ArrayOfAbsorptionLines& abs_lines, const ArrayOfAbsorptionLines& appending_lines, const Index& safe, const Verbosity&)
971{
972 if (safe) {
973 std::vector<AbsorptionLines> addedlines(0);
974
975 for (auto& alines: appending_lines) {
976 Index number_of_matching_bands = 0;
977 for (auto& tlines: abs_lines) {
978 if (tlines.Match(alines).first) {
979 number_of_matching_bands++;
980 for (auto& aline: alines.lines) {
981 Index number_of_matching_single_lines = 0;
982 for (auto& tline: tlines.lines) {
983 if (tline.localquanta.val == aline.localquanta.val) {
984 number_of_matching_single_lines++;
985 }
986 }
987 ARTS_USER_ERROR_IF (number_of_matching_single_lines not_eq 0,
988 "Error: Did match to a single single line. "
989 "This means the input data has not been understood. "
990 "This function needs exactly zero matches.");
991 tlines.AppendSingleLine(aline);
992 }
993 tlines.sort_by_frequency();
994 }
995 }
996
997 if (number_of_matching_bands == 0) addedlines.push_back(alines);
998 ARTS_USER_ERROR_IF (number_of_matching_bands not_eq 1,
999 "Error: Did not match to a single set of absorption lines. "
1000 "This means the input data has not been understood. "
1001 "This function needs exactly one or zero matches.");
1002 }
1003
1004 for (auto& lines: addedlines) {
1005 abs_lines.push_back(std::move(lines));
1006 }
1007 } else {
1008 for (auto& band: appending_lines)
1009 abs_lines.push_back(band);
1010 }
1011}
1012
1013/* Workspace method: Doxygen documentation will be auto-generated */
1014void abs_linesDeleteBadF0(ArrayOfAbsorptionLines& abs_lines, const Numeric& f0, const Index& lower, const Verbosity&)
1015{
1016 for (auto& lines: abs_lines) {
1017 std::vector<Index> hits;
1018 for (Index i=0; i<lines.NumLines(); i++) {
1019 if (lower and lines.lines[i].F0 < f0)
1020 hits.push_back(i);
1021 else if (not lower and lines.lines[i].F0 > f0)
1022 hits.push_back(i);
1023 }
1024
1025 // Remove the bad values (sort by descending firs)
1026 std::sort(hits.begin(), hits.end());
1027 while(not hits.empty()) {
1028 lines.RemoveLine(hits.back());
1029 hits.pop_back();
1030 }
1031 }
1032}
1033
1034/* Workspace method: Doxygen documentation will be auto-generated */
1036{
1037 for (auto& dlines: deleting_lines) {
1038 for (auto& tlines: abs_lines) {
1039 std::vector<Index> hits(0);
1040
1041 if (tlines.Match(dlines).first) {
1042 for (auto& dline: dlines.lines) {
1043 for (Index i=0; i<tlines.NumLines(); i++) {
1044 if (tlines.lines[i].localquanta.val == dline.localquanta.val) {
1045 hits.push_back(i);
1046 }
1047 }
1048 }
1049
1050 // Sort and test the input
1051 std::sort(hits.begin(), hits.end());
1052 auto n = hits.size();
1053 hits.erase(std::unique(hits.begin(), hits.end()), hits.end());
1054 ARTS_USER_ERROR_IF(n not_eq hits.size(),
1055 "Removing the same line more than once is not accepted");
1056
1057 // Remove the bad values
1058 while(not hits.empty()) {
1059 tlines.RemoveLine(hits.back());
1060 hits.pop_back();
1061 }
1062 }
1063 }
1064 }
1065}
1066
1067/* Workspace method: Doxygen documentation will be auto-generated */
1069{
1070 for (auto& band: abs_lines) {
1071 std::array<bool, LineShape::nVars> var_is_empty;
1072
1073 // Species by species can be empty, so loop each species by themselves
1074 for (Index ispec=0; ispec<band.NumBroadeners(); ispec++) {
1075 var_is_empty.fill(true);
1076
1077 // Check if any variable in this band for any line is non-empty
1078 for (Index iline=0; iline<band.NumLines(); iline++) {
1079 for (Index ivar=0; ivar < LineShape::nVars; ivar++) {
1080 if (not LineShape::modelparameterEmpty(band.lines[iline].lineshape.Data()[ispec].Data()[ivar])) {
1081 var_is_empty[ivar] = false;
1082 }
1083 }
1084 }
1085
1086 // Remove empty variables from the writing. This will also speed up some calculations
1087 for (Index iline=0; iline<band.NumLines(); iline++) {
1088 for (Index ivar=0; ivar < LineShape::nVars; ivar++) {
1089 if (var_is_empty[ivar]) {
1090 band.lines[iline].lineshape.Data()[ispec].Data()[ivar].type = LineShape::TemperatureModel::None;
1091 }
1092 }
1093 }
1094 }
1095 }
1096}
1097
1099{
1100 for (auto& band: abs_lines) {
1101 const Quantum::Number::StateMatch lt(qid, band.quantumidentity);
1102 while (lt not_eq Quantum::Number::StateMatchType::Full and band.NumLines()) {
1103 band.RemoveLine(0);
1104 }
1105 }
1106}
1107
1109 const Verbosity&) {
1110 const Index nb = lines.nelem();
1111 for (Index i=0; i<nb; i++) {
1112 for (Index j=i+1; j<nb; j++) {
1114 Quantum::Number::StateMatch(lines[i].quantumidentity, lines[j].quantumidentity) ==
1115 Quantum::Number::StateMatchType::Full,
1116 "Not unique, these bands match:\n", lines[i], "\nand\n", lines[j])
1117 }
1118 }
1119}
1120
1122 const ArrayOfAbsorptionLines& replacing_bands,
1123 const Verbosity&) {
1124 for (auto& replacement: replacing_bands) {
1125 struct {Index band;} pos{-1};
1126
1127 for (Index i=0; i<abs_lines.nelem(); i++) {
1128 auto& band = abs_lines[i];
1129
1130 if (Quantum::Number::StateMatch(band.quantumidentity, replacement.quantumidentity) ==
1131 Quantum::Number::StateMatchType::Full) {
1132 ARTS_USER_ERROR_IF (pos.band not_eq -1, "Duplicate band matches for replacement line:\n",
1133 replacement, "\nThese are for band indexes ", pos.band, " and ", i)
1134
1135 pos.band = i;
1136 }
1137 }
1138
1139 ARTS_USER_ERROR_IF(pos.band == -1,
1140 "There is no match for replacement band:\n", replacement,
1141 "\nYou need to append the entire band")
1142 abs_lines[pos.band] = replacement;
1143 }
1144}
1145
1147 const ArrayOfAbsorptionLines& replacing_lines,
1148 const Verbosity&) {
1149 const Index nb = abs_lines.nelem(); // Evaluate first so new bands are not counted
1150
1151 for (auto& replacement: replacing_lines) {
1152 const Index nl=replacement.NumLines();
1153
1154 struct {Index band; ArrayOfIndex lines;} pos{-1, ArrayOfIndex(nl, -1)};
1155
1156 for (Index i=0; i<nb; i++) {
1157 auto& band = abs_lines[i];
1158
1159 if (Quantum::Number::StateMatch(band.quantumidentity, replacement.quantumidentity) ==
1160 Quantum::Number::StateMatchType::Full) {
1161 ARTS_USER_ERROR_IF (pos.band not_eq -1, "Duplicate band matches for replacement line:\n",
1162 replacement, "\nThese are for band indexes ", pos.band, " and ", i, '\n',
1163 "ID ", i, ": ", band.quantumidentity, '\n',
1164 "ID ", pos.band, ": ", abs_lines[pos.band].quantumidentity, '\n',
1165 "ID target: ", replacement.quantumidentity, '\n')
1166
1167 pos.band = i;
1168
1169 for (Index k=0; k<band.NumLines(); k++) {
1170 for (Index j=0; j<nl; j++) {
1171 if (replacement.lines[j].localquanta.val == band.lines[k].localquanta.val) {
1172 // We cannot have multiple entries
1173 ARTS_USER_ERROR_IF(pos.lines[j] not_eq -1 or std::any_of(pos.lines.begin(), pos.lines.end(), [k](auto& a){return a == k;}),
1174 "Found multiple matches of lines in:\n", replacement, "\n\nin mathcing band:\n", band)
1175
1176 pos.lines[j] = k;
1177 }
1178 }
1179 }
1180 }
1181 }
1182
1183 ARTS_USER_ERROR_IF(pos.band == -1 or std::any_of(pos.lines.begin(), pos.lines.end(), [](auto& a){return a == -1;}),
1184 "There is no match for replacement line:\n", replacement,
1185 "\nYou need to append the entire band")
1186
1187 // Add or change the line catalog
1188 auto& band = abs_lines[pos.band];
1189 if (const auto [match, nullable] = band.Match(replacement); nullable) {
1190 for (Index j=0; j<nl; j++) band.lines[pos.lines[j]] = replacement.lines[j];
1191 band.MakeLineShapeModelCommon();
1192 } else {
1193 // Sort to remove from behind
1194 std::sort(pos.lines.begin(), pos.lines.end(), std::greater<>());
1195 for (auto& k: pos.lines) band.RemoveLine(k);
1196 abs_lines.push_back(replacement);
1197 }
1198 }
1199}
1200
1204
1208
1209/* Workspace method: Doxygen documentation will be auto-generated */
1211 const String& type,
1212 const Numeric& x,
1213 const Verbosity&) {
1214 auto t = Absorption::toCutoffTypeOrThrow(type);
1215 for (auto& lines : abs_lines) {
1216 lines.cutoff = t;
1217 lines.cutofffreq = x;
1218 }
1219}
1220
1221/* Workspace method: Doxygen documentation will be auto-generated */
1223 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1224 const String& type,
1225 const Numeric& x,
1226 const Verbosity& v) {
1227 for (auto& abs_lines : abs_lines_per_species)
1228 abs_linesCutoff(abs_lines, type, x, v);
1229}
1230
1231/* Workspace method: Doxygen documentation will be auto-generated */
1233 const String& type,
1234 const Numeric& x,
1235 const QuantumIdentifier& QI,
1236 const Verbosity&) {
1237 auto t = Absorption::toCutoffTypeOrThrow(type);
1238 for (auto& band : abs_lines) {
1239 const Quantum::Number::StateMatch lt(QI, band.quantumidentity);
1240 if (lt == Quantum::Number::StateMatchType::Full) {
1241 band.cutoff = t;
1242 band.cutofffreq = x;
1243 }
1244 }
1245}
1246
1247/* Workspace method: Doxygen documentation will be auto-generated */
1249 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1250 const String& type,
1251 const Numeric& x,
1252 const QuantumIdentifier& QI,
1253 const Verbosity& verbosity) {
1254 for (auto& lines : abs_lines_per_species) {
1255 abs_linesCutoffMatch(lines, type, x, QI, verbosity);
1256 }
1257}
1258
1259/* Workspace method: Doxygen documentation will be auto-generated */
1261 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1262 const ArrayOfArrayOfSpeciesTag& abs_species,
1263 const String& type,
1264 const Numeric& x,
1265 const String& species_tag,
1266 const Verbosity& verbosity) {
1267 Index t1;
1268 ArrayOfArrayOfSpeciesTag target_species;
1269 abs_speciesSet(target_species, t1, {species_tag}, verbosity);
1270 for (Index ispec = 0; ispec < abs_species.nelem(); ispec++) {
1271 if (std::equal(abs_species[ispec].begin(),
1272 abs_species[ispec].end(),
1273 target_species[0].begin())) {
1274 abs_linesCutoff(abs_lines_per_species[ispec], type, x, verbosity);
1275 }
1276 }
1277}
1278
1282
1283/* Workspace method: Doxygen documentation will be auto-generated */
1285 const String& type,
1286 const Verbosity&) {
1287 auto t = Absorption::toMirroringTypeOrThrow(type);
1288 for (auto& lines : abs_lines) lines.mirroring = t;
1289}
1290
1291/* Workspace method: Doxygen documentation will be auto-generated */
1293 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1294 const String& type,
1295 const Verbosity& v) {
1296 for (auto& abs_lines : abs_lines_per_species)
1297 abs_linesMirroring(abs_lines, type, v);
1298}
1299
1300/* Workspace method: Doxygen documentation will be auto-generated */
1302 const String& type,
1303 const QuantumIdentifier& QI,
1304 const Verbosity&) {
1305 auto t = Absorption::toMirroringTypeOrThrow(type);
1306 for (auto& band : abs_lines) {
1307 const Quantum::Number::StateMatch lt(QI, band.quantumidentity);
1308 if (lt == Quantum::Number::StateMatchType::Full) {
1309 band.mirroring = t;
1310 }
1311 }
1312}
1313
1314/* Workspace method: Doxygen documentation will be auto-generated */
1316 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1317 const String& type,
1318 const QuantumIdentifier& QI,
1319 const Verbosity& v) {
1320 for (auto& abs_lines : abs_lines_per_species)
1321 abs_linesMirroringMatch(abs_lines, type, QI, v);
1322}
1323
1324/* Workspace method: Doxygen documentation will be auto-generated */
1326 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1327 const ArrayOfArrayOfSpeciesTag& abs_species,
1328 const String& type,
1329 const String& species_tag,
1330 const Verbosity& verbosity) {
1331 Index t1;
1332 ArrayOfArrayOfSpeciesTag target_species;
1333 abs_speciesSet(target_species, t1, {species_tag}, verbosity);
1334 for (Index ispec = 0; ispec < abs_species.nelem(); ispec++) {
1335 if (std::equal(abs_species[ispec].begin(),
1336 abs_species[ispec].end(),
1337 target_species[0].begin())) {
1338 abs_linesMirroring(abs_lines_per_species[ispec], type, verbosity);
1339 }
1340 }
1341}
1342
1343/* Workspace method: Doxygen documentation will be auto-generated */
1345 const Verbosity&) {
1346 const ArrayOfAbsorptionLines abs_lines_copy = abs_lines;
1347 for (AbsorptionLines band : abs_lines_copy) {
1348 band.mirroring = Absorption::MirroringType::Manual;
1349
1352 std::find_if(abs_lines_copy.cbegin(),
1353 abs_lines_copy.cend(),
1354 [&band](const AbsorptionLines& li) {
1355 return band.Match(li).first;
1356 }) not_eq abs_lines_copy.cend(),
1357 "Dual bands with same setup is not allowed for mirroring of band:\n",
1358 band,
1359 '\n')
1360
1361 for (auto& line : band.lines) {
1362 line.F0 *= -1;
1363 }
1364
1365 abs_lines.emplace_back(std::move(band));
1366 }
1367}
1368
1369/* Workspace method: Doxygen documentation will be auto-generated */
1371 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1372 const Verbosity& verbosity) {
1373 for (auto& abs_lines : abs_lines_per_species)
1374 abs_linesManualMirroring(abs_lines, verbosity);
1375}
1376
1377/* Workspace method: Doxygen documentation will be auto-generated */
1379 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1380 const ArrayOfArrayOfSpeciesTag& abs_species,
1381 const ArrayOfSpeciesTag& species,
1382 const Verbosity& verbosity) {
1383 ARTS_USER_ERROR_IF(abs_species.size() not_eq abs_lines_per_species.size(),
1384 "Mismatch abs_species and abs_lines_per_species sizes [",
1385 abs_species.size(),
1386 " vs ",
1387 abs_lines_per_species.size(),
1388 ", respectively]")
1389
1390 if (auto ind = std::distance(
1391 abs_species.cbegin(),
1392 std::find(abs_species.cbegin(), abs_species.cend(), species));
1393 ind not_eq abs_species.nelem()) {
1394 abs_linesManualMirroring(abs_lines_per_species[ind], verbosity);
1395 } else {
1396 ARTS_USER_ERROR("Cannot find species: ",
1397 species,
1398 "\nIn abs_species: [",
1399 abs_species,
1400 ']')
1401 }
1402}
1403
1407
1408/* Workspace method: Doxygen documentation will be auto-generated */
1410 const String& type,
1411 const Verbosity&) {
1412 auto t = Absorption::toPopulationTypeOrThrow(type);
1413 for (auto& lines : abs_lines) lines.population = t;
1414}
1415
1416/* Workspace method: Doxygen documentation will be auto-generated */
1418 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1419 const String& type,
1420 const Verbosity& v) {
1421 for (auto& abs_lines : abs_lines_per_species)
1422 abs_linesPopulation(abs_lines, type, v);
1423}
1424
1425/* Workspace method: Doxygen documentation will be auto-generated */
1427 const String& type,
1428 const QuantumIdentifier& QI,
1429 const Verbosity&) {
1430 auto t = Absorption::toPopulationTypeOrThrow(type);
1431 for (auto& lines : abs_lines) {
1432 const Quantum::Number::StateMatch lt(QI, lines.quantumidentity);
1433 if (lt == Quantum::Number::StateMatchType::Full) {
1434 lines.population = t;
1435 }
1436 }
1437}
1438
1439/* Workspace method: Doxygen documentation will be auto-generated */
1441 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1442 const String& type,
1443 const QuantumIdentifier& QI,
1444 const Verbosity& v) {
1445 for (auto& abs_lines : abs_lines_per_species)
1446 abs_linesPopulationMatch(abs_lines, type, QI, v);
1447}
1448
1449/* Workspace method: Doxygen documentation will be auto-generated */
1451 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1452 const ArrayOfArrayOfSpeciesTag& abs_species,
1453 const String& type,
1454 const String& species_tag,
1455 const Verbosity& verbosity) {
1456 Index t1;
1457 ArrayOfArrayOfSpeciesTag target_species;
1458 abs_speciesSet(target_species, t1, {species_tag}, verbosity);
1459 for (Index ispec = 0; ispec < abs_species.nelem(); ispec++) {
1460 if (std::equal(abs_species[ispec].begin(),
1461 abs_species[ispec].end(),
1462 target_species[0].begin())) {
1463 abs_linesPopulation(abs_lines_per_species[ispec], type, verbosity);
1464 }
1465 }
1466}
1467
1471
1472/* Workspace method: Doxygen documentation will be auto-generated */
1474 const String& type,
1475 const Verbosity&) {
1476 auto t = Absorption::toNormalizationTypeOrThrow(type);
1477 for (auto& lines : abs_lines) lines.normalization = t;
1478}
1479
1480/* Workspace method: Doxygen documentation will be auto-generated */
1482 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1483 const String& type,
1484 const Verbosity& v) {
1485 for (auto& abs_lines : abs_lines_per_species)
1486 abs_linesNormalization(abs_lines, type, v);
1487}
1488
1489/* Workspace method: Doxygen documentation will be auto-generated */
1491 const String& type,
1492 const QuantumIdentifier& QI,
1493 const Verbosity&) {
1494 auto t = Absorption::toNormalizationTypeOrThrow(type);
1495 for (auto& lines : abs_lines) {
1496 const Quantum::Number::StateMatch lt(QI, lines.quantumidentity);
1497 if (lt == Quantum::Number::StateMatchType::Full) {
1498 lines.normalization = t;
1499 }
1500 }
1501}
1502
1503/* Workspace method: Doxygen documentation will be auto-generated */
1505 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1506 const String& type,
1507 const QuantumIdentifier& QI,
1508 const Verbosity& v) {
1509 for (auto& abs_lines : abs_lines_per_species)
1510 abs_linesNormalizationMatch(abs_lines, type, QI, v);
1511}
1512
1513/* Workspace method: Doxygen documentation will be auto-generated */
1515 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1516 const ArrayOfArrayOfSpeciesTag& abs_species,
1517 const String& type,
1518 const String& species_tag,
1519 const Verbosity& verbosity) {
1520 Index t1;
1521 ArrayOfArrayOfSpeciesTag target_species;
1522 abs_speciesSet(target_species, t1, {species_tag}, verbosity);
1523 for (Index ispec = 0; ispec < abs_species.nelem(); ispec++) {
1524 if (std::equal(abs_species[ispec].begin(),
1525 abs_species[ispec].end(),
1526 target_species[0].begin())) {
1527 abs_linesNormalization(abs_lines_per_species[ispec], type, verbosity);
1528 }
1529 }
1530}
1531
1535
1536/* Workspace method: Doxygen documentation will be auto-generated */
1538 const String& type,
1539 const Verbosity&) {
1540 auto t = LineShape::toType(type);
1541 check_enum_error(t, "Cannot understand type: ", type);
1542 for (auto& lines : abs_lines) lines.lineshapetype = t;
1543}
1544
1545/* Workspace method: Doxygen documentation will be auto-generated */
1547 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1548 const String& type,
1549 const Verbosity& v) {
1550 for (auto& abs_lines : abs_lines_per_species)
1551 abs_linesLineShapeType(abs_lines, type, v);
1552}
1553
1554/* Workspace method: Doxygen documentation will be auto-generated */
1556 const String& type,
1557 const QuantumIdentifier& QI,
1558 const Verbosity&) {
1559 auto t = LineShape::toTypeOrThrow(type);
1560 for (auto& lines : abs_lines) {
1561 const Quantum::Number::StateMatch lt(QI, lines.quantumidentity);
1562 if (lt == Quantum::Number::StateMatchType::Full) {
1563 lines.lineshapetype = t;
1564 }
1565 }
1566}
1567
1568/* Workspace method: Doxygen documentation will be auto-generated */
1570 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1571 const String& type,
1572 const QuantumIdentifier& QI,
1573 const Verbosity& v) {
1574 for (auto& abs_lines : abs_lines_per_species)
1575 abs_linesLineShapeTypeMatch(abs_lines, type, QI, v);
1576}
1577
1578/* Workspace method: Doxygen documentation will be auto-generated */
1580 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1581 const ArrayOfArrayOfSpeciesTag& abs_species,
1582 const String& type,
1583 const String& species_tag,
1584 const Verbosity& verbosity) {
1585 Index t1;
1586 ArrayOfArrayOfSpeciesTag target_species;
1587 abs_speciesSet(target_species, t1, {species_tag}, verbosity);
1588 for (Index ispec = 0; ispec < abs_species.nelem(); ispec++) {
1589 if (std::equal(abs_species[ispec].begin(),
1590 abs_species[ispec].end(),
1591 target_species[0].begin())) {
1592 abs_linesLineShapeType(abs_lines_per_species[ispec], type, verbosity);
1593 }
1594 }
1595}
1596
1600
1601/* Workspace method: Doxygen documentation will be auto-generated */
1603 const Numeric& x,
1604 const Verbosity&) {
1605 for (auto& lines : abs_lines) lines.linemixinglimit = x;
1606}
1607
1608/* Workspace method: Doxygen documentation will be auto-generated */
1610 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1611 const Numeric& x,
1612 const Verbosity& v) {
1613 for (auto& abs_lines : abs_lines_per_species)
1614 abs_linesLinemixingLimit(abs_lines, x, v);
1615}
1616
1617/* Workspace method: Doxygen documentation will be auto-generated */
1619 const Numeric& x,
1620 const QuantumIdentifier& QI,
1621 const Verbosity&) {
1622 for (auto& lines : abs_lines) {
1623 const Quantum::Number::StateMatch lt(QI, lines.quantumidentity);
1624 if (lt == Quantum::Number::StateMatchType::Full) {
1625 lines.linemixinglimit = x;
1626 }
1627 }
1628}
1629
1630/* Workspace method: Doxygen documentation will be auto-generated */
1632 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1633 const Numeric& x,
1634 const QuantumIdentifier& QI,
1635 const Verbosity& v) {
1636 for (auto& abs_lines : abs_lines_per_species)
1637 abs_linesLinemixingLimitMatch(abs_lines, x, QI, v);
1638}
1639
1640/* Workspace method: Doxygen documentation will be auto-generated */
1642 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1643 const ArrayOfArrayOfSpeciesTag& abs_species,
1644 const Numeric& x,
1645 const String& species_tag,
1646 const Verbosity& verbosity) {
1647 Index t1;
1648 ArrayOfArrayOfSpeciesTag target_species;
1649 abs_speciesSet(target_species, t1, {species_tag}, verbosity);
1650 for (Index ispec = 0; ispec < abs_species.nelem(); ispec++) {
1651 if (std::equal(abs_species[ispec].begin(),
1652 abs_species[ispec].end(),
1653 target_species[0].begin())) {
1654 abs_linesLinemixingLimit(abs_lines_per_species[ispec], x, verbosity);
1655 }
1656 }
1657}
1658
1662
1663/* Workspace method: Doxygen documentation will be auto-generated */
1665 const Numeric& x,
1666 const Verbosity&) {
1667 for (auto& lines : abs_lines) lines.T0 = x;
1668}
1669
1670/* Workspace method: Doxygen documentation will be auto-generated */
1672 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1673 const Numeric& x,
1674 const Verbosity& v) {
1675 for (auto& abs_lines : abs_lines_per_species) abs_linesT0(abs_lines, x, v);
1676}
1677
1678/* Workspace method: Doxygen documentation will be auto-generated */
1680 const Numeric& x,
1681 const QuantumIdentifier& QI,
1682 const Verbosity&) {
1683 for (auto& lines : abs_lines) {
1684 const Quantum::Number::StateMatch lt(QI, lines.quantumidentity);
1685 if (lt == Quantum::Number::StateMatchType::Full) {
1686 lines.T0 = x;
1687 }
1688 }
1689}
1690
1691/* Workspace method: Doxygen documentation will be auto-generated */
1693 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1694 const Numeric& x,
1695 const QuantumIdentifier& QI,
1696 const Verbosity& v) {
1697 for (auto& abs_lines : abs_lines_per_species)
1698 abs_linesT0Match(abs_lines, x, QI, v);
1699}
1700
1701/* Workspace method: Doxygen documentation will be auto-generated */
1703 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1704 const ArrayOfArrayOfSpeciesTag& abs_species,
1705 const Numeric& x,
1706 const String& species_tag,
1707 const Verbosity& verbosity) {
1708 Index t1;
1709 ArrayOfArrayOfSpeciesTag target_species;
1710 abs_speciesSet(target_species, t1, {species_tag}, verbosity);
1711 for (Index ispec = 0; ispec < abs_species.nelem(); ispec++) {
1712 if (std::equal(abs_species[ispec].begin(),
1713 abs_species[ispec].end(),
1714 target_species[0].begin())) {
1715 abs_linesT0(abs_lines_per_species[ispec], x, verbosity);
1716 }
1717 }
1718}
1719
1723
1724/* Workspace method: Doxygen documentation will be auto-generated */
1726 const QuantumIdentifier& QI,
1727 const String& parameter_name,
1728 const Numeric& change,
1729 const Index& relative,
1730 const Verbosity&)
1731{
1732 Index parameter_switch = -1;
1733
1734 ARTS_USER_ERROR_IF (parameter_name.nelem() == 0,
1735 "parameter_name is empty.\n");
1736
1737 if (parameter_name == "Central Frequency" or
1738 parameter_name == "Line Center")
1739 parameter_switch = 0;
1740 else if (parameter_name == "Line Strength")
1741 parameter_switch = 1;
1742 else if (parameter_name == "Lower State Energy")
1743 parameter_switch = 4;
1744 else if (parameter_name == "Einstein Coefficient")
1745 parameter_switch = 5;
1746 else if (parameter_name == "Lower Statistical Weight")
1747 parameter_switch = 6;
1748 else if (parameter_name == "Upper Statistical Weight")
1749 parameter_switch = 7;
1750 else if (parameter_name == "Lower Zeeman Coefficient")
1751 parameter_switch = 8;
1752 else if (parameter_name == "Upper Zeeman Coefficient")
1753 parameter_switch = 9;
1754
1755 for (auto& band: abs_lines) {
1756 for (Index k=0; k<band.NumLines(); k++) {
1757 const Quantum::Number::StateMatch lt(QI, band.lines[k].localquanta, band.quantumidentity);
1758 if (lt == Quantum::Number::StateMatchType::Full) {
1759 switch (parameter_switch) {
1760 case 0: // "Central Frequency":
1761 if (relative == 0)
1762 band.lines[k].F0 += change;
1763 else
1764 band.lines[k].F0 *= 1.0e0 + change;
1765 break;
1766 case 1: // "Line Strength":
1767 if (relative == 0)
1768 band.lines[k].I0 += change;
1769 else
1770 band.lines[k].I0 *= 1.0e0 + change;
1771 break;
1772 case 4: // "Lower State Energy":
1773 if (relative == 0)
1774 band.lines[k].E0 += change;
1775 else
1776 band.lines[k].E0 *= 1.0e0 + change;
1777 break;
1778 case 5: // "Einstein":
1779 if (relative == 0)
1780 band.lines[k].A += change;
1781 else
1782 band.lines[k].A *= 1.0e0 + change;
1783 break;
1784 case 6: // "Lower Statistical Weight":
1785 if (relative == 0)
1786 band.lines[k].glow += change;
1787 else
1788 band.lines[k].glow *= 1.0e0 + change;
1789 break;
1790 case 7: // "Upper Statistical Weight":
1791 if (relative == 0)
1792 band.lines[k].gupp += change;
1793 else
1794 band.lines[k].gupp *= 1.0e0 + change;
1795 break;
1796 case 8: // "Lower Zeeman Coefficient":
1797 if (relative == 0)
1798 band.lines[k].zeeman.gl() += change;
1799 else
1800 band.lines[k].zeeman.gl() *= 1.0e0 + change;
1801 break;
1802 case 9: // "Upper Zeeman Coefficient":
1803 if (relative == 0)
1804 band.lines[k].zeeman.gu() += change;
1805 else
1806 band.lines[k].zeeman.gu() *= 1.0e0 + change;
1807 break;
1808 default: {
1810 "Usupported paramter_name\n", parameter_name,
1811 "\nSee method description for supported parameter names.\n")
1812 }
1813 }
1814 }
1815 }
1816 }
1817}
1818
1819/* Workspace method: Doxygen documentation will be auto-generated */
1821 const QuantumIdentifier& QI,
1822 const String& parameter_name,
1823 const Numeric& change,
1824 const Index& relative,
1825 const Verbosity& verbosity)
1826{
1827 for (auto& lines: abs_lines_per_species)
1828 abs_linesChangeBaseParameterForMatchingLines(lines, QI, parameter_name, change, relative, verbosity);
1829}
1830
1831/* Workspace method: Doxygen documentation will be auto-generated */
1833 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1834 const ArrayOfArrayOfSpeciesTag& abs_species,
1835 const QuantumIdentifier& QI,
1836 const String& parameter_name,
1837 const Numeric& change,
1838 const Index& relative,
1839 const String& species_tag,
1840 const Verbosity& verbosity)
1841{
1842 Index t1;
1843 ArrayOfArrayOfSpeciesTag target_species;
1844 abs_speciesSet(target_species, t1, {species_tag}, verbosity);
1845 for (Index ispec=0; ispec<abs_species.nelem(); ispec++) {
1846 if (std::equal(abs_species[ispec].begin(), abs_species[ispec].end(), target_species[0].begin())) {
1847 abs_linesChangeBaseParameterForMatchingLines(abs_lines_per_species[ispec], QI, parameter_name, change, relative, verbosity);
1848 }
1849 }
1850}
1851
1852/* Workspace method: Doxygen documentation will be auto-generated */
1854 const QuantumIdentifier& QI,
1855 const String& parameter_name,
1856 const Numeric& x,
1857 const Verbosity&) {
1858 Index parameter_switch = -1;
1859
1860 ARTS_USER_ERROR_IF(parameter_name.nelem() == 0, "parameter_name is empty.\n");
1861
1862 if (parameter_name == "Central Frequency" or parameter_name == "Line Center")
1863 parameter_switch = 0;
1864 else if (parameter_name == "Line Strength")
1865 parameter_switch = 1;
1866 else if (parameter_name == "Lower State Energy")
1867 parameter_switch = 4;
1868 else if (parameter_name == "Einstein Coefficient")
1869 parameter_switch = 5;
1870 else if (parameter_name == "Lower Statistical Weight")
1871 parameter_switch = 6;
1872 else if (parameter_name == "Upper Statistical Weight")
1873 parameter_switch = 7;
1874 else if (parameter_name == "Lower Zeeman Coefficient")
1875 parameter_switch = 8;
1876 else if (parameter_name == "Upper Zeeman Coefficient")
1877 parameter_switch = 9;
1878
1879 for (auto& band : abs_lines) {
1880 for (Index k = 0; k < band.NumLines(); k++) {
1882 QI, band.lines[k].localquanta, band.quantumidentity);
1883 if (lt == Quantum::Number::StateMatchType::Full) {
1884 switch (parameter_switch) {
1885 case 0: // "Central Frequency":
1886 band.lines[k].F0 = x;
1887 break;
1888 case 1: // "Line Strength":
1889 band.lines[k].I0 = x;
1890 break;
1891 case 4: // "Lower State Energy":
1892 band.lines[k].E0 = x;
1893 break;
1894 case 5: // "Einstein":
1895 band.lines[k].A = x;
1896 break;
1897 case 6: // "Lower Statistical Weight":
1898 band.lines[k].glow = x;
1899 break;
1900 case 7: // "Upper Statistical Weight":
1901 band.lines[k].gupp = x;
1902 break;
1903 case 8:
1904 band.lines[k].zeeman.gl() = x;
1905 break;
1906 case 9:
1907 band.lines[k].zeeman.gu() = x;
1908 break;
1909 default: {
1911 "Usupported paramter_name\n",
1912 parameter_name,
1913 "\nSee method description for supported parameter names.\n")
1914 }
1915 }
1916 }
1917 }
1918 }
1919}
1920
1921/* Workspace method: Doxygen documentation will be auto-generated */
1923 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1924 const QuantumIdentifier& QI,
1925 const String& parameter_name,
1926 const Numeric& change,
1927 const Verbosity& verbosity) {
1928 for (auto& lines : abs_lines_per_species)
1930 lines, QI, parameter_name, change, verbosity);
1931}
1932
1933/* Workspace method: Doxygen documentation will be auto-generated */
1935 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1936 const ArrayOfArrayOfSpeciesTag& abs_species,
1937 const QuantumIdentifier& QI,
1938 const String& parameter_name,
1939 const Numeric& change,
1940 const String& species_tag,
1941 const Verbosity& verbosity) {
1942 Index t1;
1943 ArrayOfArrayOfSpeciesTag target_species;
1944 abs_speciesSet(target_species, t1, {species_tag}, verbosity);
1945 for (Index ispec = 0; ispec < abs_species.nelem(); ispec++) {
1946 if (std::equal(abs_species[ispec].begin(),
1947 abs_species[ispec].end(),
1948 target_species[0].begin())) {
1950 abs_lines_per_species[ispec], QI, parameter_name, change, verbosity);
1951 }
1952 }
1953}
1954
1955/* Workspace method: Doxygen documentation will be auto-generated */
1957 ArrayOfAbsorptionLines& abs_lines,
1958 const QuantumIdentifier& QI,
1959 const String& parameter,
1960 const String& species,
1961 const String& temperaturemodel,
1962 const Vector& new_values,
1963 const Verbosity& verbosity) {
1965
1966 const bool do_self = species == LineShape::self_broadening;
1967 const bool do_bath = species == LineShape::bath_broadening;
1968
1969 // Set the spec index if possible
1970 const Species::Species spec = do_self ? Species::Species::FINAL
1971 : do_bath ? Species::Species::Bath
1972 : SpeciesTag(species).Spec();
1973
1974 const LineShape::Variable var = LineShape::toVariableOrThrow(parameter);
1975
1977 "Mismatch between input and expected number of variables\n"
1978 "\tInput is: ",
1979 new_values.nelem(),
1980 " long but expects: ",
1982 " values\n")
1983
1985 newdata.type = LineShape::toTemperatureModelOrThrow(temperaturemodel);
1986 newdata.X0 = new_values[0];
1987 newdata.X1 = new_values[1];
1988 newdata.X2 = new_values[2];
1989 newdata.X3 = new_values[3];
1990
1991 for (auto& band : abs_lines) {
1992 for (Index k = 0; k < band.NumLines(); k++) {
1994 QI, band.lines[k].localquanta, band.quantumidentity);
1995 if (lt == Quantum::Number::StateMatchType::Full) {
1996 if (do_self and band.selfbroadening) {
1997 out3 << "Changing self\n";
1998 band.lines[k].lineshape.Data().front().Data()[Index(var)] = newdata;
1999 } else if (do_bath and band.bathbroadening) {
2000 out3 << "Changing bath\n";
2001 band.lines[k].lineshape.Data().back().Data()[Index(var)] = newdata;
2002 } else {
2003 for (Index i = band.selfbroadening;
2004 i < band.broadeningspecies.nelem() - band.bathbroadening;
2005 i++) {
2006 if (spec == band.broadeningspecies[i]) {
2007 out3 << "Changing species: " << Species::toShortName(spec)
2008 << '\n';
2009 band.lines[k].lineshape.Data()[i].Data()[Index(var)] = newdata;
2010 }
2011 }
2012 }
2013 }
2014 }
2015 }
2016}
2017
2018/* Workspace method: Doxygen documentation will be auto-generated */
2020 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
2021 const QuantumIdentifier& QI,
2022 const String& parameter,
2023 const String& species,
2024 const String& temperaturemodel,
2025 const Vector& new_values,
2026 const Verbosity& verbosity) {
2027 for (auto& lines : abs_lines_per_species)
2029 lines, QI, parameter, species, temperaturemodel, new_values, verbosity);
2030}
2031
2035
2036/* Workspace method: Doxygen documentation will be auto-generated */
2038 const QuantumIdentifier& QI,
2039 const String& parameter_name,
2040 const Numeric& change,
2041 const Index& relative,
2042 const Verbosity&)
2043{
2044 Index parameter_switch = -1;
2045
2046 ARTS_USER_ERROR_IF (parameter_name.nelem() == 0,
2047 "parameter_name is empty.\n");
2048 if (parameter_name == "Statistical Weight")
2049 parameter_switch = 1;
2050 else if (parameter_name == "Zeeman Coefficient")
2051 parameter_switch = 2;
2052
2053 for (auto& band: abs_lines) {
2054 for (Index k=0; k<band.NumLines(); k++) {
2055 const Quantum::Number::StateMatch lt(QI, band.lines[k].localquanta, band.quantumidentity);
2056 if (lt == Quantum::Number::StateMatchType::Level and lt.low) {
2057 switch (parameter_switch) {
2058 case 1: // "Statistical Weight":
2059 if (relative == 0)
2060 band.lines[k].glow += change;
2061 else
2062 band.lines[k].glow *= 1.0e0 + change;
2063 break;
2064 case 2: // "Zeeman Coefficient":
2065 if (relative == 0)
2066 band.lines[k].zeeman.gl() += change;
2067 else
2068 band.lines[k].zeeman.gl() *= 1.0e0 + change;
2069 break;
2070 default: {
2072 "Usupported paramter_name\n", parameter_name,
2073 "\nSee method description for supported parameter names.\n")
2074 }
2075 }
2076 }
2077
2078 if (lt == Quantum::Number::StateMatchType::Level and lt.upp) {
2079 switch (parameter_switch) {
2080 case 1: // "Statistical Weight":
2081 if (relative == 0)
2082 band.lines[k].gupp += change;
2083 else
2084 band.lines[k].gupp *= 1.0e0 + change;
2085 break;
2086 case 2: // "Zeeman Coefficient":
2087 if (relative == 0)
2088 band.lines[k].zeeman.gu() += change;
2089 else
2090 band.lines[k].zeeman.gu() *= 1.0e0 + change;
2091 break;
2092 default: {
2094 "Usupported paramter_name\n", parameter_name,
2095 "\nSee method description for supported parameter names.\n")
2096 }
2097 }
2098 }
2099 }
2100 }
2101}
2102
2103/* Workspace method: Doxygen documentation will be auto-generated */
2105 const QuantumIdentifier& QI,
2106 const String& parameter_name,
2107 const Numeric& change,
2108 const Index& relative,
2109 const Verbosity& verbosity)
2110{
2111 for (auto& lines: abs_lines_per_species)
2112 abs_linesChangeBaseParameterForMatchingLevel(lines, QI, parameter_name, change, relative, verbosity);
2113}
2114
2115/* Workspace method: Doxygen documentation will be auto-generated */
2117 const ArrayOfQuantumIdentifier& QID,
2118 const String& parameter_name,
2119 const Vector& change,
2120 const Index& relative,
2121 const Verbosity& verbosity)
2122{
2123 ARTS_USER_ERROR_IF (QID.nelem() not_eq change.nelem(),
2124 "Mismatch between QID and change input lengths not allowed");
2125
2126 for (Index iq=0; iq<QID.nelem(); iq++)
2127 abs_linesChangeBaseParameterForMatchingLevel(abs_lines, QID[iq], parameter_name, change[iq], relative, verbosity);
2128}
2129
2130/* Workspace method: Doxygen documentation will be auto-generated */
2132 const ArrayOfQuantumIdentifier& QID,
2133 const String& parameter_name,
2134 const Vector& change,
2135 const Index& relative,
2136 const Verbosity& verbosity)
2137{
2138 ARTS_USER_ERROR_IF (QID.nelem() not_eq change.nelem(),
2139 "Mismatch between QID and change input lengths not allowed");
2140
2141 for (Index iq=0; iq<QID.nelem(); iq++)
2142 for (auto& lines: abs_lines_per_species)
2143 abs_linesChangeBaseParameterForMatchingLevel(lines, QID[iq], parameter_name, change[iq], relative, verbosity);
2144}
2145
2146/* Workspace method: Doxygen documentation will be auto-generated */
2148 const QuantumIdentifier& QI,
2149 const String& parameter_name,
2150 const Numeric& x,
2151 const Verbosity&) {
2152 Index parameter_switch = -1;
2153
2154 ARTS_USER_ERROR_IF(parameter_name.nelem() == 0, "parameter_name is empty.\n");
2155 if (parameter_name == "Statistical Weight")
2156 parameter_switch = 1;
2157 else if (parameter_name == "Zeeman Coefficient")
2158 parameter_switch = 2;
2159
2160 for (auto& band : abs_lines) {
2161 for (Index k = 0; k < band.NumLines(); k++) {
2163 QI, band.lines[k].localquanta, band.quantumidentity);
2164 if (lt == Quantum::Number::StateMatchType::Level and lt.low) {
2165 switch (parameter_switch) {
2166 case 1: // "Statistical Weight":
2167 band.lines[k].glow = x;
2168 break;
2169 case 2: // "Zeeman Coefficient":
2170 band.lines[k].zeeman.gl() = x;
2171 break;
2172 default: {
2174 "Usupported paramter_name\n",
2175 parameter_name,
2176 "\nSee method description for supported parameter names.\n")
2177 }
2178 }
2179 }
2180
2181 if (lt == Quantum::Number::StateMatchType::Level and lt.upp) {
2182 switch (parameter_switch) {
2183 case 1: // "Statistical Weight":
2184 band.lines[k].gupp = x;
2185 break;
2186 case 2: // "Zeeman Coefficient":
2187 band.lines[k].zeeman.gu() = x;
2188 break;
2189 default: {
2191 "Usupported paramter_name\n",
2192 parameter_name,
2193 "\nSee method description for supported parameter names.\n")
2194 }
2195 }
2196 }
2197 }
2198 }
2199}
2200
2201/* Workspace method: Doxygen documentation will be auto-generated */
2203 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
2204 const QuantumIdentifier& QI,
2205 const String& parameter_name,
2206 const Numeric& change,
2207 const Verbosity& verbosity) {
2208 for (auto& lines : abs_lines_per_species)
2210 lines, QI, parameter_name, change, verbosity);
2211}
2212
2213/* Workspace method: Doxygen documentation will be auto-generated */
2215 const ArrayOfQuantumIdentifier& QID,
2216 const String& parameter_name,
2217 const Vector& change,
2218 const Verbosity& verbosity) {
2220 QID.nelem() not_eq change.nelem(),
2221 "Mismatch between QID and change input lengths not allowed");
2222
2223 for (Index iq = 0; iq < QID.nelem(); iq++)
2225 abs_lines, QID[iq], parameter_name, change[iq], verbosity);
2226}
2227
2228/* Workspace method: Doxygen documentation will be auto-generated */
2230 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
2231 const ArrayOfQuantumIdentifier& QID,
2232 const String& parameter_name,
2233 const Vector& change,
2234 const Verbosity& verbosity) {
2236 QID.nelem() not_eq change.nelem(),
2237 "Mismatch between QID and change input lengths not allowed");
2238
2239 for (Index iq = 0; iq < QID.nelem(); iq++)
2240 for (auto& lines : abs_lines_per_species)
2242 lines, QID[iq], parameter_name, change[iq], verbosity);
2243}
2244
2248
2249/* Workspace method: Doxygen documentation will be auto-generated */
2251 Index& nlte_do,
2252 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
2253 const EnergyLevelMap& nlte_field,
2254 const Verbosity&) {
2255 nlte_field.ThrowIfNotOK();
2256
2257 if (nlte_field.value.empty()) {
2258 nlte_do = 0;
2259 return;
2260 }
2261 nlte_do = 1;
2262
2263 const Absorption::PopulationType poptyp =
2264 nlte_field.vib_energy.empty() ? Absorption::PopulationType::NLTE
2265 : Absorption::PopulationType::VibTemps;
2266
2267 for (auto& spec_lines : abs_lines_per_species) {
2268 for (auto& band : spec_lines) {
2269 Index low = 0, upp = 0;
2270 for (auto& id : nlte_field.levels) {
2271 for (auto& line : band.lines) {
2272 const auto lt =
2273 poptyp == Absorption::PopulationType::NLTE
2275 id, line.localquanta, band.quantumidentity)
2276 : Quantum::Number::StateMatch(id, band.quantumidentity);
2277 low += lt.low;
2278 upp += lt.upp;
2279 }
2280 }
2281
2283 not(low == 0 or low == band.NumLines()) or
2284 not(upp == 0 or upp == band.NumLines()),
2285 "The band with metadata:\n",
2286 band.MetaData(),
2287 "\n\nwith lines:\n",
2288 band.lines,
2289 "\n\nThe number of levels don't add upp correctly.\n There are ",
2290 band.NumLines(),
2291 " lines but there are ",
2292 low,
2293 " lower level matches and ",
2294 upp,
2295 " upper level matches.")
2296
2297 if (upp or low) band.population = poptyp;
2298 }
2299 }
2300}
2301
2302
2306
2307/* Workspace method: Doxygen documentation will be auto-generated */
2309 ArrayOfArrayOfAbsorptionLines & abs_lines_per_species,
2310 const ArrayOfArrayOfSpeciesTag& abs_species,
2311 const Verbosity&) {
2312 abs_lines_per_species = ArrayOfArrayOfAbsorptionLines(
2313 abs_species.nelem(), ArrayOfAbsorptionLines(0));
2314}
2315
2316/* Workspace method: Doxygen documentation will be auto-generated */
2318 const Vector& f_grid,
2319 const Verbosity&) {
2320 const Numeric fmax = max(f_grid);
2321 const Numeric fmin = min(f_grid);
2322
2323 for (auto& band : abs_lines) {
2324 for (Index k = band.NumLines() - 1; k >= 0; k--) {
2325 const Numeric fcut_upp = band.CutoffFreq(k);
2326 const Numeric fcut_low = band.CutoffFreqMinus(k);
2327
2328 if (fmax < fcut_low or fmin > fcut_upp) {
2329 band.RemoveLine(k);
2330 }
2331 }
2332 }
2333}
2334
2335/* Workspace method: Doxygen documentation will be auto-generated */
2337 ArrayOfArrayOfAbsorptionLines & abs_lines_per_species,
2338 const Vector& f_grid,
2339 const Verbosity& verbosity) {
2340 for (auto& lines : abs_lines_per_species) {
2341 abs_linesCompact(lines, f_grid, verbosity);
2342 }
2343}
2344
2345/* Workspace method: Doxygen documentation will be auto-generated */
2347 const QuantumIdentifier& qid,
2348 const Verbosity&) {
2349 for (Index i = 0; i < abs_lines.nelem(); i++) {
2350 const Quantum::Number::StateMatch lt(qid, abs_lines[i].quantumidentity);
2351 if (lt == Quantum::Number::StateMatchType::Full) {
2352 abs_lines.erase(abs_lines.begin() + i);
2353 break;
2354 }
2355 }
2356}
2357
2361
2362template <class T>
2363std::vector<T> linspace(
2364 T s, T e, typename std::vector<T>::size_type count) noexcept {
2365 std::vector<T> ls(count);
2366
2367 if (count == 0) return ls;
2368 if (count == 1) {
2369 ls.front() = (e + s) / 2;
2370 return ls;
2371 }
2372 const T step = (e - s) / T(count - 1);
2373 ls.front() = s;
2374 ls.back() = e;
2375 for (typename std::vector<T>::size_type i = 1; i < count - 1; ++i)
2376 ls[i] = s + step * T(i);
2377 return ls;
2378}
2379
2380/* Workspace method: Doxygen documentation will be auto-generated */
2382 Vector & f_grid,
2383 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
2384 const Numeric& delta_f_low,
2385 const Numeric& delta_f_upp,
2386 const Index& num_freqs,
2387 const Verbosity&) {
2388 const Index n = nelem(abs_lines_per_species);
2389
2390 ARTS_USER_ERROR_IF(delta_f_low >= delta_f_upp,
2391 "The lower frequency delta has to be smaller "
2392 "than the upper frequency delta");
2393 ARTS_USER_ERROR_IF(num_freqs == 0, "Need more than zero frequency points");
2394 ARTS_USER_ERROR_IF(n < 1,
2395 "No lines found. Error? Use *VectorSet* "
2396 "to resize *f_grid*");
2397
2398 std::vector<Numeric> fout(0);
2399 for (auto& lines : abs_lines_per_species) {
2400 for (auto& band : lines) {
2401 for (Index k = 0; k < band.NumLines(); k++) {
2402 if (num_freqs > 1) {
2403 auto ftmp =
2404 linspace<Numeric>(band.lines[k].F0 + delta_f_low,
2405 band.lines[k].F0 + delta_f_upp,
2406 std::vector<Numeric>::size_type(num_freqs));
2407 for (auto& f : ftmp) {
2408 if (f > 0) fout.push_back(f);
2409 }
2410 } else {
2411 fout.push_back(band.lines[k].F0);
2412 }
2413 }
2414 }
2415 }
2416
2417 std::sort(fout.begin(), fout.end());
2418 fout.erase(std::unique(fout.begin(), fout.end()), fout.end());
2419 f_grid.resize(fout.size());
2420 for (Index i = 0; i < f_grid.nelem(); i++) f_grid[i] = fout[i];
2421}
2422
2426
2428 const ArrayOfSpeciesTag& species,
2429 const Numeric lower_frequency,
2430 const Numeric upper_frequency,
2431 const Numeric lower_intensity,
2432 const Index safe,
2433 const Index flip_flims,
2434 const Verbosity& verbosity) {
2435 const bool care_about_species = species.nelem();
2436
2437 for (auto& band : abs_lines) {
2438 if (care_about_species and
2439 species[0].Isotopologue() not_eq band.Isotopologue())
2440 continue;
2441
2442 auto& lines = band.lines;
2443
2444 if (not safe) {
2445 std::vector<std::size_t> rem;
2446 if (flip_flims) {
2447 for (std::size_t i=lines.size()-1; i<lines.size(); i--) {
2448 auto& line = lines[i];
2449 if ((line.F0 >= lower_frequency and
2450 line.F0 <= upper_frequency) or
2451 line.I0 < lower_intensity)
2452 rem.push_back(i);
2453 }
2454 } else {
2455 for (std::size_t i = lines.size() - 1; i < lines.size(); i--) {
2456 auto& line = lines[i];
2457 if (line.F0 < lower_frequency or line.F0 > upper_frequency or
2458 line.I0 < lower_intensity)
2459 rem.push_back(i);
2460 }
2461 }
2462 for (auto i : rem) band.RemoveLine(i);
2463 } else {
2464 ARTS_USER_ERROR_IF(flip_flims,
2465 "Not allowed to combine GINs flip_flims and safe.")
2466 const bool all_low = std::all_of(
2467 lines.begin(), lines.end(), [lower_frequency](auto& line) {
2468 return line.F0 < lower_frequency;
2469 });
2470 const bool all_upp = std::all_of(
2471 lines.begin(), lines.end(), [upper_frequency](auto& line) {
2472 return line.F0 > upper_frequency;
2473 });
2474 const bool low_int = std::all_of(
2475 lines.begin(), lines.end(), [lower_intensity](auto& line) {
2476 return line.I0 < lower_intensity;
2477 });
2478 if (all_low or all_upp or low_int) lines.resize(0);
2479 }
2480 }
2481
2482 // Removes empty bands
2483 abs_linesRemoveEmptyBands(abs_lines, verbosity);
2484}
2485
2487 const Numeric& lower_frequency,
2488 const Numeric& upper_frequency,
2489 const Numeric& lower_intensity,
2490 const Index& safe,
2491 const Index& flip_flims,
2492 const Verbosity& verbosity) {
2493 remove_impl(abs_lines,
2494 {},
2495 lower_frequency,
2496 upper_frequency,
2497 lower_intensity,
2498 safe,
2499 flip_flims,
2500 verbosity);
2501}
2502
2504 ArrayOfArrayOfAbsorptionLines & abs_lines_per_species,
2505 const Numeric& lower_frequency,
2506 const Numeric& upper_frequency,
2507 const Numeric& lower_intensity,
2508 const Index& safe,
2509 const Index& flip_flims,
2510 const Verbosity& verbosity) {
2511 for (auto& abs_lines : abs_lines_per_species)
2512 abs_linesRemoveLines(abs_lines,
2513 lower_frequency,
2514 upper_frequency,
2515 lower_intensity,
2516 safe,
2517 flip_flims,
2518 verbosity);
2519}
2520
2522 const ArrayOfSpeciesTag& species,
2523 const Numeric& lower_frequency,
2524 const Numeric& upper_frequency,
2525 const Numeric& lower_intensity,
2526 const Index& safe,
2527 const Index& flip_flims,
2528 const Verbosity& verbosity) {
2530 species.nelem() not_eq 1, "Must have a single species, got: ", species)
2531 ARTS_USER_ERROR_IF(species[0].Isotopologue().joker(),
2532 "Cannot give joker species, got: ",
2533 species)
2534
2535 remove_impl(abs_lines,
2536 species,
2537 lower_frequency,
2538 upper_frequency,
2539 lower_intensity,
2540 safe,
2541 flip_flims,
2542 verbosity);
2543}
2544
2546 ArrayOfArrayOfAbsorptionLines & abs_lines_per_species,
2547 const ArrayOfSpeciesTag& species,
2548 const Numeric& lower_frequency,
2549 const Numeric& upper_frequency,
2550 const Numeric& lower_intensity,
2551 const Index& safe,
2552 const Index& flip_flims,
2553 const Verbosity& verbosity) {
2554 for (auto& abs_lines : abs_lines_per_species)
2556 species,
2557 lower_frequency,
2558 upper_frequency,
2559 lower_intensity,
2560 safe,
2561 flip_flims,
2562 verbosity);
2563}
2564
2566 const String& option,
2567 const Verbosity& verbosity) {
2568 const auto opt = Options::toSortingOptionOrThrow(option);
2569
2570 abs_linesRemoveEmptyBands(abs_lines, verbosity);
2571
2572 switch (opt) {
2573 case Options::SortingOption::ByFrequency:
2574 for (auto& band : abs_lines) band.sort_by_frequency();
2575 std::sort(abs_lines.begin(), abs_lines.end(), [](auto& a, auto& b) {
2576 return a.lines[0].F0 <= b.lines[0].F0;
2577 });
2578 break;
2579 case Options::SortingOption::ByEinstein:
2580 for (auto& band : abs_lines) band.sort_by_einstein();
2581 std::sort(abs_lines.begin(), abs_lines.end(), [](auto& a, auto& b) {
2582 return a.lines[0].A <= b.lines[0].A;
2583 });
2584 break;
2585 case Options::SortingOption::FINAL: { /*leave last*/
2586 }
2587 }
2588}
Contains the absorption namespace.
Array< ArrayOfAbsorptionLines > ArrayOfArrayOfAbsorptionLines
This file contains the definition of Array.
Array< Index > ArrayOfIndex
An array of Index.
Definition: array.h:136
int arts_omp_get_max_threads()
Wrapper for omp_get_max_threads.
Definition: arts_omp.cc:46
int arts_omp_get_thread_num()
Wrapper for omp_get_thread_num.
Definition: arts_omp.cc:76
Header file for helper functions for OpenMP.
Stuff related to time in ARTS.
This can be used to make arrays out of anything.
Definition: array.h:48
Index nelem() const ARTS_NOEXCEPT
Definition: array.h:92
The ARTS XML tag class.
Definition: xml_io.h:45
void get_attribute_value(const String &aname, SpeciesTag &value)
Returns value of attribute as type SpeciesTag.
Definition: xml_io.cc:70
bool empty() const noexcept
Definition: matpackIV.h:148
Index nelem() const noexcept
Returns the number of elements.
Definition: matpackI.h:554
bool empty() const noexcept
Returns true if variable size is zero.
Definition: matpackI.h:543
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.
The Vector class.
Definition: matpackI.h:922
void resize(Index n)
Resize function.
Definition: matpackI.cc:411
void check_name(const String &expected_name)
Check tag name.
Definition: xml_io_base.cc:54
void read_from_stream(istream &is)
Reads next XML tag.
Definition: xml_io_base.cc:201
Index nelem() const
Definition: mystring.h:188
Helper macros for debugging.
#define ARTS_ASSERT(condition,...)
Definition: debug.h:83
#define ARTS_USER_ERROR(...)
Definition: debug.h:150
std::string var_string(Args &&... args)
Definition: debug.h:36
#define ARTS_USER_ERROR_IF(condition,...)
Definition: debug.h:134
constexpr void check_enum_error(EnumType type, Messages ... args)
Checks if the enum class type is good and otherwise throws an error message composed by variadic inpu...
Definition: enums.h:83
bool find_xml_file_existence(String &filename)
As find_xml_file but does not throw in the main body.
Definition: file.cc:412
void open_input_file(ifstream &file, const std::string_view name)
Open a file for reading.
Definition: file.cc:128
This file contains basic functions to handle ASCII files.
#define min(a, b)
#define max(a, b)
void abs_lines_per_speciesCreateFromLines(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfAbsorptionLines &abs_lines, const ArrayOfArrayOfSpeciesTag &abs_species, const Verbosity &)
WORKSPACE METHOD: abs_lines_per_speciesCreateFromLines.
Definition: m_abs.cc:105
void abs_speciesSet(ArrayOfArrayOfSpeciesTag &abs_species, Index &propmat_clearsky_agenda_checked, const ArrayOfString &names, const Verbosity &verbosity)
WORKSPACE METHOD: abs_speciesSet.
QuantumIdentifier global_quantumidentifier(const Array< QuantumNumberType > &qns, const QuantumIdentifier &qid)
Selects the global quantum numbers.
void abs_linesCutoff(ArrayOfAbsorptionLines &abs_lines, const String &type, const Numeric &x, const Verbosity &)
WORKSPACE METHOD: abs_linesCutoff.
void abs_lines_per_speciesManualMirroringSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfSpeciesTag &species, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesManualMirroringSpecies.
void abs_lines_per_speciesBaseParameterMatchingLines(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Verbosity &verbosity)
void abs_lines_per_speciesChangeBaseParameterForMatchingLevel(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Index &relative, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesChangeBaseParameterForMatchingLevel.
void abs_lines_per_speciesPopulation(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesPopulation.
void abs_linesPopulationMatch(ArrayOfAbsorptionLines &abs_lines, const String &type, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesPopulationMatch.
void abs_linesDeleteBadF0(ArrayOfAbsorptionLines &abs_lines, const Numeric &f0, const Index &lower, const Verbosity &)
WORKSPACE METHOD: abs_linesDeleteBadF0.
void abs_linesNormalization(ArrayOfAbsorptionLines &abs_lines, const String &type, const Verbosity &)
WORKSPACE METHOD: abs_linesNormalization.
void abs_linesT0(ArrayOfAbsorptionLines &abs_lines, const Numeric &x, const Verbosity &)
WORKSPACE METHOD: abs_linesT0.
void abs_linesFlatten(ArrayOfAbsorptionLines &abs_lines, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesFlatten.
void abs_linesNormalizationMatch(ArrayOfAbsorptionLines &abs_lines, const String &type, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesNormalizationMatch.
void abs_lines_per_speciesLinemixingLimitSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const Numeric &x, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesLinemixingLimitSpecies.
bool check_local(const Array< QuantumNumberType > &local_state)
void abs_lines_per_speciesChangeBaseParameterForMatchingLines(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Index &relative, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesChangeBaseParameterForMatchingLines.
void abs_linesKeepBand(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &qid, const Verbosity &)
WORKSPACE METHOD: abs_linesKeepBand.
void abs_lines_per_speciesLineShapeModelParametersMatchingLines(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const QuantumIdentifier &QI, const String &parameter, const String &species, const String &temperaturemodel, const Vector &new_values, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesLineShapeModelParametersMatchingLines.
void abs_lines_per_speciesLineShapeTypeMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesLineShapeTypeMatch.
void abs_linesMirroringMatch(ArrayOfAbsorptionLines &abs_lines, const String &type, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesMirroringMatch.
void ReadARTSCAT(ArrayOfAbsorptionLines &abs_lines, const String &artscat_file, const Numeric &fmin, const Numeric &fmax, const String &globalquantumnumbers, const String &localquantumnumbers, const String &normalization_option, const String &mirroring_option, const String &population_option, const String &lineshapetype_option, const String &cutoff_option, const Numeric &cutoff_value, const Numeric &linemixinglimit_value, const Verbosity &verbosity)
WORKSPACE METHOD: ReadARTSCAT.
void abs_lines_per_speciesChangeBaseParameterForSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Index &relative, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesChangeBaseParameterForSpecies.
void abs_lines_per_speciesPopulationNlteField(Index &nlte_do, ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const EnergyLevelMap &nlte_field, const Verbosity &)
WORKSPACE METHOD: abs_lines_per_speciesPopulationNlteField.
void abs_lines_per_speciesCutoff(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Numeric &x, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesCutoff.
void abs_lines_per_speciesLinemixingLimit(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Numeric &x, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesLinemixingLimit.
void abs_lines_per_speciesLinemixingLimitMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesLinemixingLimitMatch.
void abs_linesChangeBaseParameterForMatchingLines(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Index &relative, const Verbosity &)
WORKSPACE METHOD: abs_linesChangeBaseParameterForMatchingLines.
void abs_linesEmptyBroadeningParameters(ArrayOfAbsorptionLines &abs_lines, const Verbosity &)
WORKSPACE METHOD: abs_linesEmptyBroadeningParameters.
void abs_lines_per_speciesT0(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Numeric &x, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesT0.
void abs_lines_per_speciesMirroringMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesMirroringMatch.
void abs_lines_per_speciesReadSpeciesSplitCatalog(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const String &basename, const Index &robust, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesReadSpeciesSplitCatalog.
void ReadJPL(ArrayOfAbsorptionLines &abs_lines, const String &jpl_file, const Numeric &fmin, const Numeric &fmax, const String &globalquantumnumbers, const String &localquantumnumbers, const String &normalization_option, const String &mirroring_option, const String &population_option, const String &lineshapetype_option, const String &cutoff_option, const Numeric &cutoff_value, const Numeric &linemixinglimit_value, const Verbosity &verbosity)
WORKSPACE METHOD: ReadJPL.
void abs_lines_per_speciesManualMirroring(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesManualMirroring.
void abs_linesReplaceLines(ArrayOfAbsorptionLines &abs_lines, const ArrayOfAbsorptionLines &replacing_lines, const Verbosity &)
WORKSPACE METHOD: abs_linesReplaceLines.
void abs_linesLineShapeModelParametersMatchingLines(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &QI, const String &parameter, const String &species, const String &temperaturemodel, const Vector &new_values, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesLineShapeModelParametersMatchingLines.
void abs_linesCompact(ArrayOfAbsorptionLines &abs_lines, const Vector &f_grid, const Verbosity &)
WORKSPACE METHOD: abs_linesCompact.
void remove_impl(ArrayOfAbsorptionLines &abs_lines, const ArrayOfSpeciesTag &species, const Numeric lower_frequency, const Numeric upper_frequency, const Numeric lower_intensity, const Index safe, const Index flip_flims, const Verbosity &verbosity)
void abs_linesRemoveBand(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &qid, const Verbosity &)
WORKSPACE METHOD: abs_linesRemoveBand.
constexpr Index merge_local_lines_size
void abs_linesCutoffMatch(ArrayOfAbsorptionLines &abs_lines, const String &type, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesCutoffMatch.
void abs_linesChangeBaseParameterForMatchingLevels(ArrayOfAbsorptionLines &abs_lines, const ArrayOfQuantumIdentifier &QID, const String &parameter_name, const Vector &change, const Index &relative, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesChangeBaseParameterForMatchingLevels.
void abs_linesBaseParameterMatchingLevels(ArrayOfAbsorptionLines &abs_lines, const ArrayOfQuantumIdentifier &QID, const String &parameter_name, const Vector &change, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesBaseParameterMatchingLevels.
void abs_lines_per_speciesRemoveLinesFromSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfSpeciesTag &species, const Numeric &lower_frequency, const Numeric &upper_frequency, const Numeric &lower_intensity, const Index &safe, const Index &flip_flims, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesRemoveLinesFromSpecies.
void abs_linesRemoveEmptyBands(ArrayOfAbsorptionLines &abs_lines, const Verbosity &)
WORKSPACE METHOD: abs_linesRemoveEmptyBands.
void abs_lines_per_speciesNormalizationSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const String &type, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesNormalizationSpecies.
void ReadArrayOfARTSCAT(ArrayOfAbsorptionLines &abs_lines, const String &artscat_file, const Numeric &fmin, const Numeric &fmax, const String &globalquantumnumbers, const String &localquantumnumbers, const String &normalization_option, const String &mirroring_option, const String &population_option, const String &lineshapetype_option, const String &cutoff_option, const Numeric &cutoff_value, const Numeric &linemixinglimit_value, const Verbosity &verbosity)
WORKSPACE METHOD: ReadArrayOfARTSCAT.
void abs_linesDeleteWithLines(ArrayOfAbsorptionLines &abs_lines, const ArrayOfAbsorptionLines &deleting_lines, const Verbosity &)
void ReadSplitARTSCAT(ArrayOfAbsorptionLines &abs_lines, const ArrayOfArrayOfSpeciesTag &abs_species, const String &basename, const Numeric &fmin, const Numeric &fmax, const String &globalquantumnumbers, const String &localquantumnumbers, const Index &ignore_missing, const String &normalization_option, const String &mirroring_option, const String &population_option, const String &lineshapetype_option, const String &cutoff_option, const Numeric &cutoff_value, const Numeric &linemixinglimit_value, const Verbosity &verbosity)
WORKSPACE METHOD: ReadSplitARTSCAT.
void abs_linesPopulation(ArrayOfAbsorptionLines &abs_lines, const String &type, const Verbosity &)
WORKSPACE METHOD: abs_linesPopulation.
void abs_lines_per_speciesT0Species(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const Numeric &x, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesT0Species.
void abs_lines_per_speciesNormalization(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesNormalization.
void abs_lines_per_speciesRemoveLines(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Numeric &lower_frequency, const Numeric &upper_frequency, const Numeric &lower_intensity, const Index &safe, const Index &flip_flims, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesRemoveLines.
void abs_linesAppendWithLines(ArrayOfAbsorptionLines &abs_lines, const ArrayOfAbsorptionLines &appending_lines, const Index &safe, const Verbosity &)
void abs_lines_per_speciesMirroring(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesMirroring.
void abs_lines_per_speciesCutoffSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const String &type, const Numeric &x, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesCutoffSpecies.
void abs_linesLineShapeTypeMatch(ArrayOfAbsorptionLines &abs_lines, const String &type, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesLineShapeTypeMatch.
Array< QuantumNumberType > string2vecqn(std::string_view qnstr)
Get a list of quantum numbers from a string.
void abs_linesChangeBaseParameterForMatchingLevel(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Index &relative, const Verbosity &)
WORKSPACE METHOD: abs_linesChangeBaseParameterForMatchingLevel.
void abs_lines_per_speciesNormalizationMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesNormalizationMatch.
void abs_lines_per_speciesPopulationMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesPopulationMatch.
void ReadLBLRTM(ArrayOfAbsorptionLines &abs_lines, const String &lblrtm_file, const Numeric &fmin, const Numeric &fmax, const String &globalquantumnumbers, const String &localquantumnumbers, const String &normalization_option, const String &mirroring_option, const String &population_option, const String &lineshapetype_option, const String &cutoff_option, const Numeric &cutoff_value, const Numeric &linemixinglimit_value, const Verbosity &verbosity)
WORKSPACE METHOD: ReadLBLRTM.
void abs_linesBaseParameterMatchingLevel(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &x, const Verbosity &)
WORKSPACE METHOD: abs_linesBaseParameterMatchingLevel.
void abs_linesT0Match(ArrayOfAbsorptionLines &abs_lines, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesT0Match.
void abs_linesRemoveLines(ArrayOfAbsorptionLines &abs_lines, const Numeric &lower_frequency, const Numeric &upper_frequency, const Numeric &lower_intensity, const Index &safe, const Index &flip_flims, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesRemoveLines.
void abs_lines_per_speciesPopulationSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const String &type, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesPopulationSpecies.
void abs_linesManualMirroring(ArrayOfAbsorptionLines &abs_lines, const Verbosity &)
WORKSPACE METHOD: abs_linesManualMirroring.
void ReadHITRAN(ArrayOfAbsorptionLines &abs_lines, const String &hitran_file, const Numeric &fmin, const Numeric &fmax, const String &globalquantumnumbers, const String &localquantumnumbers, const String &hitran_type, const String &normalization_option, const String &mirroring_option, const String &population_option, const String &lineshapetype_option, const String &cutoff_option, const Numeric &cutoff_value, const Numeric &linemixinglimit_value, const Verbosity &verbosity)
WORKSPACE METHOD: ReadHITRAN.
void abs_lines_per_speciesT0Match(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesT0Match.
void abs_lines_per_speciesSetEmpty(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const Verbosity &)
WORKSPACE METHOD: abs_lines_per_speciesSetEmpty.
void abs_lines_per_speciesBaseParameterMatchingLevel(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesBaseParameterMatchingLevel.
void merge_external_line(ArrayOfAbsorptionLines &abs_lines, const Absorption::SingleLineExternal &sline, const QuantumIdentifier &global_qid)
Merge an external line to abs_lines.
void abs_lines_per_speciesBaseParameterMatchingLevels(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfQuantumIdentifier &QID, const String &parameter_name, const Vector &change, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesBaseParameterMatchingLevels.
void abs_linesReplaceWithLines(ArrayOfAbsorptionLines &abs_lines, const ArrayOfAbsorptionLines &replacing_lines, const Verbosity &)
void abs_lines_per_speciesMirroringSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const String &type, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesMirroringSpecies.
void abs_linesLinemixingLimit(ArrayOfAbsorptionLines &abs_lines, const Numeric &x, const Verbosity &)
WORKSPACE METHOD: abs_linesLinemixingLimit.
void abs_lines_per_speciesCutoffMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesCutoffMatch.
void abs_linesRemoveLinesFromSpecies(ArrayOfAbsorptionLines &abs_lines, const ArrayOfSpeciesTag &species, const Numeric &lower_frequency, const Numeric &upper_frequency, const Numeric &lower_intensity, const Index &safe, const Index &flip_flims, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesRemoveLinesFromSpecies.
void abs_lines_per_speciesFlatten(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesFlatten.
void abs_linesBaseParameterMatchingLines(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &x, const Verbosity &)
WORKSPACE METHOD: abs_linesBaseParameterMatchingLines.
void abs_lines_per_speciesWriteSpeciesSplitCatalog(const String &output_format, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesWriteSpeciesSplitCatalog.
void abs_linesLineShapeType(ArrayOfAbsorptionLines &abs_lines, const String &type, const Verbosity &)
WORKSPACE METHOD: abs_linesLineShapeType.
void abs_lines_per_speciesChangeBaseParameterForMatchingLevels(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfQuantumIdentifier &QID, const String &parameter_name, const Vector &change, const Index &relative, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesChangeBaseParameterForMatchingLevels.
void abs_linesReplaceBands(ArrayOfAbsorptionLines &abs_lines, const ArrayOfAbsorptionLines &replacing_bands, const Verbosity &)
WORKSPACE METHOD: abs_linesReplaceBands.
void abs_lines_per_speciesLineShapeType(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesLineShapeType.
void abs_linesMirroring(ArrayOfAbsorptionLines &abs_lines, const String &type, const Verbosity &)
WORKSPACE METHOD: abs_linesMirroring.
void f_gridFromAbsorptionLines(Vector &f_grid, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Numeric &delta_f_low, const Numeric &delta_f_upp, const Index &num_freqs, const Verbosity &)
WORKSPACE METHOD: f_gridFromAbsorptionLines.
void abs_linesLinemixingLimitMatch(ArrayOfAbsorptionLines &abs_lines, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesLinemixingLimitMatch.
void abs_linesSort(ArrayOfAbsorptionLines &abs_lines, const String &option, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesSort.
void merge_local_lines(ArrayOfAbsorptionLines &abs_lines, const ArrayOfAbsorptionLines &local_lines)
Merge lines to abs_lines.
void abs_linesWriteSpeciesSplitCatalog(const String &output_format, const ArrayOfAbsorptionLines &abs_lines, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesWriteSpeciesSplitCatalog.
void abs_lines_per_speciesLineShapeTypeSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const String &type, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesLineShapeTypeSpecies.
std::vector< T > linspace(T s, T e, typename std::vector< T >::size_type count) noexcept
void CheckUnique(const ArrayOfAbsorptionLines &lines, const Verbosity &)
WORKSPACE METHOD: CheckUnique.
void abs_linesReadSpeciesSplitCatalog(ArrayOfAbsorptionLines &abs_lines, const String &basename, const Index &robust, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesReadSpeciesSplitCatalog.
void abs_lines_per_speciesCompact(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Vector &f_grid, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesCompact.
Workspace methods and template functions for supergeneric XML IO.
void WriteXML(const String &file_format, const T &v, const String &f, const Index &no_clobber, const String &v_name, const String &f_name, const String &no_clobber_name, const Verbosity &verbosity)
WORKSPACE METHOD: WriteXML.
Definition: m_xml.h:118
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
const Joker joker
#define CREATE_OUT3
Definition: messages.h:206
#define CREATE_OUT2
Definition: messages.h:205
my_basic_string< char > String
The String type for ARTS.
Definition: mystring.h:215
SingleLineExternal ReadFromHitranOnlineStream(istream &is)
Read from HITRAN online.
SingleLineExternal ReadFromJplStream(istream &is)
Read from JPL.
SingleLineExternal ReadFromArtscat3Stream(istream &is)
Read from ARTSCAT-3.
Index nelem(const Lines &l)
Number of lines.
SingleLineExternal ReadFromLBLRTMStream(istream &is)
Read from LBLRTM.
SingleLineExternal ReadFromArtscat5Stream(istream &is)
Read from ARTSCAT-5.
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 bool modelparameterEmpty(const ModelParameters mp) noexcept
constexpr Index nVars
Current max number of line shape variables.
constexpr ValueType common_value_type(ValueType a, ValueType b) noexcept
Return a common type between a and b.
constexpr std::array global_types
A default state of global quantum numbers.
constexpr std::string_view items(std::string_view s, std::size_t i) noexcept
Get a view of a number of space-separated items from the list.
constexpr std::array local_types
A default state of local quantum numbers.
constexpr Index count_items(std::string_view s) noexcept
Count all space-separated items in s.
VectorView var(VectorView var, const Vector &y, const ArrayOfVector &ys, const Index start, const Index end_tmp)
Compute the variance of the ranged ys.
Definition: raw.cc:159
constexpr std::array Isotopologues
A list of all ARTS isotopologues, note how the species enum class input HAS to be sorted.
Definition: isotopologues.h:57
ArrayOfIsotopeRecord isotopologues(Species spec)
Definition: isotopologues.cc:6
Model GetAdvancedModel(const QuantumIdentifier &qid) ARTS_NOEXCEPT
Returns an advanced Zeeman model.
Definition: zeemandata.cc:129
constexpr Rational end(Rational Ju, Rational Jl, Polarization type) noexcept
Gives the largest M for a polarization type of this transition.
Definition: zeemandata.h:113
#define v
#define a
#define c
#define b
#define N
Definition: rng.cc:164
std::set< Species::Species > lbl_species(const ArrayOfArrayOfSpeciesTag &abs_species) noexcept
Species::Tag SpeciesTag
Definition: species_tags.h:101
Index NumLines() const noexcept
Number of 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.
Single line reading output.
Numeric F0
Central frequency.
Quantum::Number::LocalState localquanta
Local quantum numbers.
Zeeman::Model zeeman
Zeeman model.
ArrayOfQuantumIdentifier levels
void ThrowIfNotOK() const ARTS_NOEXCEPT
Coefficients and temperature model for SingleSpeciesModel.
static constexpr Index N
A logical struct for global quantum numbers with species identifiers.
Species::IsotopeRecord Isotopologue() const noexcept
StateMatchType operates so that a check less than a level should be 'better', bar None.
void xml_find_and_open_input_file(std::shared_ptr< istream > &ifs, const String &filename, const Verbosity &verbosity)
Open plain or zipped xml file.
Definition: xml_io.cc:136
void xml_read_from_file(const String &filename, T &type, const Verbosity &verbosity)
Reads data from XML file.
Definition: xml_io.h:151
void xml_read_header_from_stream(istream &is, FileType &ftype, NumericType &ntype, EndianType &etype, const Verbosity &verbosity)
Reads XML header and root tag.
Definition: xml_io_base.cc:574
@ FILE_TYPE_ASCII
Definition: xml_io_base.h:43
@ NUMERIC_TYPE_DOUBLE
Definition: xml_io_base.h:48
@ ENDIAN_TYPE_LITTLE
Definition: xml_io_base.h:49