ARTS 2.5.0 (git: 9ee3ac6c)
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 "auto_md.h"
32#include "enums.h"
33#include "file.h"
34#include "global_data.h"
35#include "m_xml.h"
36
40
42 const Verbosity&)
43{
44 abs_lines.erase(std::remove_if(abs_lines.begin(), abs_lines.end(),
45 [](auto& x){return x.NumLines() == 0;}),
46 abs_lines.end());
47}
48
50 const Verbosity& verbosity)
51{
52 const Index n = abs_lines.nelem();
53
54 for (Index i=0; i<n; i++) {
55 AbsorptionLines& band = abs_lines[i];
56 if (band.NumLines()) {
57 for (Index j=i+1; j<n; j++) {
58 if (band.Match(abs_lines[j])) {
59 for (auto& line: abs_lines[j].AllLines()) {
60 band.AppendSingleLine(line);
61 }
62 abs_lines[j].AllLines().clear();
63 }
64 }
65 }
66 }
67
68 abs_linesRemoveEmptyBands(abs_lines, verbosity);
69}
70
72 const Verbosity& verbosity)
73{
74 for (auto& abs_lines: abs_lines_per_species) abs_linesFlatten(abs_lines, verbosity);
75}
76
80
87std::vector<QuantumNumberType> string2vecqn(const String& qnstr)
88{
89 std::vector<QuantumNumberType> nums(0);
90
91 String part;
92 if (qnstr not_eq "") {
93 std::istringstream str(qnstr);
94 while (not str.eof()) {
95 str >> part;
97 "The quantum number key: \"", part, "\" is invalid.\n")
98 nums.push_back(string2quantumnumbertype(part));
99 }
100 }
101
102 return nums;
103}
104
105/* Workspace method: Doxygen documentation will be auto-generated */
107 const String& artscat_file,
108 const Numeric& fmin,
109 const Numeric& fmax,
110 const String& globalquantumnumbers,
111 const String& localquantumnumbers,
112 const String& normalization_option,
113 const String& mirroring_option,
114 const String& population_option,
115 const String& lineshapetype_option,
116 const String& cutoff_option,
117 const Numeric& cutoff_value,
118 const Numeric& linemixinglimit_value,
119 const Verbosity& verbosity)
120{
121 // Global numbers
122 const std::vector<QuantumNumberType> global_nums = string2vecqn(globalquantumnumbers);
123
124 // Local numbers
125 const std::vector<QuantumNumberType> local_nums = string2vecqn(localquantumnumbers);
126
128
129 ArtsXMLTag tag(verbosity);
130 Index nelem;
131
132 // ARTSCAT data
133 shared_ptr<istream> ifs = nullptr;
134 xml_find_and_open_input_file(ifs, artscat_file, verbosity);
135 istream& is_xml = *ifs;
136 auto a = FILE_TYPE_ASCII;
137 auto b = NUMERIC_TYPE_DOUBLE;
138 auto c = ENDIAN_TYPE_LITTLE;
140 verbosity);
141
142 tag.read_from_stream(is_xml);
143 tag.check_name("Array");
144
145 Index num_arrays;
146 tag.get_attribute_value("nelem", num_arrays);
147
148 std::vector<Absorption::SingleLineExternal> v(0);
149
150 for (Index i=0; i<num_arrays; i++) {
151 tag.read_from_stream(is_xml);
152 tag.check_name("ArrayOfLineRecord");
153
154 tag.get_attribute_value("nelem", nelem);
155
156 String version;
157 tag.get_attribute_value("version", version);
158
159 Index artscat_version;
160
161 if (version == "3") {
162 artscat_version = 3;
163 } else if (version.substr(0, 8) != "ARTSCAT-") {
165 "The ARTS line file you are trying to read does not contain a valid version tag.\n"
166 "Probably it was created with an older version of ARTS that used different units.")
167 } else {
168 istringstream is(version.substr(8));
169 is >> artscat_version;
170 }
171
172 ARTS_USER_ERROR_IF (artscat_version < 3 or artscat_version > 5,
173 "Unknown ARTS line file version: ", version)
174
175 bool go_on = true;
176 Index n = 0;
177 while (go_on and n<nelem) {
178 switch(artscat_version) {
179 case 3:
180 v.push_back(Absorption::ReadFromArtscat3Stream(is_xml));
181 break;
182 case 4:
183 v.push_back(Absorption::ReadFromArtscat4Stream(is_xml));
184 break;
185 case 5:
186 v.push_back(Absorption::ReadFromArtscat5Stream(is_xml));
187 break;
188 default:
189 ARTS_ASSERT (false, "Bad version!");
190 }
191
192 if (v.back().bad) {
193 v.pop_back();
194 go_on = false;
195 } else if (v.back().line.F0() < fmin) {
196 v.pop_back();
197 } else if (v.back().line.F0() > fmax) {
198 v.pop_back();
199 go_on = false;
200 }
201
202 n++;
203 }
204
205 tag.read_from_stream(is_xml);
206 tag.check_name("/ArrayOfLineRecord");
207 }
208
209 tag.read_from_stream(is_xml);
210 tag.check_name("/Array");
211
212 for (auto& x: v)
213 x.line.Zeeman() = Zeeman::GetAdvancedModel(x.quantumidentity);
214
215 auto x = Absorption::split_list_of_external_lines(v, local_nums, global_nums);
216 abs_lines.resize(0);
217 abs_lines.reserve(x.size());
218 while (x.size()) {
219 abs_lines.push_back(x.back());
220 abs_lines.back().sort_by_frequency();
221 x.pop_back();
222 }
223
224 abs_linesSetNormalization(abs_lines, normalization_option, verbosity);
225 abs_linesSetMirroring(abs_lines, mirroring_option, verbosity);
226 abs_linesSetPopulation(abs_lines, population_option, verbosity);
227 abs_linesSetLineShapeType(abs_lines, lineshapetype_option, verbosity);
228 abs_linesSetCutoff(abs_lines, cutoff_option, cutoff_value, verbosity);
229 abs_linesSetLinemixingLimit(abs_lines, linemixinglimit_value, verbosity);
230}
231
232/* Workspace method: Doxygen documentation will be auto-generated */
234 const String& artscat_file,
235 const Numeric& fmin,
236 const Numeric& fmax,
237 const String& globalquantumnumbers,
238 const String& localquantumnumbers,
239 const String& normalization_option,
240 const String& mirroring_option,
241 const String& population_option,
242 const String& lineshapetype_option,
243 const String& cutoff_option,
244 const Numeric& cutoff_value,
245 const Numeric& linemixinglimit_value,
246 const Verbosity& verbosity)
247{
248 // Global numbers
249 const std::vector<QuantumNumberType> global_nums = string2vecqn(globalquantumnumbers);
250
251 // Local numbers
252 const std::vector<QuantumNumberType> local_nums = string2vecqn(localquantumnumbers);
253
255
256 ArtsXMLTag tag(verbosity);
257 Index nelem;
258
259 // ARTSCAT data
260 shared_ptr<istream> ifs = nullptr;
261 xml_find_and_open_input_file(ifs, artscat_file, verbosity);
262 istream& is_xml = *ifs;
263 auto a = FILE_TYPE_ASCII;
264 auto b = NUMERIC_TYPE_DOUBLE;
265 auto c = ENDIAN_TYPE_LITTLE;
267 verbosity);
268
269 tag.read_from_stream(is_xml);
270 tag.check_name("ArrayOfLineRecord");
271
272 tag.get_attribute_value("nelem", nelem);
273
274 String version;
275 tag.get_attribute_value("version", version);
276
277 Index artscat_version;
278
279 if (version == "3") {
280 artscat_version = 3;
281 } else if (version.substr(0, 8) != "ARTSCAT-") {
283 "The ARTS line file you are trying to read does not contain a valid version tag.\n"
284 "Probably it was created with an older version of ARTS that used different units.")
285 } else {
286 istringstream is(version.substr(8));
287 is >> artscat_version;
288 }
289
290 ARTS_USER_ERROR_IF (artscat_version < 3 or artscat_version > 5,
291 "Unknown ARTS line file version: ", version)
292
293 std::vector<Absorption::SingleLineExternal> v(0);
294
295 bool go_on = true;
296 Index n = 0;
297 while (n<nelem) {
298 if (go_on) {
299 switch(artscat_version) {
300 case 3:
301 v.push_back(Absorption::ReadFromArtscat3Stream(is_xml));
302 break;
303 case 4:
304 v.push_back(Absorption::ReadFromArtscat4Stream(is_xml));
305 break;
306 case 5:
307 v.push_back(Absorption::ReadFromArtscat5Stream(is_xml));
308 break;
309 default:
310 ARTS_ASSERT(false, "Bad version!");
311 }
312
313 if (v.back().bad) {
314 v.pop_back();
315 go_on = false;
316 } else if (v.back().line.F0() < fmin) {
317 v.pop_back();
318 } else if (v.back().line.F0() > fmax) {
319 v.pop_back();
320 go_on = false;
321 }
322 } else {
323 String line;
324 getline(is_xml, line);
325 }
326
327 n++;
328 }
329
330 tag.read_from_stream(is_xml);
331 tag.check_name("/ArrayOfLineRecord");
332
333 for (auto& x: v)
334 x.line.Zeeman() = Zeeman::GetAdvancedModel(x.quantumidentity);
335
336 auto x = Absorption::split_list_of_external_lines(v, local_nums, global_nums);
337 abs_lines.resize(0);
338 abs_lines.reserve(x.size());
339 while (x.size()) {
340 abs_lines.push_back(x.back());
341 abs_lines.back().sort_by_frequency();
342 x.pop_back();
343 }
344
345 abs_linesSetNormalization(abs_lines, normalization_option, verbosity);
346 abs_linesSetMirroring(abs_lines, mirroring_option, verbosity);
347 abs_linesSetPopulation(abs_lines, population_option, verbosity);
348 abs_linesSetLineShapeType(abs_lines, lineshapetype_option, verbosity);
349 abs_linesSetCutoff(abs_lines, cutoff_option, cutoff_value, verbosity);
350 abs_linesSetLinemixingLimit(abs_lines, linemixinglimit_value, verbosity);
351}
352
353/* Workspace method: Doxygen documentation will be auto-generated */
355 const ArrayOfArrayOfSpeciesTag& abs_species,
356 const String& basename,
357 const Numeric& fmin,
358 const Numeric& fmax,
359 const String& globalquantumnumbers,
360 const String& localquantumnumbers,
361 const Index& ignore_missing,
362 const String& normalization_option,
363 const String& mirroring_option,
364 const String& population_option,
365 const String& lineshapetype_option,
366 const String& cutoff_option,
367 const Numeric& cutoff_value,
368 const Numeric& linemixinglimit_value,
369 const Verbosity& verbosity)
370{
371 // Build a set of species indices. Duplicates are ignored.
372 const std::set<Species::Species> unique_species = lbl_species(abs_species);
373
374 String tmpbasename = basename;
375 if (basename.length() && basename[basename.length() - 1] != '/') {
376 tmpbasename += '.';
377 }
378
379 // Read catalogs for each identified species and put them all into
380 // abs_lines.
381 abs_lines.resize(0);
382 for (auto& spec: unique_species) {
383 ArrayOfAbsorptionLines more_abs_lines;
384
385 try {
386 ReadARTSCAT(more_abs_lines,
387 tmpbasename + String(Species::toShortName(spec)) + ".xml",
388 fmin,
389 fmax,
390 globalquantumnumbers,
391 localquantumnumbers,
392 normalization_option,
393 mirroring_option,
394 population_option,
395 lineshapetype_option,
396 cutoff_option,
397 cutoff_value,
398 linemixinglimit_value,
399 verbosity);
400
401 // Either find a line like this in the list of lines or start a new Lines
402 for (auto& newband: more_abs_lines) {
403 bool found = false;
404 for (auto& band: abs_lines) {
405 if (band.Match(newband)) {
406 for (Index k=0; k<newband.NumLines(); k++) {
407 band.AppendSingleLine(newband.Line(k));
408 found = true;
409 }
410 }
411 }
412 if (not found) {
413 abs_lines.push_back(newband);
414 }
415 }
416 } catch (const std::exception& e) {
417 ARTS_USER_ERROR_IF (not ignore_missing,
418 "Errors in calls by *propmat_clearskyAddZeeman*:\n",
419 e.what())
420 continue;
421 }
422 }
423
424 for (auto& band: abs_lines)
425 band.sort_by_frequency();
426
427 if (normalization_option != "None")
428 abs_linesSetNormalization(abs_lines, normalization_option, verbosity);
429 if (mirroring_option != "None")
430 abs_linesSetMirroring(abs_lines, mirroring_option, verbosity);
431 if (population_option != "None")
432 abs_linesSetPopulation(abs_lines, population_option, verbosity);
433 if (lineshapetype_option != "None")
434 abs_linesSetLineShapeType(abs_lines, lineshapetype_option, verbosity);
435 if (cutoff_option != "None")
436 abs_linesSetCutoff(abs_lines, cutoff_option, cutoff_value, verbosity);
437 abs_linesSetLinemixingLimit(abs_lines, linemixinglimit_value, verbosity);
438}
439
440/* Workspace method: Doxygen documentation will be auto-generated */
442 const String& hitran_file,
443 const Numeric& fmin,
444 const Numeric& fmax,
445 const String& globalquantumnumbers,
446 const String& localquantumnumbers,
447 const String& hitran_type,
448 const String& normalization_option,
449 const String& mirroring_option,
450 const String& population_option,
451 const String& lineshapetype_option,
452 const String& cutoff_option,
453 const Numeric& cutoff_value,
454 const Numeric& linemixinglimit_value,
455 const Verbosity& verbosity)
456{
457 abs_lines.resize(0);
458
459 // Global numbers
460 const std::vector<QuantumNumberType> global_nums = string2vecqn(globalquantumnumbers);
461
462 // Local numbers
463 const std::vector<QuantumNumberType> local_nums = string2vecqn(localquantumnumbers);
464
465 // HITRAN type
466 const Options::HitranType hitran_version = Options::toHitranTypeOrThrow(hitran_type);
467
468 // Hitran data
469 ifstream is;
470 open_input_file(is, hitran_file);
471
472 bool go_on = true;
473 while (go_on) {
475 switch (hitran_version) {
476 case Options::HitranType::Post2012:
478 break;
479 case Options::HitranType::From2004To2012:
481 break;
482 case Options::HitranType::Pre2004:
484 break;
485 case Options::HitranType::Online:
487 break;
488 default:
489 ARTS_ASSERT(false, "The HitranType enum class has to be fully updated!\n");
490 }
491
492 if (sline.bad) {
493 if (is.eof())
494 break;
495 ARTS_USER_ERROR("Cannot read line ", nelem(abs_lines) + 1);
496 }
497 if (sline.line.F0() < fmin)
498 continue; // Skip this line
499 if (sline.line.F0() > fmax)
500 break; // We assume sorted so quit here
501
502 // Set Zeeman if implemented
504
505 // Get the global quantum number identifier
506 QuantumIdentifier global_qid(sline.quantumidentity.Isotopologue(), Quantum::IdentifierType::Transition);
507 for(auto qn: global_nums) {
508 global_qid.Lower()[qn] = sline.quantumidentity.Lower()[qn];
509 global_qid.Upper()[qn] = sline.quantumidentity.Upper()[qn];
510 }
511
512 // Get local quantum numbers into the line
513 sline.line.LowerQuantumNumbers().reserve(local_nums.size());
514 sline.line.UpperQuantumNumbers().reserve(local_nums.size());
515 for(auto qn: local_nums) {
516 sline.line.LowerQuantumNumbers().emplace_back(sline.quantumidentity.Lower()[qn]);
517 sline.line.UpperQuantumNumbers().emplace_back(sline.quantumidentity.Upper()[qn]);
518 }
519
520 // Either find a line like this in the list of lines or start a new Lines
521 auto band = std::find_if(abs_lines.begin(), abs_lines.end(), [&](const auto& li){return li.MatchWithExternal(sline, global_qid);});
522 if (band not_eq abs_lines.end()) {
523 band -> AppendSingleLine(sline.line);
524 } else {
525 abs_lines.emplace_back(sline.selfbroadening, sline.bathbroadening, sline.cutoff,
526 sline.mirroring, sline.population, sline.normalization,
527 sline.lineshapetype, sline.T0, sline.cutofffreq,
528 sline.linemixinglimit, global_qid, local_nums, sline.species,
529 std::vector<AbsorptionSingleLine>{sline.line});
530 }
531 }
532
533 abs_linesSetNormalization(abs_lines, normalization_option, verbosity);
534 abs_linesSetMirroring(abs_lines, mirroring_option, verbosity);
535 abs_linesSetPopulation(abs_lines, population_option, verbosity);
536 abs_linesSetLineShapeType(abs_lines, lineshapetype_option, verbosity);
537 abs_linesSetCutoff(abs_lines, cutoff_option, cutoff_value, verbosity);
538 abs_linesSetLinemixingLimit(abs_lines, linemixinglimit_value, verbosity);
539}
540
541/* Workspace method: Doxygen documentation will be auto-generated */
543 const String& lblrtm_file,
544 const Numeric& fmin,
545 const Numeric& fmax,
546 const String& globalquantumnumbers,
547 const String& localquantumnumbers,
548 const String& normalization_option,
549 const String& mirroring_option,
550 const String& population_option,
551 const String& lineshapetype_option,
552 const String& cutoff_option,
553 const Numeric& cutoff_value,
554 const Numeric& linemixinglimit_value,
555 const Verbosity& verbosity)
556{
557 // Global numbers
558 const std::vector<QuantumNumberType> global_nums = string2vecqn(globalquantumnumbers);
559
560 // Local numbers
561 const std::vector<QuantumNumberType> local_nums = string2vecqn(localquantumnumbers);
562
563 // LBLRTM data
564 ifstream is;
565 open_input_file(is, lblrtm_file);
566
567 std::vector<Absorption::SingleLineExternal> v(0);
568
569 bool go_on = true;
570 while (go_on) {
571 v.push_back(Absorption::ReadFromLBLRTMStream(is));
572
573 if (v.back().bad) {
574 v.pop_back();
575 go_on = false;
576 } else if (v.back().line.F0() < fmin) {
577 v.pop_back();
578 } else if (v.back().line.F0() > fmax) {
579 v.pop_back();
580 go_on = false;
581 }
582 }
583
584 for (auto& x: v)
585 x.line.Zeeman() = Zeeman::GetAdvancedModel(x.quantumidentity);
586
587 auto x = Absorption::split_list_of_external_lines(v, local_nums, global_nums);
588 abs_lines.resize(0);
589 abs_lines.reserve(x.size());
590 while (x.size()) {
591 abs_lines.push_back(x.back());
592 abs_lines.back().sort_by_frequency();
593 x.pop_back();
594 }
595
596 abs_linesSetNormalization(abs_lines, normalization_option, verbosity);
597 abs_linesSetMirroring(abs_lines, mirroring_option, verbosity);
598 abs_linesSetPopulation(abs_lines, population_option, verbosity);
599 abs_linesSetLineShapeType(abs_lines, lineshapetype_option, verbosity);
600 abs_linesSetCutoff(abs_lines, cutoff_option, cutoff_value, verbosity);
601 abs_linesSetLinemixingLimit(abs_lines, linemixinglimit_value, verbosity);
602}
603
604/* Workspace method: Doxygen documentation will be auto-generated */
606 const String& jpl_file,
607 const Numeric& fmin,
608 const Numeric& fmax,
609 const String& globalquantumnumbers,
610 const String& localquantumnumbers,
611 const String& normalization_option,
612 const String& mirroring_option,
613 const String& population_option,
614 const String& lineshapetype_option,
615 const String& cutoff_option,
616 const Numeric& cutoff_value,
617 const Numeric& linemixinglimit_value,
618 const Verbosity& verbosity)
619{
620 // Global numbers
621 const std::vector<QuantumNumberType> global_nums = string2vecqn(globalquantumnumbers);
622
623 // Local numbers
624 const std::vector<QuantumNumberType> local_nums = string2vecqn(localquantumnumbers);
625
626 // LBLRTM data
627 ifstream is;
628 open_input_file(is, jpl_file);
629
630 std::vector<Absorption::SingleLineExternal> v(0);
631
632 bool go_on = true;
633 while (go_on) {
634 v.push_back(Absorption::ReadFromJplStream(is));
635
636 if (v.back().bad) {
637 v.pop_back();
638 go_on = false;
639 } else if (v.back().line.F0() < fmin) {
640 v.pop_back();
641 } else if (v.back().line.F0() > fmax) {
642 v.pop_back();
643 go_on = false;
644 }
645 }
646
647 for (auto& x: v)
648 x.line.Zeeman() = Zeeman::GetAdvancedModel(x.quantumidentity);
649
650 auto x = Absorption::split_list_of_external_lines(v, local_nums, global_nums);
651 abs_lines.resize(0);
652 abs_lines.reserve(x.size());
653 while (x.size()) {
654 abs_lines.push_back(x.back());
655 abs_lines.back().sort_by_frequency();
656 x.pop_back();
657 }
658
659 abs_linesSetNormalization(abs_lines, normalization_option, verbosity);
660 abs_linesSetMirroring(abs_lines, mirroring_option, verbosity);
661 abs_linesSetPopulation(abs_lines, population_option, verbosity);
662 abs_linesSetLineShapeType(abs_lines, lineshapetype_option, verbosity);
663 abs_linesSetCutoff(abs_lines, cutoff_option, cutoff_value, verbosity);
664 abs_linesSetLinemixingLimit(abs_lines, linemixinglimit_value, verbosity);
665}
666
670
671/* Workspace method: Doxygen documentation will be auto-generated */
672void abs_linesWriteSplitXML(const String& output_format,
673 const ArrayOfAbsorptionLines& abs_lines,
674 const String& basename,
675 const Verbosity& verbosity)
676{
677 std::map<String, int> names;
678
679 String true_basename = basename;
680 if (not(true_basename.back() == '.' or true_basename.back() == '/'))
681 true_basename += '.';
682
683 for (auto& lines : abs_lines) {
684 auto name = lines.SpeciesName();
685 const String fname = true_basename + name;
686
687 WriteXML(output_format, lines,
688 fname + '.' + std::to_string(names[name]++) + ".xml",
689 0, "", "", "", verbosity);
690 }
691}
692
693/* Workspace method: Doxygen documentation will be auto-generated */
694void abs_linesWriteSpeciesSplitXML(const String& output_format,
695 const ArrayOfAbsorptionLines& abs_lines,
696 const String& basename,
697 const Verbosity& verbosity)
698{
699 // Set the true name of the saving
700 String true_basename = basename;
701 if (not(true_basename.back() == '.' or true_basename.back() == '/'))
702 true_basename += '.';
703
704 // Find all species
705 ArrayOfString specs(0);
706 for (auto& band: abs_lines) {
707 auto specname = band.SpeciesName();
708
709 bool any = false;
710 for (auto& thisname: specs) {
711 if (any) break;
712 if (thisname == specname) any = true;
713 }
714
715 if (not any)
716 specs.push_back(specname);
717 }
718
719 // Make all species into a species tag array
720 Index throwaway;
722 abs_speciesSet(as, throwaway, throwaway, specs, verbosity);
723
724 // Split lines by species
726 abs_lines_per_speciesCreateFromLines(alps, abs_lines, as, verbosity);
727
728 // Save the arrays
729 for (Index i=0; i<specs.nelem(); i++) {
730 auto& name = specs[i];
731 auto& lines = alps[i];
732
733 WriteXML(output_format, lines,
734 true_basename + name + ".xml",
735 0, "", "", "", verbosity);
736 }
737}
738
739/* Workspace method: Doxygen documentation will be auto-generated */
741 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
742 const String& basename,
743 const Verbosity& verbosity)
744{
745 std::map<String, int> names;
746
747 String true_basename = basename;
748 if (not(true_basename.back() == '.' or true_basename.back() == '/'))
749 true_basename += '.';
750
751 for (auto& spec_lines : abs_lines_per_species) {
752 for (auto& lines : spec_lines) {
753 auto name = lines.SpeciesName();
754 const String fname = true_basename + name;
755
756 WriteXML(output_format, lines,
757 fname + '.' + std::to_string(names[name]++) + ".xml",
758 0, "", "", "", verbosity);
759 }
760 }
761}
762
763/* Workspace method: Doxygen documentation will be auto-generated */
765 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
766 const String& basename,
767 const Verbosity& verbosity)
768{
769 // Compact to abs_lines
770 ArrayOfAbsorptionLines abs_lines(0);
771 for (auto& lines: abs_lines_per_species) {
772 for (auto& band: lines) {
773 abs_lines.push_back(band);
774 }
775 }
776
777 // Save using the other function
778 abs_linesWriteSpeciesSplitXML(output_format, abs_lines, basename, verbosity);
779}
780
781/* Workspace method: Doxygen documentation will be auto-generated */
783 const ArrayOfArrayOfSpeciesTag& abs_species,
784 const String& basename,
785 const Verbosity& verbosity)
786{
787 // Build a set of species indices. Duplicates are ignored.
788 const std::set<Species::Species> unique_species = lbl_species(abs_species);
789
790 String tmpbasename = basename;
791 if (basename.length() && basename[basename.length() - 1] != '/') {
792 tmpbasename += '.';
793 }
794
795 // Read catalogs for each identified species and put them all into
796 // abs_lines
797 ArrayOfAbsorptionLines abs_lines(0);
798 for (auto& spec: unique_species) {
800 for (auto& ir: specs) {
801 Index i=0;
802 do {
803 String filename = tmpbasename + ir.FullName() + '.' + std::to_string(i) + ".xml";
804 if (find_xml_file_existence(filename)) {
805 abs_lines.push_back(AbsorptionLines());
806 xml_read_from_file(filename, abs_lines.back(), verbosity);
807 i++;
808 } else {
809 break;
810 }
811 } while (true);
812 }
813 }
814
815 abs_lines_per_speciesCreateFromLines(abs_lines_per_species, abs_lines, abs_species, verbosity);
816}
817
818/* Workspace method: Doxygen documentation will be auto-generated */
820 const String& basename,
821 const Index& robust,
822 const Verbosity& verbosity)
823{
825 std::size_t bands_found{0};
826
827 String tmpbasename = basename;
828 if (basename.length() && basename[basename.length() - 1] != '/') {
829 tmpbasename += '.';
830 }
831
832 // Read catalogs for each identified species and put them all into
833 // abs_lines
834 abs_lines.resize(0);
835 for (auto& ir: Species::Isotopologues) {
836 String filename = tmpbasename + ir.FullName() + ".xml";
837 if (find_xml_file_existence(filename)) {
838 ArrayOfAbsorptionLines speclines;
839 xml_read_from_file(filename, speclines, verbosity);
840 for (auto& band: speclines) {
841 abs_lines.push_back(band);
842 bands_found++;
843 }
844 }
845 }
846
847 ARTS_USER_ERROR_IF (not bands_found and not robust,
848 "Cannot find any bands in the directory you are reading");
849 out3 << "Found " << bands_found << " bands\n";
850}
851
852/* Workspace method: Doxygen documentation will be auto-generated */
854 const ArrayOfArrayOfSpeciesTag& abs_species,
855 const String& basename,
856 const Index& robust,
857 const Verbosity& verbosity)
858{
860 std::size_t bands_found{0};
861
862 // Build a set of species indices. Duplicates are ignored.
863 const std::set<Species::Species> unique_species = lbl_species(abs_species);
864
865 String tmpbasename = basename;
866 if (basename.length() && basename[basename.length() - 1] != '/') {
867 tmpbasename += '.';
868 }
869
870 // Read catalogs for each identified species and put them all into
871 // abs_lines
872 ArrayOfAbsorptionLines abs_lines(0);
873 for (auto spec: unique_species) {
874 auto isots = Species::isotopologues(spec);
875 for (auto& isot: isots) {
876 String filename = tmpbasename + isot.FullName() + ".xml";
877 if (find_xml_file_existence(filename)) {
878 ArrayOfAbsorptionLines speclines;
879 xml_read_from_file(filename, speclines, verbosity);
880 for (auto& band: speclines) {
881 abs_lines.push_back(band);
882 bands_found++;
883 }
884 }
885 }
886 }
887
888 ARTS_USER_ERROR_IF (not bands_found and not robust,
889 "Cannot find any bands in the directory you are reading");
890 out3 << "Found " << bands_found << " bands\n";
891
892 abs_lines_per_speciesCreateFromLines(abs_lines_per_species, abs_lines, abs_species, verbosity);
893}
894
898
899/* Workspace method: Doxygen documentation will be auto-generated */
901 const String& qn,
902 const Rational& x,
903 const QuantumIdentifier& QI,
904 const Verbosity&)
905{
906 auto QN = string2quantumnumbertype(qn);
907 ARTS_USER_ERROR_IF (QN == QuantumNumberType::FINAL,
908 "Usupported quantum number key: ", qn, '\n')
909
910 for (auto& band: abs_lines) {
911 for (Index k=0; k<band.NumLines(); k++) {
912 const Absorption::QuantumIdentifierLineTarget lt = Absorption::QuantumIdentifierLineTarget(QI, band, k);
913 if (lt == Absorption::QuantumIdentifierLineTargetType::Level and lt.lower) {
914 band.LowerQuantumNumber(k, QN) = x;
915 } else if (lt == Absorption::QuantumIdentifierLineTargetType::Level and lt.upper) {
916 band.UpperQuantumNumber(k, QN) = x;
917 }
918 }
919 }
920}
921
922/* Workspace method: Doxygen documentation will be auto-generated */
924 const String& qn,
925 const Rational& x,
926 const QuantumIdentifier& QI,
927 const Verbosity& v)
928{
929 for (auto& band: abs_lines_per_species)
930 abs_linesSetQuantumNumberForMatch(band, qn, x, QI, v);
931}
932
933/* Workspace method: Doxygen documentation will be auto-generated */
935 const Verbosity&)
936{
938
939 for (auto& lines: abs_lines) {
940 lines.truncate_global_quantum_numbers();
941
942 Index match = -1;
943 for (Index ind=0; ind<x.nelem(); ind++) {
944 if (x[ind].Match(lines)) {
945 match = ind;
946 break;
947 }
948 }
949
950 if (match < 0) {
951 x.push_back(lines);
952 } else {
953 for (auto& line: lines.AllLines()) {
954 x[match].AppendSingleLine(line);
955 }
956 }
957 }
958
959 abs_lines = std::move(x);
960 for (auto& lines: abs_lines) {
961 lines.sort_by_frequency();
962 }
963}
964
965/* Workspace method: Doxygen documentation will be auto-generated */
967 const Verbosity& verbosity) {
968 for (auto& band: abs_lines) {
969 band.truncate_local_quantum_numbers();
970 }
971
972 abs_linesTruncateGlobalQuantumNumbers(abs_lines, verbosity);
973}
974
975/* Workspace method: Doxygen documentation will be auto-generated */
977 const Index& pos,
978 const Verbosity& verbosity) {
979 if (pos == -1) {
980 for (auto& abs_lines: abs_lines_per_species) {
981 abs_linesTruncateQuantumNumbers(abs_lines, verbosity);
982 }
983 } else {
984 ARTS_USER_ERROR_IF(pos < 0 or pos >= abs_lines_per_species.nelem(),
985 "Not a valid position for current line lists")
986 abs_linesTruncateQuantumNumbers(abs_lines_per_species[pos], verbosity);
987 }
988}
989
990/* Workspace method: Doxygen documentation will be auto-generated */
992 const Verbosity&)
993{
994 for (auto& lines: abs_lines) {
995 lines.RemoveUnusedLocalQuantums();
996 }
997}
998
1002
1003/* Workspace method: Doxygen documentation will be auto-generated */
1005{
1006 for (auto& rlines: replacing_lines) {
1007 Index number_of_matching_bands = 0;
1008 for (auto& tlines: abs_lines) {
1009 if (tlines.Match(rlines)) {
1010 number_of_matching_bands++;
1011 for (auto& rline: rlines.AllLines()) {
1012 Index number_of_matching_single_lines = 0;
1013 for (auto& tline: tlines.AllLines()) {
1014 if (tline.SameQuantumNumbers(rline)) {
1015 number_of_matching_single_lines++;
1016 tline = rline;
1017 }
1018 }
1019
1020 ARTS_USER_ERROR_IF (number_of_matching_single_lines not_eq 1,
1021 "Error: Did not match to a single single line. "
1022 "This means the input data has not been understood. "
1023 "This function needs exactly one match.");
1024 }
1025 tlines.sort_by_frequency();
1026 }
1027 }
1028
1029 ARTS_USER_ERROR_IF (number_of_matching_bands not_eq 1,
1030 "Error: Did not match to a single set of absorption lines. "
1031 "This means the input data has not been understood. "
1032 "This function needs exactly one match.");
1033 }
1034}
1035
1036/* Workspace method: Doxygen documentation will be auto-generated */
1037void abs_linesAppendWithLines(ArrayOfAbsorptionLines& abs_lines, const ArrayOfAbsorptionLines& appending_lines, const Index& safe, const Verbosity&)
1038{
1039 if (safe) {
1040 std::vector<AbsorptionLines> addedlines(0);
1041
1042 for (auto& alines: appending_lines) {
1043 Index number_of_matching_bands = 0;
1044 for (auto& tlines: abs_lines) {
1045 if (tlines.Match(alines)) {
1046 number_of_matching_bands++;
1047 for (auto& aline: alines.AllLines()) {
1048 Index number_of_matching_single_lines = 0;
1049 for (auto& tline: tlines.AllLines()) {
1050 if (tline.SameQuantumNumbers(aline)) {
1051 number_of_matching_single_lines++;
1052 }
1053 }
1054 ARTS_USER_ERROR_IF (number_of_matching_single_lines not_eq 0,
1055 "Error: Did match to a single single line. "
1056 "This means the input data has not been understood. "
1057 "This function needs exactly zero matches.");
1058 tlines.AppendSingleLine(aline);
1059 }
1060 tlines.sort_by_frequency();
1061 }
1062 }
1063
1064 if (number_of_matching_bands == 0) addedlines.push_back(alines);
1065 ARTS_USER_ERROR_IF (number_of_matching_bands not_eq 1,
1066 "Error: Did not match to a single set of absorption lines. "
1067 "This means the input data has not been understood. "
1068 "This function needs exactly one or zero matches.");
1069 }
1070
1071 for (auto& lines: addedlines) {
1072 abs_lines.push_back(std::move(lines));
1073 }
1074 } else {
1075 for (auto& band: appending_lines)
1076 abs_lines.push_back(band);
1077 }
1078}
1079
1080/* Workspace method: Doxygen documentation will be auto-generated */
1081void abs_linesDeleteBadF0(ArrayOfAbsorptionLines& abs_lines, const Numeric& f0, const Index& lower, const Verbosity&)
1082{
1083 for (auto& lines: abs_lines) {
1084 std::vector<Index> hits;
1085 for (Index i=0; i<lines.NumLines(); i++) {
1086 if (lower and lines.F0(i) < f0)
1087 hits.push_back(i);
1088 else if (not lower and lines.F0(i) > f0)
1089 hits.push_back(i);
1090 }
1091
1092 // Remove the bad values (sort by descending firs)
1093 std::sort(hits.begin(), hits.end());
1094 while(not hits.empty()) {
1095 lines.RemoveLine(hits.back());
1096 hits.pop_back();
1097 }
1098 }
1099}
1100
1101/* Workspace method: Doxygen documentation will be auto-generated */
1103{
1104 for (auto& dlines: deleting_lines) {
1105 for (auto& tlines: abs_lines) {
1106 std::vector<Index> hits(0);
1107
1108 if (tlines.Match(dlines)) {
1109 for (auto& dline: dlines.AllLines()) {
1110 for (Index i=0; i<tlines.NumLines(); i++) {
1111 if (tlines.AllLines()[i].SameQuantumNumbers(dline)) {
1112 hits.push_back(i);
1113 }
1114 }
1115 }
1116
1117 // Sort and test the input
1118 std::sort(hits.begin(), hits.end());
1119 auto n = hits.size();
1120 hits.erase(std::unique(hits.begin(), hits.end()), hits.end());
1121 ARTS_USER_ERROR_IF(n not_eq hits.size(),
1122 "Removing the same line more than once is not accepted");
1123
1124 // Remove the bad values
1125 while(not hits.empty()) {
1126 tlines.RemoveLine(hits.back());
1127 hits.pop_back();
1128 }
1129 }
1130 }
1131 }
1132}
1133
1134/* Workspace method: Doxygen documentation will be auto-generated */
1136{
1138 Index i = 0;
1139
1140 for (auto& band: abs_lines) {
1141 std::vector<Index> deleters(0);
1142
1143 for (Index iline=0; iline<band.NumLines(); iline++) {
1144 if (std::any_of(band.Line(iline).LowerQuantumNumbers().cbegin(),
1145 band.Line(iline).LowerQuantumNumbers().cend(),
1146 [](auto x) -> bool {return x.isUndefined();}) or
1147 std::any_of(band.Line(iline).UpperQuantumNumbers().cbegin(),
1148 band.Line(iline).UpperQuantumNumbers().cend(),
1149 [](auto x) -> bool {return x.isUndefined();})) {
1150 deleters.push_back(iline);
1151 }
1152 }
1153
1154 while (deleters.size()) {
1155 band.RemoveLine(deleters.back());
1156 deleters.pop_back();
1157 i++;
1158 }
1159 }
1160
1161 out2 << "Deleted " << i << " lines.\n";
1162}
1163
1164/* Workspace method: Doxygen documentation will be auto-generated */
1166{
1168 Index i = 0;
1169
1170 for (auto& band: abs_lines) {
1171 std::vector<Index> deleters(0);
1172
1173 for (Index iline=0; iline<band.NumLines(); iline++) {
1174 auto Jlo = band.LowerQuantumNumber(iline, QuantumNumberType::J);
1175 auto Jup = band.UpperQuantumNumber(iline, QuantumNumberType::J);
1176
1177 if (Jlo.isUndefined() or Jup.isUndefined() or 1 < abs(Jup - Jlo)) {
1178 deleters.push_back(iline);
1179 }
1180 }
1181
1182 while (deleters.size()) {
1183 band.RemoveLine(deleters.back());
1184 deleters.pop_back();
1185 i++;
1186 }
1187 }
1188
1189 out2 << "Deleted " << i << " lines.\n";
1190}
1191
1192/* Workspace method: Doxygen documentation will be auto-generated */
1194{
1195 const auto qn = string2quantumnumbertype(qn_id);
1196
1197 for (auto& band: abs_lines) {
1198 for (Index i=band.NumLines() - 1; i>=0; i--) {
1199 if (band.UpperQuantumNumber(i, qn) > qn_val or band.LowerQuantumNumber(i, qn) > qn_val) {
1200 band.RemoveLine(i);
1201 }
1202 }
1203 }
1204}
1205
1206/* Workspace method: Doxygen documentation will be auto-generated */
1208{
1209 for (auto& band: abs_lines) {
1210 std::array<bool, LineShape::nVars> var_is_empty;
1211
1212 // Species by species can be empty, so loop each species by themselves
1213 for (Index ispec=0; ispec<band.NumBroadeners(); ispec++) {
1214 var_is_empty.fill(true);
1215
1216 // Check if any variable in this band for any line is non-empty
1217 for (Index iline=0; iline<band.NumLines(); iline++) {
1218 for (Index ivar=0; ivar < LineShape::nVars; ivar++) {
1219 if (not LineShape::modelparameterEmpty(band.Line(iline).LineShape().Data()[ispec].Data()[ivar])) {
1220 var_is_empty[ivar] = false;
1221 }
1222 }
1223 }
1224
1225 // Remove empty variables from the writing. This will also speed up some calculations
1226 for (Index iline=0; iline<band.NumLines(); iline++) {
1227 for (Index ivar=0; ivar < LineShape::nVars; ivar++) {
1228 if (var_is_empty[ivar]) {
1229 band.Line(iline).LineShape().Data()[ispec].Data()[ivar].type = LineShape::TemperatureModel::None;
1230 }
1231 }
1232 }
1233 }
1234 }
1235}
1236
1238{
1239 for (auto& band: abs_lines) {
1240 const Absorption::QuantumIdentifierLineTarget lt(qid, band);
1241 while (lt not_eq Absorption::QuantumIdentifierLineTargetType::Band and band.NumLines()) {
1242 band.RemoveLine(0);
1243 }
1244 }
1245}
1246
1250
1254
1255/* Workspace method: Doxygen documentation will be auto-generated */
1257 const String& type,
1258 const Numeric& x,
1259 const Verbosity&)
1260{
1261 auto t = Absorption::toCutoffTypeOrThrow(type);
1262 for (auto& lines: abs_lines) {
1263 lines.Cutoff(t);
1264 lines.CutoffFreqValue(x);
1265 }
1266}
1267
1268/* Workspace method: Doxygen documentation will be auto-generated */
1270 const String& type,
1271 const Numeric& x,
1272 const Verbosity& v)
1273{
1274 for (auto& abs_lines: abs_lines_per_species)
1275 abs_linesSetCutoff(abs_lines, type, x, v);
1276}
1277
1278/* Workspace method: Doxygen documentation will be auto-generated */
1280 ArrayOfAbsorptionLines& abs_lines,
1281 const String& type,
1282 const Numeric& x,
1283 const QuantumIdentifier& QI,
1284 const Verbosity&)
1285{
1286 auto t = Absorption::toCutoffTypeOrThrow(type);
1287 for (auto& band: abs_lines) {
1288 const Absorption::QuantumIdentifierLineTarget lt(QI, band);
1289 if (lt == Absorption::QuantumIdentifierLineTargetType::Band) {
1290 band.Cutoff(t);
1291 band.CutoffFreqValue(x);
1292 }
1293 }
1294}
1295
1296/* Workspace method: Doxygen documentation will be auto-generated */
1298 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1299 const String& type,
1300 const Numeric& x,
1301 const QuantumIdentifier& QI,
1302 const Verbosity& verbosity)
1303{
1304 for (auto& lines: abs_lines_per_species) {
1305 abs_linesSetCutoffForMatch(lines, type, x, QI, verbosity);
1306 }
1307}
1308
1309/* Workspace method: Doxygen documentation will be auto-generated */
1311 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1312 const ArrayOfArrayOfSpeciesTag& abs_species,
1313 const String& type,
1314 const Numeric& x,
1315 const String& species_tag,
1316 const Verbosity& verbosity)
1317{
1318 Index t1, t2;
1319 ArrayOfArrayOfSpeciesTag target_species;
1320 abs_speciesSet(target_species, t1, t2, {species_tag}, verbosity);
1321 for (Index ispec=0; ispec<abs_species.nelem(); ispec++) {
1322 if (std::equal(abs_species[ispec].begin(), abs_species[ispec].end(), target_species[0].begin())) {
1323 abs_linesSetCutoff(abs_lines_per_species[ispec], type, x, verbosity);
1324 }
1325 }
1326}
1327
1331
1332/* Workspace method: Doxygen documentation will be auto-generated */
1334 const String& type,
1335 const Verbosity&)
1336{
1337 auto t = Absorption::toMirroringTypeOrThrow(type);
1338 for (auto& lines: abs_lines)
1339 lines.Mirroring(t);
1340}
1341
1342/* Workspace method: Doxygen documentation will be auto-generated */
1344 const String& type,
1345 const Verbosity& v)
1346{
1347 for (auto& abs_lines: abs_lines_per_species)
1348 abs_linesSetMirroring(abs_lines, type, v);
1349}
1350
1351/* Workspace method: Doxygen documentation will be auto-generated */
1353 const String& type,
1354 const QuantumIdentifier& QI,
1355 const Verbosity&)
1356{
1357 auto t = Absorption::toMirroringTypeOrThrow(type);
1358 for (auto& band: abs_lines) {
1359 const Absorption::QuantumIdentifierLineTarget lt(QI, band);
1360 if (lt == Absorption::QuantumIdentifierLineTargetType::Band) {
1361 band.Mirroring(t);
1362 }
1363 }
1364}
1365
1366/* Workspace method: Doxygen documentation will be auto-generated */
1368 const String& type,
1369 const QuantumIdentifier& QI,
1370 const Verbosity& v)
1371{
1372 for (auto& abs_lines: abs_lines_per_species)
1373 abs_linesSetMirroringForMatch(abs_lines, type, QI, v);
1374}
1375
1376/* Workspace method: Doxygen documentation will be auto-generated */
1378 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1379 const ArrayOfArrayOfSpeciesTag& abs_species,
1380 const String& type,
1381 const String& species_tag,
1382 const Verbosity& verbosity)
1383{
1384 Index t1, t2;
1385 ArrayOfArrayOfSpeciesTag target_species;
1386 abs_speciesSet(target_species, t1, t2, {species_tag}, verbosity);
1387 for (Index ispec=0; ispec<abs_species.nelem(); ispec++) {
1388 if (std::equal(abs_species[ispec].begin(), abs_species[ispec].end(), target_species[0].begin())) {
1389 abs_linesSetMirroring(abs_lines_per_species[ispec], type, verbosity);
1390 }
1391 }
1392}
1393
1394/* Workspace method: Doxygen documentation will be auto-generated */
1396 const Verbosity&)
1397{
1398 const ArrayOfAbsorptionLines abs_lines_copy = abs_lines;
1399 for (AbsorptionLines band: abs_lines_copy) {
1400 band.Mirroring(Absorption::MirroringType::Manual);
1401
1404 std::find_if(abs_lines_copy.cbegin(), abs_lines_copy.cend(),
1405 [&band](const AbsorptionLines& li){
1406 return band.Match(li);
1407 }) not_eq abs_lines_copy.cend(),
1408 "Dual bands with same setup is not allowed for mirroring of band:\n",
1409 band, '\n')
1410
1411 for (auto& line: band.AllLines()) {
1412 line.F0() *= -1;
1413 }
1414
1415 abs_lines.emplace_back(std::move(band));
1416 }
1417}
1418
1419/* Workspace method: Doxygen documentation will be auto-generated */
1421 const Verbosity& verbosity)
1422{
1423 for (auto& abs_lines: abs_lines_per_species) abs_linesMakeManualMirroring(abs_lines, verbosity);
1424}
1425
1426/* Workspace method: Doxygen documentation will be auto-generated */
1428 const ArrayOfArrayOfSpeciesTag& abs_species,
1429 const ArrayOfSpeciesTag& species,
1430 const Verbosity& verbosity)
1431{
1432 ARTS_USER_ERROR_IF(abs_species.size() not_eq abs_lines_per_species.size(),
1433 "Mismatch abs_species and abs_lines_per_species sizes [",
1434 abs_species.size(), " vs ", abs_lines_per_species.size(),
1435 ", respectively]")
1436
1437 if (auto ind = std::distance(abs_species.cbegin(), std::find(abs_species.cbegin(), abs_species.cend(), species)); ind not_eq abs_species.nelem()) {
1438 abs_linesMakeManualMirroring(abs_lines_per_species[ind], verbosity);
1439 } else {
1440 ARTS_USER_ERROR("Cannot find species: ", species, "\nIn abs_species: [", abs_species, ']')
1441 }
1442}
1443
1444
1448
1449/* Workspace method: Doxygen documentation will be auto-generated */
1451 const String& type,
1452 const Verbosity&)
1453{
1454 auto t = Absorption::toPopulationTypeOrThrow(type);
1455 for (auto& lines: abs_lines)
1456 lines.Population(t);
1457}
1458
1459/* Workspace method: Doxygen documentation will be auto-generated */
1461 const String& type,
1462 const Verbosity& v)
1463{
1464 for (auto& abs_lines: abs_lines_per_species)
1465 abs_linesSetPopulation(abs_lines, type, v);
1466}
1467
1468/* Workspace method: Doxygen documentation will be auto-generated */
1470 const String& type,
1471 const QuantumIdentifier& QI,
1472 const Verbosity&)
1473{
1474 auto t = Absorption::toPopulationTypeOrThrow(type);
1475 for (auto& lines: abs_lines) {
1476 const Absorption::QuantumIdentifierLineTarget lt(QI, lines);
1477 if (lt == Absorption::QuantumIdentifierLineTargetType::Band) {
1478 lines.Population(t);
1479 }
1480 }
1481}
1482
1483/* Workspace method: Doxygen documentation will be auto-generated */
1485 const String& type,
1486 const QuantumIdentifier& QI,
1487 const Verbosity& v)
1488{
1489 for (auto& abs_lines: abs_lines_per_species)
1490 abs_linesSetPopulationForMatch(abs_lines, type, QI, v);
1491}
1492
1493/* Workspace method: Doxygen documentation will be auto-generated */
1495 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1496 const ArrayOfArrayOfSpeciesTag& abs_species,
1497 const String& type,
1498 const String& species_tag,
1499 const Verbosity& verbosity)
1500{
1501 Index t1, t2;
1502 ArrayOfArrayOfSpeciesTag target_species;
1503 abs_speciesSet(target_species, t1, t2, {species_tag}, verbosity);
1504 for (Index ispec=0; ispec<abs_species.nelem(); ispec++) {
1505 if (std::equal(abs_species[ispec].begin(), abs_species[ispec].end(), target_species[0].begin())) {
1506 abs_linesSetPopulation(abs_lines_per_species[ispec], type, verbosity);
1507 }
1508 }
1509}
1510
1514
1515/* Workspace method: Doxygen documentation will be auto-generated */
1517 const String& type,
1518 const Verbosity&)
1519{
1520 auto t = Absorption::toNormalizationTypeOrThrow(type);
1521 for (auto& lines: abs_lines)
1522 lines.Normalization(t);
1523}
1524
1525/* Workspace method: Doxygen documentation will be auto-generated */
1527 const String& type,
1528 const Verbosity& v)
1529{
1530 for (auto& abs_lines: abs_lines_per_species)
1531 abs_linesSetNormalization(abs_lines, type, v);
1532}
1533
1534/* Workspace method: Doxygen documentation will be auto-generated */
1536 const String& type,
1537 const QuantumIdentifier& QI,
1538 const Verbosity&)
1539{
1540 auto t = Absorption::toNormalizationTypeOrThrow(type);
1541 for (auto& lines: abs_lines) {
1542 const Absorption::QuantumIdentifierLineTarget lt(QI, lines);
1543 if (lt == Absorption::QuantumIdentifierLineTargetType::Band) {
1544 lines.Normalization(t);
1545 }
1546 }
1547}
1548
1549/* Workspace method: Doxygen documentation will be auto-generated */
1551 const String& type,
1552 const QuantumIdentifier& QI,
1553 const Verbosity& v)
1554{
1555 for (auto& abs_lines: abs_lines_per_species)
1556 abs_linesSetNormalizationForMatch(abs_lines, type, QI, v);
1557}
1558
1559/* Workspace method: Doxygen documentation will be auto-generated */
1561 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1562 const ArrayOfArrayOfSpeciesTag& abs_species,
1563 const String& type,
1564 const String& species_tag,
1565 const Verbosity& verbosity)
1566{
1567 Index t1, t2;
1568 ArrayOfArrayOfSpeciesTag target_species;
1569 abs_speciesSet(target_species, t1, t2, {species_tag}, verbosity);
1570 for (Index ispec=0; ispec<abs_species.nelem(); ispec++) {
1571 if (std::equal(abs_species[ispec].begin(), abs_species[ispec].end(), target_species[0].begin())) {
1572 abs_linesSetNormalization(abs_lines_per_species[ispec], type, verbosity);
1573 }
1574 }
1575}
1576
1580
1581/* Workspace method: Doxygen documentation will be auto-generated */
1583 const String& type,
1584 const Verbosity&)
1585{
1586 auto t = LineShape::toType(type);
1587 check_enum_error(t, "Cannot understand type: ", type);
1588 for (auto& lines: abs_lines)
1589 lines.LineShapeType(t);
1590}
1591
1592/* Workspace method: Doxygen documentation will be auto-generated */
1594 const String& type,
1595 const Verbosity& v)
1596{
1597 for (auto& abs_lines: abs_lines_per_species)
1598 abs_linesSetLineShapeType(abs_lines, type, v);
1599}
1600
1601/* Workspace method: Doxygen documentation will be auto-generated */
1603 const String& type,
1604 const QuantumIdentifier& QI,
1605 const Verbosity&)
1606{
1607 auto t = LineShape::toTypeOrThrow(type);
1608 for (auto& lines: abs_lines) {
1609 const Absorption::QuantumIdentifierLineTarget lt(QI, lines);
1610 if (lt == Absorption::QuantumIdentifierLineTargetType::Band) {
1611 lines.LineShapeType(t);
1612 }
1613 }
1614}
1615
1616/* Workspace method: Doxygen documentation will be auto-generated */
1618 const String& type,
1619 const QuantumIdentifier& QI,
1620 const Verbosity& v)
1621{
1622 for (auto& abs_lines: abs_lines_per_species)
1623 abs_linesSetLineShapeTypeForMatch(abs_lines, type, QI, v);
1624}
1625
1626/* Workspace method: Doxygen documentation will be auto-generated */
1628 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1629 const ArrayOfArrayOfSpeciesTag& abs_species,
1630 const String& type,
1631 const String& species_tag,
1632 const Verbosity& verbosity)
1633{
1634 Index t1, t2;
1635 ArrayOfArrayOfSpeciesTag target_species;
1636 abs_speciesSet(target_species, t1, t2, {species_tag}, verbosity);
1637 for (Index ispec=0; ispec<abs_species.nelem(); ispec++) {
1638 if (std::equal(abs_species[ispec].begin(), abs_species[ispec].end(), target_species[0].begin())) {
1639 abs_linesSetLineShapeType(abs_lines_per_species[ispec], type, verbosity);
1640 }
1641 }
1642}
1643
1647
1648/* Workspace method: Doxygen documentation will be auto-generated */
1650 const Numeric& x,
1651 const Verbosity&)
1652{
1653 for (auto& lines: abs_lines)
1654 lines.LinemixingLimit(x);
1655}
1656
1657/* Workspace method: Doxygen documentation will be auto-generated */
1659 const Numeric& x,
1660 const Verbosity& v)
1661{
1662 for (auto& abs_lines: abs_lines_per_species)
1663 abs_linesSetLinemixingLimit(abs_lines, x, v);
1664}
1665
1666/* Workspace method: Doxygen documentation will be auto-generated */
1668 const Numeric& x,
1669 const QuantumIdentifier& QI,
1670 const Verbosity&)
1671{
1672 for (auto& lines: abs_lines) {
1673 const Absorption::QuantumIdentifierLineTarget lt(QI, lines);
1674 if (lt == Absorption::QuantumIdentifierLineTargetType::Band) {
1675 lines.LinemixingLimit(x);
1676 }
1677 }
1678}
1679
1680/* Workspace method: Doxygen documentation will be auto-generated */
1682 const Numeric& x,
1683 const QuantumIdentifier& QI,
1684 const Verbosity& v)
1685{
1686 for (auto& abs_lines: abs_lines_per_species)
1687 abs_linesSetLinemixingLimitForMatch(abs_lines, x, QI, v);
1688}
1689
1690/* Workspace method: Doxygen documentation will be auto-generated */
1692 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1693 const ArrayOfArrayOfSpeciesTag& abs_species,
1694 const Numeric& x,
1695 const String& species_tag,
1696 const Verbosity& verbosity)
1697{
1698 Index t1, t2;
1699 ArrayOfArrayOfSpeciesTag target_species;
1700 abs_speciesSet(target_species, t1, t2, {species_tag}, verbosity);
1701 for (Index ispec=0; ispec<abs_species.nelem(); ispec++) {
1702 if (std::equal(abs_species[ispec].begin(), abs_species[ispec].end(), target_species[0].begin())) {
1703 abs_linesSetLinemixingLimit(abs_lines_per_species[ispec], x, verbosity);
1704 }
1705 }
1706}
1707
1711
1712/* Workspace method: Doxygen documentation will be auto-generated */
1714 const Numeric& x,
1715 const Verbosity&)
1716{
1717 for (auto& lines: abs_lines)
1718 lines.T0(x);
1719}
1720
1721/* Workspace method: Doxygen documentation will be auto-generated */
1723 const Numeric& x,
1724 const Verbosity& v)
1725{
1726 for (auto& abs_lines: abs_lines_per_species)
1727 abs_linesSetT0(abs_lines, x, v);
1728}
1729
1730/* Workspace method: Doxygen documentation will be auto-generated */
1732 const Numeric& x,
1733 const QuantumIdentifier& QI,
1734 const Verbosity&)
1735{
1736 for (auto& lines: abs_lines) {
1737 const Absorption::QuantumIdentifierLineTarget lt(QI, lines);
1738 if (lt == Absorption::QuantumIdentifierLineTargetType::Band) {
1739 lines.T0(x);
1740 }
1741 }
1742}
1743
1744/* Workspace method: Doxygen documentation will be auto-generated */
1746 const Numeric& x,
1747 const QuantumIdentifier& QI,
1748 const Verbosity& v)
1749{
1750 for (auto& abs_lines: abs_lines_per_species)
1751 abs_linesSetT0ForMatch(abs_lines, x, QI, v);
1752}
1753
1754/* Workspace method: Doxygen documentation will be auto-generated */
1756 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1757 const ArrayOfArrayOfSpeciesTag& abs_species,
1758 const Numeric& x,
1759 const String& species_tag,
1760 const Verbosity& verbosity)
1761{
1762 Index t1, t2;
1763 ArrayOfArrayOfSpeciesTag target_species;
1764 abs_speciesSet(target_species, t1, t2, {species_tag}, verbosity);
1765 for (Index ispec=0; ispec<abs_species.nelem(); ispec++) {
1766 if (std::equal(abs_species[ispec].begin(), abs_species[ispec].end(), target_species[0].begin())) {
1767 abs_linesSetT0(abs_lines_per_species[ispec], x, verbosity);
1768 }
1769 }
1770}
1771
1775
1776/* Workspace method: Doxygen documentation will be auto-generated */
1778 const QuantumIdentifier& QI,
1779 const String& parameter_name,
1780 const Numeric& change,
1781 const Index& relative,
1782 const Index& loose_matching,
1783 const Verbosity&)
1784{
1785 Index parameter_switch = -1;
1786
1787 ARTS_USER_ERROR_IF (parameter_name.nelem() == 0,
1788 "parameter_name is empty.\n");
1789
1790 if (parameter_name == "Central Frequency" or
1791 parameter_name == "Line Center")
1792 parameter_switch = 0;
1793 else if (parameter_name == "Line Strength")
1794 parameter_switch = 1;
1795 else if (parameter_name == "Lower State Energy")
1796 parameter_switch = 4;
1797 else if (parameter_name == "Einstein Coefficient")
1798 parameter_switch = 5;
1799 else if (parameter_name == "Lower Statistical Weight")
1800 parameter_switch = 6;
1801 else if (parameter_name == "Upper Statistical Weight")
1802 parameter_switch = 7;
1803 else if (parameter_name == "Lower Zeeman Coefficient")
1804 parameter_switch = 8;
1805 else if (parameter_name == "Upper Zeeman Coefficient")
1806 parameter_switch = 9;
1807
1808 for (auto& band: abs_lines) {
1809 for (Index k=0; k<band.NumLines(); k++) {
1810 const Absorption::QuantumIdentifierLineTarget lt = Absorption::QuantumIdentifierLineTarget(QI, band, k);
1811 if (loose_matching ? lt == Absorption::QuantumIdentifierLineTargetType::Band
1812 : lt == Absorption::QuantumIdentifierLineTargetType::Line) {
1813 switch (parameter_switch) {
1814 case 0: // "Central Frequency":
1815 if (relative == 0)
1816 band.F0(k) += change;
1817 else
1818 band.F0(k) *= 1.0e0 + change;
1819 break;
1820 case 1: // "Line Strength":
1821 if (relative == 0)
1822 band.I0(k) += change;
1823 else
1824 band.I0(k) *= 1.0e0 + change;
1825 break;
1826 case 4: // "Lower State Energy":
1827 if (relative == 0)
1828 band.E0(k) += change;
1829 else
1830 band.E0(k) *= 1.0e0 + change;
1831 break;
1832 case 5: // "Einstein":
1833 if (relative == 0)
1834 band.A(k) += change;
1835 else
1836 band.A(k) *= 1.0e0 + change;
1837 break;
1838 case 6: // "Lower Statistical Weight":
1839 if (relative == 0)
1840 band.g_low(k) += change;
1841 else
1842 band.g_low(k) *= 1.0e0 + change;
1843 break;
1844 case 7: // "Upper Statistical Weight":
1845 if (relative == 0)
1846 band.g_upp(k) += change;
1847 else
1848 band.g_upp(k) *= 1.0e0 + change;
1849 break;
1850 case 8: // "Lower Zeeman Coefficient":
1851 if (relative == 0)
1852 band.Line(k).Zeeman().gl() += change;
1853 else
1854 band.Line(k).Zeeman().gl() *= 1.0e0 + change;
1855 break;
1856 case 9: // "Upper Zeeman Coefficient":
1857 if (relative == 0)
1858 band.Line(k).Zeeman().gu() += change;
1859 else
1860 band.Line(k).Zeeman().gu() *= 1.0e0 + change;
1861 break;
1862 default: {
1864 "Usupported paramter_name\n", parameter_name,
1865 "\nSee method description for supported parameter names.\n")
1866 }
1867 }
1868 }
1869 }
1870 }
1871}
1872
1873/* Workspace method: Doxygen documentation will be auto-generated */
1875 const QuantumIdentifier& QI,
1876 const String& parameter_name,
1877 const Numeric& change,
1878 const Index& relative,
1879 const Index& loose_matching,
1880 const Verbosity& verbosity)
1881{
1882 for (auto& lines: abs_lines_per_species)
1883 abs_linesChangeBaseParameterForMatchingLines(lines, QI, parameter_name, change, relative, loose_matching, verbosity);
1884}
1885
1886/* Workspace method: Doxygen documentation will be auto-generated */
1888 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1889 const ArrayOfArrayOfSpeciesTag& abs_species,
1890 const QuantumIdentifier& QI,
1891 const String& parameter_name,
1892 const Numeric& change,
1893 const Index& relative,
1894 const Index& loose_matching,
1895 const String& species_tag,
1896 const Verbosity& verbosity)
1897{
1898 Index t1, t2;
1899 ArrayOfArrayOfSpeciesTag target_species;
1900 abs_speciesSet(target_species, t1, t2, {species_tag}, verbosity);
1901 for (Index ispec=0; ispec<abs_species.nelem(); ispec++) {
1902 if (std::equal(abs_species[ispec].begin(), abs_species[ispec].end(), target_species[0].begin())) {
1903 abs_linesChangeBaseParameterForMatchingLines(abs_lines_per_species[ispec], QI, parameter_name, change, relative, loose_matching, verbosity);
1904 }
1905 }
1906}
1907
1908/* Workspace method: Doxygen documentation will be auto-generated */
1910 const QuantumIdentifier& QI,
1911 const String& parameter_name,
1912 const Numeric& x,
1913 const Index& loose_matching,
1914 const Verbosity&)
1915{
1916 Index parameter_switch = -1;
1917
1918 ARTS_USER_ERROR_IF (parameter_name.nelem() == 0,
1919 "parameter_name is empty.\n");
1920
1921 if (parameter_name == "Central Frequency" or
1922 parameter_name == "Line Center")
1923 parameter_switch = 0;
1924 else if (parameter_name == "Line Strength")
1925 parameter_switch = 1;
1926 else if (parameter_name == "Lower State Energy")
1927 parameter_switch = 4;
1928 else if (parameter_name == "Einstein Coefficient")
1929 parameter_switch = 5;
1930 else if (parameter_name == "Lower Statistical Weight")
1931 parameter_switch = 6;
1932 else if (parameter_name == "Upper Statistical Weight")
1933 parameter_switch = 7;
1934 else if (parameter_name == "Lower Zeeman Coefficient")
1935 parameter_switch = 8;
1936 else if (parameter_name == "Upper Zeeman Coefficient")
1937 parameter_switch = 9;
1938
1939 for (auto& band: abs_lines) {
1940 for (Index k=0; k<band.NumLines(); k++) {
1941 const Absorption::QuantumIdentifierLineTarget lt(QI, band, k);
1942 if (loose_matching ? lt == Absorption::QuantumIdentifierLineTargetType::Line
1943 : lt == Absorption::QuantumIdentifierLineTargetType::Band) {
1944 switch (parameter_switch) {
1945 case 0: // "Central Frequency":
1946 band.F0(k) = x;
1947 break;
1948 case 1: // "Line Strength":
1949 band.I0(k) = x;
1950 break;
1951 case 4: // "Lower State Energy":
1952 band.E0(k) = x;
1953 break;
1954 case 5: // "Einstein":
1955 band.A(k) = x;
1956 break;
1957 case 6: // "Lower Statistical Weight":
1958 band.g_low(k) = x;
1959 break;
1960 case 7: // "Upper Statistical Weight":
1961 band.g_upp(k) = x;
1962 break;
1963 case 8:
1964 band.Line(k).Zeeman().gl() = x;
1965 break;
1966 case 9:
1967 band.Line(k).Zeeman().gu() = x;
1968 break;
1969 default: {
1971 "Usupported paramter_name\n", parameter_name,
1972 "\nSee method description for supported parameter names.\n")
1973 }
1974 }
1975 }
1976 }
1977 }
1978}
1979
1980/* Workspace method: Doxygen documentation will be auto-generated */
1982 const QuantumIdentifier& QI,
1983 const String& parameter_name,
1984 const Numeric& change,
1985 const Index& loose_matching,
1986 const Verbosity& verbosity)
1987{
1988 for (auto& lines: abs_lines_per_species)
1989 abs_linesSetBaseParameterForMatchingLines(lines, QI, parameter_name, change, loose_matching, verbosity);
1990}
1991
1992/* Workspace method: Doxygen documentation will be auto-generated */
1994 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
1995 const ArrayOfArrayOfSpeciesTag& abs_species,
1996 const QuantumIdentifier& QI,
1997 const String& parameter_name,
1998 const Numeric& change,
1999 const Index& loose_matching,
2000 const String& species_tag,
2001 const Verbosity& verbosity)
2002{
2003 Index t1, t2;
2004 ArrayOfArrayOfSpeciesTag target_species;
2005 abs_speciesSet(target_species, t1, t2, {species_tag}, verbosity);
2006 for (Index ispec=0; ispec<abs_species.nelem(); ispec++) {
2007 if (std::equal(abs_species[ispec].begin(), abs_species[ispec].end(), target_species[0].begin())) {
2008 abs_linesSetBaseParameterForMatchingLines(abs_lines_per_species[ispec], QI, parameter_name, change, loose_matching, verbosity);
2009 }
2010 }
2011}
2012
2013/* Workspace method: Doxygen documentation will be auto-generated */
2015 ArrayOfAbsorptionLines& abs_lines,
2016 const QuantumIdentifier& QI,
2017 const String& parameter,
2018 const String& species,
2019 const String& temperaturemodel,
2020 const Vector& new_values,
2021 const Verbosity& verbosity)
2022{
2024
2025 const bool do_self = species == LineShape::self_broadening;
2026 const bool do_bath = species == LineShape::bath_broadening;
2027
2028 // Set the spec index if possible
2029 const Species::Species spec = do_self ? Species::Species::FINAL : do_bath ? Species::Species::Bath : SpeciesTag(species).Spec();
2030
2031 const LineShape::Variable var = LineShape::toVariableOrThrow(parameter);
2032
2034 "Mismatch between input and expected number of variables\n"
2035 "\tInput is: ", new_values.nelem(), " long but expects: ", LineShape::ModelParameters::N, " values\n")
2036
2038 newdata.type = LineShape::toTemperatureModelOrThrow(temperaturemodel);
2039 newdata.X0 = new_values[0];
2040 newdata.X1 = new_values[1];
2041 newdata.X2 = new_values[2];
2042 newdata.X3 = new_values[3];
2043
2044 for (auto& band: abs_lines) {
2045 for (Index k=0; k<band.NumLines(); k++) {
2046 const Absorption::QuantumIdentifierLineTarget lt = Absorption::QuantumIdentifierLineTarget(QI, band, k);
2047 if (lt == Absorption::QuantumIdentifierLineTargetType::Line) {
2048 if (do_self and band.Self()) {
2049 out3 << "Changing self\n";
2050 band.Line(k).LineShape().Data().front().Data()[Index(var)] = newdata;
2051 } else if (do_bath and band.Bath()) {
2052 out3 << "Changing bath\n";
2053 band.Line(k).LineShape().Data().back().Data()[Index(var)] = newdata;
2054 } else {
2055 for (Index i=band.Self(); i<band.BroadeningSpecies().nelem()-band.Bath(); i++) {
2056 if (spec == band.BroadeningSpecies()[i]) {
2057 out3 << "Changing species: " << Species::toShortName(spec) << '\n';
2058 band.Line(k).LineShape().Data()[i].Data()[Index(var)] = newdata;
2059 }
2060 }
2061 }
2062 }
2063 }
2064 }
2065}
2066
2067/* Workspace method: Doxygen documentation will be auto-generated */
2069 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
2070 const QuantumIdentifier& QI,
2071 const String& parameter,
2072 const String& species,
2073 const String& temperaturemodel,
2074 const Vector& new_values,
2075 const Verbosity& verbosity)
2076{
2077 for (auto& lines: abs_lines_per_species)
2079 lines, QI, parameter, species, temperaturemodel, new_values, verbosity);
2080}
2081
2085
2086/* Workspace method: Doxygen documentation will be auto-generated */
2088 const QuantumIdentifier& QI,
2089 const String& parameter_name,
2090 const Numeric& change,
2091 const Index& relative,
2092 const Verbosity&)
2093{
2094 ARTS_USER_ERROR_IF (QI.type not_eq Quantum::IdentifierType::EnergyLevel,
2095 "Bad input. Must be energy level. Is: ", QI, '\n')
2096
2097 Index parameter_switch = -1;
2098
2099 ARTS_USER_ERROR_IF (parameter_name.nelem() == 0,
2100 "parameter_name is empty.\n");
2101 if (parameter_name == "Statistical Weight")
2102 parameter_switch = 1;
2103 else if (parameter_name == "Zeeman Coefficient")
2104 parameter_switch = 2;
2105
2106 for (auto& band: abs_lines) {
2107 for (Index k=0; k<band.NumLines(); k++) {
2108 const Absorption::QuantumIdentifierLineTarget lt = Absorption::QuantumIdentifierLineTarget(QI, band, k);
2109 if (lt == Absorption::QuantumIdentifierLineTargetType::Level and lt.lower) {
2110 switch (parameter_switch) {
2111 case 1: // "Statistical Weight":
2112 if (relative == 0)
2113 band.g_low(k) += change;
2114 else
2115 band.g_low(k) *= 1.0e0 + change;
2116 break;
2117 case 2: // "Zeeman Coefficient":
2118 if (relative == 0)
2119 band.Line(k).Zeeman().gl() += change;
2120 else
2121 band.Line(k).Zeeman().gl() *= 1.0e0 + change;
2122 break;
2123 default: {
2125 "Usupported paramter_name\n", parameter_name,
2126 "\nSee method description for supported parameter names.\n")
2127 }
2128 }
2129 }
2130
2131 if (lt == Absorption::QuantumIdentifierLineTargetType::Level and lt.upper) {
2132 switch (parameter_switch) {
2133 case 1: // "Statistical Weight":
2134 if (relative == 0)
2135 band.g_upp(k) += change;
2136 else
2137 band.g_upp(k) *= 1.0e0 + change;
2138 break;
2139 case 2: // "Zeeman Coefficient":
2140 if (relative == 0)
2141 band.Line(k).Zeeman().gu() += change;
2142 else
2143 band.Line(k).Zeeman().gu() *= 1.0e0 + change;
2144 break;
2145 default: {
2147 "Usupported paramter_name\n", parameter_name,
2148 "\nSee method description for supported parameter names.\n")
2149 }
2150 }
2151 }
2152 }
2153 }
2154}
2155
2156/* Workspace method: Doxygen documentation will be auto-generated */
2158 const QuantumIdentifier& QI,
2159 const String& parameter_name,
2160 const Numeric& change,
2161 const Index& relative,
2162 const Verbosity& verbosity)
2163{
2164 for (auto& lines: abs_lines_per_species)
2165 abs_linesChangeBaseParameterForMatchingLevel(lines, QI, parameter_name, change, relative, verbosity);
2166}
2167
2168/* Workspace method: Doxygen documentation will be auto-generated */
2170 const ArrayOfQuantumIdentifier& QID,
2171 const String& parameter_name,
2172 const Vector& change,
2173 const Index& relative,
2174 const Verbosity& verbosity)
2175{
2176 ARTS_USER_ERROR_IF (QID.nelem() not_eq change.nelem(),
2177 "Mismatch between QID and change input lengths not allowed");
2178
2179 for (Index iq=0; iq<QID.nelem(); iq++)
2180 abs_linesChangeBaseParameterForMatchingLevel(abs_lines, QID[iq], parameter_name, change[iq], relative, verbosity);
2181}
2182
2183/* Workspace method: Doxygen documentation will be auto-generated */
2185 const ArrayOfQuantumIdentifier& QID,
2186 const String& parameter_name,
2187 const Vector& change,
2188 const Index& relative,
2189 const Verbosity& verbosity)
2190{
2191 ARTS_USER_ERROR_IF (QID.nelem() not_eq change.nelem(),
2192 "Mismatch between QID and change input lengths not allowed");
2193
2194 for (Index iq=0; iq<QID.nelem(); iq++)
2195 for (auto& lines: abs_lines_per_species)
2196 abs_linesChangeBaseParameterForMatchingLevel(lines, QID[iq], parameter_name, change[iq], relative, verbosity);
2197}
2198
2199/* Workspace method: Doxygen documentation will be auto-generated */
2201 const QuantumIdentifier& QI,
2202 const String& parameter_name,
2203 const Numeric& x,
2204 const Verbosity&)
2205{
2206 ARTS_USER_ERROR_IF (QI.type not_eq Quantum::IdentifierType::EnergyLevel,
2207 "Bad input. Must be energy level. Is: ", QI, '\n')
2208
2209 Index parameter_switch = -1;
2210
2211 ARTS_USER_ERROR_IF (parameter_name.nelem() == 0,
2212 "parameter_name is empty.\n");
2213 if (parameter_name == "Statistical Weight")
2214 parameter_switch = 1;
2215 else if (parameter_name == "Zeeman Coefficient")
2216 parameter_switch = 2;
2217
2218 for (auto& band: abs_lines) {
2219 for (Index k=0; k<band.NumLines(); k++) {
2220 const Absorption::QuantumIdentifierLineTarget lt = Absorption::QuantumIdentifierLineTarget(QI, band, k);
2221 if (lt == Absorption::QuantumIdentifierLineTargetType::Level and lt.lower) {
2222 switch (parameter_switch) {
2223 case 1: // "Statistical Weight":
2224 band.g_low(k) = x;
2225 break;
2226 case 2: // "Zeeman Coefficient":
2227 band.Line(k).Zeeman().gl() = x;
2228 break;
2229 default: {
2231 "Usupported paramter_name\n", parameter_name,
2232 "\nSee method description for supported parameter names.\n")
2233 }
2234 }
2235 }
2236
2237 if (lt == Absorption::QuantumIdentifierLineTargetType::Level and lt.upper) {
2238 switch (parameter_switch) {
2239 case 1: // "Statistical Weight":
2240 band.g_upp(k) = x;
2241 break;
2242 case 2: // "Zeeman Coefficient":
2243 band.Line(k).Zeeman().gu() = x;
2244 break;
2245 default: {
2247 "Usupported paramter_name\n", parameter_name,
2248 "\nSee method description for supported parameter names.\n")
2249 }
2250 }
2251 }
2252 }
2253 }
2254}
2255
2256/* Workspace method: Doxygen documentation will be auto-generated */
2258 const QuantumIdentifier& QI,
2259 const String& parameter_name,
2260 const Numeric& change,
2261 const Verbosity& verbosity)
2262{
2263 for (auto& lines: abs_lines_per_species)
2264 abs_linesSetBaseParameterForMatchingLevel(lines, QI, parameter_name, change, verbosity);
2265}
2266
2267/* Workspace method: Doxygen documentation will be auto-generated */
2269 const ArrayOfQuantumIdentifier& QID,
2270 const String& parameter_name,
2271 const Vector& change,
2272 const Verbosity& verbosity)
2273{
2274 ARTS_USER_ERROR_IF (QID.nelem() not_eq change.nelem(),
2275 "Mismatch between QID and change input lengths not allowed");
2276
2277 for (Index iq=0; iq<QID.nelem(); iq++)
2278 abs_linesSetBaseParameterForMatchingLevel(abs_lines, QID[iq], parameter_name, change[iq], verbosity);
2279}
2280
2281/* Workspace method: Doxygen documentation will be auto-generated */
2283 const ArrayOfQuantumIdentifier& QID,
2284 const String& parameter_name,
2285 const Vector& change,
2286 const Verbosity& verbosity)
2287{
2288 ARTS_USER_ERROR_IF (QID.nelem() not_eq change.nelem(),
2289 "Mismatch between QID and change input lengths not allowed");
2290
2291 for (Index iq=0; iq<QID.nelem(); iq++)
2292 for (auto& lines: abs_lines_per_species)
2293 abs_linesSetBaseParameterForMatchingLevel(lines, QID[iq], parameter_name, change[iq], verbosity);
2294}
2295
2299
2300/* Workspace method: Doxygen documentation will be auto-generated */
2302 Index& nlte_do,
2303 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
2304 const EnergyLevelMap& nlte_field,
2305 const Verbosity&)
2306{
2307 nlte_field.ThrowIfNotOK();
2308
2309 if (nlte_field.Data().empty()) {
2310 nlte_do = 0;
2311 return;
2312 }
2313 nlte_do = 1;
2314
2315 const Absorption::PopulationType poptyp = nlte_field.Energies().empty() ?
2316 Absorption::PopulationType::NLTE :
2317 Absorption::PopulationType::VibTemps;
2318
2319 for (auto& spec_lines: abs_lines_per_species) {
2320 for (auto& band: spec_lines) {
2321 for (auto& id: nlte_field.Levels()) {
2322 const Absorption::QuantumIdentifierLineTarget lt(id, band);
2323 if (lt == Absorption::QuantumIdentifierLineTargetType::Level) {
2324 for (Index k=0; k<band.NumLines(); k++) {
2325 ARTS_USER_ERROR_IF (poptyp==Absorption::PopulationType::NLTE and
2326 (not std::isnormal(band.A(k)) or band.A(k) < 0),
2327 "Error in band deemed for NLTE calculations by population distribution\n"
2328 "some of the lines in the band below have a bad Einstein coefficient:\n",
2329 band, '\n')
2330 }
2331 band.Population(poptyp);
2332 }
2333 }
2334 }
2335 }
2336}
2337
2341
2342/* Workspace method: Doxygen documentation will be auto-generated */
2344 const ArrayOfArrayOfSpeciesTag& abs_species,
2345 const Verbosity&) {
2346 abs_lines_per_species = ArrayOfArrayOfAbsorptionLines(abs_species.nelem(), ArrayOfAbsorptionLines(0));
2347}
2348
2349/* Workspace method: Doxygen documentation will be auto-generated */
2351 const Vector& f_grid,
2352 const Verbosity&)
2353{
2354 const Numeric fmax = max(f_grid);
2355 const Numeric fmin = min(f_grid);
2356
2357 for (auto& band: abs_lines) {
2358 for (Index k=band.NumLines()-1; k>=0; k--) {
2359 const Numeric fcut_upp = band.CutoffFreq(k);
2360 const Numeric fcut_low = band.CutoffFreqMinus(k);
2361
2362 if (fmax < fcut_low or fmin > fcut_upp) {
2363 band.RemoveLine(k);
2364 }
2365 }
2366 }
2367}
2368
2369/* Workspace method: Doxygen documentation will be auto-generated */
2371 const Vector& f_grid,
2372 const Verbosity& verbosity)
2373{
2374 for (auto& lines: abs_lines_per_species) {
2375 abs_linesCompact(lines, f_grid, verbosity);
2376 }
2377}
2378
2379/* Workspace method: Doxygen documentation will be auto-generated */
2381 const QuantumIdentifier& qid,
2382 const Verbosity&)
2383{
2384 for (Index i=0; i<abs_lines.nelem(); i++) {
2385 const Absorption::QuantumIdentifierLineTarget lt(qid, abs_lines[i]);
2386 if (lt == Absorption::QuantumIdentifierLineTargetType::Band) {
2387 abs_lines.erase(abs_lines.begin()+i);
2388 break;
2389 }
2390 }
2391}
2392
2396
2397template <class T>
2398std::vector<T> linspace(T s, T e, typename std::vector<T>::size_type count) noexcept {
2399 std::vector<T> ls(count);
2400
2401 if (count == 0)
2402 return ls;
2403 if (count == 1) {
2404 ls.front() = (e + s) / 2;
2405 return ls;
2406 }
2407 const T step = (e - s) / T(count - 1);
2408 ls.front() = s;
2409 ls.back() = e;
2410 for (typename std::vector<T>::size_type i = 1; i < count - 1; ++i)
2411 ls[i] = s + step * T(i);
2412 return ls;
2413}
2414
2415/* Workspace method: Doxygen documentation will be auto-generated */
2417 const ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
2418 const Numeric& delta_f_low,
2419 const Numeric& delta_f_upp,
2420 const Index& num_freqs,
2421 const Verbosity&)
2422{
2423 const Index n=nelem(abs_lines_per_species);
2424
2425 ARTS_USER_ERROR_IF (delta_f_low >= delta_f_upp,
2426 "The lower frequency delta has to be smaller "
2427 "than the upper frequency delta");
2428 ARTS_USER_ERROR_IF (num_freqs == 0,
2429 "Need more than zero frequency points");
2430 ARTS_USER_ERROR_IF(n < 1,
2431 "No lines found. Error? Use *VectorSet* "
2432 "to resize *f_grid*");
2433
2434 std::vector<Numeric> fout(0);
2435 for (auto& lines: abs_lines_per_species) {
2436 for (auto& band: lines) {
2437 for (Index k=0; k<band.NumLines(); k++) {
2438 if (num_freqs > 1) {
2439 auto ftmp = linspace<Numeric>(band.F0(k)+delta_f_low, band.F0(k)+delta_f_upp, std::vector<Numeric>::size_type(num_freqs));
2440 for (auto& f: ftmp) {
2441 if (f > 0) fout.push_back(f);
2442 }
2443 } else {
2444 fout.push_back(band.F0(k));
2445 }
2446 }
2447 }
2448 }
2449
2450 std::sort(fout.begin(), fout.end());
2451 fout.erase(std::unique(fout.begin(), fout.end()), fout.end());
2452 f_grid.resize(fout.size());
2453 for (Index i=0; i<f_grid.nelem(); i++)
2454 f_grid[i] = fout[i];
2455}
2456
2460
2461/* Workspace method: Doxygen documentation will be auto-generated */
2463 const Verbosity& verbosity)
2464{
2466
2467 std::map<Index, Index> qns;
2468
2469 for (auto& band: abs_lines) {
2470 for (Index iline=0; iline<band.NumLines(); iline++) {
2471 for (Index iqn=0; iqn<Index(QuantumNumberType::FINAL); iqn++) {
2472 if (band.LowerQuantumNumber(iline, QuantumNumberType(iqn)).isDefined() or
2473 band.UpperQuantumNumber(iline, QuantumNumberType(iqn)).isDefined()) {
2474 qns[iqn]++;
2475 }
2476 }
2477 }
2478 }
2479
2480 for (auto& qn: qns) {
2481 out0 << QuantumNumberType(qn.first) << ':' << ' ' << qn.second << '\n';
2482 }
2483}
2484
2485
2489
2491 const ArrayOfSpeciesTag& species,
2492 const Numeric lower_frequency,
2493 const Numeric upper_frequency,
2494 const Numeric lower_intensity,
2495 const Index safe,
2496 const Verbosity& verbosity) {
2497 const bool care_about_species = species.nelem();
2498
2499 for (auto& band: abs_lines) {
2500 if (care_about_species and species[0].Isotopologue() not_eq band.Isotopologue()) continue;
2501
2502 auto& lines = band.AllLines();
2503
2504 if (not safe) {
2505 std::vector<std::size_t> rem;
2506 for (std::size_t i=lines.size()-1; i<lines.size(); i--) {
2507 auto& line = lines[i];
2508 if (line.F0() < lower_frequency or
2509 line.F0() > upper_frequency or
2510 line.I0() < lower_intensity)
2511 rem.push_back(i);
2512 }
2513
2514 for (auto i: rem) band.RemoveLine(i);
2515 } else {
2516 const bool all_low = std::all_of(lines.begin(), lines.end(),
2517 [lower_frequency] (auto& line) {return line.F0() < lower_frequency;});
2518 const bool all_upp = std::all_of(lines.begin(), lines.end(),
2519 [upper_frequency] (auto& line) {return line.F0() > upper_frequency;});
2520 const bool low_int = std::all_of(lines.begin(), lines.end(),
2521 [lower_intensity] (auto& line) {return line.I0() < lower_intensity;});
2522 if (all_low or all_upp or low_int) lines.resize(0);
2523 }
2524 }
2525
2526 // Removes empty bands
2527 abs_linesRemoveEmptyBands(abs_lines, verbosity);
2528}
2529
2531 const Numeric& lower_frequency,
2532 const Numeric& upper_frequency,
2533 const Numeric& lower_intensity,
2534 const Index& safe,
2535 const Verbosity& verbosity) {
2536 remove_impl(abs_lines, {}, lower_frequency, upper_frequency, lower_intensity, safe, verbosity);
2537}
2538
2540 const Numeric& lower_frequency,
2541 const Numeric& upper_frequency,
2542 const Numeric& lower_intensity,
2543 const Index& safe,
2544 const Verbosity& verbosity) {
2545 for (auto& abs_lines: abs_lines_per_species)
2546 abs_linesRemoveLines(abs_lines, lower_frequency, upper_frequency, lower_intensity, safe, verbosity);
2547}
2548
2550 const ArrayOfSpeciesTag& species,
2551 const Numeric& lower_frequency,
2552 const Numeric& upper_frequency,
2553 const Numeric& lower_intensity,
2554 const Index& safe,
2555 const Verbosity& verbosity) {
2556 ARTS_USER_ERROR_IF(species.nelem() not_eq 1, "Must have a single species, got: ", species)
2557 ARTS_USER_ERROR_IF(species[0].Isotopologue().joker(), "Cannot give joker species, got: ", species)
2558
2559 remove_impl(abs_lines, species, lower_frequency, upper_frequency, lower_intensity, safe, verbosity);
2560}
2561
2563 const ArrayOfSpeciesTag& species,
2564 const Numeric& lower_frequency,
2565 const Numeric& upper_frequency,
2566 const Numeric& lower_intensity,
2567 const Index& safe,
2568 const Verbosity& verbosity) {
2569 for (auto& abs_lines: abs_lines_per_species)
2570 abs_linesRemoveLinesFromSpecies(abs_lines, species, lower_frequency, upper_frequency, lower_intensity, safe, verbosity);
2571}
Contains the absorption namespace.
Array< ArrayOfAbsorptionLines > ArrayOfArrayOfAbsorptionLines
Absorption::Lines AbsorptionLines
char * isot
char * newdata
type upp char * str
void abs_lines_per_speciesCreateFromLines(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfAbsorptionLines &abs_lines, const ArrayOfArrayOfSpeciesTag &abs_species, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesCreateFromLines.
Definition: m_abs.cc:93
void WriteXML(const String &output_file_format, const T &in, const String &filename, const Index &no_clobber, const String &in_wsvname, const String &filename_wsvname, const String &no_clobber_wsvname, const Verbosity &verbosity)
WORKSPACE METHOD: WriteXML.
Definition: m_xml.h:118
void abs_speciesSet(ArrayOfArrayOfSpeciesTag &abs_species, Index &abs_xsec_agenda_checked, Index &propmat_clearsky_agenda_checked, const ArrayOfString &species, const Verbosity &verbosity)
WORKSPACE METHOD: abs_speciesSet.
Index NumLines() const noexcept
Number of lines.
void AppendSingleLine(SingleLine &&sl)
Appends a single line to the absorption lines.
bool Match(const Lines &l) const noexcept
Checks if another line list matches this structure.
const std::vector< Rational > & UpperQuantumNumbers() const noexcept
Upper level quantum numbers.
const std::vector< Rational > & LowerQuantumNumbers() const noexcept
Lower level quantum numbers.
Numeric F0() const noexcept
Central frequency.
Zeeman::Model Zeeman() const noexcept
Zeeman model.
This can be used to make arrays out of anything.
Definition: array.h:107
Index nelem() const ARTS_NOEXCEPT
Number of elements.
Definition: array.h:195
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
Check if variable is empty.
Definition: matpackIV.cc:50
bool empty() const ARTS_NOEXCEPT
Returns true if variable size is zero.
Definition: matpackI.cc:46
Index nelem() const ARTS_NOEXCEPT
Returns the number of elements.
Definition: matpackI.cc:48
const ArrayOfQuantumIdentifier & Levels() const noexcept
Energy level type.
const Vector & Energies() const noexcept
Energy level type.
const Tensor4 & Data() const noexcept
Energy level type.
void ThrowIfNotOK() const ARTS_NOEXCEPT
Implements rational numbers to work with other ARTS types.
Definition: rational.h:52
The Vector class.
Definition: matpackI.h:876
void resize(Index n)
Resize function.
Definition: matpackI.cc:408
void check_name(const String &expected_name)
Check tag name.
Definition: xml_io_base.cc:48
void read_from_stream(istream &is)
Reads next XML tag.
Definition: xml_io_base.cc:179
Index nelem() const
Number of elements.
Definition: mystring.h:253
#define ARTS_ASSERT(condition,...)
Definition: debug.h:83
#define ARTS_USER_ERROR(...)
Definition: debug.h:150
#define ARTS_USER_ERROR_IF(condition,...)
Definition: debug.h:134
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:96
bool find_xml_file_existence(String &filename)
As find_xml_file but does not throw in the main body.
Definition: file.cc:452
void open_input_file(ifstream &file, const String &name)
Open a file for reading.
Definition: file.cc:147
This file contains basic functions to handle ASCII files.
#define abs(x)
#define min(a, b)
#define max(a, b)
void abs_lines_per_speciesSetPopulationForSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const String &type, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesSetPopulationForSpecies.
void abs_linesSetCutoffForMatch(ArrayOfAbsorptionLines &abs_lines, const String &type, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesSetCutoffForMatch.
void abs_lines_per_speciesSetCutoffForMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesSetCutoffForMatch.
void abs_lines_per_speciesChangeBaseParameterForMatchingLines(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Index &relative, const Index &loose_matching, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesChangeBaseParameterForMatchingLines.
void abs_lines_per_speciesSetNormalization(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetNormalization.
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_linesDeleteBadF0(ArrayOfAbsorptionLines &abs_lines, const Numeric &f0, const Index &lower, const Verbosity &)
WORKSPACE METHOD: abs_linesDeleteBadF0.
void abs_linesDeleteLinesWithBadOrHighChangingJs(ArrayOfAbsorptionLines &abs_lines, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesDeleteLinesWithBadOrHighChangingJs.
void abs_linesChangeBaseParameterForMatchingLines(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Index &relative, const Index &loose_matching, const Verbosity &)
WORKSPACE METHOD: abs_linesChangeBaseParameterForMatchingLines.
void abs_lines_per_speciesReadSplitCatalog(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesReadSplitCatalog.
void abs_linesFlatten(ArrayOfAbsorptionLines &abs_lines, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesFlatten.
void abs_lines_per_speciesSetT0ForMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetT0ForMatch.
void abs_lines_per_speciesSetLinemixingLimit(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Numeric &x, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetLinemixingLimit.
void abs_lines_per_speciesSetNormalizationForMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetNormalizationForMatch.
void abs_linesKeepBand(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &qid, const Verbosity &)
WORKSPACE METHOD: abs_linesKeepBand.
void abs_linesSetLineShapeTypeForMatch(ArrayOfAbsorptionLines &abs_lines, const String &type, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesSetLineShapeTypeForMatch.
void abs_lines_per_speciesWriteSpeciesSplitXML(const String &output_format, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesWriteSpeciesSplitXML.
void abs_linesTruncateGlobalQuantumNumbers(ArrayOfAbsorptionLines &abs_lines, const Verbosity &)
WORKSPACE METHOD: abs_linesTruncateGlobalQuantumNumbers.
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_speciesSetLineShapeModelParametersForMatchingLines(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_speciesSetLineShapeModelParametersForMatchingLines.
void abs_lines_per_speciesSetT0(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Numeric &x, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetT0.
void abs_linesSetNormalization(ArrayOfAbsorptionLines &abs_lines, const String &type, const Verbosity &)
WORKSPACE METHOD: abs_linesSetNormalization.
void abs_lines_per_speciesSetLineShapeType(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetLineShapeType.
void abs_linesSetNormalizationForMatch(ArrayOfAbsorptionLines &abs_lines, const String &type, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesSetNormalizationForMatch.
void abs_lines_per_speciesSetNormalizationForSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const String &type, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesSetNormalizationForSpecies.
void abs_lines_per_speciesWriteSplitXML(const String &output_format, const ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesWriteSplitXML.
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 abs_linesDeleteLinesWithUndefinedLocalQuanta(ArrayOfAbsorptionLines &abs_lines, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesDeleteLinesWithUndefinedLocalQuanta.
void abs_lines_per_speciesSetQuantumNumberForMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &qn, const Rational &x, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetQuantumNumberForMatch.
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 Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesRemoveLines.
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_speciesSetMirroringForSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const String &type, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesSetMirroringForSpecies.
void abs_linesCompact(ArrayOfAbsorptionLines &abs_lines, const Vector &f_grid, const Verbosity &)
WORKSPACE METHOD: abs_linesCompact.
void abs_linesMakeManualMirroring(ArrayOfAbsorptionLines &abs_lines, const Verbosity &)
WORKSPACE METHOD: abs_linesMakeManualMirroring.
void abs_linesRemoveBand(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &qid, const Verbosity &)
WORKSPACE METHOD: abs_linesRemoveBand.
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 Verbosity &verbosity)
WORKSPACE METHOD: abs_linesRemoveLinesFromSpecies.
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_lines_per_speciesSetCutoff(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Numeric &x, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetCutoff.
void abs_linesRemoveEmptyBands(ArrayOfAbsorptionLines &abs_lines, const Verbosity &)
WORKSPACE METHOD: abs_linesRemoveEmptyBands.
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 &)
WORKSPACE METHOD: abs_linesDeleteWithLines.
void abs_linesSetLineShapeType(ArrayOfAbsorptionLines &abs_lines, const String &type, const Verbosity &)
WORKSPACE METHOD: abs_linesSetLineShapeType.
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_lines_per_speciesSetLineShapeTypeForSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const String &type, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesSetLineShapeTypeForSpecies.
void abs_linesWriteSplitXML(const String &output_format, const ArrayOfAbsorptionLines &abs_lines, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesWriteSplitXML.
void abs_linesSetQuantumNumberForMatch(ArrayOfAbsorptionLines &abs_lines, const String &qn, const Rational &x, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesSetQuantumNumberForMatch.
void abs_lines_per_speciesSetLineShapeTypeForMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetLineShapeTypeForMatch.
void abs_lines_per_speciesSetPopulation(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetPopulation.
void abs_linesSetCutoff(ArrayOfAbsorptionLines &abs_lines, const String &type, const Numeric &x, const Verbosity &)
WORKSPACE METHOD: abs_linesSetCutoff.
void abs_lines_per_speciesSetBaseParameterForMatchingLevel(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesSetBaseParameterForMatchingLevel.
void abs_linesAppendWithLines(ArrayOfAbsorptionLines &abs_lines, const ArrayOfAbsorptionLines &appending_lines, const Index &safe, const Verbosity &)
WORKSPACE METHOD: abs_linesAppendWithLines.
void abs_linesRemoveLines(ArrayOfAbsorptionLines &abs_lines, const Numeric &lower_frequency, const Numeric &upper_frequency, const Numeric &lower_intensity, const Index &safe, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesRemoveLines.
void abs_lines_per_speciesMakeManualMirroring(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesMakeManualMirroring.
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_linesSetT0ForMatch(ArrayOfAbsorptionLines &abs_lines, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesSetT0ForMatch.
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_lines_per_speciesSetBaseParameterForSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Index &loose_matching, const String &species_tag, const Verbosity &verbosity)
void abs_linesSetBaseParameterForMatchingLines(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &x, const Index &loose_matching, const Verbosity &)
WORKSPACE METHOD: abs_linesSetBaseParameterForMatchingLines.
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 Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesRemoveLinesFromSpecies.
void abs_lines_per_speciesSetLinemixingLimitForSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const Numeric &x, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesSetLinemixingLimitForSpecies.
void abs_linesSetLineShapeModelParametersForMatchingLines(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_linesSetLineShapeModelParametersForMatchingLines.
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_speciesSetEmpty(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const Verbosity &)
WORKSPACE METHOD: abs_lines_per_speciesSetEmpty.
void abs_lines_per_speciesTruncateQuantumNumbers(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Index &pos, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesTruncateQuantumNumbers.
void abs_linesSetBaseParameterForMatchingLevel(ArrayOfAbsorptionLines &abs_lines, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &x, const Verbosity &)
WORKSPACE METHOD: abs_linesSetBaseParameterForMatchingLevel.
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 Verbosity &verbosity)
void abs_linesReplaceWithLines(ArrayOfAbsorptionLines &abs_lines, const ArrayOfAbsorptionLines &replacing_lines, const Verbosity &)
WORKSPACE METHOD: abs_linesReplaceWithLines.
void abs_linesSetLinemixingLimit(ArrayOfAbsorptionLines &abs_lines, const Numeric &x, const Verbosity &)
WORKSPACE METHOD: abs_linesSetLinemixingLimit.
void abs_linesSetPopulationForMatch(ArrayOfAbsorptionLines &abs_lines, const String &type, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesSetPopulationForMatch.
void abs_linesSetLinemixingLimitForMatch(ArrayOfAbsorptionLines &abs_lines, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesSetLinemixingLimitForMatch.
void abs_linesRemoveUnusedLocalQuantumNumbers(ArrayOfAbsorptionLines &abs_lines, const Verbosity &)
WORKSPACE METHOD: abs_linesRemoveUnusedLocalQuantumNumbers.
void abs_lines_per_speciesFlatten(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesFlatten.
std::vector< QuantumNumberType > string2vecqn(const String &qnstr)
Get a list of quantum numbers from a string.
void abs_linesSetMirroring(ArrayOfAbsorptionLines &abs_lines, const String &type, const Verbosity &)
WORKSPACE METHOD: abs_linesSetMirroring.
void abs_linesSetEmptyBroadeningParametersToEmpty(ArrayOfAbsorptionLines &abs_lines, const Verbosity &)
WORKSPACE METHOD: abs_linesSetEmptyBroadeningParametersToEmpty.
void abs_lines_per_speciesMakeManualMirroringSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const ArrayOfSpeciesTag &species, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesMakeManualMirroringSpecies.
void abs_linesSetMirroringForMatch(ArrayOfAbsorptionLines &abs_lines, const String &type, const QuantumIdentifier &QI, const Verbosity &)
WORKSPACE METHOD: abs_linesSetMirroringForMatch.
void abs_lines_per_speciesSetT0ForSpecies(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfArrayOfSpeciesTag &abs_species, const Numeric &x, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesSetT0ForSpecies.
void abs_linesTruncateQuantumNumbers(ArrayOfAbsorptionLines &abs_lines, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesTruncateQuantumNumbers.
void abs_linesDeleteLinesWithQuantumNumberAbove(ArrayOfAbsorptionLines &abs_lines, const String &qn_id, const Index &qn_val, const Verbosity &)
WORKSPACE METHOD: abs_linesDeleteLinesWithQuantumNumberAbove.
void abs_lines_per_speciesSetMirroringForMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetMirroringForMatch.
void abs_lines_per_speciesSetBaseParameterForMatchingLines(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const QuantumIdentifier &QI, const String &parameter_name, const Numeric &change, const Index &loose_matching, const Verbosity &verbosity)
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 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_linesPrintDefinedQuantumNumbers(const ArrayOfAbsorptionLines &abs_lines, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesPrintDefinedQuantumNumbers.
void abs_linesSetT0(ArrayOfAbsorptionLines &abs_lines, const Numeric &x, const Verbosity &)
WORKSPACE METHOD: abs_linesSetT0.
void abs_linesSetPopulation(ArrayOfAbsorptionLines &abs_lines, const String &type, const Verbosity &)
WORKSPACE METHOD: abs_linesSetPopulation.
void abs_lines_per_speciesSetMirroring(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetMirroring.
void nlteSetByQuantumIdentifiers(Index &nlte_do, ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const EnergyLevelMap &nlte_field, const Verbosity &)
WORKSPACE METHOD: nlteSetByQuantumIdentifiers.
void abs_lines_per_speciesSetPopulationForMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const String &type, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetPopulationForMatch.
void abs_linesWriteSpeciesSplitXML(const String &output_format, const ArrayOfAbsorptionLines &abs_lines, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesWriteSpeciesSplitXML.
void abs_linesSetBaseParameterForMatchingLevels(ArrayOfAbsorptionLines &abs_lines, const ArrayOfQuantumIdentifier &QID, const String &parameter_name, const Vector &change, const Verbosity &verbosity)
WORKSPACE METHOD: abs_linesSetBaseParameterForMatchingLevels.
void abs_lines_per_speciesSetCutoffForSpecies(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_speciesSetCutoffForSpecies.
void abs_lines_per_speciesSetBaseParameterForMatchingLevels(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const ArrayOfQuantumIdentifier &QID, const String &parameter_name, const Vector &change, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesSetBaseParameterForMatchingLevels.
std::vector< T > linspace(T s, T e, typename std::vector< T >::size_type count) noexcept
void abs_lines_per_speciesSetLinemixingLimitForMatch(ArrayOfArrayOfAbsorptionLines &abs_lines_per_species, const Numeric &x, const QuantumIdentifier &QI, const Verbosity &v)
WORKSPACE METHOD: abs_lines_per_speciesSetLinemixingLimitForMatch.
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 Index &loose_matching, const String &species_tag, const Verbosity &verbosity)
WORKSPACE METHOD: abs_lines_per_speciesChangeBaseParameterForSpecies.
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.
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
const Joker joker
#define CREATE_OUT3
Definition: messages.h:207
#define CREATE_OUT2
Definition: messages.h:206
#define CREATE_OUT0
Definition: messages.h:204
my_basic_string< char > String
The String type for ARTS.
Definition: mystring.h:287
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 ReadFromHitran2012Stream(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.
char Type type
constexpr bool modelparameterEmpty(const ModelParameters mp) noexcept
constexpr Index nVars
Current max number of line shape variables.
X3 LineCenter None
Definition: constants.h:576
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
Definition: isotopologues.h:50
ArrayOfIsotopeRecord isotopologues(Species spec)
Definition: isotopologues.cc:6
Model GetAdvancedModel(const QuantumIdentifier &qid) ARTS_NOEXCEPT
Returns an advanced Zeeman model.
Definition: zeemandata.cc:103
constexpr Rational end(Rational Ju, Rational Jl, Polarization type) noexcept
Gives the largest M for a polarization type of this transition.
Definition: zeemandata.h:109
#define v
#define a
#define c
#define b
Quantum::Identifier QuantumIdentifier
Definition: quantum.h:471
constexpr bool IsValidQuantumNumberName(const std::string_view name)
Check for valid quantum number name.
Definition: quantum.h:482
std::set< Species::Species > lbl_species(const ArrayOfArrayOfSpeciesTag &abs_species) noexcept
Species::Tag SpeciesTag
Definition: species_tags.h:99
Single line reading output.
Coefficients and temperature model for SingleSpeciesModel.
static constexpr Index N
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:154
void xml_read_from_file(const String &filename, T &type, const Verbosity &verbosity)
Reads data from XML file.
Definition: xml_io.h:160
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:552
@ FILE_TYPE_ASCII
Definition: xml_io_base.h:41
@ NUMERIC_TYPE_DOUBLE
Definition: xml_io_base.h:46
@ ENDIAN_TYPE_LITTLE
Definition: xml_io_base.h:47