| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | /// @begin_doc_example{cpp} | ||
| 2 | /// @file iguana_ex_cpp_00_run_functions.cc | ||
| 3 | /// @brief Example using **full HIPO banks** with Iguana algorithms' `Run` functions, using `hipo::banklist` | ||
| 4 | /// | ||
| 5 | /// This example requires the user to have the C++ `hipo::banklist` objects, which are lists of `hipo::bank` objects. | ||
| 6 | /// | ||
| 7 | /// - see `iguana_ex_cpp_00_run_functions_with_banks.cc` if you prefer just `hipo::bank` objects, rather than `hipo::banklist` | ||
| 8 | /// - see other examples if you do not have `hipo::bank` objects | ||
| 9 | /// | ||
| 10 | /// @par Usage | ||
| 11 | /// ```bash | ||
| 12 | /// iguana_ex_cpp_00_run_functions [HIPO_FILE] [NUM_EVENTS] | ||
| 13 | /// | ||
| 14 | /// HIPO_FILE the HIPO file to analyze | ||
| 15 | /// | ||
| 16 | /// NUM_EVENTS the number of events to analyze; | ||
| 17 | /// set to zero to analyze all events | ||
| 18 | /// ``` | ||
| 19 | /// @end_doc_example | ||
| 20 | |||
| 21 | #include <hipo4/reader.h> | ||
| 22 | #include <iguana/algorithms/AlgorithmSequence.h> | ||
| 23 | |||
| 24 | /// main function | ||
| 25 | 1 | int main(int argc, char** argv) | |
| 26 | { | ||
| 27 | |||
| 28 | // parse arguments | ||
| 29 |
1/2✓ Branch 2 → 3 taken 1 time.
✗ Branch 2 → 12 not taken.
|
1 | char const* inFileName = argc > 1 ? argv[1] : "data.hipo"; |
| 30 |
1/2✓ Branch 3 → 4 taken 1 time.
✗ Branch 3 → 12 not taken.
|
1 | int const numEvents = argc > 2 ? std::stoi(argv[2]) : 3; |
| 31 | |||
| 32 | // read input file | ||
| 33 |
1/2✓ Branch 13 → 14 taken 1 time.
✗ Branch 13 → 186 not taken.
|
1 | hipo::reader reader(inFileName, {0}); |
| 34 | |||
| 35 | // set list of banks to be read | ||
| 36 |
1/2✓ Branch 17 → 18 taken 1 time.
✗ Branch 17 → 302 not taken.
|
2 | hipo::banklist banks = reader.getBanks({"RUN::config", |
| 37 | "REC::Particle", | ||
| 38 | "REC::Calorimeter", | ||
| 39 | "REC::Track", | ||
| 40 |
1/2✓ Branch 18 → 19 taken 1 time.
✗ Branch 18 → 190 not taken.
|
1 | "REC::Scintillator"}); |
| 41 | |||
| 42 | // iguana algorithm sequence | ||
| 43 | // NOTE: the order that they are added to the sequence here will be the same order in which they will be run | ||
| 44 |
2/4✓ Branch 20 → 21 taken 1 time.
✗ Branch 20 → 300 not taken.
✓ Branch 21 → 22 taken 1 time.
✗ Branch 21 → 298 not taken.
|
1 | iguana::AlgorithmSequence seq; |
| 45 |
3/8✓ Branch 22 → 23 taken 1 time.
✗ Branch 22 → 198 not taken.
✓ Branch 23 → 24 taken 1 time.
✗ Branch 23 → 192 not taken.
✓ Branch 29 → 30 taken 1 time.
✗ Branch 29 → 32 not taken.
✗ Branch 198 → 199 not taken.
✗ Branch 198 → 201 not taken.
|
2 | seq.Add("clas12::EventBuilderFilter"); // filter by Event Builder PID (a filter algorithm) |
| 46 |
3/8✓ Branch 35 → 36 taken 1 time.
✗ Branch 35 → 210 not taken.
✓ Branch 36 → 37 taken 1 time.
✗ Branch 36 → 204 not taken.
✓ Branch 42 → 43 taken 1 time.
✗ Branch 42 → 45 not taken.
✗ Branch 210 → 211 not taken.
✗ Branch 210 → 213 not taken.
|
2 | seq.Add("clas12::SectorFinder"); // get the sector for each particle (a creator algorithm) |
| 47 |
3/8✓ Branch 48 → 49 taken 1 time.
✗ Branch 48 → 222 not taken.
✓ Branch 49 → 50 taken 1 time.
✗ Branch 49 → 216 not taken.
✓ Branch 55 → 56 taken 1 time.
✗ Branch 55 → 58 not taken.
✗ Branch 222 → 223 not taken.
✗ Branch 222 → 225 not taken.
|
2 | seq.Add("clas12::rga::MomentumCorrection"); // momentum corrections (a transformer algorithm) |
| 48 | // seq.PrintSequence(); | ||
| 49 | |||
| 50 | // set log levels | ||
| 51 | // NOTE: this can also be done in a config file | ||
| 52 |
3/8✓ Branch 60 → 61 taken 1 time.
✗ Branch 60 → 298 not taken.
✓ Branch 61 → 62 taken 1 time.
✗ Branch 61 → 234 not taken.
✓ Branch 70 → 71 taken 1 time.
✗ Branch 70 → 73 not taken.
✗ Branch 234 → 235 not taken.
✗ Branch 234 → 237 not taken.
|
2 | seq.SetOption("clas12::EventBuilderFilter", "log", "info"); |
| 53 |
3/8✓ Branch 75 → 76 taken 1 time.
✗ Branch 75 → 298 not taken.
✓ Branch 76 → 77 taken 1 time.
✗ Branch 76 → 246 not taken.
✓ Branch 85 → 86 taken 1 time.
✗ Branch 85 → 88 not taken.
✗ Branch 246 → 247 not taken.
✗ Branch 246 → 249 not taken.
|
2 | seq.SetOption("clas12::SectorFinder", "log", "info"); |
| 54 |
3/8✓ Branch 90 → 91 taken 1 time.
✗ Branch 90 → 298 not taken.
✓ Branch 91 → 92 taken 1 time.
✗ Branch 91 → 258 not taken.
✓ Branch 100 → 101 taken 1 time.
✗ Branch 100 → 103 not taken.
✗ Branch 258 → 259 not taken.
✗ Branch 258 → 261 not taken.
|
2 | seq.SetOption("clas12::rga::MomentumCorrection", "log", "info"); |
| 55 | |||
| 56 | // set algorithm options | ||
| 57 | // NOTE: this can also be done in a config file, but setting options here OVERRIDES config file settings | ||
| 58 |
6/16✓ Branch 105 → 106 taken 1 time.
✗ Branch 105 → 298 not taken.
✓ Branch 106 → 107 taken 1 time.
✗ Branch 106 → 276 not taken.
✓ Branch 107 → 108 taken 1 time.
✗ Branch 107 → 270 not taken.
✓ Branch 108 → 109 taken 1 time.
✗ Branch 108 → 264 not taken.
✓ Branch 114 → 115 taken 1 time.
✗ Branch 114 → 117 not taken.
✓ Branch 119 → 120 taken 1 time.
✗ Branch 119 → 122 not taken.
✗ Branch 270 → 271 not taken.
✗ Branch 270 → 273 not taken.
✗ Branch 276 → 277 not taken.
✗ Branch 276 → 279 not taken.
|
3 | seq.SetOption<std::vector<int>>("clas12::EventBuilderFilter", "pids", {11, 211, -211}); |
| 59 | |||
| 60 | // start the algorithms | ||
| 61 |
1/2✓ Branch 122 → 123 taken 1 time.
✗ Branch 122 → 298 not taken.
|
1 | seq.Start(banks); |
| 62 | |||
| 63 | // get bank index, for each bank we want to use after Iguana algorithms run | ||
| 64 | // NOTE: new banks from creator algorithms are initialized by `Start` | ||
| 65 |
2/4✓ Branch 123 → 124 taken 1 time.
✗ Branch 123 → 298 not taken.
✓ Branch 124 → 125 taken 1 time.
✗ Branch 124 → 280 not taken.
|
1 | auto b_config = iguana::tools::GetBankIndex(banks, "RUN::config"); |
| 66 |
2/4✓ Branch 130 → 131 taken 1 time.
✗ Branch 130 → 298 not taken.
✓ Branch 131 → 132 taken 1 time.
✗ Branch 131 → 286 not taken.
|
1 | auto b_particle = iguana::tools::GetBankIndex(banks, "REC::Particle"); |
| 67 |
2/4✓ Branch 137 → 138 taken 1 time.
✗ Branch 137 → 298 not taken.
✓ Branch 138 → 139 taken 1 time.
✗ Branch 138 → 292 not taken.
|
1 | auto b_sector = seq.GetCreatedBankIndex(banks, "clas12::SectorFinder"); // newly created bank; string parameter is algorithm name, not bank name |
| 68 | |||
| 69 | // run the algorithm sequence on each event | ||
| 70 | int iEvent = 0; | ||
| 71 |
5/8✓ Branch 170 → 171 taken 101 times.
✗ Branch 170 → 298 not taken.
✓ Branch 171 → 172 taken 101 times.
✗ Branch 171 → 175 not taken.
✓ Branch 172 → 173 taken 101 times.
✗ Branch 172 → 174 not taken.
✓ Branch 173 → 174 taken 100 times.
✓ Branch 173 → 175 taken 1 time.
|
101 | while(reader.next(banks) && (numEvents == 0 || iEvent++ < numEvents)) { |
| 72 | |||
| 73 | // references to this event's banks | ||
| 74 |
1/2✓ Branch 174 → 145 taken 100 times.
✗ Branch 174 → 298 not taken.
|
100 | auto& bank_config = banks.at(b_config); |
| 75 |
1/2✓ Branch 145 → 146 taken 100 times.
✗ Branch 145 → 298 not taken.
|
100 | auto& bank_particle = banks.at(b_particle); |
| 76 |
1/2✓ Branch 146 → 147 taken 100 times.
✗ Branch 146 → 298 not taken.
|
100 | auto& bank_sector = banks.at(b_sector); |
| 77 | |||
| 78 | // print the event number | ||
| 79 |
2/4✓ Branch 148 → 149 taken 100 times.
✗ Branch 148 → 298 not taken.
✓ Branch 149 → 150 taken 100 times.
✗ Branch 149 → 298 not taken.
|
100 | fmt::println("===== EVENT {} =====", bank_config.getInt("event", 0)); |
| 80 | |||
| 81 | // print the particle bank before Iguana algorithms | ||
| 82 |
1/2✓ Branch 149 → 150 taken 100 times.
✗ Branch 149 → 298 not taken.
|
100 | fmt::println("----- BEFORE IGUANA -----"); |
| 83 |
1/2✓ Branch 150 → 151 taken 100 times.
✗ Branch 150 → 298 not taken.
|
100 | bank_particle.show(); // the original particle bank |
| 84 | |||
| 85 | // run the sequence of Iguana algorithms | ||
| 86 |
1/2✓ Branch 151 → 152 taken 100 times.
✗ Branch 151 → 298 not taken.
|
100 | seq.Run(banks); |
| 87 | |||
| 88 | // print the banks after Iguana algorithms | ||
| 89 |
1/2✓ Branch 152 → 153 taken 100 times.
✗ Branch 152 → 298 not taken.
|
100 | fmt::println("----- AFTER IGUANA -----"); |
| 90 |
1/2✓ Branch 153 → 154 taken 100 times.
✗ Branch 153 → 298 not taken.
|
100 | bank_particle.show(); // the filtered particle bank, with corrected momenta |
| 91 |
1/2✓ Branch 154 → 155 taken 100 times.
✗ Branch 154 → 298 not taken.
|
100 | bank_sector.show(); // the new sector bank |
| 92 | |||
| 93 | // print a table; first the header | ||
| 94 |
1/2✓ Branch 156 → 157 taken 100 times.
✗ Branch 156 → 298 not taken.
|
100 | fmt::print("----- Analysis Particles -----\n"); |
| 95 | 100 | fmt::print(" {:<20} {:<20} {:<20} {:<20}\n", "row == pindex", "PDG", "|p|", "sector"); | |
| 96 | // then print a row for each particle | ||
| 97 | // - use the `hipo::bank::getRowList()` method to loop over the bank rows that PASS the filter | ||
| 98 | // - if you'd rather loop over ALL bank rows, iterate from `i=0` up to `i < hipo::bank::getRows()` instead | ||
| 99 |
3/4✓ Branch 157 → 158 taken 100 times.
✗ Branch 157 → 298 not taken.
✓ Branch 167 → 159 taken 201 times.
✓ Branch 167 → 168 taken 100 times.
|
301 | for(auto const& row : bank_particle.getRowList()) { |
| 100 | 201 | auto p = std::hypot( | |
| 101 | bank_particle.getFloat("px", row), | ||
| 102 | bank_particle.getFloat("py", row), | ||
| 103 | bank_particle.getFloat("pz", row)); | ||
| 104 | 201 | auto pdg = bank_particle.getInt("pid", row); | |
| 105 | 201 | auto sector = bank_sector.getInt("sector", row); | |
| 106 | 201 | fmt::print(" {:<20} {:<20} {:<20.3f} {:<20}\n", row, pdg, p, sector); | |
| 107 | } | ||
| 108 | 100 | fmt::print("\n"); | |
| 109 | } | ||
| 110 | |||
| 111 | // stop algorithms | ||
| 112 |
1/2✓ Branch 175 → 176 taken 1 time.
✗ Branch 175 → 298 not taken.
|
1 | seq.Stop(); |
| 113 | return 0; | ||
| 114 | 1 | } | |
| 115 |