Hall-D Software  alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
JEventProcessor_CDC_amp.cc
Go to the documentation of this file.
1 // $Id$
2 //
3 // File: JEventProcessor_CDC_amp.cc
4 // Created: Tue Sep 6 10:13:02 EDT 2016
5 // Creator: njarvis (on Linux egbert 2.6.32-642.3.1.el6.x86_64 x86_64)
6 //
7 
9 
10 
11 
12 // Routine used to create our JEventProcessor
13 #include <JANA/JApplication.h>
14 #include <JANA/JFactory.h>
15 extern "C"{
16 void InitPlugin(JApplication *app){
17  InitJANAPlugin(app);
18  app->AddProcessor(new JEventProcessor_CDC_amp());
19 }
20 } // "C"
21 
22 
23 //------------------
24 // JEventProcessor_CDC_amp (Constructor)
25 //------------------
27 {
28 
29 }
30 
31 //------------------
32 // ~JEventProcessor_CDC_amp (Destructor)
33 //------------------
35 {
36 
37 }
38 
39 //------------------
40 // init
41 //------------------
43 {
44  // This is called once at program startup.
45 
46  TDirectory *main = gDirectory;
47  gDirectory->mkdir("CDC_amp")->cd();
48 
49  time = new TH1I("time","CDC time (hits on tracks); raw time",200,0,2000);
50 
51  a = new TH1I("a","CDC amplitude (hits on tracks);amplitude - pedestal",4096,0,4096);
52 
53  an = new TH2I("an","CDC amplitude vs n (hits on tracks); n; amplitude - pedestal",3522,0,3522,4096,0,4096);
54 
55  atime = new TH2D("atime","CDC amplitude vs time (hits on tracks); raw time; amplitude - pedestal",200,0,2000,4096,0,4096);
56 
57  atheta = new TH2D("atheta","CDC amplitude vs theta (hits on tracks); theta; amplitude - pedestal",180,0,180,4096,0,4096);
58 
59 
60  a30 = new TH1I("a30","CDC amplitude (tracks, theta 28-32 deg, z 52-78cm);amplitude - pedestal",4096,0,4096);
61 
62  a30_100ns = new TH1I("a30_100ns","CDC amplitude, hits with drift time < 100ns (tracks, theta 28-32 deg, z 52-78cm);amplitude - pedestal",4096,0,4096);
63 
64  an30_100ns = new TH2I("an30_100ns","CDC amplitude vs n, hits with drift time < 100ns (tracks, theta 28-32 deg, z 52-78cm); n; amplitude - pedestal",3522,0,3522,4096,0,4096);
65 
66  atime30 = new TH2D("atime30","CDC amplitude vs time (tracks, theta 28-32 deg, z 52-78cm); raw time; amplitude - pedestal",200,0,2000,4096,0,4096);
67 
68  adoca30 = new TH2D("adoca30","CDC amplitude vs DOCA (tracks, theta 28-32 deg, z 52-78cm); DOCA (cm); amplitude - pedestal",200,0,1.0,4096,0,4096);
69 
70 
71 
72  a45 = new TH1I("a45","CDC amplitude (tracks, theta 43-47 deg, z 52-78cm);amplitude - pedestal",4096,0,4096);
73 
74  a45_100ns = new TH1I("a45_100ns","CDC amplitude, hits with drift time < 100ns (tracks, theta 43-47 deg, z 52-78cm); amplitude - pedestal",4096,0,4096);
75 
76  an45_100ns = new TH2I("an45_100ns","CDC amplitude vs n, hits with drift time < 100ns (tracks, theta 43-47 deg, z 52-78cm); n; amplitude - pedestal",3522,0,3522,4096,0,4096);
77 
78  atime45 = new TH2D("atime45","CDC amplitude vs time (tracks, theta 43-47 deg, z 52-78cm); raw time; amplitude - pedestal",200,0,2000,4096,0,4096);
79 
80  adoca45 = new TH2D("adoca45","CDC amplitude vs DOCA (tracks, theta 43-47 deg, z 52-78cm); DOCA (cm); amplitude - pedestal",200,0,1.0,4096,0,4096);
81 
82 
83 
84 
85 
86 
87  a90 = new TH1I("a90","CDC amplitude (tracks, theta 85-95 deg, z 52-78cm);amplitude - pedestal",4096,0,4096);
88 
89  a90_100ns = new TH1I("a90_100ns","CDC amplitude, hits with drift time < 100ns (tracks, theta 85-95 deg, z 52-78cm);amplitude - pedestal",4096,0,4096);
90 
91  an90_100ns = new TH2I("an90_100ns","CDC amplitude vs n, hits with drift time < 100ns (tracks, theta 85-95 deg, z 52-78cm); n; amplitude - pedestal",3522,0,3522,4096,0,4096);
92 
93  atime90 = new TH2D("atime90","CDC amplitude vs time (tracks, theta 85-95 deg, z 52-78cm); raw time; amplitude - pedestal",200,0,2000,4096,0,4096);
94 
95  adoca90 = new TH2D("adoca90","CDC amplitude vs DOCA (tracks, theta 85-95 deg, z 52-78cm); DOCA (cm); amplitude - pedestal",200,0,1.0,4096,0,4096);
96 
97  time90 = new TH1I("time90","CDC time (tracks, theta 85-95 deg, z 52-78cm); raw time",200,0,2000);
98 
99  an90 = new TH2I("an90","CDC amplitude vs n (tracks, theta 85-95 deg, z 52-78cm); n; amplitude - pedestal",3522,0,3522,4096,0,4096);
100 
101  xt90 = new TH2D("xt90","CDC DOCA vs time (tracks, theta 85-95 deg, z 52-78cm); time (ns); DOCA (cm)",400,0,1200,400,0,1.0);
102 
103 
104  main->cd();
105 
106  return NOERROR;
107 }
108 
109 //------------------
110 // brun
111 //------------------
112 jerror_t JEventProcessor_CDC_amp::brun(JEventLoop *eventLoop, int32_t runnumber)
113 {
114  // This is called whenever the run number changes
115 
116  if (runnumber<40000) ASCALE = 8; // default for ASCALE before run 40,000 to be used if Df125config is not present
117 
118  if (runnumber>41497) ASCALE = 8; // default for ASCALE before run 40,000 to be used if Df125config is not present
119 
120  return NOERROR;
121 }
122 
123 //------------------
124 // evnt
125 //------------------
126 jerror_t JEventProcessor_CDC_amp::evnt(JEventLoop *loop, uint64_t eventnumber)
127 {
128  // This is called for every event. Use of common resources like writing
129  // to a file or filling a histogram should be mutex protected. Using
130  // loop->Get(...) to get reconstructed objects (and thereby activating the
131  // reconstruction algorithm) should be done outside of any mutex lock
132  // since multiple threads may call this method at the same time.
133  // Here's an example:
134  //
135  // vector<const MyDataClass*> mydataclasses;
136  // loop->Get(mydataclasses);
137  //
138  // japp->RootFillLock(this);
139  // ... fill histograms or trees ...
140  // japp->RootFillUnLock(this);
141 
142  int ring, straw, n; // ring number, straw number within ring, straw number overall (1 to 3522)
143 
144  uint32_t amp,ped,t; // dcdcdigihits raw quantities: time, pedestal, amplitude, quality factor, overflow count
145 
146  //scaling factors will be overridden by Df125Config if presqnt
147  // uint16_t ASCALE = 8; //amplitude
148  // uint16_t PSCALE = 1; //ped
149 
150  //add extra 0 at front to use offset[1] for ring 1 //used to calculate straw number n
151  int straw_offset[29] = {0,0,42,84,138,192,258,324,404,484,577,670,776,882,1005,1128,1263,1398,1544,1690,1848,2006,2176,2346,2528,2710,2907,3104,3313};
152 
153 
154  // select events with physics events, i.e., not LED and other front panel triggers
155  const DTrigger* locTrigger = NULL;
156  loop->GetSingle(locTrigger);
157  if(locTrigger->Get_L1FrontPanelTriggerBits() != 0)
158  return NOERROR;
159 
160  // use only events with track vertex in the target region
161  const DVertex* locVertex = NULL;
162  loop->GetSingle(locVertex);
163  double vertexz = locVertex->dSpacetimeVertex.Z();
164  if (vertexz < 52 || vertexz > 78) return NOERROR;
165 
166 
167  // // test whether this is simulated or real data (skip digihits for sim data)
168  // int SIMULATION;
169  // vector<const DMCThrown*> MCThrowns;
170  // loop->Get(MCThrowns);
171  // if (MCThrowns.empty()) SIMULATION = 0;
172  // if (!MCThrowns.empty()) SIMULATION = 1;
173 
174 
175  const DCDCHit *hit = NULL;
176  const DCDCDigiHit *digihit = NULL;
177  const Df125CDCPulse *cp = NULL;
178 
179  int netamp = 0;
180  float scaledped;
181 
182  int used[3523] = {0};
183 
184  //--------tracks---------------------------
185 
186 
187  vector<const DTrackTimeBased*> tracks;
188 
189  loop->Get(tracks);
190 
191 
192  for (uint32_t i=0; i<tracks.size(); i++) {
193 
194  DVector3 mom = tracks[i]->momentum();
195  double theta = mom.Theta();
196  theta = 180.0*theta/3.14159;
197 
198  vector<DTrackFitter::pull_t> pulls = tracks[i]->pulls;
199 
200  for (uint32_t j=0; j<pulls.size(); j++) {
201 
202  if (pulls[j].cdc_hit == NULL) continue;
203 
204  hit = NULL;
205  pulls[j].cdc_hit->GetSingle(hit);
206 
207  double doca = pulls[j].d;
208  double tdrift = pulls[j].tdrift;
209 
210  netamp = 0;
211 
212  // if (!SIMULATION) {
213 
214  digihit = NULL;
215  hit->GetSingle(digihit);
216  if (!digihit) continue;
217 
218  cp = NULL;
219  digihit->GetSingle(cp);
220  if (!cp) continue; //no CDCPulseData (happens occasionally)
221 
222  if (cp->time_quality_bit) continue;
223  if (cp->overflow_count) continue;
224 
225  amp = cp->first_max_amp;
226  ped = cp->pedestal;
227  t = cp->le_time;
228 
229  scaledped = ped*(float)PSCALE/(float)ASCALE;
230 
231  netamp = (int)amp - (int)scaledped;
232 
233  // }
234 
235  if (netamp ==0) continue;
236 
237  ring = (int)hit->ring;
238  straw = (int)hit->straw;
239 
240  n = straw_offset[ring] + straw;
241 
242 
243  japp->RootFillLock(this); //ACQUIRE ROOT LOCK!!
244 
245  if (!used[n]) {
246 
247  used[n] = 1;
248 
249  a->Fill(netamp);
250  an->Fill(n,netamp);
251  atheta->Fill(theta,netamp);
252  atime->Fill((int)t,netamp);
253  time->Fill((int)t);
254 
255  if ((theta>85) && (theta<95)) {
256 
257  a90->Fill(netamp);
258  atime90->Fill((int)t,netamp);
259  adoca90->Fill(doca,netamp);
260 
261  if (hit->t < 100.0) a90_100ns->Fill(netamp);
262  if (hit->t < 100.0) an90_100ns->Fill(n,netamp);
263 
264  an90->Fill(n,netamp);
265  time90->Fill((int)t);
266  xt90->Fill(tdrift,doca);
267 
268  } else if ((theta > 28.05) && (theta < 32)) {
269 
270  a30->Fill(netamp);
271  atime30->Fill((int)t,netamp);
272  adoca30->Fill(doca,netamp);
273 
274  if (hit->t < 100.0) a30_100ns->Fill(netamp);
275  if (hit->t < 100.0) an30_100ns->Fill(n,netamp);
276 
277  } else if ((theta > 43.07) && (theta < 47)) {
278 
279  a45->Fill(netamp);
280  atime45->Fill((int)t,netamp);
281  adoca45->Fill(doca,netamp);
282 
283  if (hit->t < 100.0) a45_100ns->Fill(netamp);
284  if (hit->t < 100.0) an45_100ns->Fill(n,netamp);
285 
286 
287  }
288 
289  } //if !used
290 
291 
292  japp->RootFillUnLock(this);
293 
294  }
295  }
296 
297  return NOERROR;
298 }
299 
300 //------------------
301 // erun
302 //------------------
304 {
305  // This is called whenever the run number changes, before it is
306  // changed to give you a chance to clean up before processing
307  // events from the next run number.
308  return NOERROR;
309 }
310 
311 //------------------
312 // fini
313 //------------------
315 {
316  // Called before program exit after event processing is finished.
317 
318 
319  return NOERROR;
320 }
uint32_t first_max_amp
from second word
Definition: Df125CDCPulse.h:68
TVector3 DVector3
Definition: DVector3.h:14
jerror_t fini(void)
Called after last event of last event source has been processed.
uint32_t pedestal
from second word
Definition: Df125CDCPulse.h:66
uint32_t Get_L1FrontPanelTriggerBits(void) const
uint32_t overflow_count
from first word
Definition: Df125CDCPulse.h:65
JApplication * japp
uint32_t le_time
from first word
Definition: Df125CDCPulse.h:63
float t
Definition: DCDCHit.h:22
jerror_t erun(void)
Called everytime run number changes, provided brun has been called.
InitPlugin_t InitPlugin
jerror_t init(void)
Called once at program start.
int ring
Definition: DCDCHit.h:18
uint32_t time_quality_bit
from first word
Definition: Df125CDCPulse.h:64
DLorentzVector dSpacetimeVertex
Definition: DVertex.h:28
void cdc_hit(Int_t &, Int_t &, Int_t &, Int_t[], Int_t, Int_t, Int_t, Int_t, Int_t, Int_t)
Definition: fa125fns.h:55
jerror_t brun(jana::JEventLoop *eventLoop, int32_t runnumber)
Called everytime a new run number is detected.
jerror_t evnt(jana::JEventLoop *eventLoop, uint64_t eventnumber)
Called every event.
TCanvas * cp
int straw
Definition: DCDCHit.h:19
int main(int argc, char *argv[])
Definition: gendoc.cc:6