63 for (
Index i = 0; i < n_new_grid; ++i) {
68 while (
abs(new_grid[i] - old_grid[j]) >
69 max(
abs(new_grid[i]),
abs(old_grid[j])) * DBL_EPSILON) {
71 if (j >= n_old_grid) {
73 os <<
"Cannot find new frequency " << i <<
" (" << new_grid[i]
74 <<
"Hz) in the lookup table frequency grid.";
75 throw runtime_error(os.str());
80 out3 <<
" " << new_grid[i] <<
" found, index = " << pos[i] <<
".\n";
124 const Index n_current_species = current_species.
nelem();
125 const Index n_current_f_grid = current_f_grid.
nelem();
133 out2 <<
" Original table: " << n_species <<
" species, " << n_f_grid
135 <<
" Adapt to: " << n_current_species <<
" species, "
136 << n_current_f_grid <<
" frequencies.\n";
139 out2 <<
" Table contains no nonlinear species.\n";
144 for (
Index s = 0; s < n_nls; ++s) {
149 out2 <<
" Table contains no temperature perturbations.\n";
161 if (0 == n_species) {
163 os <<
"The lookup table should have at least one species.";
164 throw runtime_error(os.str());
171 os <<
"The table must not have duplicate nonlinear species.\n"
173 throw runtime_error(os.str());
177 for (
Index i = 0; i < n_nls; ++i) {
179 os <<
"nonlinear_species[" << i <<
"]";
205 if (0 == n_nls_pert) {
207 os <<
"The vector nls_pert should contain the perturbations\n"
208 <<
"for the nonlinear species, but it is empty.";
209 throw runtime_error(os.str());
228 chk_size(
"xsec",
xsec, 1, n_species, n_f_grid, n_p_grid);
246 Index b = n_species + n_nls * (n_nls_pert - 1);
257 for (
Index i = 0, sp = 0; i < n_species; ++i) {
258 original_spec_pos_in_xsec[i] = sp;
268 if (0 == n_current_species) {
270 os <<
"The list of current species should not be empty.";
271 throw runtime_error(os.str());
281 out3 <<
" Looking for species in lookup table:\n";
282 for (
Index i = 0; i < n_current_species; ++i) {
283 out3 <<
" " << current_species[i].Name() <<
": ";
286 i_current_species[i] =
291 const auto spec_type = current_species[i].Type();
292 if (spec_type == Species::TagType::Zeeman ||
293 spec_type == Species::TagType::FreeElectrons ||
294 spec_type == Species::TagType::Particles) {
296 i_current_species[i] = -1;
299 os <<
"Species " << current_species[i].Name()
300 <<
" is missing in absorption lookup table.";
301 throw runtime_error(os.str());
305 out3 <<
"found (or trivial).\n";
309 Index n_current_nonlinear_species = 0;
315 out3 <<
" Finding out which of the current species are nonlinear:\n";
316 for (
Index i = 0; i < n_current_species; ++i) {
317 if (i_current_species[i] >= 0)
320 if (non_linear[i_current_species[i]]) {
321 out3 <<
" " << current_species[i] <<
"\n";
323 current_non_linear[i] = 1;
324 ++n_current_nonlinear_species;
337 out3 <<
" Looking for Frequencies in lookup table:\n";
343 i_current_f_grid,
f_grid, current_f_grid, verbosity);
349 new_table.
species.resize(n_current_species);
350 for (
Index i = 0; i < n_current_species; ++i) {
351 if (i_current_species[i] >= 0) {
360 new_table.
species[i] = current_species[i];
366 for (
Index i = 0; i < n_current_f_grid; ++i) {
376 for (
Index i = 0; i < n_current_species; ++i) {
377 if (i_current_species[i] >= 0) {
405 n_current_species + n_current_nonlinear_species * (n_nls_pert - 1),
414 for (
Index i_s = 0, sp = 0; i_s < n_current_species; ++i_s) {
417 if (current_non_linear[i_s])
426 for (
Index i_f = 0; i_f < n_current_f_grid; ++i_f) {
427 if (i_current_species[i_s] >= 0) {
430 Range(original_spec_pos_in_xsec[i_current_species[i_s]], n_v),
431 i_current_f_grid[i_f],
515 const Index& p_interp_order,
516 const Index& t_interp_order,
517 const Index& h2o_interp_order,
518 const Index& f_interp_order,
523 const Numeric& extpolfac)
const {
546 const Index n_new_f_grid = new_f_grid.
nelem();
557 Index h2o_index = -1;
565 if (h2o_index == -1) {
567 os <<
"With nonlinear species, at least one species must be a H2O species.";
568 throw runtime_error(os.str());
585 b = n_species + n_nls * (n_nls_pert - 1);
604 os <<
"The lookup table internal variable log_p_grid is not initialized.\n"
605 <<
"Use the abs_lookupAdapt method!";
606 throw runtime_error(os.str());
615 if ((n_p_grid < p_interp_order + 1)) {
617 os <<
"The number of pressure grid points in the table (" << n_p_grid
618 <<
") is not enough for the desired order of interpolation ("
619 << p_interp_order <<
").";
620 throw runtime_error(os.str());
623 if ((n_nls_pert != 0) && (n_nls_pert < h2o_interp_order + 1)) {
625 os <<
"The number of humidity perturbation grid points in the table ("
627 <<
") is not enough for the desired order of interpolation ("
628 << h2o_interp_order <<
").";
629 throw runtime_error(os.str());
632 if ((n_t_pert != 0) && (n_t_pert < t_interp_order + 1)) {
634 os <<
"The number of temperature perturbation grid points in the table ("
635 << n_t_pert <<
") is not enough for the desired order of interpolation ("
636 << t_interp_order <<
").";
637 throw runtime_error(os.str());
640 if ((n_f_grid < f_interp_order + 1)) {
642 os <<
"The number of frequency grid points in the table (" << n_f_grid
643 <<
") is not enough for the desired order of interpolation ("
644 << f_interp_order <<
").";
645 throw runtime_error(os.str());
651 if (!
is_size(abs_vmrs, n_species)) {
653 os <<
"Number of species in lookup table does not match number\n"
654 <<
"of species for which you want to extract absorption.\n"
655 <<
"Have you used abs_lookupAdapt? Or did you miss to add\n"
656 <<
"some VRM fields (e.g. for free electrons or particles)?\n";
657 throw runtime_error(os.str());
673 if (f_interp_order == 0) {
680 const Numeric allowed_f_margin = 0.09;
682 if (n_new_f_grid == n_f_grid) {
688 if (
abs(
f_grid[0] - new_f_grid[0]) > allowed_f_margin)
691 os <<
"First frequency in f_grid inconsistent with lookup table.\n"
692 <<
"f_grid[0] = " <<
f_grid[0] <<
"\n"
693 <<
"new_f_grid[0] = " << new_f_grid[0] <<
".";
694 throw runtime_error(os.str());
698 if (
abs(
f_grid[n_f_grid - 1] - new_f_grid[n_new_f_grid - 1]) > allowed_f_margin)
701 os <<
"Last frequency in f_grid inconsistent with lookup table.\n"
702 <<
"f_grid[n_f_grid-1] = " <<
f_grid[n_f_grid - 1]
704 <<
"new_f_grid[n_new_f_grid-1] = " << new_f_grid[n_new_f_grid - 1]
706 throw runtime_error(os.str());
708 }
else if (n_new_f_grid == 1) {
713 if (
abs(
f_grid[flag_local[0].pos] - new_f_grid[0]) > allowed_f_margin)
716 os <<
"Cannot find a matching lookup table frequency for frequency "
717 << new_f_grid[0] <<
".\n"
718 <<
"(This check has not been properly tested, so perhaps this is\n"
719 <<
"a false alarm. Check for this in file gas_abs_lookup.cc.)";
720 throw runtime_error(os.str());
724 os <<
"With f_interp_order 0 the frequency grid has to have the same\n"
725 <<
"size as in the lookup table, or exactly one element.";
726 throw runtime_error(os.str());
732 if (new_f_grid[0] < f_min) {
734 os <<
"Problem with gas absorption lookup table.\n"
735 <<
"At least one frequency is outside the range covered by the lookup table.\n"
736 <<
"Your new frequency value is " << new_f_grid[0] <<
" Hz.\n"
737 <<
"The allowed range is " << f_min <<
" to " << f_max <<
" Hz.";
738 throw runtime_error(os.str());
740 if (new_f_grid[n_new_f_grid - 1] > f_max) {
742 os <<
"Problem with gas absorption lookup table.\n"
743 <<
"At least one frequency is outside the range covered by the lookup table.\n"
744 <<
"Your new frequency value is " << new_f_grid[n_new_f_grid - 1]
746 <<
"The allowed range is " << f_min <<
" to " << f_max <<
" Hz.";
747 throw runtime_error(os.str());
759 const Index do_T = n_t_pert;
763 for (
Index s = 0; s < n_nls; ++s) {
779 if ((p > p_max) || (p < p_min)) {
781 os <<
"Problem with gas absorption lookup table.\n"
782 <<
"Pressure p is outside the range covered by the lookup table.\n"
783 <<
"Your p value is " << p <<
" Pa.\n"
784 <<
"The allowed range is " << p_min <<
" to " << p_max <<
".\n"
785 <<
"The pressure grid range in the table is " <<
p_grid[n_p_grid - 1]
786 <<
" to " <<
p_grid[0] <<
".\n"
787 <<
"We allow a bit of extrapolation, but NOT SO MUCH!";
788 throw runtime_error(os.str());
843 xsec_pre_interpolated.
resize(
844 p_interp_order + 1, n_species, 1, 1, n_new_f_grid);
852 for (
Index pi = 0; pi < p_interp_order + 1; ++pi) {
871 const Index this_p_grid_index = plag[0].pos + pi;
903 const Numeric effective_T_ref =
t_ref[this_p_grid_index];
906 const Numeric T_offset = T - effective_T_ref;
915 extpolfac * (
t_pert[n_t_pert - 1] -
t_pert[n_t_pert - 2]);
916 if ((T_offset > t_max) || (T_offset < t_min)) {
918 os <<
"Problem with gas absorption lookup table.\n"
919 <<
"Temperature T is outside the range covered by the lookup table.\n"
920 <<
"Your temperature was " << T <<
" K at a pressure of " << p
922 <<
"The temperature offset value is " << T_offset <<
".\n"
923 <<
"The allowed range is " << t_min <<
" to " << t_max <<
".\n"
924 <<
"The temperature perturbation grid range in the table is "
926 <<
"We allow a bit of extrapolation, but NOT SO MUCH!";
927 throw runtime_error(os.str());
957 const Numeric effective_vmr_ref =
vmrs_ref(h2o_index, this_p_grid_index);
960 const Numeric VMR_frac = abs_vmrs[h2o_index] / effective_vmr_ref;
971 if ((VMR_frac > x_max) || (VMR_frac < x_min)) {
973 os <<
"Problem with gas absorption lookup table.\n"
974 <<
"VMR for H2O (species " << h2o_index
975 <<
") is outside the range covered by the lookup table.\n"
976 <<
"Your VMR was " << abs_vmrs[h2o_index] <<
" at a pressure of "
978 <<
"The reference VMR value there is " << effective_vmr_ref <<
"\n"
979 <<
"The fractional VMR relative to the reference value is "
981 <<
"The allowed range is " << x_min <<
" to " << x_max <<
".\n"
982 <<
"The fractional VMR perturbation grid range in the table is "
984 <<
"We allow a bit of extrapolation, but NOT SO MUCH!";
985 throw runtime_error(os.str());
994 if (n_nls < n_species) {
1008 for (
Index si = 0; si < n_species; ++si) {
1011 const Index do_VMR = non_linear[si];
1024 os <<
"Problem with gas absorption lookup table.\n"
1025 <<
"VMR interpolation is not allowed for species \""
1026 <<
species[si][0].Name() <<
"\"";
1027 throw runtime_error(os.str());
1035 Index this_h2o_extent;
1038 this_h2o_extent = n_nls_pert;
1041 vlag = &lag_trivial;
1042 this_h2o_extent = 1;
1049 Range(fpi, this_h2o_extent),
1083 sga.
resize(n_species, n_new_f_grid);
1085 for (
Index pi = 0; pi < p_interp_order + 1; ++pi) {
1093 xsec_pre_interpolated(
1105 for (
Index si = 0; si < n_species; ++si) {
1106 if (select_abs_species.
nelem()) {
1107 if (
species[si] == select_abs_species)
1125 os <<
"GasAbsLookup: Output operator not implemented";
base max(const Array< base > &x)
Max function.
Index nelem() const ARTS_NOEXCEPT
A constant view of a Tensor3.
Index ncols() const noexcept
Index nbooks() const noexcept
Index npages() const noexcept
A constant view of a Vector.
Index nelem() const noexcept
Returns the number of elements.
An absorption lookup table.
Vector p_grid
The pressure grid for the table [Pa].
Tensor4 xsec
Absorption cross sections.
Vector log_p_grid
The natural log of the pressure grid.
Vector t_pert
The vector of temperature perturbations [K].
const Vector & GetPgrid() const
Vector f_grid
The frequency grid [Hz].
Matrix vmrs_ref
The reference VMR profiles.
const Vector & GetFgrid() const
ArrayOfLagrangeInterpolation flag_default
Frequency grid positions.
Vector nls_pert
The vector of perturbations for the VMRs of the nonlinear species.
void Adapt(const ArrayOfArrayOfSpeciesTag ¤t_species, ConstVectorView current_f_grid, const Verbosity &verbosity)
Adapt lookup table to current calculation.
Vector t_ref
The reference temperature profile [K].
void Extract(Matrix &sga, const ArrayOfSpeciesTag &select_abs_species, const Index &p_interp_order, const Index &t_interp_order, const Index &h2o_interp_order, const Index &f_interp_order, const Numeric &p, const Numeric &T, ConstVectorView abs_vmrs, ConstVectorView new_f_grid, const Numeric &extpolfac) const
Extract scalar gas absorption coefficients from the lookup table.
ArrayOfIndex nonlinear_species
The species tags with non-linear treatment.
ArrayOfArrayOfSpeciesTag species
The species tags for which the table is valid.
void resize(Index r, Index c)
Resize function.
void resize(Index b, Index p, Index r, Index c)
Resize function.
void resize(Index s, Index b, Index p, Index r, Index c)
Resize function.
void resize(Index n)
Resize function.
Subclasses of runtime_error.
#define ARTS_ASSERT(condition,...)
void find_new_grid_in_old_grid(ArrayOfIndex &pos, ConstVectorView old_grid, ConstVectorView new_grid, const Verbosity &verbosity)
Find positions of new grid points in old grid.
ostream & operator<<(ostream &os, const GasAbsLookup &)
Output operatior for GasAbsLookup.
Declarations for the gas absorption lookup table.
void interpweights(VectorView itw, const GridPos &tc)
Red 1D interpolation weights.
Header file for interpolation.cc.
Interpolation::Lagrange LagrangeInterpolation
bool is_size(ConstVectorView x, const Index &n)
Verifies that the size of x is l.
bool is_unique(const ArrayOfIndex &x)
Checks if an ArrayOfIndex is unique, i.e., has no duplicate values.
Header file for logic.cc.
void abs(Sparse &A, const Sparse &B)
Absolute value of sparse matrix elements.
void transform(VectorView y, double(&my_func)(double), ConstVectorView x)
A generic transform function for vectors, which can be used to implement mathematical functions opera...
NUMERIC Numeric
The type to use for all floating point numbers.
INDEX Index
The type to use for all integer numbers and indices.
Declarations having to do with the four output streams.
Array< Lagrange > LagrangeVector(const ConstVectorView &xs, const ConstVectorView &xi, const Index polyorder, const Numeric extrapol, const bool do_derivs, const GridType type, const std::pair< Numeric, Numeric > cycle)
Implements Zeeman modeling.
This file contains declerations of functions of physical character.
constexpr Numeric number_density(Numeric p, Numeric t) noexcept
number_density