Back to home page

EIC code displayed by LXR

 
 

    


Warning, /EICrecon/docs/tutorial/04-factory.md is written in an unsupported language. File is not indexed.

0001 
0002 | title                                                                                  | teaching | exercises | questions                                            | objectives                                                                                                                                                                                                                                                                           | keypoints                                                                                 |
0003 |----------------------------------------------------------------------------------------|----------|-----------|------------------------------------------------------|--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------|-------------------------------------------------------------------------------------------|
0004 | Creating or modifying a JANA factory in order to implement a reconstruction algorithm  | 15       | 20        | How to write a reconstruction algorithm in EICrecon? | Learn how to create a new factory in EICrecon that supplies a reconstruction algorithm for all to use.<br/><br/>Understand the directory structure for where the factory should be placed in the source tree.<br/><br/>Understand how to use a generic algorithm in a JANA factory.  | Create a factory for reconstructing single subdetector data or for global reconstruction. |
0005 
0006 ## Notice
0007 
0008 **This section of the tutorial is outdated due to switch from `JFactory` to `JOmniFactory`. Refer to a [tutorial on reconstruction algorithms](https://eic.github.io/tutorial-reconstruction-algorithms/) for up-to-date information.**
0009 
0010 ## Introduction
0011 
0012 
0013 Now that you've learned about JANA plugins and JEventProcessors, let's talk about JFactories. JFactories are another essential JANA component just like JEventProcessors and JEventSources. While JEventProcessors are used for _aggregating_ results from each event into a structured output such as a histogram or a file, JFactories are used for computing those results in an organized way.
0014 
0015 
0016 ### When do I use a JFactory?
0017 
0018 - If you have an input file and need to read data model objects from it, use a JEventSource.
0019 - If you have an output file (or histogram) and wish to write data model objects to it, use a JEventProcessor.
0020 - If you have some data model objects and wish to produce a new data model object, use a JFactory.
0021 
0022 
0023 ### Why should I prefer writing a JFactory?
0024 
0025 1. They make your code reusable. Different people can use your results later without having to understand the specifics of what you did.
0026 
0027 2. If you are consuming some data which doesn't look right to you, JFactories make it extremely easy to pinpoint exactly which code produced this data.
0028 
0029 3. EICrecon needs to run multithreaded, and using JFactories can help steer you away from introducing thorny parallelism bugs.
0030 
0031 4. You can simply ask for the results you need and the JFactory will provide it. If nobody needs the results from the JFactory, it won't be run. If the results were already in the input file, it won't be run. If there are multiple consumers, the results are only computed once and then cached. If the JFactory relies on results from other JFactories, it will call them transparently and recursively.
0032 
0033 
0034 ### When do I create my own plugin?
0035 
0036 - If you are doing a one-off prototype, it's fine to just use a ROOT macro.
0037 - If you are writing code you'll probably return to, we recommend putting the code in a standalone (i.e. outside of the EICrecon source tree) plugin.
0038 - If you are writing code other people will probably want to run, we recommend adding your plugin to the EICrecon source tree.
0039 - If you are writing a JFactory, we recommend adding it to the EICrecon source tree, either to an existing plugin or to a new one.
0040 
0041 
0042 ## Algorithms vs Factories
0043 
0044 In general, a Factory is a programming pattern for constructing objects in an abstract way. Oftentimes, the Factory is calling an algorithm under the hood.
0045 This algorithm may be very generic. For instance, we may have a Factory that produces Cluster objects for a barrel calorimeter, and it calls a clustering algorithm
0046 that doesn't care at all about barrel calorimeters, just the position and energy of energy of each CalorimeterHit object. Perhaps multiple factories for creating clusters
0047 for completely different detectors are all using the same algorithm.
0048 
0049 Note that Gaudi provides an abstraction called "Algorithm" which is essentially its own version of a JFactory. In EICrecon, we have been separating out _generic algorithms_ from the old Gaudi and new JANA code so that these can be developed and tested independently. To see an example of how a generic algorithm is being implemented, look at these examples:
0050 
0051 src/detectors/EEMC/RawCalorimeterHit_factory_EcalEndcapNRawHits.h
0052 src/algorithms/calorimetry/CalorimeterHitDigi.h
0053 src/algorithms/calorimetry/CalorimeterHitDigi.cc
0054 
0055 Using generic algorithms makes things slightly more complex. However, the generic algorithms can be recycled for use in multiple detector systems which adds some simplification.
0056 
0057 
0058 ## Parallelism considerations
0059 
0060 
0061 JEventProcessors observe the entire _event stream_, and require a _critical section_ where only one thread is allowed to modify a shared resource (such as a histogram) at any time.
0062 JFactories, on the other hand, only observe a single event at a time, and work on each event independently. Each worker thread is given an independent event with its own set of factories. This means that for a given JFactory instance, there will be only one thread working on one event at any time. You get the benefits of multithreading _without_ having to make each JFactory thread-safe.
0063 
0064 
0065 You can write JFactories in an almost-functional style, but you can also cache some data on the JFactory that will stick around from event-to-event. This is useful for things like conditions and geometry data, where for performance reasons you don't want to be doing a deep lookup on every event. Instead, you can write callbacks such as `BeginRun()`, where you can update your cached values when the run number changes.
0066 
0067 
0068 Note that just because the JFactory _can_ be called in parallel doesn't mean it always will. If you call event->Get() from inside `JEventProcessor::ProcessSequential`, in particular, the factory will run single-threaded and slow everything down. However, if you call it using `Prefetch` instead, it will run in parallel and you may get a speed boost.
0069 
0070 
0071 ## How do I use an existing JFactory? ##
0072 
0073 Using an existing JFactory is extremely easy! Any time you are someplace where you have access to a `JEvent` object, do this:
0074 
0075 
0076 ~~~ c++
0077 
0078 auto clusters = event->Get<edm4eic::Cluster>("EcalEndcapNIslandClusters");
0079 
0080 for (auto c : clusters) {
0081   // ... do something with a cluster
0082 }
0083 
0084 ~~~
0085 
0086 As you can see, it doesn't matter whether the `Cluster` objects were calculated from some simpler objects, or were simply loaded from a file. This is a very powerful concept.
0087 
0088 One thing we might want to do is to swap one factory for another, possibly even at runtime. This is easy to do if you just make the factory tag be a parameter:
0089 
0090 
0091 ~~~ c++
0092 
0093 std::string my_cluster_source = "EcalEndcapNIslandClusters";  // Make this be a parameter
0094 app->SetDefaultParameter("MyPlugin:MyAnalysis:my_cluster_source", my_cluster_source, "Cluster source for MyAnalysis");
0095 auto clusters = event->Get<edm4eic::Cluster>(my_cluster_source);
0096 ~~~
0097 
0098 ## How do I create a new JFactory?
0099 
0100 We are going to add a new JFactory inside EICrecon.
0101 
0102 `src/detectors/EEMC/Cluster_factory_EcalEndcapNIslandClusters.h`:
0103 ~~~ c++
0104 #pragma once
0105 
0106 #include <edm4eic/Cluster.h>
0107 #include <JANA/JFactoryT.h>
0108 #include <services/log/Log_service.h>
0109 
0110 class Cluster_factory_EcalEndcapNIslandClusters : public JFactoryT<edm4eic::Cluster> {
0111 public:
0112 
0113     Cluster_factory_EcalEndcapNIslandClusters(); // Constructor
0114 
0115     void Init() override;
0116     // Gets called exactly once at the beginning of the JFactory's life
0117 
0118     void ChangeRun(const std::shared_ptr<const JEvent> &event) override {};
0119     // Gets called on events where the run number has changed (before Process())
0120 
0121     void Process(const std::shared_ptr<const JEvent> &event) override;
0122     // Gets called on every event
0123 
0124     void Finish() override {};
0125     // Gets called exactly once at the end of the JFactory's life
0126 
0127 private:
0128     float m_scaleFactor;
0129     std::shared_ptr<spdlog::logger> m_log;
0130 
0131 };
0132 
0133 ~~~
0134 
0135 `src/detectors/EEMC/Cluster_factory_EcalEndcapNIslandClusters.cc`:
0136 ~~~ c++
0137 #include "Cluster_factory_EcalEndcapNIslandClusters.h"
0138 
0139 #include <edm4eic/ProtoCluster.h>
0140 #include <JANA/JEvent.h>
0141 
0142 
0143 Cluster_factory_EcalEndcapNIslandClusters::Cluster_factory_EcalEndcapNIslandClusters() {
0144 
0145     SetTag("EcalEndcapNIslandClusters");
0146 }
0147 
0148 
0149 void Cluster_factory_EcalEndcapNIslandClusters::Init() {
0150     auto app = GetApplication();
0151 
0152     // This is an example of how to declare a configuration parameter that
0153     // can be set at run time. e.g. with -PEEMC:EcalEndcapNIslandClusters:scaleFactor=0.97
0154     m_scaleFactor =0.98;
0155     app->SetDefaultParameter("EEMC:EcalEndcapNIslandClusters:scaleFactor", m_scaleFactor, "Energy scale factor");
0156 
0157     // This is how you access shared resources using the JService interface
0158     m_log = app->GetService<Log_service>()->logger("EcalEndcapNIslandClusters");
0159 }
0160 
0161 
0162 void Cluster_factory_EcalEndcapNIslandClusters::Process(const std::shared_ptr<const JEvent> &event) {
0163 
0164     m_log->info("Processing event {}", event->GetEventNumber());
0165 
0166     // Grab inputs
0167     auto protoclusters = event->Get<edm4eic::ProtoCluster>("EcalEndcapNIslandProtoClusters");
0168 
0169     // Loop over protoclusters and turn each into a cluster
0170     std::vector<edm4eic::Cluster*> outputClusters;
0171     for( auto proto : protoclusters ) {
0172 
0173         // ======================
0174         // Algorithm goes here!
0175         // ======================
0176 
0177         auto cluster = new edm4eic::Cluster(
0178             0, // type
0179             energy * m_scaleFactor,
0180             sqrt(energyError_squared),
0181             time,
0182             timeError,
0183             proto->hits_size(),
0184             position,
0185             edm4eic::Cov3f(), // positionError,
0186             0.0,              // intrinsicTheta,
0187             0.0,              // intrinsicPhi,
0188             edm4eic::Cov2f()  // intrinsicDirectionError
0189             );
0190 
0191         outputClusters.push_back( cluster );
0192     }
0193 
0194     // Hand ownership of algorithm objects over to JANA
0195     Set(outputClusters);
0196 }
0197 
0198 
0199 ~~~
0200 
0201 We can now fill in the algorithm with anything we like!
0202 
0203 
0204 ~~~ c++
0205         // Grab inputs
0206         auto protoclusters = event->Get<edm4eic::ProtoCluster>("EcalEndcapNIslandProtoClusters");
0207 
0208         // Loop over protoclusters and turn each into a cluster
0209         std::vector<edm4eic::Cluster*> outputClusters;
0210         for( auto proto : protoclusters ) {
0211 
0212             // Fill cumulative values by looping over all hits in proto cluster
0213             float energy = 0;
0214             double energyError_squared = 0.0;
0215             float time = 1.0E8;
0216             float timeError;
0217             edm4hep::Vector3f position;
0218             double sum_weights = 0.0;
0219             for( uint32_t ihit=0; ihit<proto->hits_size() ; ihit++){
0220                 auto const &hit = proto->getHits(ihit);
0221                 auto weight = proto->getWeights(ihit);
0222                 energy += hit.getEnergy();
0223                 energyError_squared += std::pow(hit.getEnergyError(), 2.0);
0224                 if( hit.getTime() < time ){
0225                     time = hit.getTime();            // use earliest time
0226                     timeError = hit.getTimeError();  // use error of earliest time
0227                 }
0228                 auto &p = hit.getPosition();
0229                 position.x += p.x*weight;
0230                 position.y += p.y*weight;
0231                 position.z += p.z*weight;
0232                 sum_weights += weight;
0233             }
0234 
0235             // Normalize position
0236             position.x /= sum_weights;
0237             position.y /= sum_weights;
0238             position.z /= sum_weights;
0239 
0240             // Create a cluster object from values accumulated from hits above
0241             auto cluster = new edm4eic::Cluster(
0242                 0, // type (?))
0243                 energy * m_scaleFactor,
0244                 sqrt(energyError_squared),
0245                 time,
0246                 timeError,
0247                 proto->hits_size(),
0248                 position,
0249 
0250                 // Not sure how to calculate these last few
0251                 edm4eic::Cov3f(), // positionError,
0252                 0.0, // intrinsicTheta,
0253                 0.0, // intrinsicPhi,
0254                 edm4eic::Cov2f() // intrinsicDirectionError
0255                 );
0256 
0257             outputClusters.push_back( cluster );
0258         }
0259 
0260         // Hand ownership of algorithm objects over to JANA
0261         Set(outputClusters);
0262 
0263 ~~~
0264 
0265 
0266 You can't pass JANA a JFactory directly (because it needs to create an arbitrary number of them on the fly). Instead you register a `JFactoryGenerator` object:
0267 
0268 `src/detectors/EEMC/EEMC.cc`
0269 ~~~ c++
0270 
0271 // In your plugin's init
0272 
0273 #include <JANA/JFactoryGenerator.h>
0274 // ...
0275 #include "Cluster_factory_EcalEndcapNIslandClusters.h"
0276 
0277 extern "C" {
0278     void InitPlugin(JApplication *app) {
0279         InitJANAPlugin(app);
0280         // ...
0281 
0282         app->Add(new JFactoryGeneratorT<Cluster_factory_EcalEndcapNIslandClusters>());
0283      }
0284 
0285 ~~~
0286 
0287 Finally, we go ahead and trigger the factory (remember, factories won't do anything unless activated by a JEventProcessor). You can open the
0288 
0289 ~~~ bash
0290 eicrecon in.root -Ppodio:output_file=out.root -Ppodio:output_collections=EcalEndcapNIslandClusters -Pjana:nevents=10
0291 ~~~
0292 
0293 
0294 Your exercise is to get this JFactory working! You can tweak the algorithm, add log messages, add additional config parameters, etc.
0295 
0296 
0297 
0298 {% include links.md %}