Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

fix: fix #238 – wrong handling of offsets #274

Draft
wants to merge 6 commits into
base: main
Choose a base branch
from
Draft

fix: fix #238 – wrong handling of offsets #274

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
42ed4e6
add test for masking in exclusive_constituents
chrispap95 Feb 16, 2024
549040b
Merge branch 'main' into fix-238-offsets
lgray Feb 17, 2024
8ce0a84
Merge branch 'main' into fix-238-offsets
chrispap95 Mar 13, 2024
ece217a
First attempt to fix the issue. Implemented starts & stops.
chrispap95 Mar 13, 2024
64b3777
should use the pointer
chrispap95 Mar 13, 2024
9500865
start/stop are smaller than offsets array
chrispap95 Mar 13, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
38 changes: 18 additions & 20 deletions src/_ext.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -33,9 +33,9 @@ typedef struct {
PyObject *next;
} SwigPyObject;

template <typename T>
T swigtocpp(py::object obj) { // unwraps python object to get the cpp pointer
// from the swig bindings
template <typename T> T swigtocpp(py::object obj) {
// unwraps python object to get the cpp pointer
// from the swig bindings
auto upointer = obj.attr("this").ptr();
auto swigpointer = reinterpret_cast<SwigPyObject *>(upointer);
auto objpointervoid = swigpointer->ptr;
Expand All @@ -59,47 +59,45 @@ output_wrapper interfacemulti(
py::array_t<double, py::array::c_style | py::array::forcecast> pyi,
py::array_t<double, py::array::c_style | py::array::forcecast> pzi,
py::array_t<double, py::array::c_style | py::array::forcecast> Ei,
py::array_t<int, py::array::c_style | py::array::forcecast> offsets,
py::array_t<int, py::array::c_style | py::array::forcecast> starts,
py::array_t<int, py::array::c_style | py::array::forcecast> stops,
py::object jetdef) {
py::buffer_info infooff = offsets.request();
// requesting buffer information of the input
py::buffer_info infostarts = starts.request();
py::buffer_info infostops = stops.request();
py::buffer_info infopx = pxi.request();
py::buffer_info infopy =
pyi.request(); // requesting buffer information of the input
py::buffer_info infopy = pyi.request();
py::buffer_info infopz = pzi.request();
py::buffer_info infoE = Ei.request();

auto offptr = static_cast<int *>(infooff.ptr);
// pointers to the initial values
auto startsptr = static_cast<int *>(infostarts.ptr);
auto stopsptr = static_cast<int *>(infostops.ptr);
auto pxptr = static_cast<double *>(infopx.ptr);
auto pyptr =
static_cast<double *>(infopy.ptr); // pointer to the initial value
auto pyptr = static_cast<double *>(infopy.ptr);
auto pzptr = static_cast<double *>(infopz.ptr);
auto Eptr = static_cast<double *>(infoE.ptr);

int dimoff = infooff.shape[0];
int dimoff = infostarts.shape[0];
output_wrapper ow;
std::vector<double> nevents;
std::vector<double> offidx;
std::vector<double> constphi;
std::vector<double> idx;
std::vector<double> idxo;
for (int i = 0; i < dimoff - 1; i++) {
for (int i = 0; i < dimoff; i++) {
std::vector<fj::PseudoJet> particles;
for (int j = *offptr; j < *(offptr + 1); j++) {
for (int j = *startsptr; j < *stopsptr; j++) {
particles.push_back(fj::PseudoJet(*pxptr, *pyptr, *pzptr, *Eptr));
pxptr++;
pyptr++;
pzptr++;
Eptr++;
}

std::vector<fj::PseudoJet> jets;
auto jet_def = swigtocpp<fj::JetDefinition *>(jetdef);
std::shared_ptr<std::vector<fj::PseudoJet>> pj =
std::make_shared<std::vector<fj::PseudoJet>>(particles);
std::shared_ptr<fastjet::ClusterSequence> cs =
std::make_shared<fastjet::ClusterSequence>(*pj, *jet_def);
auto j = cs->inclusive_jets();
offptr++;
startsptr++;
stopsptr++;
ow.cse.push_back(cs);
ow.parts.push_back(pj);
}
Expand Down
15 changes: 9 additions & 6 deletions src/fastjet/_generalevent.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -21,14 +21,15 @@ def __init__(self, data, jetdef):
self._clusterable_level[i] = ak.Array(
self._clusterable_level[i].layout.to_ListOffsetArray64(True)
)
px, py, pz, E, offsets = self.extract_cons(self._clusterable_level[i])
px, py, pz, E, starts, stops = self.extract_cons(self._clusterable_level[i])
px = self.correct_byteorder(px)
py = self.correct_byteorder(py)
pz = self.correct_byteorder(pz)
E = self.correct_byteorder(E)
offsets = self.correct_byteorder(offsets)
starts = self.correct_byteorder(starts)
stops = self.correct_byteorder(stops)
self._results.append(
fastjet._ext.interfacemulti(px, py, pz, E, offsets, jetdef)
fastjet._ext.interfacemulti(px, py, pz, E, starts, stops, jetdef)
)

def _check_listoffset_subtree(self, data):
Expand Down Expand Up @@ -215,9 +216,11 @@ def extract_cons(self, array):
py = np.asarray(ak.Array(array.layout.content, behavior=array.behavior).py)
pz = np.asarray(ak.Array(array.layout.content, behavior=array.behavior).pz)
E = np.asarray(ak.Array(array.layout.content, behavior=array.behavior).E)
off = np.asarray(array.layout.stops)
off = np.insert(off, 0, 0)
return px, py, pz, E, off
starts = np.asarray(array.layout.starts)
stops = np.asarray(array.layout.stops)
starts = np.insert(starts, 0, 0)
stops = np.insert(stops, 0, 0)
return px, py, pz, E, starts, stops

def _replace_multi(self):
self._mod_data = self.data
Expand Down
18 changes: 11 additions & 7 deletions src/fastjet/_multievent.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -9,15 +9,17 @@
class _classmultievent:
def __init__(self, data, jetdef):
self.jetdef = jetdef
# data = ak.Array(data.layout.to_ListOffsetArray64(True)) If I do this vector can't convert the coordinates
self.data = data
px, py, pz, E, offsets = self.extract_cons(self.data)
px, py, pz, E, starts, stops = self.extract_cons(self.data)
px = self.correct_byteorder(px)
py = self.correct_byteorder(py)
pz = self.correct_byteorder(pz)
E = self.correct_byteorder(E)
offsets = self.correct_byteorder(offsets)
self._results = fastjet._ext.interfacemulti(px, py, pz, E, offsets, jetdef)
starts = self.correct_byteorder(starts)
stops = self.correct_byteorder(stops)
self._results = fastjet._ext.interfacemulti(
px, py, pz, E, starts, stops, jetdef
)

def _check_record(self, data):
return data.layout.is_record or data.layout.is_numpy
Expand All @@ -34,9 +36,11 @@ def extract_cons(self, array):
py = np.asarray(ak.Array(array.layout.content, behavior=array.behavior).py)
pz = np.asarray(ak.Array(array.layout.content, behavior=array.behavior).pz)
E = np.asarray(ak.Array(array.layout.content, behavior=array.behavior).E)
off = np.asarray(array.layout.stops)
off = np.insert(off, 0, 0)
return px, py, pz, E, off
starts = np.asarray(array.layout.starts)
stops = np.asarray(array.layout.stops)
starts = np.insert(starts, 0, 0)
stops = np.insert(stops, 0, 0)
return px, py, pz, E, starts, stops

def single_to_jagged(self, array):
single = ak.Array(
Expand Down
17 changes: 11 additions & 6 deletions src/fastjet/_singleevent.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -10,13 +10,16 @@ class _classsingleevent:
def __init__(self, data, jetdef):
self.jetdef = jetdef
self.data = self.single_to_jagged(data)
px, py, pz, E, offsets = self.extract_cons(self.data)
px, py, pz, E, starts, stops = self.extract_cons(self.data)
px = self.correct_byteorder(px)
py = self.correct_byteorder(py)
pz = self.correct_byteorder(pz)
E = self.correct_byteorder(E)
offsets = self.correct_byteorder(offsets)
self._results = fastjet._ext.interfacemulti(px, py, pz, E, offsets, jetdef)
starts = self.correct_byteorder(starts)
stops = self.correct_byteorder(stops)
self._results = fastjet._ext.interfacemulti(
px, py, pz, E, starts, stops, jetdef
)

def correct_byteorder(self, data):
if data.dtype.byteorder == "=":
Expand All @@ -36,9 +39,11 @@ def extract_cons(self, array):
py = np.asarray(ak.Array(array.layout.content, behavior=array.behavior).py)
pz = np.asarray(ak.Array(array.layout.content, behavior=array.behavior).pz)
E = np.asarray(ak.Array(array.layout.content, behavior=array.behavior).E)
off = np.asarray(array.layout.stops)
off = np.insert(off, 0, 0)
return px, py, pz, E, off
starts = np.asarray(array.layout.starts)
stops = np.asarray(array.layout.stops)
starts = np.insert(starts, 0, 0)
stops = np.insert(stops, 0, 0)
return px, py, pz, E, starts, stops

def _check_record(self, data):
return data.layout.is_record or data.layout.is_numpy
Expand Down
35 changes: 35 additions & 0 deletions tests/test_002-exclusive_jets.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -534,6 +534,41 @@ def test_exclusive_constituents_multi():
)


def test_exclusive_constituents_masking_multi():
array = ak.Array(
[
[
{"px": -0.181, "py": -0.798, "pz": -0.652, "E": 1.06},
{"px": 0.377, "py": 0.116, "pz": -0.0749, "E": 0.425},
],
[
{"px": 0.54, "py": -0.65, "pz": -0.00527, "E": 0.857},
{"px": 0.253, "py": -0.0342, "pz": -0.0731, "E": 0.3},
],
[
{"px": 0.294, "py": 0.254, "pz": -0.259, "E": 0.467},
{"px": 4.65, "py": -1.88, "pz": -3.29, "E": 6},
],
[
{"px": 1.45, "py": -0.179, "pz": -0.876, "E": 1.71},
{"px": 12.8, "py": -1.8, "pz": -7.18, "E": 14.8},
],
[
{"px": -3.55, "py": -1.64, "pz": -0.0941, "E": 3.91},
{"px": -1.33, "py": -1.03, "pz": 0.147, "E": 1.7},
],
],
with_name="Momentum4D",
)
mask = ak.Array([True, True, False, True, True])
jetdef = fastjet.JetDefinition(fastjet.kt_algorithm, 1)
result = fastjet.ClusterSequence(
array[mask], jetdef
).exclusive_jets_constituent_index(njets=2)
output = ak.Array([[[1], [0]], [[0], [1]], [[1], [0]], [[0], [1]]])
assert ak.all(result == output)


def test_exclusive_ycut():
array = ak.Array(
[
Expand Down
Loading