ARTS 2.5.0 (git: 9ee3ac6c)
lineshapemodel.cc
Go to the documentation of this file.
1/* Copyright (C) 2018
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
35#include "constants.h"
36#include "lineshapemodel.h"
37#include "matpackI.h"
38
40 const String& coeff) {
41 // Test viability of model variables
42 const auto var_type = LineShape::toVariableOrThrow(var);
43
44 // Test viability of model coefficients
45 const auto coeff_type = Options::toLineShapeCoeffOrThrow(coeff);
46
47// Define a repetitive pattern. Update if/when there are more coefficients in the future
48#define ReturnJacPropMatType(ID) \
49 case LineShape::Variable::ID: \
50 switch (coeff_type) { \
51 case Options::LineShapeCoeff::X0: \
52 return Jacobian::Line::Shape##ID##X0; \
53 case Options::LineShapeCoeff::X1: \
54 return Jacobian::Line::Shape##ID##X1; \
55 case Options::LineShapeCoeff::X2: \
56 return Jacobian::Line::Shape##ID##X2; \
57 case Options::LineShapeCoeff::X3: \
58 return Jacobian::Line::Shape##ID##X3; \
59 case Options::LineShapeCoeff::FINAL: \
60 return Jacobian::Line::FINAL; \
61 } break
62
63 switch (var_type) {
73 case LineShape::Variable::FINAL: {
74 /* Leave last */ }
75 }
76
77#undef ReturnJacPropMatType
78
79 return Jacobian::Line::FINAL;
80}
81
82std::istream& LineShape::from_artscat4(std::istream& is,
83 Type& mtype,
84 bool& self,
85 bool& bath,
86 Model& m,
87 ArrayOfSpecies& species,
88 const QuantumIdentifier& qid) {
89 // Set or reset variables
90 mtype = Type::VP;
91 self = true;
92 bath = false;
93 m.mdata = std::vector<SingleSpeciesModel>(7);
94 species = ArrayOfSpecies(7);
95
96 // Set species
97 species[1] = Species::fromShortName("N2");
98 species[2] = Species::fromShortName("O2");
99 species[3] = Species::fromShortName("H2O");
100 species[4] = Species::fromShortName("CO2");
101 species[5] = Species::fromShortName("H2");
102 species[6] = Species::fromShortName("He");
103
104 // Temperature types
105 for (auto& v : m.mdata) {
106 v.G0().type = TemperatureModel::T1;
107 v.D0().type = TemperatureModel::T5;
108 }
109
110 // G0 main coefficient
111 for (auto& v : m.mdata) is >> v.G0().X0;
112
113 // G0 exponent is same as D0 exponent
114 for (auto& v : m.mdata) {
115 is >> v.G0().X1;
116 v.D0().X1 = v.G0().X1;
117 }
118
119 // D0 coefficient
120 m.mdata.front().D0().X0 = 0;
121 for (int k = 1; k < 7; k++)
122 is >> m.mdata[k].D0().X0;
123
124 // Special case when self is part of this list, it needs to be removed
125 for (int k = 1; k < 7; k++) {
126 if (qid.Species() == species[k]) {
127 ARTS_USER_ERROR_IF(m.mdata.front().G0().X0 not_eq m.mdata[k].G0().X0 or
128 m.mdata.front().G0().X1 not_eq m.mdata[k].G0().X1 or
129 m.mdata.front().D0().X1 not_eq m.mdata[k].D0().X1,
130 "Species is ", qid.Isotopologue(), " and this is a broadening species in ARTSCAT-4.\n"
131 "Despite this, values representing self and ", qid.Isotopologue(), " does not match "
132 "in input string\n")
133 m.mdata.front().D0().X0 = m.mdata[k].D0().X0;
134 m.Remove(k, species);
135 break;
136 }
137 }
138
139 return is;
140}
141
142std::istream& LineShape::from_linefunctiondata(std::istream& data,
143 Type& mtype,
144 bool& self,
145 bool& bath,
146 Model& m,
147 ArrayOfSpecies& species) {
148 self = bath = false;
149 Index specs;
150 String s;
151
152 // The first tag should give the line shape scheme
153 data >> mtype;
154 check_enum_error(mtype, "Bad Data");
155
156 // Order of elements for line shape
157 const auto shapeparams =
159
160 // The second tag should give the line mixing scheme
161 data >> s;
162
163 // Order of elements for line mixing
164 const auto mixingparams =
166
167 // The third tag should contain the number of species
168 data >> specs;
169 species.resize(specs);
170 m.mdata.resize(specs);
171
172 ARTS_USER_ERROR_IF (not specs and mtype not_eq Type::DP,
173 "Need at least one species for non-Doppler line shapes");
174
175 // For all species, we need to set the methods to compute them
176 for (Index i = 0; i < specs; i++) {
177 // This should be a species tag or one of the specials, SELF or BATH
178 data >> s;
179 if (s == self_broadening) {
180 // If the species is self, then we need to flag this
181 self = true;
182 ARTS_USER_ERROR_IF (i not_eq 0,
183 "Self broadening must be first, it is not\n");
184 } else if (s == bath_broadening) {
185 // If the species is air, then we need to flag this
186 bath = true;
187 species[i] = Species::Species::Bath;
188 ARTS_USER_ERROR_IF (i not_eq specs - 1,
189 "Air/bath broadening must be last, it is not\n");
190 } else {
191 // Otherwise, we hope we find a species
192 species[i] = Species::fromShortName(s);
193 ARTS_USER_ERROR_IF (not good_enum(species[i]),
194 "Encountered ", s, " in a position where a species should have been "
195 "defined.\nPlease check your pressure broadening data structure and ensure "
196 "that it follows the correct conventions.\n")
197 }
198
199 // For all parameters
200 for (auto& params : {shapeparams, mixingparams}) {
201 for (auto& param : params) {
202 data >> s; // Should contain a temperature tag
203
204 const auto type = toTemperatureModel(s);
205 const Index ntemp =
207
208 m.mdata[i].Data()[Index(param)].type = type;
209 if (ntemp <= ModelParameters::N) {
210 switch (ntemp) {
211 case 1:
212 data >> m.mdata[i].Data()[Index(param)].X0;
213 break;
214 case 2:
215 data >> m.mdata[i].Data()[Index(param)].X0 >>
216 m.mdata[i].Data()[Index(param)].X1;
217 break;
218 case 3:
219 data >> m.mdata[i].Data()[Index(param)].X0 >>
220 m.mdata[i].Data()[Index(param)].X1 >>
221 m.mdata[i].Data()[Index(param)].X2;
222 break;
223 case 0:
224 break;
225 default:
227 "Unknown number of input parameters in Legacy mode.");
228 }
229 } else { // Has to be the only allowed interpolation case
230 ARTS_USER_ERROR_IF (ntemp > 12,
231 "Too many input parameters in interpolation results Legacy mode.");
233 data >> temp; // should be 200
234 data >> temp; // should be 250
235 data >> temp; // should be 296
236 data >> temp; // should be 340
237 data >> m.mdata[i].Y().X0 >> m.mdata[i].Y().X1 >> m.mdata[i].Y().X2 >> m.mdata[i].Y().X3
238 >> m.mdata[i].G().X0 >> m.mdata[i].G().X1 >> m.mdata[i].G().X2 >> m.mdata[i].G().X3;
239 }
240 }
241 }
242 }
243
244 return data;
245}
246
248 std::istream& data,
249 LineShape::Type& mtype,
250 bool& self,
251 bool& bath,
252 Model& m,
253 ArrayOfSpecies& species,
254 const QuantumIdentifier& qid) {
255 String s;
256 data >> s;
257
260 const auto self_in_list = LegacyPressureBroadeningData::self_listed(qid, type);
261
262 Vector x(n);
263 for (auto& num : x) data >> num;
264
266 self,
267 bath,
268 m,
269 species,
270 x,
271 type,
272 self_in_list,
273 qid.Species());
274
275 return data;
276}
277
278std::istream& LineShape::from_linemixingdata(std::istream& data,
279 LineShape::Model& lsc) {
280 String s;
281 data >> s;
282
285
286 Vector x(n);
287 for (auto& num : x) data >> num;
288
290
291 return data;
292}
293
294#pragma GCC diagnostic push
295#pragma GCC diagnostic ignored "-Wreturn-type"
297 LineShape::Type& mtype,
298 bool& self,
299 bool& bath,
300 Model& m,
301 ArrayOfSpecies& species,
302 Vector x,
304 bool self_in_list,
305 Species::Species self_spec) {
306 switch (type) {
307 case TypePB::PB_NONE:
308 mtype = LineShape::Type::DP;
309 self = bath = false;
310 m = Model();
311 species.resize(0);
312 return;
313 case TypePB::PB_AIR_BROADENING:
314 mtype = LineShape::Type::VP;
315 self = bath = true;
316 m = Model(x[0], x[1], x[2], x[3], x[4]);
317 species.resize(2);
318 species[0] = self_spec;
319 species[1] = Species::Species::Bath;
320 return;
321 case TypePB::PB_AIR_AND_WATER_BROADENING:
322 if (self_in_list) {
323 mtype = LineShape::Type::VP;
324 self = false;
325 bath = true;
326 m.Data().resize(2);
327 m.Data()[0].G0() = {TemperatureModel::T1, x[0], x[1], 0, 0};
328 m.Data()[0].D0() = {TemperatureModel::T5, x[2], x[1], 0, 0};
329 m.Data()[1].G0() = {TemperatureModel::T1, x[3], x[4], 0, 0};
330 m.Data()[1].D0() = {TemperatureModel::T5, x[5], x[4], 0, 0};
331 species.resize(2);
332 species[0] = Species::fromShortName("H2O");
333 species[1] = Species::Species::Bath;
334 return;
335 } else {
336 mtype = LineShape::Type::VP;
337 self = bath = true;
338 m.Data().resize(2);
339 m.Data()[0].G0() = {TemperatureModel::T1, x[0], x[1], 0, 0};
340 m.Data()[0].D0() = {TemperatureModel::T5, x[2], x[1], 0, 0};
341 m.Data()[2].G0() = {TemperatureModel::T1, x[3], x[4], 0, 0};
342 m.Data()[2].D0() = {TemperatureModel::T5, x[5], x[4], 0, 0};
343 m.Data()[1].G0() = {TemperatureModel::T1, x[6], x[7], 0, 0};
344 m.Data()[1].D0() = {TemperatureModel::T5, x[8], x[7], 0, 0};
345 species.resize(3);
346 species[0] = self_spec;
347 species[1] = Species::fromShortName("H2O");
348 species[2] = Species::Species::Bath;
349 return;
350 }
351 case TypePB::PB_PLANETARY_BROADENING:
352 if (self_in_list) {
353 mtype = LineShape::Type::VP;
354 self = bath = false;
355 m.Data().resize(6);
356 m.Data()[0].G0() = {TemperatureModel::T1, x[1], x[8], 0, 0};
357 m.Data()[0].D0() = {TemperatureModel::T5, x[14], x[8], 0, 0};
358 m.Data()[1].G0() = {TemperatureModel::T1, x[2], x[9], 0, 0};
359 m.Data()[1].D0() = {TemperatureModel::T5, x[15], x[9], 0, 0};
360 m.Data()[2].G0() = {TemperatureModel::T1, x[3], x[10], 0, 0};
361 m.Data()[2].D0() = {TemperatureModel::T5, x[16], x[10], 0, 0};
362 m.Data()[3].G0() = {TemperatureModel::T1, x[4], x[11], 0, 0};
363 m.Data()[3].D0() = {TemperatureModel::T5, x[17], x[11], 0, 0};
364 m.Data()[4].G0() = {TemperatureModel::T1, x[5], x[12], 0, 0};
365 m.Data()[4].D0() = {TemperatureModel::T5, x[18], x[12], 0, 0};
366 m.Data()[5].G0() = {TemperatureModel::T1, x[6], x[13], 0, 0};
367 m.Data()[5].D0() = {TemperatureModel::T5, x[19], x[13], 0, 0};
368 species = {Species::fromShortName("N2"),
369 Species::fromShortName("O2"),
370 Species::fromShortName("H2O"),
371 Species::fromShortName("CO2"),
372 Species::fromShortName("H2"),
373 Species::fromShortName("He")};
374 return;
375 } else {
376 mtype = LineShape::Type::VP;
377 self = true;
378 bath = false;
379 m.Data().resize(7);
380 m.Data()[0].G0() = {TemperatureModel::T1, x[0], x[7], 0, 0};
381 // ssm[0].D0() = ...
382 m.Data()[1].G0() = {TemperatureModel::T1, x[1], x[8], 0, 0};
383 m.Data()[1].D0() = {TemperatureModel::T5, x[14], x[8], 0, 0};
384 m.Data()[2].G0() = {TemperatureModel::T1, x[2], x[9], 0, 0};
385 m.Data()[2].D0() = {TemperatureModel::T5, x[15], x[9], 0, 0};
386 m.Data()[3].G0() = {TemperatureModel::T1, x[3], x[10], 0, 0};
387 m.Data()[3].D0() = {TemperatureModel::T5, x[16], x[10], 0, 0};
388 m.Data()[4].G0() = {TemperatureModel::T1, x[4], x[11], 0, 0};
389 m.Data()[4].D0() = {TemperatureModel::T5, x[17], x[11], 0, 0};
390 m.Data()[5].G0() = {TemperatureModel::T1, x[5], x[12], 0, 0};
391 m.Data()[5].D0() = {TemperatureModel::T5, x[18], x[12], 0, 0};
392 m.Data()[6].G0() = {TemperatureModel::T1, x[6], x[13], 0, 0};
393 m.Data()[6].D0() = {TemperatureModel::T5, x[19], x[13], 0, 0};
394 species.resize(7);
395 species[0] = self_spec;
396 species[1] = Species::fromShortName("N2");
397 species[2] = Species::fromShortName("O2");
398 species[3] = Species::fromShortName("H2O");
399 species[4] = Species::fromShortName("CO2");
400 species[5] = Species::fromShortName("H2");
401 species[6] = Species::fromShortName("He");
402 return;
403 }
404 }
405}
406#pragma GCC diagnostic pop
407
410 Model y(1);
411 switch (type) {
412 case TypeLM::LM_NONE:
413 break;
414 case TypeLM::LM_LBLRTM:
415 y.Data().front().Y().type = LineShape::TemperatureModel::LM_AER;
416 y.Data().front().G().type = LineShape::TemperatureModel::LM_AER;
417 y.Data().front().Y().X0 = x[4];
418 y.Data().front().Y().X1 = x[5];
419 y.Data().front().Y().X2 = x[6];
420 y.Data().front().Y().X3 = x[7];
421 y.Data().front().G().X0 = x[8];
422 y.Data().front().G().X1 = x[9];
423 y.Data().front().G().X2 = x[10];
424 y.Data().front().G().X3 = x[11];
425 break;
426 case TypeLM::LM_LBLRTM_O2NonResonant:
427 y.Data().front().G().type = LineShape::TemperatureModel::T0;
428 y.Data().front().G().X0 = x[0];
429 break;
430 case TypeLM::LM_2NDORDER:
431 y.Data().front().Y().type = LineShape::TemperatureModel::T4;
432 y.Data().front().Y().X0 = x[0];
433 y.Data().front().Y().X1 = x[1];
434 y.Data().front().Y().X2 = x[7];
435 y.Data().front().G().type = LineShape::TemperatureModel::T4;
436 y.Data().front().G().X0 = x[2];
437 y.Data().front().G().X1 = x[3];
438 y.Data().front().G().X2 = x[8];
439 y.Data().front().DV().type = LineShape::TemperatureModel::T4;
440 y.Data().front().DV().X0 = x[4];
441 y.Data().front().DV().X1 = x[5];
442 y.Data().front().DV().X2 = x[9];
443 break;
444 case TypeLM::LM_1STORDER:
445 y.Data().front().Y().type = LineShape::TemperatureModel::T1;
446 y.Data().front().Y().X0 = x[1];
447 y.Data().front().Y().X1 = x[2];
448 break;
449 case TypeLM::LM_BYBAND:
450 break;
451 }
452 return y;
453}
454
455
456Vector LineShape::vmrs(const ConstVectorView& atmospheric_vmrs,
457 const ArrayOfArrayOfSpeciesTag& atmospheric_species,
458 const ArrayOfSpecies& lineshape_species) ARTS_NOEXCEPT {
459 ARTS_ASSERT (atmospheric_species.nelem() == atmospheric_vmrs.nelem(), "Bad atmospheric inputs");
460
461 const Index n = lineshape_species.nelem();
462 const Index back = n - 1; // Last index
463
464 // Initialize list of VMRS to 0
465 Vector line_vmrs(n, 0);
466
467 // We need to know if bath is an actual species
468 const bool bath = lineshape_species.back() == Species::Species::Bath;
469
470 // Loop species
471 for (Index i = 0; i < n - bath; i++) {
472 const Species::Species target = lineshape_species[i];
473
474 // Find species in list or do nothing at all
475 Index this_species_index = -1;
476 for (Index j = 0; j < atmospheric_species.nelem(); j++) {
477 if (atmospheric_species[j].Species() == target) {
478 this_species_index = j;
479 }
480 }
481
482 // Set to non-zero in-case species exists
483 if (this_species_index not_eq -1) {
484 line_vmrs[i] = atmospheric_vmrs[this_species_index];
485 }
486 }
487
488 // Renormalize, if bath-species exist this is automatic.
489 if (bath) {
490 line_vmrs[back] = 1.0 - line_vmrs.sum();
491 } else if(line_vmrs.sum() == 0) { // Special case
492 } else {
493 line_vmrs /= line_vmrs.sum();
494 }
495
496 return line_vmrs;
497}
498
499Vector LineShape::mass(const ConstVectorView& atmospheric_vmrs,
500 const ArrayOfArrayOfSpeciesTag& atmospheric_species,
501 const ArrayOfSpecies& lineshape_species,
503 ARTS_ASSERT (atmospheric_species.nelem() == atmospheric_vmrs.nelem(),
504 "Bad atmospheric inputs");
505
506 const Index n = lineshape_species.nelem();
507 const Index back = n - 1; // Last index
508
509 // Initialize list of VMRS to 0
510 Vector line_vmrs(lineshape_species.nelem(), 0);
511 Vector line_mass(lineshape_species.nelem(), 0);
512
513 // We need to know if bath is an actual species
514 const bool bath = lineshape_species.back() == Species::Species::Bath;
515
516 // Loop species ignoring self and bath
517 for (Index i = 0; i < n - bath; i++) {
518 // Select target in-case this is self-broadening
519 const Species::Species target = lineshape_species[i];
520
521 // Find species in list or do nothing at all
522 Index this_species_index = -1;
523 for (Index j = 0; j < atmospheric_species.nelem(); j++) {
524 if (atmospheric_species[j].Species() == target) {
525 this_species_index = j;
526 }
527 }
528
529 // Set to non-zero in-case species exists
530 if (this_species_index not_eq -1) {
531 line_vmrs[i] = atmospheric_vmrs[this_species_index];
532 line_mass[i] = Species::mean_mass(target, ir);
533 }
534 }
535
536 // Renormalize, if bath-species exist this is automatic.
537 if(line_vmrs.sum() == 0) { // Special case
538 } else if (bath) {
539 line_mass[back] = (line_vmrs * line_mass) / line_vmrs.sum();
540 }
541
542 return line_mass;
543}
544
545std::ostream& LineShape::operator<<(std::ostream& os, const LineShape::Model& m)
546{
547 for(auto& data: m.Data())
548 os << data;
549 return os;
550}
551
552std::istream& LineShape::operator>>(std::istream& is, Model& m)
553{
554 for(auto& data: m.Data())
555 is >> data;
556 return is;
557}
558
559
561{
562 String out = "";
563 const auto& vars = enumtyps::VariableTypes;
564
565 for (auto& var: vars) {
566 if (std::any_of(m.Data().cbegin(), m.Data().cend(),
567 [var](auto& x){return x.Get(var).type not_eq TemperatureModel::None;})) {
568 out += String(toString(var)) + ' ';
569 for (auto& ssm: m.Data())
570 out += String(toString(ssm.Get(var).type)) + ' ';
571 }
572 }
573
574 if(out.size())
575 out.pop_back();
576
577 return out;
578}
579
580
582{
583 if (s.nelem() == 0)
584 return LineShape::Model();
585
586 const auto& names = enumstrs::VariableNames;
587
588 std::istringstream str(s);
589 String part;
590 Variable var=Variable::ETA;
591 TemperatureModel tm=TemperatureModel::None;
592 Index i=-100000;
593
594 std::vector<SingleSpeciesModel> ssms(0);
595 while (not str.eof()) {
596 str >> part;
597 if(std::any_of(names.cbegin(), names.cend(),
598 [part](auto x){return part == x;})) {
599 i=-1;
600 var = toVariable(part);
601 }
602 else {
603 i++;
604 tm = toTemperatureModel(part);
605 }
606
607 if (i < 0)
608 continue;
609 if (i < Index(ssms.size()))
610 goto add_var;
611 else {
612 ssms.emplace_back();
613 add_var:
614 auto mp = ssms[i].Get(var);
615 mp.type = tm;
616 ssms[i].Set(var, mp);
617 }
618 }
619
620 return Model(ssms);
621}
622
624{
625 std::ostringstream os;
626 switch (mp.type) {
628 os << 0;
629 break;
630 case TemperatureModel::T0:
631 os << mp.X0;
632 break;
633 case TemperatureModel::T1:
634 os << mp.X0 << " * (" << T0 << "/T)^" << mp.X1;
635 break;
636 case TemperatureModel::T2:
637 os << mp.X0 << " * (" << T0 << "/T)^" << mp.X1 << " / (1 + " << mp.X2 << " * log(T/" << T0 << "))";
638 break;
639 case TemperatureModel::T3:
640 os << mp.X0 << " + " << mp.X1 << " * (" << T0 << " - T)";
641 break;
642 case TemperatureModel::T4:
643 os << "(" << mp.X0 << " + " << mp.X1 << " * (" << T0 << "/T - 1)) * (" << T0 << "/T)^" << mp.X2;
644 break;
645 case TemperatureModel::T5:
646 os << mp.X0 << " * (" << T0 << "/T)^(0.25 + 1.5 * " << mp.X1 << ")";
647 break;
648 case TemperatureModel::LM_AER:
649 os << '(' << "Linear interpolation to y(x) from x-ref = [200, 250, 296, 340] and y-ref = [" << mp.X0 << ", " << mp.X1 << ", " << mp.X2 << ", " << mp.X3 << ']' << ')';
650 break;
651 case TemperatureModel::DPL:
652 os << '(' << mp.X0 << " * (" << T0 << "/T)^" << mp.X1 << " + " << mp.X2 << " * (" << T0 << "/T)^" << mp.X3 << ')';
653 break;
654 case TemperatureModel::POLY:
655 os << '(' << mp.X0 << " + " << mp.X1 << " * T + " << mp.X2 << " * T * T + " << mp.X3 << " * T * T * T)";
656 case TemperatureModel::FINAL: break;
657 }
658
659 return os.str();
660}
661
663 const bool self,
664 const ArrayOfSpecies& sts,
665 const Numeric T0)
666{
667 const auto& vars = enumtyps::VariableTypes;
668 ArrayOfString as(0);
669
670 for (Index i=0; i<Index(Variable::FINAL); i++) {
671 Variable var = vars[i];
672
673 if (std::any_of(m.Data().cbegin(), m.Data().cend(),
674 [var](auto& x){return x.Get(var).type not_eq TemperatureModel::None;})) {
675
676 std::ostringstream os;
677 os << var << " ~ ";
678 for (Index j=0; j<sts.nelem(); j++) {
679 if (j == 0 and self)
680 os << "VMR(" << self_broadening << ") * " << modelparameters2metadata(m.Data().front().Get(var), T0);
681 else
682 os << "VMR(" << toShortName(sts[j]) << ") * " << modelparameters2metadata(m.Data()[j].Get(var), T0);
683
684 if (sts[j] not_eq sts.back()) os << " + ";
685 }
686 as.push_back(os.str());
687 }
688 }
689
690 return as;
691}
692
693
694namespace LineShape {
695#pragma GCC diagnostic push
696#pragma GCC diagnostic ignored "-Wreturn-type"
698 if (type == "X0")
699 return mp.X0;
700 if (type == "X1")
701 return mp.X1;
702 if (type == "X2")
703 return mp.X2;
704 if (type == "X3")
705 return mp.X3;
707 "Type: ", type, ", is not accepted. "
708 "See documentation for accepted types\n")
709}
710#pragma GCC diagnostic pop
711
712std::ostream& operator<<(std::ostream& os, const ModelParameters& mp) {
713 os << mp.type << ' ' << mp.X0 << ' ' << mp.X1 << ' '
714 << mp.X2 << ' ' << mp.X3 << ' ';
715 return os;
716}
717
718std::istream& operator>>(std::istream& is, ModelParameters& mp) {
719 is >> mp.type >> double_imanip() >> mp.X0 >> mp.X1 >> mp.X2 >> mp.X3;
720 return is;
721}
722
724 using std::log;
725 using std::pow;
726
727 switch (type) {
729 return 0;
730 case TemperatureModel::T0:
731 return X0;
732 case TemperatureModel::T1:
733 return X0 * pow(T0 / T, X1);
734 case TemperatureModel::T2:
735 return X0 * pow(T0 / T, X1) * (1 + X2 * log(T / T0));
736 case TemperatureModel::T3:
737 return X0 + X1 * (T - T0);
738 case TemperatureModel::T4:
739 return (X0 + X1 * (T0 / T - 1.)) * pow(T0 / T, X2);
740 case TemperatureModel::T5:
741 return X0 * pow(T0 / T, 0.25 + 1.5 * X1);
742 case TemperatureModel::LM_AER:
743 return special_linemixing_aer(T);
744 case TemperatureModel::DPL:
745 return X0 * pow(T0 / T, X1) + X2 * pow(T0 / T, X3);
746 case TemperatureModel::POLY:
747 return X0 + X1 * T + X2 * T * T + X3 * T * T * T;
748 case TemperatureModel::FINAL: { /* Leave last */ }
749 }
750 return std::numeric_limits<Numeric>::quiet_NaN();
751}
752
754 using std::log;
755 using std::pow;
756
757 switch (type) {
759 return 0;
760 case TemperatureModel::T0:
761 return 1;
762 case TemperatureModel::T1:
763 return pow(T0 / T, X1);
764 case TemperatureModel::T2:
765 return pow(T0 / T, X1) * (1 + X2 * log(T / T0));
766 case TemperatureModel::T3:
767 return 1;
768 case TemperatureModel::T4:
769 return pow(T0 / T, X2);
770 case TemperatureModel::T5:
771 return pow(T0 / T, 1.5 * X1 + 0.25);
772 case TemperatureModel::LM_AER:
773 return special_linemixing_aer_dX0(T);
774 case TemperatureModel::DPL:
775 return pow(T0 / T, X1);
776 case TemperatureModel::POLY:
777 return 1;
778 case TemperatureModel::FINAL: { /* Leave last */ }
779 }
780 return std::numeric_limits<Numeric>::quiet_NaN();
781}
782
784 using std::log;
785 using std::pow;
786
787 switch (type) {
789 return 0;
790 case TemperatureModel::T0:
791 return 0;
792 case TemperatureModel::T1:
793 return X0 * pow(T0 / T, X1) * log(T0 / T);
794 case TemperatureModel::T2:
795 return X0 * pow(T0 / T, X1) * (X2 * log(T / T0) + 1.) * log(T0 / T);
796 case TemperatureModel::T3:
797 return (T - T0);
798 case TemperatureModel::T4:
799 return pow(T0 / T, X2) * (T0 / T - 1.);
800 case TemperatureModel::T5:
801 return 1.5 * X0 * pow(T0 / T, 1.5 * X1 + 0.25) * log(T0 / T);
802 case TemperatureModel::LM_AER:
803 return special_linemixing_aer_dX1(T);
804 case TemperatureModel::DPL:
805 return X0 * pow(T0 / T, X1) * log(T0 / T);
806 case TemperatureModel::POLY:
807 return T;
808 case TemperatureModel::FINAL: {/* Leave last */ }
809 }
810 return std::numeric_limits<Numeric>::quiet_NaN();
811}
812
814 using std::log;
815 using std::pow;
816
817 switch (type) {
819 return 0;
820 case TemperatureModel::T0:
821 return 0;
822 case TemperatureModel::T1:
823 return 0;
824 case TemperatureModel::T2:
825 return X0 * pow(T0 / T, X1) * log(T / T0);
826 case TemperatureModel::T3:
827 return 0;
828 case TemperatureModel::T4:
829 return pow(T0 / T, X2) * (X0 + X1 * (T0 / T - 1)) * log(T0 / T);
830 case TemperatureModel::T5:
831 return 0;
832 case TemperatureModel::LM_AER:
833 return special_linemixing_aer_dX2(T);
834 case TemperatureModel::DPL:
835 return pow(T0 / T, X3);
836 case TemperatureModel::POLY:
837 return T * T;
838 case TemperatureModel::FINAL: {/* Leave last */ }
839 }
840 return std::numeric_limits<Numeric>::quiet_NaN();
841}
842
844 using std::log;
845 using std::pow;
846
847 switch (type) {
849 return 0;
850 case TemperatureModel::T0:
851 return 0;
852 case TemperatureModel::T1:
853 return 0;
854 case TemperatureModel::T2:
855 return 0;
856 case TemperatureModel::T3:
857 return 0;
858 case TemperatureModel::T4:
859 return 0;
860 case TemperatureModel::T5:
861 return 0;
862 case TemperatureModel::LM_AER:
863 return special_linemixing_aer_dX3(T);
864 case TemperatureModel::DPL:
865 return X2 * pow(T0 / T, X3) * log(T0 / T);
866 case TemperatureModel::POLY:
867 return T * T * T;
868 case TemperatureModel::FINAL: {/* Leave last */ }
869 }
870 return std::numeric_limits<Numeric>::quiet_NaN();
871}
872
874 using std::log;
875 using std::pow;
876
877 switch (type) {
879 return 0;
880 case TemperatureModel::T0:
881 return 0;
882 case TemperatureModel::T1:
883 return -X0 * X1 * pow(T0 / T, X1) / T;
884 case TemperatureModel::T2:
885 return -X0 * X1 * pow(T0 / T, X1) * (X2 * log(T / T0) + 1.) / T +
886 X0 * X2 * pow(T0 / T, X1) / T;
887 case TemperatureModel::T3:
888 return X1;
889 case TemperatureModel::T4:
890 return -X2 * pow(T0 / T, X2) * (X0 + X1 * (T0 / T - 1.)) / T -
891 T0 * X1 * pow(T0 / T, X2) / pow(T, 2);
892 case TemperatureModel::T5:
893 return -X0 * pow(T0 / T, 1.5 * X1 + 0.25) * (1.5 * X1 + 0.25) / T;
894 case TemperatureModel::LM_AER:
895 return special_linemixing_aer_dT(T);
896 case TemperatureModel::DPL:
897 return -X0 * X1 * pow(T0 / T, X1) / T + -X2 * X3 * pow(T0 / T, X3) / T;
898 case TemperatureModel::POLY:
899 return X1 + 2 * X2 * T + 3 * X3 * T * T;
900 case TemperatureModel::FINAL: {/* Leave last */ }
901 }
902 return std::numeric_limits<Numeric>::quiet_NaN();
903}
904
906 using std::log;
907 using std::pow;
908
909 switch (type) {
911 return 0;
912 case TemperatureModel::T0:
913 return 0;
914 case TemperatureModel::T1:
915 return X0 * X1 * pow(T0 / T, X1) / T0;
916 case TemperatureModel::T2:
917 return X0 * X1 * pow(T0 / T, X1) * (X2 * log(T / T0) + 1.) / T0 -
918 X0 * X2 * pow(T0 / T, X1) / T0;
919 case TemperatureModel::T3:
920 return -X1;
921 case TemperatureModel::T4:
922 return X2 * pow(T0 / T, X2) * (X0 + X1 * (T0 / T - 1.)) / T0 +
923 X1 * pow(T0 / T, X2) / T;
924 case TemperatureModel::T5:
925 return X0 * pow(T0 / T, 1.5 * X1 + 0.25) * (1.5 * X1 + 0.25) / T0;
926 case TemperatureModel::LM_AER:
927 return 0;
928 case TemperatureModel::DPL:
929 return X0 * X1 * pow(T0 / T, X1) / T0 + X2 * X3 * pow(T0 / T, X3) / T0;
930 case TemperatureModel::POLY:
931 return 0;
932 case TemperatureModel::FINAL: {/* Leave last */ }
933 }
934 return std::numeric_limits<Numeric>::quiet_NaN();
935}
936
938 for (auto& data: X) {
939 Index x;
940 bif >> x >> data.X0 >> data.X1 >> data.X2 >> data.X3;
941 data.type = TemperatureModel(x);
942 }
943 return bif;
944}
945
947 for (auto& data: X) {
948 bof << Index(data.type) << data.X0 << data.X1 << data.X2 << data.X3;
949 }
950 return bof;
951}
952
953bool SingleSpeciesModel::MatchTypes(const SingleSpeciesModel& other) const noexcept {
954 return std::equal (X.cbegin(), X.cend(), other.X.cbegin(), other.X.cend(), [](const auto& a, const auto& b){return a.type == b.type;});
955}
956
957std::ostream& operator<<(std::ostream& os, const SingleSpeciesModel& ssm) {
958 for (const auto& mp : ssm.Data())
959 if (mp.type not_eq TemperatureModel::None)
960 os << mp.X0 << ' ' << mp.X1 << ' ' << mp.X2 << ' ' << mp.X3 << ' ';
961 return os;
962}
963
964std::istream& operator>>(std::istream& is, SingleSpeciesModel& ssm) {
965 for (auto& mp : ssm.Data())
966 if(mp.type not_eq TemperatureModel::None)
967 is >> double_imanip() >> mp.X0 >> mp.X1 >> mp.X2 >> mp.X3;
968 return is;
969}
970
971std::ostream& operator<<(std::ostream& os, Output x) {
972 return os << "G0: " << x.G0 << " D0: " << x.D0 << " G2: " << x.G2
973 << " D2: " << x.D2 << " FVC: " << x.FVC << " ETA: " << x.ETA
974 << " Y: " << x.Y << " G: " << x.G << " DV: " << x.DV;
975}
976
978 Numeric nself,
979 Numeric agam,
980 Numeric nair,
981 Numeric psf,
982 std::array<Numeric, 12> aer_interp) noexcept : mdata(2) {
983 mdata.front().G0() = {TemperatureModel::T1, sgam, nself, 0, 0};
984 mdata.front().D0() = {TemperatureModel::T5, psf, nair, 0, 0};
985
986 mdata.back().G0() = {TemperatureModel::T1, agam, nair, 0, 0};
987 mdata.back().D0() = {TemperatureModel::T5, psf, nair, 0, 0};
988
989 if (std::any_of(aer_interp.cbegin(), aer_interp.cend(), [](auto x){return x not_eq 0;})) {
990 mdata.front().Y().type = TemperatureModel::LM_AER;
991 mdata.front().Y().X0 = aer_interp[4];
992 mdata.front().Y().X1 = aer_interp[5];
993 mdata.front().Y().X2 = aer_interp[6];
994 mdata.front().Y().X3 = aer_interp[7];
995 mdata.front().G().type = TemperatureModel::LM_AER;
996 mdata.front().G().X0 = aer_interp[8];
997 mdata.front().G().X1 = aer_interp[9];
998 mdata.front().G().X2 = aer_interp[10];
999 mdata.front().G().X3 = aer_interp[11];
1000
1001 mdata.back().Y().type = TemperatureModel::LM_AER;
1002 mdata.back().Y().X0 = aer_interp[4];
1003 mdata.back().Y().X1 = aer_interp[5];
1004 mdata.back().Y().X2 = aer_interp[6];
1005 mdata.back().Y().X3 = aer_interp[7];
1006 mdata.back().G().type = TemperatureModel::LM_AER;
1007 mdata.back().G().X0 = aer_interp[8];
1008 mdata.back().G().X1 = aer_interp[9];
1009 mdata.back().G().X2 = aer_interp[10];
1010 mdata.back().G().X3 = aer_interp[11];
1011 }
1012}
1013
1015 Numeric nself,
1016 Numeric agam,
1017 Numeric nair,
1018 Numeric psf) {
1019 Model m(2);
1020
1021 m.Data().front().G0() = {TemperatureModel::T1, sgam, nself, 0, 0};
1022 m.Data().front().D0() = {TemperatureModel::T0, psf, 0, 0, 0};
1023
1024 m.Data().back().G0() = {TemperatureModel::T1, agam, nair, 0, 0};
1025 m.Data().back().D0() = {TemperatureModel::T0, psf, 0, 0, 0};
1026
1027 return m;
1028}
1029
1031 Numeric nself,
1032 Numeric agam,
1033 Numeric nair,
1034 Numeric psf,
1035 std::array<Numeric, 12> aer_interp) {
1036 Model m(2);
1037
1038 m.Data().front().G0() = {TemperatureModel::T1, sgam, nself, 0, 0};
1039 m.Data().front().D0() = {TemperatureModel::T0, psf, 0, 0, 0};
1040
1041 m.Data().back().G0() = {TemperatureModel::T1, agam, nair, 0, 0};
1042 m.Data().back().D0() = {TemperatureModel::T0, psf, 0, 0, 0};
1043
1044 if (std::any_of(aer_interp.cbegin(), aer_interp.cend(), [](auto x){return x not_eq 0;})) {
1045 m.Data().front().Y().type = TemperatureModel::LM_AER;
1046 m.Data().front().Y().X0 = aer_interp[4];
1047 m.Data().front().Y().X1 = aer_interp[5];
1048 m.Data().front().Y().X2 = aer_interp[6];
1049 m.Data().front().Y().X3 = aer_interp[7];
1050 m.Data().front().G().type = TemperatureModel::LM_AER;
1051 m.Data().front().G().X0 = aer_interp[8];
1052 m.Data().front().G().X1 = aer_interp[9];
1053 m.Data().front().G().X2 = aer_interp[10];
1054 m.Data().front().G().X3 = aer_interp[11];
1055
1056 m.Data().back().Y().type = TemperatureModel::LM_AER;
1057 m.Data().back().Y().X0 = aer_interp[4];
1058 m.Data().back().Y().X1 = aer_interp[5];
1059 m.Data().back().Y().X2 = aer_interp[6];
1060 m.Data().back().Y().X3 = aer_interp[7];
1061 m.Data().back().G().type = TemperatureModel::LM_AER;
1062 m.Data().back().G().X0 = aer_interp[8];
1063 m.Data().back().G().X1 = aer_interp[9];
1064 m.Data().back().G().X2 = aer_interp[10];
1065 m.Data().back().G().X3 = aer_interp[11];
1066 }
1067
1068 return m;
1069}
1070
1071bool Model::OK(Type type, bool self, bool bath,
1072 const std::size_t nspecies) const noexcept {
1073 Index n = mdata.size();
1074 Index m = Index(self) + Index(bath);
1075 bool needs_any = type not_eq Type::DP;
1076 return not (n not_eq Index(nspecies) or m > n or (needs_any and not n));
1077}
1078
1079#define LSPC(XVAR, PVAR) \
1080 Numeric Model::XVAR( \
1081 Numeric T, Numeric T0, Numeric P [[maybe_unused]], \
1082 const ConstVectorView& vmrs) \
1083 const noexcept { \
1084 return PVAR * \
1085 std::inner_product( \
1086 mdata.cbegin(), \
1087 mdata.cend(), \
1088 vmrs.begin(), \
1089 0.0, \
1090 std::plus<>(), \
1091 [&](auto& x, auto vmr) -> Numeric { \
1092 return vmr * x.XVAR().at(T, T0); \
1093 }); \
1094 }
1095 LSPC(G0, P)
1096 LSPC(D0, P)
1097 LSPC(G2, P)
1098 LSPC(D2, P)
1099 LSPC(FVC, P)
1100 LSPC(ETA, 1)
1101 LSPC(Y, P)
1102 LSPC(G, P* P)
1103 LSPC(DV, P* P)
1104#undef LSPC
1105
1106#define LSPCV(XVAR, PVAR) \
1107 Numeric Model::d##XVAR##_dVMR(Numeric T, \
1108 Numeric T0, \
1109 Numeric P [[maybe_unused]], \
1110 const Index deriv_pos) const noexcept { \
1111 if (deriv_pos not_eq -1) \
1112 return PVAR * mdata[deriv_pos].XVAR().at(T, T0); \
1113 return 0; \
1114 }
1115 LSPCV(G0, P)
1116 LSPCV(D0, P)
1117 LSPCV(G2, P)
1118 LSPCV(D2, P)
1119 LSPCV(FVC, P)
1120 LSPCV(ETA, 1)
1121 LSPCV(Y, P)
1122 LSPCV(G, P* P)
1123 LSPCV(DV, P* P)
1124#undef LSPCV
1125
1126#define LSPCT(XVAR, PVAR) \
1127 Numeric Model::d##XVAR##_dT( \
1128 Numeric T, Numeric T0, Numeric P [[maybe_unused]], \
1129 const ConstVectorView& vmrs) \
1130 const noexcept { \
1131 return PVAR * \
1132 std::inner_product( \
1133 mdata.cbegin(), \
1134 mdata.cend(), \
1135 vmrs.begin(), \
1136 0.0, \
1137 std::plus<>(), \
1138 [&](auto& x, auto vmr) -> Numeric { \
1139 return vmr * x.XVAR().dT(T, T0); \
1140 }); \
1141 }
1142 LSPCT(G0, P)
1143 LSPCT(D0, P)
1144 LSPCT(G2, P)
1145 LSPCT(D2, P)
1146 LSPCT(FVC, P)
1147 LSPCT(ETA, 1)
1148 LSPCT(Y, P)
1149 LSPCT(G, P* P)
1150 LSPCT(DV, P* P)
1151#undef LSPCT
1152
1153#define LSPDC(XVAR, DERIV, PVAR) \
1154 Numeric Model::d##XVAR##_##DERIV(Numeric T, \
1155 Numeric T0, \
1156 Numeric P [[maybe_unused]], \
1157 Index deriv_pos, \
1158 const ConstVectorView& vmrs) const noexcept { \
1159 if (deriv_pos not_eq -1) \
1160 return vmrs[deriv_pos] * PVAR * \
1161 mdata[deriv_pos].XVAR().DERIV(T, T0); \
1162 return 0; \
1163 }
1164 LSPDC(G0, dT0, P)
1165 LSPDC(G0, dX0, P)
1166 LSPDC(G0, dX1, P)
1167 LSPDC(G0, dX2, P)
1168 LSPDC(G0, dX3, P)
1169 LSPDC(D0, dT0, P)
1170 LSPDC(D0, dX0, P)
1171 LSPDC(D0, dX1, P)
1172 LSPDC(D0, dX2, P)
1173 LSPDC(D0, dX3, P)
1174 LSPDC(G2, dT0, P)
1175 LSPDC(G2, dX0, P)
1176 LSPDC(G2, dX1, P)
1177 LSPDC(G2, dX2, P)
1178 LSPDC(G2, dX3, P)
1179 LSPDC(D2, dT0, P)
1180 LSPDC(D2, dX0, P)
1181 LSPDC(D2, dX1, P)
1182 LSPDC(D2, dX2, P)
1183 LSPDC(D2, dX3, P)
1184 LSPDC(FVC, dT0, P)
1185 LSPDC(FVC, dX0, P)
1186 LSPDC(FVC, dX1, P)
1187 LSPDC(FVC, dX2, P)
1188 LSPDC(FVC, dX3, P)
1189 LSPDC(ETA, dT0, 1)
1190 LSPDC(ETA, dX0, 1)
1191 LSPDC(ETA, dX1, 1)
1192 LSPDC(ETA, dX2, 1)
1193 LSPDC(ETA, dX3, 1)
1194 LSPDC(Y, dT0, P)
1195 LSPDC(Y, dX0, P)
1196 LSPDC(Y, dX1, P)
1197 LSPDC(Y, dX2, P)
1198 LSPDC(Y, dX3, P)
1199 LSPDC(G, dT0, P* P)
1200 LSPDC(G, dX0, P* P)
1201 LSPDC(G, dX1, P* P)
1202 LSPDC(G, dX2, P* P)
1203 LSPDC(G, dX3, P* P)
1204 LSPDC(DV, dT0, P* P)
1205 LSPDC(DV, dX0, P* P)
1206 LSPDC(DV, dX1, P* P)
1207 LSPDC(DV, dX2, P* P)
1208 LSPDC(DV, dX3, P* P)
1209#undef LSPDC
1210
1212 Numeric T0,
1213 Numeric P,
1214 const ConstVectorView& vmrs) const noexcept {
1215 return {G0(T, T0, P, vmrs),
1216 D0(T, T0, P, vmrs),
1217 G2(T, T0, P, vmrs),
1218 D2(T, T0, P, vmrs),
1219 FVC(T, T0, P, vmrs),
1220 ETA(T, T0, P, vmrs),
1221 Y(T, T0, P, vmrs),
1222 G(T, T0, P, vmrs),
1223 DV(T, T0, P, vmrs)};
1224}
1225
1227 Numeric T0,
1228 Numeric P,
1229 size_t k) const noexcept {
1230 return {P * mdata[k].G0().at(T, T0),
1231 P * mdata[k].D0().at(T, T0),
1232 P * mdata[k].G2().at(T, T0),
1233 P * mdata[k].D2().at(T, T0),
1234 P * mdata[k].FVC().at(T, T0),
1235 mdata[k].ETA().at(T, T0),
1236 P * mdata[k].Y().at(T, T0),
1237 P * P * mdata[k].G().at(T, T0),
1238 P * P * mdata[k].DV().at(T, T0)};
1239}
1240
1242 Numeric T0,
1243 Numeric P,
1244 ConstVectorView vmrs) const noexcept {
1245 return {dG0_dT(T, T0, P, vmrs),
1246 dD0_dT(T, T0, P, vmrs),
1247 dG2_dT(T, T0, P, vmrs),
1248 dD2_dT(T, T0, P, vmrs),
1249 dFVC_dT(T, T0, P, vmrs),
1250 dETA_dT(T, T0, P, vmrs),
1251 dY_dT(T, T0, P, vmrs),
1252 dG_dT(T, T0, P, vmrs),
1253 dDV_dT(T, T0, P, vmrs)};
1254}
1255
1256Output Model::GetVMRDerivs(Numeric T, Numeric T0, Numeric P, const Index pos) const noexcept {
1257 return {dG0_dVMR(T, T0, P, pos),
1258 dD0_dVMR(T, T0, P, pos),
1259 dG2_dVMR(T, T0, P, pos),
1260 dD2_dVMR(T, T0, P, pos),
1261 dFVC_dVMR(T, T0, P, pos),
1262 dETA_dVMR(T, T0, P, pos),
1263 dY_dVMR(T, T0, P, pos),
1264 dG_dVMR(T, T0, P, pos),
1265 dDV_dVMR(T, T0, P, pos)};
1266}
1267
1269 Numeric T0,
1270 Numeric P,
1271 Index pos,
1272 const Vector& vmrs,
1273 Jacobian::Line deriv) const noexcept {
1274 if (pos < 0) return 0;
1275
1276#define RETURNINTERNALDERIVATIVE(TYPE) \
1277 case Jacobian::Line::Shape##TYPE##X0: \
1278 return d##TYPE##_dX0(T, T0, P, pos, vmrs); \
1279 case Jacobian::Line::Shape##TYPE##X1: \
1280 return d##TYPE##_dX1(T, T0, P, pos, vmrs); \
1281 case Jacobian::Line::Shape##TYPE##X2: \
1282 return d##TYPE##_dX2(T, T0, P, pos, vmrs); \
1283 case Jacobian::Line::Shape##TYPE##X3: \
1284 return d##TYPE##_dX3(T, T0, P, pos, vmrs)
1285 switch (deriv) {
1295 default:
1296 return std::numeric_limits<Numeric>::quiet_NaN();
1297 }
1298#undef RETURNINTERNALDERIVATIVE
1299}
1300
1302 mdata.erase(mdata.begin() + i);
1303 specs.erase(specs.begin() + i);
1304}
1305
1307 mdata.erase(mdata.begin() + i);
1308 specs.erase(specs.begin() + i);
1309}
1310
1312 for (auto& ssm : mdata) {
1313 ssm.Y() = x.Y();
1314 ssm.G() = x.G();
1315 ssm.DV() = x.DV();
1316 }
1317}
1318
1319bool Model::Match(const Model& other) const noexcept {
1320 return std::equal (mdata.cbegin(), mdata.cend(), other.mdata.cbegin(), other.mdata.cend(), [](auto& a, auto& b) {return a.MatchTypes(b);});
1321}
1322
1324 for (auto& data: mdata)
1325 data.read(bif);
1326 return bif;
1327}
1328
1330 for (auto& data: mdata)
1331 data.write(bof);
1332 return bof;
1333}
1334
1335namespace LegacyLineFunctionData {
1337 if (type == String("DP"))
1338 return {};
1339 if (type == String("LP"))
1340 return {Variable::G0, Variable::D0};
1341 if (type == String("VP"))
1342 return {Variable::G0, Variable::D0};
1343 if (type == String("SDVP"))
1345 if (type == String("HTP"))
1346 return {Variable::G0,
1353 "Type: ", type, ", is not accepted. "
1354 "See documentation for accepted types\n")
1355}
1356
1358 if (type == "#")
1359 return {};
1360 if (type == "LM1")
1361 return {Variable::Y};
1362 if (type == "LM2")
1363 return {Variable::Y, Variable::G, Variable::DV};
1364 if (type == "INT")
1365 return {};
1366 if (type == "ConstG")
1367 return {Variable::G};
1369 "Type: ", type, ", is not accepted. "
1370 "See documentation for accepted types\n")
1371}
1372} // namespace LegacyLineFunctionData
1373
1374namespace LegacyLineMixingData {
1376 if (type == "NA") // The standard case
1377 return TypeLM::LM_NONE;
1378 if (type == "LL") // The LBLRTM case
1379 return TypeLM::LM_LBLRTM;
1380 if (type == "NR") // The LBLRTM O2 non-resonant case
1382 if (type == "L2") // The 2nd order case
1383 return TypeLM::LM_2NDORDER;
1384 if (type == "L1") // The 2nd order case
1385 return TypeLM::LM_1STORDER;
1386 if (type == "BB") // The band class
1387 return TypeLM::LM_BYBAND;
1389 "Type: ", type, ", is not accepted. "
1390 "See documentation for accepted types\n")
1391}
1392} // namespace LegacyLineMixingData
1393
1394
1395namespace LegacyPressureBroadeningData {
1397 if (type == "NA") // The none case
1398 return TypePB::PB_NONE;
1399 if (type == "N2") // Air Broadening is N2 broadening mostly...
1401 if (type == "WA") // Water and Air Broadening
1403 if (type == "AP") // Planetary broadening
1406 "Type: ", type, ", is not accepted. "
1407 "See documentation for accepted types\n")
1408}
1409
1413 (qid.Species() == Species::fromShortName("N2") or
1414 qid.Species() == Species::fromShortName("O2") or
1415 qid.Species() == Species::fromShortName("H2O") or
1416 qid.Species() == Species::fromShortName("CO2") or
1417 qid.Species() == Species::fromShortName("H2") or
1418 qid.Species() == Species::fromShortName("He")))
1419 return true;
1420 if (t == TypePB::PB_AIR_AND_WATER_BROADENING and qid.Species() == Species::fromShortName("H2O"))
1421 return true;
1422 return false;
1423}
1424} // namespace LegacyPressureBroadeningData
1425} // namespace LineShape
void * data
type upp char * str
Index nelem() const ARTS_NOEXCEPT
Number of elements.
Definition: array.h:195
A constant view of a Vector.
Definition: matpackI.h:489
Numeric sum() const ARTS_NOEXCEPT
The sum of all elements of a Vector.
Definition: matpackI.cc:52
Main line shape model class.
void SetLineMixingModel(SingleSpeciesModel x)
Sets the same line mixing model to all species.
Output GetVMRDerivs(Numeric T, Numeric T0, Numeric P, const Index pos) const noexcept
Derivative of GetParams(...) wrt VMR.
void Remove(Index i, ArrayOfSpeciesTag &specs)
Remove species and data at position.
bool Match(const Model &other) const noexcept
std::vector< SingleSpeciesModel > mdata
Model(Index n=0) noexcept
Default init just sets the size.
Numeric GetInternalDeriv(Numeric T, Numeric T0, Numeric P, Index pos, const Vector &vmrs, Jacobian::Line deriv) const noexcept
Derivative of GetParams(...) wrt Coefficient.
Output GetTemperatureDerivs(Numeric T, Numeric T0, Numeric P, ConstVectorView vmrs) const noexcept
Derivative of GetParams(...) wrt T.
bool OK(Type type, bool self, bool bath, const std::size_t nspecies) const noexcept
The Model is good to use.
bifstream & read(bifstream &bif)
Binary read for Model.
Output GetParams(Numeric T, Numeric T0, Numeric P, const ConstVectorView &vmrs) const noexcept
Compute all shape parameters.
bofstream & write(bofstream &bof) const
Binary write for Model.
const std::vector< SingleSpeciesModel > & Data() const noexcept
The line shape model data.
Compute the line shape parameters for a single broadening species.
constexpr std::array< ModelParameters, nVars > & Data() noexcept
Get internal Data reference.
std::array< ModelParameters, nVars > X
bofstream & write(bofstream &bof) const
Binary write for SingleSpeciesModel.
bifstream & read(bifstream &bif)
Binary read for SingleSpeciesModel.
bool MatchTypes(const SingleSpeciesModel &other) const noexcept
The Vector class.
Definition: matpackI.h:876
Binary output file stream class.
Definition: bifstream.h:42
Binary output file stream class.
Definition: bofstream.h:42
Input manipulator class for doubles to enable nan and inf parsing.
Definition: double_imanip.h:64
Index nelem() const
Number of elements.
Definition: mystring.h:253
Constants of physical expressions as constexpr.
#define ARTS_NOEXCEPT
Definition: debug.h:80
#define ARTS_ASSERT(condition,...)
Definition: debug.h:83
#define ARTS_USER_ERROR(...)
Definition: debug.h:150
#define ARTS_USER_ERROR_IF(condition,...)
Definition: debug.h:134
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
constexpr bool good_enum(EnumType x) noexcept
Checks if the enum number is good.
Definition: enums.h:21
Array< Species::Species > ArrayOfSpecies
#define temp
#define LSPCV(XVAR, PVAR)
#define ReturnJacPropMatType(ID)
#define RETURNINTERNALDERIVATIVE(TYPE)
#define LSPCT(XVAR, PVAR)
#define LSPC(XVAR, PVAR)
#define LSPDC(XVAR, DERIV, PVAR)
Jacobian::Line select_derivativeLineShape(const String &var, const String &coeff)
Return the derivative type based on string input.
Contains the line shape namespace.
Implementation of Matrix, Vector, and such stuff.
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
my_basic_string< char > String
The String type for ARTS.
Definition: mystring.h:287
char Type type
std::vector< Variable > linemixingtag2variablesvector(String type)
Line mixing models from string.
std::vector< Variable > lineshapetag2variablesvector(String type)
Line shape models from string.
constexpr Index temperaturemodel2legacynelem(TemperatureModel type) noexcept
Length per variable for temperature model.
LegacyLineMixingData::TypeLM string2typelm(String type)
Line mixing types from string.
constexpr Index typelm2nelem(LegacyLineMixingData::TypeLM type)
Line mixing types to number.
Model vector2modellm(Vector x, LegacyLineMixingData::TypeLM type)
LineShape::Model from legacy input vector.
TypeLM
Line mixing types that used to exist.
Index self_listed(const QuantumIdentifier &qid, LegacyPressureBroadeningData::TypePB t)
Pressure broadening if self exist.
void vector2modelpb(LineShape::Type &mtype, bool &self, bool &bath, Model &m, ArrayOfSpecies &species, Vector x, LegacyPressureBroadeningData::TypePB type, bool self_in_list, Species::Species self_spec)
LineShape::Model from legacy input vector.
constexpr Index typepb2nelem(LegacyPressureBroadeningData::TypePB type)
Pressure broadening types to number of elements.
LegacyPressureBroadeningData::TypePB string2typepb(String type)
Pressure broadening types from string.
TypePB
Pressure broadening types that used to exist.
Computations of line shape derived parameters.
Definition: lineshape.cc:23
ArrayOfString ModelMetaDataArray(const Model &m, const bool self, const ArrayOfSpecies &sts, const Numeric T0)
Model lblrtm_model(Numeric sgam, Numeric nself, Numeric agam, Numeric nair, Numeric psf, std::array< Numeric, 12 > aer_interp)
Numeric & SingleModelParameter(ModelParameters &mp, const String &type)
Get a coefficient from ModelParameters by name.
std::istream & from_linefunctiondata(std::istream &data, Type &type, bool &self, bool &bath, Model &m, ArrayOfSpecies &species)
Model hitran_model(Numeric sgam, Numeric nself, Numeric agam, Numeric nair, Numeric psf)
Vector mass(const ConstVectorView &atmospheric_vmrs, const ArrayOfArrayOfSpeciesTag &atmospheric_species, const ArrayOfSpecies &lineshape_species, const SpeciesIsotopologueRatios &ir) ARTS_NOEXCEPT
Returns a mass vector for this model's main calculations.
std::istream & from_pressurebroadeningdata(std::istream &data, LineShape::Type &type, bool &self, bool &bath, Model &m, ArrayOfSpecies &species, const QuantumIdentifier &qid)
Legacy reading of old deprecated PressureBroadeningData class.
std::istream & operator>>(std::istream &is, ModelParameters &mp)
Input operator for ModelParameters.
String modelparameters2metadata(const ModelParameters mp, const Numeric T0)
String ModelShape2MetaData(const Model &m)
std::ostream & operator<<(std::ostream &os, const ModelParameters &mp)
Output operator for ModelParameters.
Model MetaData2ModelShape(const String &s)
std::istream & from_artscat4(std::istream &is, Type &type, bool &self, bool &bath, Model &m, ArrayOfSpecies &species, const QuantumIdentifier &qid)
std::istream & from_linemixingdata(std::istream &data, Model &lsc)
Legacy reading of old deprecated LineMixingData class.
Vector vmrs(const ConstVectorView &atmospheric_vmrs, const ArrayOfArrayOfSpeciesTag &atmospheric_species, const ArrayOfSpecies &lineshape_species) ARTS_NOEXCEPT
Returns a VMR vector for this model's main calculations.
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 Numeric mean_mass(Species spec, const IsotopologueRatios &ir) noexcept
#define v
#define a
#define b
Quantum::Identifier QuantumIdentifier
Definition: quantum.h:471
Numeric pow(const Rational base, Numeric exp)
Power of.
Definition: rational.h:747
#define N
Definition: rng.cc:164
Coefficients and temperature model for SingleSpeciesModel.
Numeric dT0(Numeric T, Numeric T0) const noexcept
Numeric dX0(Numeric T, Numeric T0) const noexcept
Numeric dT(Numeric T, Numeric T0) const noexcept
Numeric at(Numeric T, Numeric T0) const noexcept
Numeric dX1(Numeric T, Numeric T0) const noexcept
Numeric dX3(Numeric T, Numeric T0) const noexcept
Numeric dX2(Numeric T, Numeric T0) const noexcept