Lmtest

doc/arts/references.bib

@@ -2326,3 +2326,50 @@ @article{Ellison2007
     url = {https://doi.org/10.1063/1.2360986},
     eprint = {https://pubs.aip.org/aip/jpr/article-pdf/36/1/1/14719718/1_1_online.pdf},
 }
+
+@article{NIRO2004483,
+title = {Spectra calculations in central and wing regions of CO2 IR bands between 10 and 20μm. I: model and laboratory measurements},
+journal = {Journal of Quantitative Spectroscopy and Radiative Transfer},
+volume = {88},
+number = {4},
+pages = {483-498},
+year = {2004},
+issn = {0022-4073},
+doi = {https://doi.org/10.1016/j.jqsrt.2004.04.003},
+url = {https://www.sciencedirect.com/science/article/pii/S0022407304001049},
+author = {F Niro and C Boulet and J.-M Hartmann},
+keywords = {CO, Infrared, Absorption, Shape, Model, Laboratory, Spectra},
+}
+
+@article{Pieroni1999,
+    author = {Pieroni, D. and Nguyen-Van-Thanh and Brodbeck, C. and
+    Claveau, C. and Valentin, A. and Hartmann, J. M. and Gabard, T. and
+    Champion, J.-P. and Bermejo, D. and Domenech, J.-L.},
+    title = {Experimental and theoretical study of line mixing in methane spectra. I.
+    The N$_2$-broadened $\nu_3$ band at room temperature},
+    journal = {The Journal of Chemical Physics},
+    volume = {110},
+    number = {16},
+    pages = {7717-7732},
+    year = {1999},
+    month = {04},
+    issn = {0021-9606},
+    doi = {10.1063/1.478724},
+    url = {https://doi.org/10.1063/1.478724},
+    eprint = {https://pubs.aip.org/aip/jcp/article-pdf/110/16/7717/19328457/7717_1_online.pdf},
+}
+
+@article{Hadded2002,
+    author = {Hadded, S. and Thibault, F. and Flaud, P.-M. and Aroui, H. and Hartmann, J.-M.},
+    title = {Experimental and theoretical study of line mixing in NH$_3$ spectra. I. Scaling analysis of parallel bands perturbed by He},
+    journal = {The Journal of Chemical Physics},
+    volume = {116},
+    number = {17},
+    pages = {7544-7557},
+    year = {2002},
+    month = {05},
+    issn = {0021-9606},
+    doi = {10.1063/1.1463442},
+    url = {https://doi.org/10.1063/1.1463442},
+    eprint = {https://pubs.aip.org/aip/jcp/article-pdf/116/17/7544/19021600/7544_1_online.pdf},
+}

src/core/lbl/CMakeLists.txt

@@ -8,6 +8,8 @@ add_library(lbl STATIC
   lbl_lineshape_voigt_ecs.cpp
   lbl_lineshape_voigt_ecs_hartmann.cpp
   lbl_lineshape_voigt_ecs_makarov.cpp
+  lbl_lineshape_voigt_ecs_stotop.cpp
+  lbl_lineshape_voigt_ecs_sphtop.cpp
   lbl_lineshape_voigt_lte.cpp
   lbl_lineshape_voigt_lte_mirrored.cpp
   lbl_lineshape_voigt_lte_matrix.cpp

src/core/lbl/lbl_lineshape.cpp

@@ -62,7 +62,9 @@ std::unique_ptr<voigt::ecs::ComputeData> init_voigt_abs_ecs_data(
     const Vector2 los) {
   if (stdr::any_of(bnds | stdv::values, [](auto& bnd) {
         return bnd.lineshape == LineByLineLineshape::VP_ECS_MAKAROV or
-               bnd.lineshape == LineByLineLineshape::VP_ECS_HARTMANN;
+               bnd.lineshape == LineByLineLineshape::VP_ECS_HARTMANN or
+               bnd.lineshape == LineByLineLineshape::VP_ECS_STOTOP or
+               bnd.lineshape == LineByLineLineshape::VP_ECS_SPHTOP;
       }))
     return std::make_unique<voigt::ecs::ComputeData>(
         f_grid, atm, los, ZeemanPolarization::no);
@@ -177,7 +179,9 @@ void calculate(PropmatVectorView pm,
         calc_voigt_line_nlte(bnd_key, bnd, pol);
         break;
       case LineByLineLineshape::VP_ECS_MAKAROV: [[fallthrough]];
-      case LineByLineLineshape::VP_ECS_HARTMANN:
+      case LineByLineLineshape::VP_ECS_HARTMANN: [[fallthrough]];
+      case LineByLineLineshape::VP_ECS_STOTOP: [[fallthrough]];
+      case LineByLineLineshape::VP_ECS_SPHTOP:
         calc_voigt_ecs_linemixing(bnd_key, bnd, pol);
         break;
     }

src/core/lbl/lbl_lineshape_voigt_ecs.cpp

@@ -21,6 +21,8 @@
 #include "lbl_lineshape_model.h"
 #include "lbl_lineshape_voigt_ecs_hartmann.h"
 #include "lbl_lineshape_voigt_ecs_makarov.h"
+#include "lbl_lineshape_voigt_ecs_sphtop.h"
+#include "lbl_lineshape_voigt_ecs_stotop.h"
 
 #undef WIGNER3
 #undef WIGNER6
@@ -187,6 +189,15 @@ void ComputeData::adapt_multi(const QuantumIdentifier& bnd_qid,
       auto& l2 = bnd_qid.state.at(QuantumNumberType::l2);
       dipr[i]  = hartmann::reduced_dipole(J.upper, J.lower, l2.upper, l2.lower);
       dip[i] *= std::signbit(dipr[i]) ? -1 : 1;
+    } else if (bnd.lineshape == LineByLineLineshape::VP_ECS_STOTOP) {
+      auto& J  = bnd.lines[i].qn.at(QuantumNumberType::J);
+      auto& Kq = bnd.lines[i].qn.at(QuantumNumberType::K);
+      dipr[i]  = stotop::reduced_dipole(J.upper, J.lower, Kq.lower);
+      dip[i]  *= std::signbit(dipr[i]) ? -1 : 1;
+    } else if (bnd.lineshape == LineByLineLineshape::VP_ECS_SPHTOP) {
+      auto& J  = bnd.lines[i].qn.at(QuantumNumberType::J);
+      dipr[i]  = sphtop::reduced_dipole(J.upper, J.lower);
+      dip[i]  *= std::signbit(dipr[i]) ? -1 : 1;
     }
   }
 
@@ -254,6 +265,12 @@ void ComputeData::adapt_multi(const QuantumIdentifier& bnd_qid,
     } else if (bnd.lineshape == LineByLineLineshape::VP_ECS_HARTMANN) {
       hartmann::relaxation_matrix_offdiagonal(
           Wimag, bnd_qid, bnd, sort, spec, rovib_data_it->second, dipr, atm);
+    } else if (bnd.lineshape == LineByLineLineshape::VP_ECS_STOTOP) {
+      stotop::relaxation_matrix_offdiagonal(
+          Wimag, bnd_qid, bnd, sort, spec, rovib_data_it->second, dipr, atm);
+    } else if (bnd.lineshape == LineByLineLineshape::VP_ECS_SPHTOP) {
+      sphtop::relaxation_matrix_offdiagonal(
+          Wimag, bnd_qid, bnd, sort, spec, rovib_data_it->second, dipr, atm);
     } else {
       ARTS_USER_ERROR("UNKNOWN ECS LINE SHAPE {}", bnd.lineshape)
     }
@@ -314,6 +331,15 @@ void ComputeData::adapt_single(const QuantumIdentifier& bnd_qid,
       auto& l2 = bnd_qid.state.at(QuantumNumberType::l2);
       dipr[i]  = hartmann::reduced_dipole(J.upper, J.lower, l2.upper, l2.lower);
       dip[i] *= std::signbit(dipr[i]) ? -1 : 1;
+    } else if (bnd.lineshape == LineByLineLineshape::VP_ECS_STOTOP) {
+      auto& J  = bnd.lines[i].qn.at(QuantumNumberType::J);
+      auto& Kq = bnd.lines[i].qn.at(QuantumNumberType::K);
+      dipr[i]  = stotop::reduced_dipole(J.upper, J.lower, Kq.lower);
+      dip[i]  *= std::signbit(dipr[i]) ? -1 : 1;
+    } else if (bnd.lineshape == LineByLineLineshape::VP_ECS_SPHTOP) {
+      auto& J  = bnd.lines[i].qn.at(QuantumNumberType::J);
+      dipr[i]  = sphtop::reduced_dipole(J.upper, J.lower);
+      dip[i]  *= std::signbit(dipr[i]) ? -1 : 1;
     }
   }
 
@@ -382,6 +408,12 @@ void ComputeData::adapt_single(const QuantumIdentifier& bnd_qid,
     } else if (bnd.lineshape == LineByLineLineshape::VP_ECS_HARTMANN) {
       hartmann::relaxation_matrix_offdiagonal(
           Wimag, bnd_qid, bnd, sort, spec, rovib_data_it->second, dipr, atm);
+    } else if (bnd.lineshape == LineByLineLineshape::VP_ECS_STOTOP) {
+      stotop::relaxation_matrix_offdiagonal(
+          Wimag, bnd_qid, bnd, sort, spec, rovib_data_it->second, dipr, atm);
+    } else if (bnd.lineshape == LineByLineLineshape::VP_ECS_SPHTOP) {
+      sphtop::relaxation_matrix_offdiagonal(
+          Wimag, bnd_qid, bnd, sort, spec, rovib_data_it->second, dipr, atm);
     } else {
       ARTS_USER_ERROR("UNKNOWN ECS LINE SHAPE {}", bnd.lineshape)
     }

src/core/lbl/lbl_lineshape_voigt_ecs_hartmann.cpp

@@ -35,7 +35,7 @@ Numeric wig6(const Rational& a,
 }
 
 std::function<Numeric(Rational)> erot_selection(const SpeciesIsotope& isot) {
-  if (isot.spec == SpeciesEnum::CarbonDioxide and isot.isotname == "626") {
+  if (isot == "CO2-626"_isot) {
     return [](const Rational J) -> Numeric {
       return Conversion::kaycm2joule(0.39021) * Numeric(J * (J + 1));
     };
@@ -79,7 +79,6 @@ void relaxation_matrix_offdiagonal(MatrixView& W,
   using std::swap;
   const bool swap_order = li > lf;
   if (swap_order) swap(li, lf);
-  const int sgn = iseven(li + lf + 1) ? -1 : 1;
   if (abs(li - lf) > 1) return;
 
   const Numeric T = atm.temperature;
@@ -129,20 +128,20 @@ void relaxation_matrix_offdiagonal(MatrixView& W,
       if (Jf_p > Jf) continue;
 
       Index L         = std::max(std::abs((Ji - Ji_p).toIndex()),
-                         std::abs((Jf - Jf_p).toIndex()));
+                                 std::abs((Jf - Jf_p).toIndex()));
       L              += L % 2;
       const Index Lf  = std::min((Ji + Ji_p).toIndex(), (Jf + Jf_p).toIndex());
 
       Numeric sum = 0;
       for (; L <= Lf; L += 2) {
-        const Numeric a  = wig3(Ji_p, Rational{L}, Ji, li, Rational{0}, -li);
-        const Numeric b  = wig3(Jf_p, Rational{L}, Jf, lf, Rational{0}, -lf);
+        const Numeric a  = wig3(Ji, Ji_p, Rational{L}, li, -li, Rational{0});
+        const Numeric b  = wig3(Jf, Jf_p, Rational{L}, lf, -lf, Rational{0});
         const Numeric c  = wig6(Ji, Jf, Rational{1}, Jf_p, Ji_p, Rational{L});
         sum             += a * b * c * Numeric(2 * L + 1) * Q[L] / Om[L];
       }
       const Numeric ECS = Om[Ji.toIndex()];
-      const Numeric scl = sgn * ECS * Numeric(2 * Ji_p + 1) *
-                          sqrtr((2 * Jf + 1) * (2 * Jf_p + 1));
+      const Numeric scl =
+          ECS * Numeric(2 * Ji_p + 1) * sqrtr((2 * Jf + 1) * (2 * Jf_p + 1));
       sum *= scl;
 
       // Add to W and rescale to upwards element by the populations
@@ -154,21 +153,16 @@ void relaxation_matrix_offdiagonal(MatrixView& W,
 
   ARTS_USER_ERROR_IF(errno == EDOM, "Cannot compute the wigner symbols")
 
-  // Undocumented negative absolute sign
-  for (Size i = 0; i < n; i++)
-    for (Size j = 0; j < n; j++)
-      if (j not_eq i and W[i, j] > 0) W[i, j] *= -1;
-
   // Sum rule correction
   for (Size i = 0; i < n; i++) {
     Numeric sumlw = 0.0;
     Numeric sumup = 0.0;
 
     for (Size j = 0; j < n; j++) {
       if (j > i) {
-        sumlw += std::abs(dipr[j]) * W[j, i];  // Undocumented abs-sign
+        sumlw += dipr[j] * W[j, i];
       } else {
-        sumup += std::abs(dipr[j]) * W[j, i];  // Undocumented abs-sign
+        sumup += dipr[j] * W[j, i];
       }
     }
 
@@ -182,7 +176,7 @@ void relaxation_matrix_offdiagonal(MatrixView& W,
       } else {
         W[j, i] *= -sumup / sumlw;
         W[i, j]  = W[j, i] * std::exp((erot(Ji) - erot(Jj)) /
-                                     kelvin2joule(T));  // This gives LTE
+                                      kelvin2joule(T));  // This gives LTE
       }
     }
   }