Back to home page

EIC code displayed by LXR

 
 

    


Warning, /EICrecon/docs/tutorial/03-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 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 ## Algorithms vs Factories
0035 
0036 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.
0037 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
0038 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
0039 for completely different detectors are all using the same algorithm.
0040 
0041 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:
0042 
0043 src/detectors/EEMC/RawCalorimeterHit_factory_EcalEndcapNRawHits.h
0044 src/algorithms/calorimetry/CalorimeterHitDigi.h
0045 src/algorithms/calorimetry/CalorimeterHitDigi.cc
0046 
0047 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.
0048 
0049 
0050 ## Parallelism considerations
0051 
0052 
0053 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.
0054 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.
0055 
0056 
0057 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.
0058 
0059 
0060 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.
0061 
0062 
0063 ## How do I use an existing JFactory? ##
0064 
0065 Using an existing JFactory is extremely easy! Any time you are someplace where you have access to a `JEvent` object, do this:
0066 
0067 
0068 ~~~ c++
0069 
0070 auto clusters = event->Get<edm4eic::Cluster>("EcalEndcapNIslandClusters");
0071 
0072 for (auto c : clusters) {
0073   // ... do something with a cluster
0074 }
0075 
0076 ~~~
0077 
0078 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.
0079 
0080 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:
0081 
0082 
0083 ~~~ c++
0084 
0085 std::string my_cluster_source = "EcalEndcapNIslandClusters";  // Make this be a parameter
0086 app->SetDefaultParameter("MyPlugin:MyAnalysis:my_cluster_source", my_cluster_source, "Cluster source for MyAnalysis");
0087 auto clusters = event->Get<edm4eic::Cluster>(my_cluster_source);
0088 ~~~
0089 
0090 ## How do I create a new JFactory?
0091 
0092 We are going to add a new JFactory inside EICrecon.
0093 
0094 `src/detectors/EEMC/Cluster_factory_EcalEndcapNIslandClusters.h`:
0095 ~~~ c++
0096 #pragma once
0097 
0098 #include <edm4eic/Cluster.h>
0099 #include <JANA/JFactoryT.h>
0100 #include <services/log/Log_service.h>
0101 
0102 class Cluster_factory_EcalEndcapNIslandClusters : public JFactoryT<edm4eic::Cluster> {
0103 public:
0104 
0105     Cluster_factory_EcalEndcapNIslandClusters(); // Constructor
0106 
0107     void Init() override;
0108     // Gets called exactly once at the beginning of the JFactory's life
0109 
0110     void ChangeRun(const std::shared_ptr<const JEvent> &event) override {};
0111     // Gets called on events where the run number has changed (before Process())
0112 
0113     void Process(const std::shared_ptr<const JEvent> &event) override;
0114     // Gets called on every event
0115 
0116     void Finish() override {};
0117     // Gets called exactly once at the end of the JFactory's life
0118 
0119 private:
0120     float m_scaleFactor;
0121     std::shared_ptr<spdlog::logger> m_log;
0122 
0123 };
0124 
0125 ~~~
0126 
0127 `src/detectors/EEMC/Cluster_factory_EcalEndcapNIslandClusters.cc`:
0128 ~~~ c++
0129 #include "Cluster_factory_EcalEndcapNIslandClusters.h"
0130 
0131 #include <edm4eic/ProtoCluster.h>
0132 #include <JANA/JEvent.h>
0133 
0134 
0135 Cluster_factory_EcalEndcapNIslandClusters::Cluster_factory_EcalEndcapNIslandClusters() {
0136 
0137     SetTag("EcalEndcapNIslandClusters");
0138 }
0139 
0140 
0141 void Cluster_factory_EcalEndcapNIslandClusters::Init() {
0142     auto app = GetApplication();
0143 
0144     // This is an example of how to declare a configuration parameter that
0145     // can be set at run time. e.g. with -PEEMC:EcalEndcapNIslandClusters:scaleFactor=0.97
0146     m_scaleFactor =0.98;
0147     app->SetDefaultParameter("EEMC:EcalEndcapNIslandClusters:scaleFactor", m_scaleFactor, "Energy scale factor");
0148 
0149     // This is how you access shared resources using the JService interface
0150     m_log = app->GetService<Log_service>()->logger("EcalEndcapNIslandClusters");
0151 }
0152 
0153 
0154 void Cluster_factory_EcalEndcapNIslandClusters::Process(const std::shared_ptr<const JEvent> &event) {
0155 
0156     m_log->info("Processing event {}", event->GetEventNumber());
0157 
0158     // Grab inputs
0159     auto protoclusters = event->Get<edm4eic::ProtoCluster>("EcalEndcapNIslandProtoClusters");
0160 
0161     // Loop over protoclusters and turn each into a cluster
0162     std::vector<edm4eic::Cluster*> outputClusters;
0163     for( auto proto : protoclusters ) {
0164 
0165         // ======================
0166         // Algorithm goes here!
0167         // ======================
0168 
0169         auto cluster = new edm4eic::Cluster(
0170             0, // type
0171             energy * m_scaleFactor,
0172             sqrt(energyError_squared),
0173             time,
0174             timeError,
0175             proto->hits_size(),
0176             position,
0177             edm4eic::Cov3f(), // positionError,
0178             0.0,              // intrinsicTheta,
0179             0.0,              // intrinsicPhi,
0180             edm4eic::Cov2f()  // intrinsicDirectionError
0181             );
0182 
0183         outputClusters.push_back( cluster );
0184     }
0185 
0186     // Hand ownership of algorithm objects over to JANA
0187     Set(outputClusters);
0188 }
0189 
0190 
0191 ~~~
0192 
0193 We can now fill in the algorithm with anything we like!
0194 
0195 
0196 ~~~ c++
0197         // Grab inputs
0198         auto protoclusters = event->Get<edm4eic::ProtoCluster>("EcalEndcapNIslandProtoClusters");
0199 
0200         // Loop over protoclusters and turn each into a cluster
0201         std::vector<edm4eic::Cluster*> outputClusters;
0202         for( auto proto : protoclusters ) {
0203 
0204             // Fill cumulative values by looping over all hits in proto cluster
0205             float energy = 0;
0206             double energyError_squared = 0.0;
0207             float time = 1.0E8;
0208             float timeError;
0209             edm4hep::Vector3f position;
0210             double sum_weights = 0.0;
0211             for( uint32_t ihit=0; ihit<proto->hits_size() ; ihit++){
0212                 auto const &hit = proto->getHits(ihit);
0213                 auto weight = proto->getWeights(ihit);
0214                 energy += hit.getEnergy();
0215                 energyError_squared += std::pow(hit.getEnergyError(), 2.0);
0216                 if( hit.getTime() < time ){
0217                     time = hit.getTime();            // use earliest time
0218                     timeError = hit.getTimeError();  // use error of earliest time
0219                 }
0220                 auto &p = hit.getPosition();
0221                 position.x += p.x*weight;
0222                 position.y += p.y*weight;
0223                 position.z += p.z*weight;
0224                 sum_weights += weight;
0225             }
0226 
0227             // Normalize position
0228             position.x /= sum_weights;
0229             position.y /= sum_weights;
0230             position.z /= sum_weights;
0231 
0232             // Create a cluster object from values accumulated from hits above
0233             auto cluster = new edm4eic::Cluster(
0234                 0, // type (?))
0235                 energy * m_scaleFactor,
0236                 sqrt(energyError_squared),
0237                 time,
0238                 timeError,
0239                 proto->hits_size(),
0240                 position,
0241 
0242                 // Not sure how to calculate these last few
0243                 edm4eic::Cov3f(), // positionError,
0244                 0.0, // intrinsicTheta,
0245                 0.0, // intrinsicPhi,
0246                 edm4eic::Cov2f() // intrinsicDirectionError
0247                 );
0248 
0249             outputClusters.push_back( cluster );
0250         }
0251 
0252         // Hand ownership of algorithm objects over to JANA
0253         Set(outputClusters);
0254 
0255 ~~~
0256 
0257 
0258 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:
0259 
0260 `src/detectors/EEMC/EEMC.cc`
0261 ~~~ c++
0262 
0263 // In your plugin's init
0264 
0265 #include <JANA/JFactoryGenerator.h>
0266 // ...
0267 #include "Cluster_factory_EcalEndcapNIslandClusters.h"
0268 
0269 extern "C" {
0270     void InitPlugin(JApplication *app) {
0271         InitJANAPlugin(app);
0272         // ...
0273 
0274         app->Add(new JFactoryGeneratorT<Cluster_factory_EcalEndcapNIslandClusters>());
0275      }
0276 
0277 ~~~
0278 
0279 Finally, we go ahead and trigger the factory (remember, factories won't do anything unless activated by a JEventProcessor). You can open the
0280 
0281 ~~~ bash
0282 eicrecon in.root -Ppodio:output_file=out.root -Ppodio:output_collections=EcalEndcapNIslandClusters -Pjana:nevents=10
0283 ~~~
0284 
0285 
0286 Your exercise is to get this JFactory working! You can tweak the algorithm, add log messages, add additional config parameters, etc.
0287 
0288 
0289 
0290 {% include links.md %}