Fixed binary convert for ELL1H

CHANGELOG-unreleased.md

@@ -28,6 +28,7 @@ the released changes.
 - Fixed bug in changing epoch for ELL1k model
 - Fixed `gridutils` behavior for 1 CPU
 - Fixed bug in `GaussianRV_gen`, where the probability distribution function was not normalized correctly. Changed to use `scipy.stats.truncnorm` instead of the custom `GaussianRV_gen`.
+- Fixed `convert_binary()` for ELL1H models to run `setup()` and not use H4 when not desired
 - Fixed bug in `model.compare()` where it failed for `PosixPath` objects
 - Fixed bug in printing of parameter correlation/covariance matrices
 - `make_fake_toas_fromMJDs` now does not assume `PLANET_SHAPIRO` is in the model - it checks.

src/pint/binaryconvert.py

@@ -553,7 +553,7 @@ def _transfer_params(
 def convert_binary(
     model: pint.models.TimingModel,
     output: str,
-    NHARMS: int = 3,
+    NHARMS: int = 7,
     useSTIGMA: bool = False,
     KOM: u.Quantity = 0 * u.deg,
 ) -> pint.models.TimingModel:
@@ -645,10 +645,13 @@ def convert_binary(
                 outmodel.STIGMA.uncertainty = stigma_unc
                 outmodel.STIGMA.frozen = outmodel.H3.frozen
             else:
-                # use H4 and H3
-                outmodel.H4.quantity = h4
-                outmodel.H4.uncertainty = h4_unc
-                outmodel.H4.frozen = outmodel.H3.frozen
+                # use H4?
+                if NHARMS > 3:
+                    outmodel.H4.quantity = h4
+                    outmodel.H4.uncertainty = h4_unc
+                    outmodel.H4.frozen = outmodel.H3.frozen
+                else:
+                    outmodel.H4._quantity = None
         elif output in ["ELL1"]:
             if model.BINARY.value == "ELL1H":
                 # ELL1H -> ELL1
@@ -1201,7 +1204,9 @@ def convert_binary(
                         outmodel.SINI.frozen = model.KIN.frozen
                 else:
                     tempmodel = convert_binary(model, "DD")
-                    outmodel = convert_binary(tempmodel, output)
+                    outmodel = convert_binary(
+                        tempmodel, output, NHARMS=NHARMS, useSTIGMA=useSTIGMA
+                    )
             if output == "ELL1H":
                 if binary_component.binary_model_name in ["DDGR", "DDH", "DDK"]:
                     model = convert_binary(model, "DD")
@@ -1218,11 +1223,13 @@ def convert_binary(
                     outmodel.STIGMA.uncertainty = stigma_unc
                     outmodel.STIGMA.frozen = outmodel.H3.frozen
                 else:
-                    # use H4 and H3
-                    outmodel.H4.quantity = h4
-                    outmodel.H4.uncertainty = h4_unc
-                    outmodel.H4.frozen = outmodel.H3.frozen
-
+                    # use H4?
+                    if NHARMS > 3:
+                        outmodel.H4.quantity = h4
+                        outmodel.H4.uncertainty = h4_unc
+                        outmodel.H4.frozen = outmodel.H3.frozen
+                    else:
+                        outmodel.H4._quantity = None
     if (
         output == "DDS"
         and binary_component.binary_model_name != "DDGR"
@@ -1265,5 +1272,6 @@ def convert_binary(
                 )
                 outmodel.KIN.frozen = model.SINI.frozen
     outmodel.validate()
+    outmodel.setup()
 
     return outmodel

src/pint/models/binary_ell1.py

@@ -410,7 +410,11 @@ def setup(self):
                 )
             self.binary_instance.ds_func = self.binary_instance.delayS_H3_STIGMA_exact
             if self.STIGMA.quantity <= 0:
-                raise ValueError("STIGMA must be greater than zero.")
+                if not self.H3.value == 0:
+                    raise ValueError("STIGMA must be greater than zero.")
+                log.warning(
+                    f"ELL1H model has STIGMA=0 and H3=0: this may require adjusting STIGMA to be >0 before fitting"
+                )
         self.update_binary_object(None)
 
     def validate(self):

tests/test_binconvert.py

@@ -349,3 +349,17 @@ def test_DDFB0_roundtrip(output):
                 )
         else:
             assert getattr(m, p).value == getattr(mback, p).value
+
+
+def test_ELL1_ELL1H():
+    m = get_model(io.StringIO(parELL1))
+    mout = pint.binaryconvert.convert_binary(m, "ELL1H", useSTIGMA=True)
+    assert "STIGMA" in mout.components["BinaryELL1H"].binary_instance.fit_params
+
+    mout = pint.binaryconvert.convert_binary(m, "ELL1H", useSTIGMA=False)
+    assert "H4" in mout.components["BinaryELL1H"].binary_instance.fit_params
+    assert mout.NHARMS.value == 7
+
+    mout = pint.binaryconvert.convert_binary(m, "ELL1H", useSTIGMA=False, NHARMS=3)
+    assert "H4" not in mout.components["BinaryELL1H"].binary_instance.fit_params
+    assert mout.NHARMS.value == 3