Added tests for spincut and D0->nld(Sn)->D0

tests/test_normalizer_nld.py

@@ -0,0 +1,40 @@
+import pytest
+from numpy.testing import assert_allclose
+from ompy.normalizer_nld import NormalizerNLD
+
+
+# expected values from d2rho
+testdata = [(2.2, 0., {"sigma": 8.43}, 6.506E+07),
+            (2.2, 0.5, {"sigma": 8.43}, 3.27e7),
+            (5, 1, {"sigma": 7.}, 6.736E+06),
+            (5, 5, {"sigma": 7.}, 2.441E+06),
+            ]
+@pytest.mark.parametrize("D0, Jtarget, spincutPars, expected", testdata)
+def test_nldSn_from_D0(D0, Jtarget, spincutPars, expected):
+    spincutModel = "const"
+    Sn = 321  # dummy value
+    nld = NormalizerNLD.nldSn_from_D0(D0, Sn, Jtarget, spincutModel,
+                                      spincutPars)
+    assert_allclose(nld[1], expected, rtol=0.01)
+
+
+@pytest.mark.parametrize("expected, Jtarget, spincutPars, nldSn", testdata)
+def test_D0_from_nldSn(nldSn, Jtarget, spincutPars, expected) -> float:
+    """ Wrapper for ompy.NormalizerNLD.D0_from_nldSn with given nld(Sn)
+
+    Args:
+        nldSn: Level density at Sn
+        pars: Normalization parameters (spin cut model ...).
+            Can be provided through ompy.NormalizationParameters.asdict().
+
+    Returns:
+        D0: Average s-wave resonance spacing
+    """
+
+    def nld_model_dummy(x):
+        return nldSn
+
+    D0 = NormalizerNLD.D0_from_nldSn(nld_model_dummy, Sn=321, Jtarget=Jtarget,
+                                     spincutModel="const",
+                                     spincutPars=spincutPars)
+    assert_allclose(D0, expected, rtol=0.01)

tests/test_spindists.py

@@ -1,6 +1,7 @@
 import pytest
 import ompy as om
 import numpy as np
+from numpy.testing import assert_allclose
 
 @pytest.mark.parametrize(
                 "spinpars",
@@ -14,3 +15,34 @@ def test_missing_parameter(spinpars):
         spindist = om.SpinFunctions(E, Js,
                                     model="EB05", pars=spinpars).distibution()
         spindist /= spindist.sum(axis=1)[:, np.newaxis]
+
+
+# tests comparing to rhobin / and or derived values
+def y_interpol(x, x0, x1, y0, y1):
+    return y0 + (x - x0) * (y1 - y0) / (x1 - x0)
+testdata = [("EB05", [4.2, 6.534],
+             {"mass": 240, "NLDa": 25.506, "Eshift": 0.162, "Sn": 6.534},
+             [7.521, 8.387]),
+            ("EB09_emp", [4.2, 6.534],
+             {"mass": 240, "Pa_prime": 1.336},
+             [4.841, 5.239]),
+            ("EB09_CT", [4.2, 6.534],
+             {"mass": 240},
+             [4.803, 4.803]),
+            ("const", [4.2, 6.534],  # not in robin, but just sigma=const
+             {"sigma": 645.456},
+             [645.456, 645.456]),
+            ("Disc_and_EB05", [1.5, 4.2, 6.534],   # combination of the above
+             {"mass": 240, "NLDa": 25.506, "Eshift": 0.162, "Sn": 6.534,
+              "sigma2_disc": [1.5, 3.6]},
+             [np.sqrt(3.6),
+              np.sqrt(y_interpol(4.2, 1.5, 6.534, 3.6, 8.387**2)), 8.387])
+            ]
+
+
+@pytest.mark.parametrize("model, Ex, pars, expected", testdata)
+def test_sigma2(model, Ex, pars, expected):
+    J = None
+    spinfunc = om.SpinFunctions(Ex, J, model=model, pars=pars)
+    spincut = np.sqrt(spinfunc.get_sigma2())
+    assert_allclose(spincut, expected, atol=0.01)