Hall-D Software  alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DEventProcessor_fcal_hists.cc
Go to the documentation of this file.
1 // $Id$
2 //
3 // File: DEventProcessor_fcal_hists.cc
4 // Created: Mon Apr 3 11:38:03 EDT 2006
5 // Creator: davidl (on Darwin swire-b241.jlab.org 8.4.0 powerpc)
6 //
7 
8 #include <map>
9 using namespace std;
10 
12 
13 #include <TLorentzVector.h>
14 
15 #include <DANA/DApplication.h>
16 #include <FCAL/DFCALCluster.h>
17 #include <TRACKING/DMCThrown.h>
18 #include <FCAL/DFCALHit.h>
19 #include <FCAL/DFCALGeometry.h>
20 
21 
22 // Routine used to create our DEventProcessor
23 extern "C"{
24 void InitPlugin(JApplication *app){
25  InitJANAPlugin(app);
26  app->AddProcessor(new DEventProcessor_fcal_hists());
27 }
28 } // "C"
29 
30 
31 #define FCAL_Z_OFFSET 640.0-65.0 // I don't know what this value is ???
32 //#define FCAL_Z_OFFSET 170.0 // I don't know what this value is ???
33 #define PI_ZERO_MASS 0.13497
34 
35 //------------------
36 // init
37 //------------------
39 {
40 
41  dE_over_E_vs_E = new TH2D("dE_over_E_vs_E","Smeared-unsmeared energy diff for single FCAL blocks", 50, 0.0, 1.0, 50, -0.1, 0.1);
42 
43  return NOERROR;
44 }
45 
46 //------------------
47 // brun
48 //------------------
49 jerror_t DEventProcessor_fcal_hists::brun(JEventLoop *eventLoop, int32_t runnumber)
50 {
51  return NOERROR;
52 }
53 
54 //------------------
55 // evnt
56 //------------------
57 jerror_t DEventProcessor_fcal_hists::evnt(JEventLoop *loop, uint64_t eventnumber)
58 {
59  // extract the FCAL Geometry (for isBlockActive() and positionOnFace())
60  vector<const DFCALGeometry*> fcalGeomVect;
61  vector<const DFCALHit*> hits;
62  vector<const DFCALHit*> truthhits;
63  loop->Get(fcalGeomVect);
64  loop->Get(hits);
65  loop->Get(truthhits,"TRUTH");
66 
67  const DFCALGeometry& fcalGeom = *(fcalGeomVect[0]);
68  if(fcalGeomVect.size()<1)return OBJECT_NOT_AVAILABLE;
69 
70  // Create STL map of "real" hits that can be used below to match
71  // with MC hits.
72  map<int,const DFCALHit*> hit_map;
73  for(unsigned int i=0; i<hits.size(); i++){
74  const DFCALHit *hit = hits[i];
75 
76  int channel = fcalGeom.channel(hit->row, hit->column);
77  hit_map[channel] = hit;
78  }
79 
80  // Since we are modifying histograms local to this plugin, it will not interfere with other ROOT operations: can use plugin-wide ROOT fill lock
81  japp->RootFillLock(this); //ACQUIRE ROOT FILL LOCK
82 
83  // Loop over "truth" hits and match them to real hits by using assuming
84  // one hit per channel
85  for(unsigned int i=0; i<truthhits.size(); i++){
86  const DFCALHit *truthhit = truthhits[i];
87  map<int,const DFCALHit*>::iterator iter = hit_map.find(fcalGeom.channel(truthhit->row, truthhit->column));
88  if(iter==hit_map.end())continue;
89 
90  const DFCALHit *hit = iter->second;
91  double deltaE = hit->E - truthhit->E;
92  dE_over_E_vs_E->Fill(truthhit->E, deltaE/truthhit->E);
93  }
94 
95  japp->RootFillUnLock(this); //RELEASE ROOT FILL LOCK
96 
97  return NOERROR;
98 }
99 
100 //------------------
101 // erun
102 //------------------
104 {
105  // Any final calculations on histograms (like dividing them)
106  // should be done here. This may get called more than once.
107  return NOERROR;
108 }
109 
110 //------------------
111 // fini
112 //------------------
114 {
115  return NOERROR;
116 }
117 
jerror_t fini(void)
Called after last event of last event source has been processed.
JApplication * japp
InitPlugin_t InitPlugin
int row
Definition: DFCALHit.h:23
jerror_t evnt(JEventLoop *eventLoop, uint64_t eventnumber)
Called every event.
jerror_t init(void)
Called once at program start.
int channel(int row, int column) const
jerror_t erun(void)
Called everytime run number changes, provided brun has been called.
float E
Definition: DFCALHit.h:27
jerror_t brun(JEventLoop *eventLoop, int32_t runnumber)
Called everytime a new run number is detected.
int column
Definition: DFCALHit.h:24