Add a base algorithm to create NTuples

[jmcarcell]

Anyone can create an algorithm that creates a TFile, a TTree, creates all the member variables that will hold the values and sets the addresses of the branches to them and then they are filled in the event loop. Therefore, if we provide something here, it has to be different. The following proposal allows

• Automatic selection of output type based on a list of choices specified only once (instead of having to specify it once for every variable). Compilation fails if the existing type can not be casted to any in the list, and it is possible to explicitly request a type by using static_cast.
• Very easy filling, one line per variable: NTupleMap["NameOfTheBranch"] = value;
• Easily switch between TTree and RNTuple with a Gaudi property
• Look like other functional algorithms; inputs are given the same way, arbitrary number of collections are supported for a given type and the output types are where typically we put the output type of the algorithm
• Performance: addresses are set only once and values are written to those always. I measured that resetting the address for a single branch every time a value is written makes writing 30% slower.
• Multiple TTrees or RNTuples in the same file in case we want to have one entry per event and one entry per track, for example.

Runs single-threaded, to have a solution running multithreaded the options provided by ROOT to write in a multithreaded way should also be studied in addition to an algorithm that can write multithreaded TTrees or RNTuples the same way as in this PR but locking.

Why not use the NTuple writers from Gaudi? None of the existing algorithms will support reading an arbitrary number of collections like we do in this repository. The GenericNTupleWriter takes inputs directly from the store (see https://gitlab.cern.ch/gaudi/Gaudi/-/blob/master/GaudiTestSuite/tests/pytest/NTuple/test_GenericNTupleWriter_st.py#L69) so it won't understand our collections, doesn't support RNTuples and sets a new address for each event.

BEGINRELEASENOTES

• Add a base algorithm to create NTuples, based on std::variant that allows writing TTrees and RNTuples with simple map[key] = value assignments and has some protections at compile time related to output types.
• Add tests to check that the algorithm is working as intended, both when writing TTrees and RNTuples.
• Add some documentation on how to use this NTuple writer.

ENDRELEASENOTES

[tmadlener]

I have played only a little bit with this now and I generally like the concept. There are a few things that are not entirely clear to me yet, and a I stumbled over a few others, where I think we might need to improve a bit.

• IIUC, each of these algorithms will write to its own TTree / RNTuple model, right? I.e. for example for an analysis tuple creation there would have to be one algorithm that does everything in one go, rather than potentially several smaller ones that each handle e.g. one edm4hep type. (I am not sure at the moment whether that is a feature that would need to be covered by this base algorithm, just pointing it out for now).
• IIUC, each instance of such an algorithm will also create its own output file. I.e. it's not possible to run, e.g. two of these algorithms in a chain and simply make them write to different TTrees in the same file? (Again not sure if this would even be considered a valid use case).
• I am not sure if the map<string, variant<...>> is the best approach as it allows quite some potentially unwanted flexibility that is neither caught at compile time nor at run time (see below for a few examples).
  • (This is mainly a guess at this point) I think the repeated string key lookups and the assignment through a variant might be quite a performance hit in tight loops and I currently don't see a way to simply hold on to a handle to the map entry such that it could be repeatedly assigned between filling without having to repeatedly go through the map lookup.
• As already described in NTuple service example #332 (comment) this doesn't necessarily solve the issue for having an "ntuple service for monitoring purposes"

Some issues with flexibility

In general the approach has some type-safety built in, e.g. it's not possible to assign a value that is not part of the Out... variadic types (since they will not be part of the internal variant). However, as soon as a type is in the Out... types, quite a bit of flexibility is possible, including things that can make for some very confusing errors (or rather outcomes, as they are not considered as errors at the moment), I think. Especially, for algorithms that produce ntuples that are useful for analysis (i.e. probably rather large ones).

Runtime type changes

Consider a slightly changed example, where the assignment to the PDG field changes type at runtime.

// ... all the rest
    if (m_counter > 5) {
      NTupleMap["PDG"] = particle.getPDG();
    } else {
      NTupleMap["PDG"] = particle.getMomentum();
    }
// ... all the rest

This will pass assignment, since edm4hep::Vector3d is in the output tuple, so it's also part of the variant that we can assign to inside the map. Hence, compilation is fine. Naively, I would have expected that this at least crashes at runtime, given that ROOT is usually rather picky about fields / branches changing types at runtime. However, the test still runs fine and it produces a file that can be opened and read. I don't think this should be possible. Ideally we catch this at compile time somehow.

The type of the branch will (quite expectedly) be the first type that is used.

Conditionally filling slots

Again consider the slightly changed example, now the following way

// ... all the rest
  if (m_counter < 5) {
    NTupleMap["AConditionallyFilledSlot"] = 42;
  } else {
    NTuplesMap["AnotherConditionallyFilledSlot"] = 54;
  }
// ... all the rest

This will result in an output that only has the AConditionallyFilledSlot but not the AnotherConditionallyFilledSlot. It will still compile and run without any problems. However, the AConditionallyFilledSlot field will still contain the same number of entries as all other slots, for the ones where no new value has been assigned, it will simply re-use the one that has been previously assigned.

Given, that we have quite a few things in place that mandate that all events should look the same in other places (e.g. for regular EDM4hep output), IMO this breaks that consistency and can lead to quite unexpected results.

Some random thoughts

I haven't yet fully though many of these through, but I would still like to write them down, while the impressions are fresh. Maybe they can serve as a starting point for further discussion.

• I am wondering whether we should provide this as a base algorithm at all. It might be easier to have this as some sort of "tool" (not necessarily in the Gaudi sense), that handles all the ROOT parts of dynamically creating branches, etc. I think pretty much all of the current functionality (e.g. configurability via Gaudi::Propertys) should also be possible with that, but we might gain a bit of freedom by not having to integrate it as a proper algorithm. On the other hand, I don't know for sure whether we would gain anything by doing that.
• I am fairly certain this doesn't scale (and I am also not sure if it's easily possible to implement this); It would be nice if one had to define a branch / field name and a type at compile time. One could than replace the map<string, variant<...>> with an array<tuple<string, variant<...>>> and create field indices (maybe even at compile time, but definitely at initialize). Those could probably be wrapped into some form of handle, and assignment could work through those. In this way one would skip all the map lookups during the event loop. (The wrapping into handles for skipping the lookup is almost certainly also possible with the current map, I think). The main advantage of this approach would be that it could probably be used to eliminate both of the issues described above, by simply not allowing them to even happen (ideally at compile time). However, for large ntuples the number of handles one has to define is fairly large. It would probably effectively double the current number of lines, as currently one simply has to index into the map, without having to declare a field first.