Repository: scikit-hep/pyjet (language: python)

Issue:
Recombination scheme
Dear authors,

thanks a lot for developing this code.  I find it very useful.  I have a user question, I have been trying for a long time.  I want to use the Antikt algorithm R=0.4 but with the Et recombination scheme.  And I tried many things and none of them worked.  In fact if I do:
this_jetdefinition = JetDefinition('genkt',jetsize,-1) and try to add more arguments it says it does not admit any.  If I do: 

this_jetdefinition = JetDefinition('antikt',jetsize,-1) it says that only one parameters is needed, and if I do 

this_jetdefinition.set_recombination_scheme() it says it does not recognise the function.

So, that's was it, thanks a lot again and please have my best regards,

Esteban



Discussion / hints:

Fix patch:
diff --git a/pyjet/src/_libpyjet.pyx b/pyjet/src/_libpyjet.pyx
index 14a5bed..0339777 100644
--- a/pyjet/src/_libpyjet.pyx
+++ b/pyjet/src/_libpyjet.pyx
@@ -59,6 +59,17 @@ JET_AREA = {
     'voronoi': fastjet.voronoi_area,
 }
 
+RECOMB_SCHEME = {
+  'E_scheme': fastjet.E_scheme,
+  'pt_scheme': fastjet.pt_scheme,
+  'pt2_scheme': fastjet.pt2_scheme,
+  'Et_scheme': fastjet.Et_scheme,
+  'Et2_scheme': fastjet.Et2_scheme,
+  'BIpt_scheme': fastjet.BIpt_scheme,
+  'BIpt2_scheme': fastjet.BIpt2_scheme,
+  'WTA_pt_scheme': fastjet.WTA_pt_scheme,
+  'WTA_modp_scheme': fastjet.WTA_modp_scheme
+}
 
 cdef class JetDefinition:
     cdef fastjet.JetDefinition* jdef
@@ -66,10 +77,11 @@ cdef class JetDefinition:
     def __cinit__(self):
         self.jdef = NULL
 
-    def __init__(self, algo='undefined', R=None, p=None):
+    def __init__(self, algo='undefined', R=None, p=None, recomb_scheme=None):
         if self.jdef != NULL:
             del self.jdef
         cdef fastjet.JetAlgorithm _algo
+        cdef fastjet.RecombinationScheme _recomb_scheme
         try:
             _algo = JET_ALGORITHM[algo]
         except KeyError:
@@ -82,6 +94,14 @@ cdef class JetDefinition:
         else:
             self.jdef = new fastjet.JetDefinition(_algo)
 
+        if recomb_scheme is not None:
+            try:
+                _recomb_scheme = RECOMB_SCHEME[recomb_scheme]
+            except KeyError:
+                raise ValueError("{} is not a valid recombination scheme".format(recomb_scheme))
+                
+            self.jdef.set_recombination_scheme(_recomb_scheme)
+            
     def __dealloc__(self):
         del self.jdef
 
diff --git a/pyjet/src/fastjet.pxd b/pyjet/src/fastjet.pxd
index 7ce56de..4634d0c 100644
--- a/pyjet/src/fastjet.pxd
+++ b/pyjet/src/fastjet.pxd
@@ -66,11 +66,23 @@ cdef extern from "fastjet.h" namespace "fastjet":
         ee_genkt_algorithm,
         plugin_algorithm,
         undefined_jet_algorithm
+           
+    cdef enum RecombinationScheme "fastjet::RecombinationScheme":
+         E_scheme, 
+         pt_scheme,        
+         pt2_scheme,
+         Et_scheme,
+         Et2_scheme, 
+         BIpt_scheme,      
+         BIpt2_scheme, 
+         WTA_pt_scheme, 
+         WTA_modp_scheme
 
     cdef cppclass JetDefinition:
         JetDefinition(JetAlgorithm) except +raise_py_error
         JetDefinition(JetAlgorithm, double R) except +raise_py_error
         JetDefinition(JetAlgorithm, double R, double extra) except +raise_py_error
+        void set_recombination_scheme(RecombinationScheme recomb_scheme)
 
     cdef enum AreaType "fastjet::AreaType":
         invalid_area,


Test patch:
diff --git a/tests/test_cluster.py b/tests/test_cluster.py
index 05930a0..5e09b27 100644
--- a/tests/test_cluster.py
+++ b/tests/test_cluster.py
@@ -46,6 +46,29 @@ def test_cluster_vectors_wrong_type():
         vectors = np.zeros(10, dtype=[('a', 'f8'), ('b', 'f8'), ('c', 'f4'), ('d', 'f8')])
         cluster(vectors, R=0.6, p=-1)
 
+def test_cluster_recombination_schemes():
+    recomb_schemes = {
+        'E_scheme': (983.28, -0.8676, 36.46),
+        'pt_scheme': (983.52, -0.8672, 0.00),
+        'pt2_scheme': (983.52, -0.8679, 0.00),
+        'Et_scheme': (983.55, -0.8672, 0.00),
+        'Et2_scheme': (983.55, -0.8679, 0.00),
+        'BIpt_scheme': (983.52, -0.8671, 0.00),
+        'BIpt2_scheme': (983.52, -0.8679, 0.00),
+        'WTA_pt_scheme': (983.52, -0.8684, 0.14),
+        'WTA_modp_scheme': (983.03, -0.8684, 0.14),
+    }
+
+    for recomb_scheme, jet in recomb_schemes.items():
+        sequence = cluster(get_event(), R=0.6, p=-1, recomb_scheme=recomb_scheme)
+        jets = sequence.inclusive_jets()
+        
+        assert jets[0].pt == approx(jet[0], abs=1e-2)
+        assert jets[0].eta == approx(jet[1], abs=1e-4)
+        assert jets[0].mass == approx(jet[2], abs=1e-2)
+
+    with pytest.raises(ValueError):
+        sequence = cluster(get_event(), R=0.6, p=-1, recomb_scheme='invalid_scheme')
 
 def test_userinfo():
     event = get_event()