Enable calorimeter hit merging by functions

[ruse-traveler]

Briefly, what does this PR introduce?

This PR extends the functionality of the CalorimeterHitsMerger algorithm. Previously, reconstructed hits could only be merged across a given field of the readout of a calorimeter. This presents a challenge for calorimeters such as the Barrel HCal where

1. our readout has no segmentation beyond just eta & phi, and
2. it could be useful to study the response of the detector as a function of how many readout channels we gang together after reconstruction.

This PR addresses this point by utilizing the EvaluatorSvc in a manner similar to the adjacencyMatrix, peakNeighbourhoodMatrix of the CalorimeterIslandClustering and the sampFrac of CalorimeterHitReco. Now the user has an option to specify a mapping for a specific field of the readout via the mappings option, which defines a function to map the indices of a field onto the desired reference indices.

For example:

app->Add(new JOmniFactoryGeneratorT<CalorimeterHitsMerger_factory>(
  "HcalBarrelMergedHits", {"HcalBarrelRecHits"}, {"HcalBarrelMergedHits"},
  {
    .readout = "HcalBarrelHits",
    .fields = {"phi"},
    .mappings = {"phi-fmod(phi,5)"}
  },
  app   // TODO: Remove me once fixed
));

The HcalBarrelMergedHits collection will merge 5 hits (i.e. scintillator tiles for the BHCal) adjacent in phi into a one with the position and cellID of the 1st of the 5. However, no hits will be merged along eta in this example.

If a mapping is provided, the algorithm will construct the map of which hits to merge according the provided functions. If not, it will attempt to merge hits across the provided fields into the indices defined by refs.

An example script of to change a mapping (and how to update the adjacency matrix accordingly) from the command-line is provided in the snippets repo here.

What kind of change does this PR introduce?

• Bug fix (issue #__)
• New feature (issue Generalize caloriemter hit merging #1669 )
• Documentation update
• Other: __

Please check if this PR fulfills the following:

• Tests for the changes have been added
• Documentation has been added / updated
• Changes have been communicated to collaborators

Does this PR introduce breaking changes? What changes might users need to make to their code?

No.

Does this PR change default behavior?

No.

[github-actions]

Capybara summary for PR 1668

• rec_dis_10x100_minQ2=0_craterlake
• rec_dis_10x100_minQ2=1000_craterlake_tracking_only
• rec_dis_18x275_minQ2=0_craterlake_18x275
• rec_dis_18x275_minQ2=1000_craterlake_18x275
• rec_dis_5x41_minQ2=0_craterlake_5x41
• rec_e_1GeV_20GeV_craterlake
• rec_pi_1GeV_20GeV_craterlake
  ^{_{Last updated 2024-12-19T13:59-05:00 3b6e445}}

[ruse-traveler]

I think this is ready to go. I'll post some plots here showing the test case in the BHCal at work shortly...

[ruse-traveler]

Here are some (eta, phi) maps of the tiles before and after merging. I tried out two cases for merging (using the example map from the PR description): once with the no. of tiles to merge set to 1 (so in other words, don't merge the tiles, and once with the no. set to 5 (merge into an sPHENIX-style tower).

Broadly, things look like I would expect them to. Except for the giant hole in the hit map...

[ruse-traveler]

Then another odd feature I noticed is that even when the no. to merge is 1, there are still slight changes in the eta, phi distributions of the tiles even though it should just be copying one collection onto another...

[veprbl]

I'm wondering if we should just have only one function that can map any cellID to any cellID, and to that function we'd provide the precomputed bit masks. That would be most general and simple implementation, at the cost of complexitry. E.g. your example could look like: (cellID & ~_phi_mask) | ((((cellID & _phi_mask) >> _phi_mask_shift) / 5) << _phi_mask_shift) and we could introduce some helper functions to not have to deal with most raw C++. Just throwing this in as an idea.

[veprbl]

Mappings sounds a bit generic to me, and it doesn't self-describe that there should be as many mappings as fields. Should we call it something like field transformations or rewrites?

Suggested change

	ParameterRef<std::vector<std::string>> m_mappings {this, "mappings", config().mappings};
	ParameterRef<std::vector<std::string>> m_field_transformations {this, "fieldTransformations", config().field_transformations};

[ruse-traveler]

Good point! I also like that name!

[veprbl]

Can we put mask operations as closures to the same ref_maps, so that the application of transofrmation is uniform?

[ruse-traveler]

I don't quite follow... So by closure do you mean if a user provides both a bitmask and a mapping, they both get applied?

[veprbl]

Allow me to clarify. I meant that this code can be replaced with something like

for (size_t ix = 0; std::string &field : m_cfg.fields) {
  uint64_t ref = m_cfg.refs.at(ix);
  ref_maps[field] = [ref](const edm4eic::CalorimeterHit&) -> int {
    return ref;
  };
  ix++;
}

Then we can remove all the different logic for handling ref_mask and id_mask.

[veprbl]

I think this it's more like:

for (size_t ix = 0; std::string &field : m_cfg.fields) {
  uint64_t ref_mask = m_cfg.refs.at(ix);
  auto field_element = id_desc.field(m_cfg.fields[i]);
  ref_maps[field] = [ref_mask,field_element](const edm4eic::CalorimeterHit& hit) -> int {
    return (field_element->value(hit.cellID()) & ~ref_mask) | ref_mask;
  };
  ix++;
}

[veprbl]

This looks a bit complicated, and I doubt that it actually works for mixing of "mappings" transformations and "refs" bitmasks. Could it be simplified assuming 1-to-1 loop?

[ruse-traveler]

Yes, as-is it's set up so that if a user provides a mapping, it defaults to using that and will ignore any provided bitmasks. But it can definitely be just turned into a 1-to-1 loop where it checks for the mapping or the bitmask.

[veprbl]

This also breaks my brain a bit. A check for matching length and a 1-to-1 loop would be nicer.

[ruse-traveler]

Sure! This was trying to work in the above design choice that providing a map makes the algorithm ignores any provided bitmasks.

[ruse-traveler]

I'm wondering if we should just have only one function that can map any cellID to any cellID, and to that function we'd provide the precomputed bit masks. That would be most general and simple implementation, at the cost of complexitry. E.g. your example could look like: (cellID & ~_phi_mask) | ((((cellID & _phi_mask) >> _phi_mask_shift) / 5) << _phi_mask_shift) and we could introduce some helper functions to not have to deal with most raw C++. Just throwing this in as an idea.

Definitely an option! And it might be nice to have that functionality available even if we keep the field-specific mappings...

To be frank, I still find bitwise operations a bit mystifying and opaque, which is part of why I shied away from using them here 😅
One thing I was hoping to accomplish with this was to provide a more "user friendly" a way to express mappings, which is another reason why I opted for the more "functional" expressions...

[Chao1009]

To be frank, I still find bitwise operations a bit mystifying and opaque, which is part of why I shied away from using them here 😅 One thing I was hoping to accomplish with this was to provide a more "user friendly" a way to express mappings, which is another reason why I opted for the more "functional" expressions...

Agreed.
It should be better to use functional mapping and remove all bitwise operations. To do this, can we merge the input "fields" and "mappings", so it supports both a string for the mapping operation and simple masking from field names?

[veprbl]

@Chao1009 Is that not the suggestion from #1668 (comment) ?

[Chao1009]

@Chao1009 Is that not the suggestion from #1668 (comment) ?

It's similar. The only difference is that I suggest using a single input for the fields and mappings (or call it field_transformations). For example, assuming the input is a vector of strings, like
mapping={"phi:0", "theta:<transformation_string>"}
The former does a simple merge of the phi field (transform it to 0), while the latter does some sophisticated transformation for theta.
In this way, we also avoid one-to-one matching of the fields and mappings in the code.

[veprbl]

I think the current implementation in this PR already allows to provide arbitrary functions of arbitrary other fields to re-define fields. The old merging behavior was to just mask out bits of the original values and set them to all binary 1.