Allow different matching scale for the photon

n3fit/src/evolven3fit_new/eko_utils.py

@@ -22,7 +22,7 @@
 }
 
 EVOLVEN3FIT_CONFIGS_DEFAULTS_EXA = {
-    "ev_op_iterations": 10,
+    "ev_op_iterations": 30,
     "ev_op_max_order": (1, 0),
     "evolution_method": "iterate-exact",
     "inversion_method": "exact",
@@ -77,6 +77,11 @@ def construct_eko_cards(
     legacy_class = runcards.Legacy(theory, {})
     theory_card = legacy_class.new_theory
 
+    # if Qedref = Qref alphaem is running, otherwise it's fixed
+    if theory["QED"] > 0:
+        if np.isclose(theory["Qref"], theory["Qedref"]):
+            theory_card.couplings.em_running = True
+
     # construct operator card
     q2_grid = utils.generate_q2grid(
         theory["Q0"],

validphys2/src/validphys/photon/compute.py

@@ -17,7 +17,6 @@
 
 log = logging.getLogger(__name__)
 
-Q_IN = 100
 FIATLUX_DEFAULT = {
     "apfel": False,
     "eps_base": 1e-5,  # precision on final integration of double integral.
@@ -71,6 +70,9 @@ def __init__(self, theoryid, lux_params, replicas):
         path_to_F2 = theoryid.path / "fastkernel/fiatlux_dis_F2.pineappl.lz4"
         path_to_FL = theoryid.path / "fastkernel/fiatlux_dis_FL.pineappl.lz4"
         self.path_to_eko_photon = theoryid.path / "eko_photon.tar"
+        with EKO.read(self.path_to_eko_photon) as eko:
+            self.q_in = np.sqrt(eko.mu20)
+
 
         # set fiatlux
         self.lux = {}
@@ -129,7 +131,7 @@ def compute_photon_array(self, replica):
         """
         # Compute photon PDF
         log.info(f"Computing photon")
-        photon_qin = np.array([self.lux[replica].EvaluatePhoton(x, Q_IN**2).total for x in XGRID])
+        photon_qin = np.array([self.lux[replica].EvaluatePhoton(x, self.q_in**2).total for x in XGRID])
         photon_qin += self.generate_errors(replica)
         # fiatlux computes x * gamma(x)
         photon_qin /= XGRID
@@ -141,21 +143,20 @@ def compute_photon_array(self, replica):
             # it has to be done inside vp-setupfit
 
             # construct PDFs
-            pdfs_init = np.zeros((len(eko.rotations.inputpids), len(XGRID)))
-            for j, pid in enumerate(eko.rotations.inputpids):
+            pdfs_init = np.zeros((len(eko.bases.inputpids), len(XGRID)))
+            for j, pid in enumerate(eko.bases.inputpids):
                 if pid == 22:
                     pdfs_init[j] = photon_qin
                     ph_id = j
                 else:
                     if pid not in self.luxpdfset.flavors:
                         continue
                     pdfs_init[j] = np.array(
-                        [self.luxpdfset.xfxQ(x, Q_IN, replica, pid) / x for x in XGRID]
+                        [self.luxpdfset.xfxQ(x, self.q_in, replica, pid) / x for x in XGRID]
                     )
 
             # Apply EKO to PDFs
-            q2 = eko.mu2grid[0]
-            with eko.operator(q2) as elem:
+            for _, elem in eko.items():
                 pdfs_final = np.einsum("ajbk,bk", elem.operator, pdfs_init)
 
         photon_Q0 = pdfs_final[ph_id]
@@ -188,7 +189,7 @@ def error_matrix(self):
         if not self.additional_errors:
             return None
         extra_set = self.additional_errors.load()
-        qs = [Q_IN] * len(XGRID)
+        qs = [self.q_in] * len(XGRID)
         res_central = np.array(extra_set.central_member.xfxQ(22, XGRID, qs))
         res = []
         for idx_member in range(101, 107 + 1):

validphys2/src/validphys/tests/photon/test_compute.py

@@ -8,6 +8,7 @@
 from validphys.core import PDF as PDFset
 
 from ..conftest import PDF
+from eko.io import EKO
 
 
 class FakeTheory:
@@ -68,6 +69,17 @@ def PlugStructureFunctions(self, f2, fl, f2lo):
     def EvaluatePhoton(self, x, q):
         return self.res
 
+class FakeEKO:
+    def __init__(self, path):
+        self.path = path
+        self.mu20 = 100**2
+
+    def __enter__(self):
+        return self
+
+    def __exit__(self, exc_type: type, _exc_value, _traceback):
+        pass
+
 
 class FakeStructureFunction:
     def __init__(self, path, pdfs):
@@ -95,6 +107,7 @@ def test_parameters_init(monkeypatch):
     monkeypatch.setattr(structure_functions, "F2LO", FakeF2LO)
     monkeypatch.setattr(fiatlux, "FiatLux", FakeFiatlux)
     monkeypatch.setattr(Photon, "compute_photon_array", lambda *args: np.zeros(196))
+    monkeypatch.setattr(EKO, "read", FakeEKO)
 
     photon = Photon(FakeTheory(), fiatlux_runcard, [1, 2, 3])
     alpha = Alpha(FakeTheory().get_description())