1 |
#include <CaloEnergy.h> |
2 |
#include <PamLevel2.h> |
3 |
|
4 |
//-------------------------------------- |
5 |
/** |
6 |
* Default constructor |
7 |
*/ |
8 |
CaloEnergy::CaloEnergy(){ |
9 |
Clear(); |
10 |
} |
11 |
|
12 |
CaloEnergy::CaloEnergy(PamLevel2 *l2p){ |
13 |
// |
14 |
Clear(); |
15 |
// |
16 |
L2 = l2p; |
17 |
// |
18 |
if ( !L2->IsORB() ) printf(" WARNING: OrbitalInfo Tree is needed, the plugin could not work properly without it \n"); |
19 |
// |
20 |
fSimu = false; |
21 |
this->Set(); |
22 |
// |
23 |
} |
24 |
|
25 |
CaloEnergy::CaloEnergy(PamLevel2 *l2p, Bool_t simulation){ |
26 |
// |
27 |
Clear(); |
28 |
// |
29 |
L2 = l2p; |
30 |
// |
31 |
if ( !L2->IsORB() ) printf(" WARNING: OrbitalInfo Tree is needed, the plugin could not work properly without it \n"); |
32 |
// |
33 |
fSimu = simulation; |
34 |
this->Set(); |
35 |
// |
36 |
} |
37 |
|
38 |
void CaloEnergy::Delete(){ |
39 |
Clear(); |
40 |
delete this; |
41 |
} |
42 |
|
43 |
void CaloEnergy::UseLevel2(){ |
44 |
if ( clong ){ |
45 |
delete clong; |
46 |
clong = new CaloLong(L2); |
47 |
clong->SplitInto(0,22); |
48 |
// clong->HeavyTail(true); |
49 |
clong->HeavyTail(false); |
50 |
}; |
51 |
if ( cp ) delete cp; |
52 |
cp = NULL; |
53 |
} |
54 |
|
55 |
void CaloEnergy::UseCaloPreSampler(){ |
56 |
// |
57 |
// use the presampler setting forcefitmode to 1000 means to force the DV routine to find the track inside the calorimeter using the "shower" approach developed for electrons |
58 |
// |
59 |
if ( !cp ) cp = new CaloPreSampler(L2); |
60 |
cp->SplitInto(0,22); |
61 |
cp->SetForceFitMode(1000); |
62 |
// cp->UseTracker(false); |
63 |
// cp->ForceCaloFit(); |
64 |
// cp->SetDebug(true); |
65 |
// cp->Process(); |
66 |
if ( clong ) clong->SetCaloLevel2Pointer(cp->GetCaloLevel2Pointer()); |
67 |
} |
68 |
|
69 |
|
70 |
void CaloEnergy::UseLongFit(){ |
71 |
fPl = 0; |
72 |
fLong = true; |
73 |
if ( !clong ){ |
74 |
clong = new CaloLong(L2); |
75 |
if ( cp ) clong->SetCaloLevel2Pointer(cp->GetCaloLevel2Pointer()); |
76 |
clong->SplitInto(0,22); |
77 |
// clong->HeavyTail(true); |
78 |
clong->HeavyTail(false); |
79 |
}; |
80 |
// |
81 |
} |
82 |
|
83 |
void CaloEnergy::Set(){ |
84 |
// |
85 |
// set default values, NB default conversion factor for energy is just very approximated! |
86 |
// |
87 |
OBT = 0; |
88 |
PKT = 0; |
89 |
atime = 0; |
90 |
sntr = "start"; |
91 |
// |
92 |
AOBT = 0; |
93 |
APKT = 0; |
94 |
aatime = 0; |
95 |
asntr = "start"; |
96 |
// |
97 |
debug = false; |
98 |
usepl18x = false; |
99 |
// |
100 |
indep = false; |
101 |
// |
102 |
fAllpl = true; |
103 |
fLong = false; |
104 |
fPl = 1; |
105 |
fRad = -1; |
106 |
cp = NULL; |
107 |
clong = NULL; |
108 |
x0max = -1.; |
109 |
// |
110 |
multicol = false; |
111 |
// |
112 |
this->DefineGeometry(); |
113 |
fXosel =true; |
114 |
fXesel = true; |
115 |
fYosel = true; |
116 |
fYesel = true; |
117 |
fConv_rxe = 44.4; |
118 |
fConv_rxo = 44.4; |
119 |
fConv_ryo = 44.4; |
120 |
fConv_rye = 44.4; |
121 |
fXomin = 1000; |
122 |
fXemin = 1000; |
123 |
fYomin = 1000; |
124 |
fYemin = 1000; |
125 |
// |
126 |
this->UseCaloPreSampler(); // use it by default, to go back to "standard" mode use CaloEnergy::UseLevel2(). |
127 |
// |
128 |
} |
129 |
void CaloEnergy::DefineGeometry(){ |
130 |
// |
131 |
// Use CaloStrip to determine once the position of border strips for each section |
132 |
// |
133 |
// fM = 2. + 0.096; // real position from cbar BUG the 0.096 is already taken into account in the border calculation made by Giovanna |
134 |
fM = 2. ; // real position from cbar |
135 |
// fM1 = 2. - 0.122 - 0.096; // due to calculation of xe1 etc. BUG! this way we count from the silicon border not from the silicon sensitive area |
136 |
fM1 = 2. - 0.122 - 0.096 + 0.096; // due to calculation of xe1 etc. |
137 |
if ( fM1 < 0. ) fM1 = 0.; |
138 |
// |
139 |
CaloStrip *cs = new CaloStrip(fSimu); |
140 |
// |
141 |
// view y plane 0 strip 0 |
142 |
cs->Set(1,0,0); |
143 |
xe1= cs->GetY(); |
144 |
// view y plane 0 strip 31 |
145 |
cs->Set(1,0,31); |
146 |
xe2= cs->GetY(); |
147 |
// view y plane 0 strip 32 |
148 |
cs->Set(1,0,32); |
149 |
xe3= cs->GetY(); |
150 |
// view y plane 0 strip 63 |
151 |
cs->Set(1,0,63); |
152 |
xe4= cs->GetY(); |
153 |
// view y plane 0 strip 64 |
154 |
cs->Set(1,0,64); |
155 |
xe5= cs->GetY(); |
156 |
// view y plane 0 strip 95 |
157 |
cs->Set(1,0,95); |
158 |
xe6= cs->GetY(); |
159 |
// view x plane 0 strip 0 |
160 |
cs->Set(0,0,0); |
161 |
yo1= cs->GetX(); |
162 |
// view x plane 0 strip 31 |
163 |
cs->Set(0,0,31); |
164 |
yo2= cs->GetX(); |
165 |
// view x plane 0 strip 32 |
166 |
cs->Set(0,0,32); |
167 |
yo3= cs->GetX(); |
168 |
// view x plane 0 strip 63 |
169 |
cs->Set(0,0,63); |
170 |
yo4= cs->GetX(); |
171 |
// view x plane 0 strip 64 |
172 |
cs->Set(0,0,64); |
173 |
yo5= cs->GetX(); |
174 |
// view x plane 0 strip 95 |
175 |
cs->Set(0,0,95); |
176 |
yo6= cs->GetX(); |
177 |
// view y plane 1 strip 0 |
178 |
cs->Set(1,1,0); |
179 |
xo1= cs->GetY(); |
180 |
// view y plane 1 strip 31 |
181 |
cs->Set(1,1,31); |
182 |
xo2= cs->GetY(); |
183 |
// view y plane 1 strip 32 |
184 |
cs->Set(1,1,32); |
185 |
xo3= cs->GetY(); |
186 |
// view y plane 1 strip 63 |
187 |
cs->Set(1,1,63); |
188 |
xo4= cs->GetY(); |
189 |
// view y plane 1 strip 64 |
190 |
cs->Set(1,1,64); |
191 |
xo5= cs->GetY(); |
192 |
// view y plane 1 strip 95 |
193 |
cs->Set(1,1,95); |
194 |
xo6= cs->GetY(); |
195 |
// view x plane 1 strip 0 |
196 |
cs->Set(0,1,0); |
197 |
ye1= cs->GetX(); |
198 |
// view x plane 1 strip 31 |
199 |
cs->Set(0,1,31); |
200 |
ye2= cs->GetX(); |
201 |
// view x plane 1 strip 32 |
202 |
cs->Set(0,1,32); |
203 |
ye3= cs->GetX(); |
204 |
// view x plane 1 strip 63 |
205 |
cs->Set(0,1,63); |
206 |
ye4= cs->GetX(); |
207 |
// view x plane 1 strip 64 |
208 |
cs->Set(0,1,64); |
209 |
ye5= cs->GetX(); |
210 |
// view x plane 1 strip 95 |
211 |
cs->Set(0,1,95); |
212 |
ye6= cs->GetX(); |
213 |
// |
214 |
for (Int_t p = 0; p<22; p ++){ |
215 |
for (Int_t v = 0; v<2; v++ ){ |
216 |
cs->Set(v,p,0); |
217 |
trk_z[p][v]= cs->GetZ(); // Z coord for each plane |
218 |
}; |
219 |
}; |
220 |
// |
221 |
delete cs; |
222 |
// |
223 |
} |
224 |
|
225 |
void CaloEnergy::Clear(){ |
226 |
// |
227 |
// clear variables |
228 |
// |
229 |
fPartsel = false; |
230 |
fSel = false; |
231 |
fXosel = false; |
232 |
fXesel = false; |
233 |
fYosel = false; |
234 |
fYesel = false; |
235 |
fCount = 0.; |
236 |
fEnergy = 0.; |
237 |
fEnergyxe = 0.; |
238 |
fEnergyxo = 0.; |
239 |
fEnergyye = 0.; |
240 |
fEnergyyo = 0.; |
241 |
fMax_plane = 0; |
242 |
fMax_planexo = 0; |
243 |
fMax_planexe = 0; |
244 |
fMax_planeyo = 0; |
245 |
fMax_planeye = 0; |
246 |
xomax_en= 0.; |
247 |
xemax_en= 0.; |
248 |
yomax_en= 0.; |
249 |
yemax_en= 0.; |
250 |
// |
251 |
memset(enstrip,0,2*22*96*(sizeof(Float_t))); |
252 |
en = 0.; |
253 |
view = 0; |
254 |
plane = 0; |
255 |
strip = 0; |
256 |
energyxe = 0.; |
257 |
energyyo = 0.; |
258 |
energyxo = 0.; |
259 |
energyye = 0.; |
260 |
fYoout = 0; |
261 |
fYeout = 0; |
262 |
fXoout = 0; |
263 |
fXeout = 0; |
264 |
fXEen_maxplane = 0.; |
265 |
fXOen_maxplane = 0.; |
266 |
fYEen_maxplane = 0.; |
267 |
fYOen_maxplane = 0.; |
268 |
memset(en_xep,0,11*sizeof(Float_t)); |
269 |
memset(en_yep,0,11*sizeof(Float_t)); |
270 |
memset(en_xop,0,11*sizeof(Float_t)); |
271 |
memset(en_yop,0,11*sizeof(Float_t)); |
272 |
// |
273 |
fColumn = -1; |
274 |
fColXE = -1; |
275 |
fColXO = -1; |
276 |
fColYE = -1; |
277 |
fColYO = -1; |
278 |
memset(encol,0,2*3*sizeof(Float_t)); |
279 |
entot[0] = 0.; |
280 |
entot[1] = 0.; |
281 |
// |
282 |
X0pl = 0.76; |
283 |
// |
284 |
} |
285 |
|
286 |
void CaloEnergy::Print(){ |
287 |
// |
288 |
printf("========================================================================\n"); |
289 |
printf(" OBT: %u PKT: %u ATIME: %u \n",OBT,PKT,atime); |
290 |
printf(" fEnergy :.............. %f \n",fEnergy); |
291 |
printf(" fMax_plane :........... %f \n",fMax_plane); |
292 |
printf(" fMax_planexo :......... %i \n",fMax_planexo); |
293 |
printf(" fMax_planexe :......... %i \n",fMax_planexe); |
294 |
printf(" fMax_planeyo :......... %i \n",fMax_planeyo); |
295 |
printf(" fMax_planeye :......... %i \n",fMax_planeye); |
296 |
printf(" fCount :.............. %f \n",fCount); |
297 |
printf(" fSel :.............. %i \n",fSel); |
298 |
printf(" fPartsel:.............. %i \n",fPartsel); |
299 |
printf(" fXesel :.............. %i \n",fXesel); |
300 |
printf(" fXosel :.............. %i \n",fXosel); |
301 |
printf(" fYesel :.............. %i \n",fYesel); |
302 |
printf(" fYosel :.............. %i \n",fYosel); |
303 |
printf(" fXemin :.............. %i \n",fXemin); |
304 |
printf(" fXomin :.............. %i \n",fXomin); |
305 |
printf(" fYemin :.............. %i \n",fYemin); |
306 |
printf(" fYomin :.............. %i \n",fYomin); |
307 |
printf(" fXeout :.............. %i \n",fXeout); |
308 |
printf(" fXoout :.............. %i \n",fXoout); |
309 |
printf(" fYeout :.............. %i \n",fYeout); |
310 |
printf(" fYoout :.............. %i \n",fYoout); |
311 |
printf(" fSimu :.............. %i \n",fSimu); |
312 |
printf(" fM :.............. %f \n",fM); |
313 |
printf(" fM1 :.............. %f \n",fM1); |
314 |
printf(" fRad :.............. %i \n",fRad); |
315 |
printf(" fPl :.............. %i \n",fPl); |
316 |
printf(" fColumn :.............. %i \n",fColumn); |
317 |
printf(" multicol:.............. %i \n",multicol); |
318 |
printf(" encol x :.............. %f \n",this->GetEncol(0)); |
319 |
printf(" encol y :.............. %f \n",this->GetEncol(1)); |
320 |
printf(" entot x :.............. %f \n",this->GetEntot(0)); |
321 |
printf(" entot y :.............. %f \n",this->GetEntot(1)); |
322 |
printf(" fColXE :.............. %i \n",fColXE); |
323 |
printf(" fColXO :.............. %i \n",fColXO); |
324 |
printf(" fColYE :.............. %i \n",fColYE); |
325 |
printf(" fColYO :.............. %i \n",fColYO); |
326 |
printf(" fConv_rxe ............. %f \n",fConv_rxe); |
327 |
printf(" fConv_rxo ............. %f \n",fConv_rxo); |
328 |
printf(" fConv_ryo ............. %f \n",fConv_ryo); |
329 |
printf(" fConv_rye ............. %f \n",fConv_rye); |
330 |
printf(" fLong :.............. %i \n",fLong); |
331 |
printf(" energyxe:.............. %f \n",energyxe); |
332 |
printf(" energyxo:.............. %f \n",energyxo); |
333 |
printf(" energyye:.............. %f \n",energyye); |
334 |
printf(" energyyo:.............. %f \n",energyyo); |
335 |
printf(" fXEen_maxplane:........ %f \n",fXEen_maxplane); |
336 |
printf(" fXOen_maxplane:........ %f \n",fXOen_maxplane); |
337 |
printf(" fYEen_maxplane:........ %f \n",fYEen_maxplane); |
338 |
printf(" fYOen_maxplane:........ %f \n",fYOen_maxplane); |
339 |
printf(" x0max :.............. %f \n",x0max); |
340 |
printf(" X0pl :.............. %f \n",X0pl); |
341 |
printf(" debug :.............. %i \n",debug); |
342 |
|
343 |
printf("========================================================================\n"); |
344 |
// |
345 |
} |
346 |
|
347 |
void CaloEnergy::SetMinimumContainment(TString section, Int_t plane){ |
348 |
section.ToUpper(); |
349 |
if ( section.Contains("XO") ) fXomin = plane; |
350 |
if ( section.Contains("XE") ) fXemin = plane; |
351 |
if ( section.Contains("YO") ) fYomin = plane; |
352 |
if ( section.Contains("YE") ) fYemin = plane; |
353 |
} |
354 |
|
355 |
void CaloEnergy::SetMinimumContainment(Int_t plane){ |
356 |
this->SetMinimumContainment("XEXOYEYO",plane); |
357 |
} |
358 |
|
359 |
void CaloEnergy::SetConversionFactor(TString section, Float_t conv){ |
360 |
section.ToUpper(); |
361 |
if ( section.Contains("XO") ) fConv_rxo = conv; |
362 |
if ( section.Contains("XE") ) fConv_rxe = conv; |
363 |
if ( section.Contains("YO") ) fConv_ryo = conv; |
364 |
if ( section.Contains("YE") ) fConv_rye = conv; |
365 |
} |
366 |
|
367 |
void CaloEnergy::SetConversionFactor(Float_t conv){ |
368 |
this->SetConversionFactor("XEXOYEYO",conv); |
369 |
} |
370 |
|
371 |
Int_t CaloEnergy::GetMinimumContainment(TString section){ |
372 |
section.ToUpper(); |
373 |
if ( section.Contains("XO") ) return(fXomin); |
374 |
if ( section.Contains("XE") ) return(fXemin); |
375 |
if ( section.Contains("YE") ) return(fYemin); |
376 |
if ( section.Contains("YO") ) return(fYomin); |
377 |
printf(" ERROR: section not recognized \n"); |
378 |
return(-1000); |
379 |
} |
380 |
|
381 |
Float_t CaloEnergy::GetConversionFactor(TString section){ |
382 |
section.ToUpper(); |
383 |
if ( section.Contains("XO") ) return(fConv_rxo); |
384 |
if ( section.Contains("XE") ) return(fConv_rxe); |
385 |
if ( section.Contains("YO") ) return(fConv_ryo); |
386 |
if ( section.Contains("YE") ) return(fConv_rye); |
387 |
printf(" ERROR: section not recognized \n"); |
388 |
return(-1000.); |
389 |
} |
390 |
|
391 |
Int_t CaloEnergy::GetMaxplane(TString section){ |
392 |
section.ToUpper(); |
393 |
if ( section.Contains("XO") ) return fMax_planexo; |
394 |
if ( section.Contains("XE") ) return fMax_planexe; |
395 |
if ( section.Contains("YO") ) return fMax_planeyo; |
396 |
if ( section.Contains("YE") ) return fMax_planeye; |
397 |
return(-1); |
398 |
} |
399 |
|
400 |
Int_t CaloEnergy::GetColumn(TString section){ |
401 |
section.ToUpper(); |
402 |
if ( section.Contains("XO") ) return fColXO; |
403 |
if ( section.Contains("XE") ) return fColXE; |
404 |
if ( section.Contains("YO") ) return fColYO; |
405 |
if ( section.Contains("YE") ) return fColYE; |
406 |
return(-1); |
407 |
} |
408 |
|
409 |
Float_t CaloEnergy::GetMipEnergyAtMaxplane(TString section){ |
410 |
printf(" WARNING: OBSOLETE METHOD, use GetMipEnergyAtMaxplane(TString) instead! \n"); |
411 |
return (this->GetEnergyAtMaxplane(section)); |
412 |
} |
413 |
|
414 |
Float_t CaloEnergy::GetEnergyAtMaxplane(TString section){ |
415 |
section.ToUpper(); |
416 |
if ( section.Contains("XO") ) return xomax_en; |
417 |
if ( section.Contains("XE") ) return xemax_en; |
418 |
if ( section.Contains("YO") ) return yomax_en; |
419 |
if ( section.Contains("YE") ) return yemax_en; |
420 |
return(-1); |
421 |
} |
422 |
|
423 |
Float_t CaloEnergy::GetMaxEnergy(TString section){ |
424 |
printf(" WARNING: OBSOLETE METHOD, use GetMipEnergy(TString) instead! \n"); |
425 |
return (this->GetMipEnergy(section)); |
426 |
} |
427 |
|
428 |
Float_t CaloEnergy::GetMipEnergy(TString section){ |
429 |
section.ToUpper(); |
430 |
if ( fLong ){ |
431 |
this->Process(section); |
432 |
return fXOen_maxplane; |
433 |
} else { |
434 |
if ( section.Contains("XO") ) return fXOen_maxplane; |
435 |
if ( section.Contains("XE") ) return fXEen_maxplane; |
436 |
if ( section.Contains("YO") ) return fYOen_maxplane; |
437 |
if ( section.Contains("YE") ) return fYEen_maxplane; |
438 |
}; |
439 |
return(-1); |
440 |
} |
441 |
|
442 |
Float_t CaloEnergy::GetEncol(Int_t i){ |
443 |
if ( fColumn > -1 && (((fXesel || fXosel)&&i==1) || ((fYesel || fYosel)&&i==0)) ){ |
444 |
Int_t t = -1; |
445 |
if ( i == 0 ){ |
446 |
if ( fColumn == 0 || fColumn == 3 || fColumn == 6 ) t = 0; |
447 |
if ( fColumn == 1 || fColumn == 4 || fColumn == 7 ) t = 1; |
448 |
if ( fColumn == 2 || fColumn == 5 || fColumn == 8 ) t = 2; |
449 |
} else { |
450 |
if ( fColumn == 0 || fColumn == 1 || fColumn == 2 ) t = 0; |
451 |
if ( fColumn == 3 || fColumn == 4 || fColumn == 5 ) t = 1; |
452 |
if ( fColumn == 6 || fColumn == 7 || fColumn == 8 ) t = 2; |
453 |
}; |
454 |
if ( debug ) printf(" encol: i %i t %i encol %f \n",i,t,encol[i][t]); |
455 |
return encol[i][t]; |
456 |
}; |
457 |
return(-1.); |
458 |
} |
459 |
|
460 |
Float_t CaloEnergy::GetMaxEnergy(){ |
461 |
printf(" WARNING: OBSOLETE METHOD, use GetMipEnergy() instead! \n"); |
462 |
return (this->GetMipEnergy()); |
463 |
} |
464 |
|
465 |
Float_t CaloEnergy::GetMipEnergy(){ |
466 |
if ( fLong ){ |
467 |
if ( debug ) printf(" oh! call process! with asntr %s and sntr %s \n",asntr.Data(),sntr.Data()); |
468 |
this->Process(asntr); |
469 |
}; |
470 |
return((fXEen_maxplane+fYOen_maxplane+fYEen_maxplane+fXOen_maxplane)); |
471 |
} |
472 |
|
473 |
|
474 |
Bool_t CaloEnergy::IsInsideAcceptance(TString section){ |
475 |
// |
476 |
// check if the event is inside the acceptance of the given section(s) |
477 |
// |
478 |
TString ntr = section; |
479 |
if ( !L2 ){ |
480 |
printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n"); |
481 |
printf(" ERROR: CaloEnergy variables not filled \n"); |
482 |
return false; |
483 |
}; |
484 |
// |
485 |
Bool_t newentry = false; |
486 |
// |
487 |
if ( L2->IsORB() ){ |
488 |
if ( L2->GetOrbitalInfo()->pkt_num != APKT || L2->GetOrbitalInfo()->OBT != AOBT || L2->GetOrbitalInfo()->absTime != aatime || strcmp(ntr.Data(),asntr.Data()) ){ |
489 |
newentry = true; |
490 |
AOBT = L2->GetOrbitalInfo()->OBT; |
491 |
APKT = L2->GetOrbitalInfo()->pkt_num; |
492 |
aatime = L2->GetOrbitalInfo()->absTime; |
493 |
asntr = ntr; |
494 |
}; |
495 |
} else { |
496 |
newentry = true; |
497 |
}; |
498 |
// |
499 |
// if we have already called this method for this event and no input changed then return fSel and exit |
500 |
// |
501 |
if ( !newentry ) return fSel; |
502 |
// |
503 |
// process the event |
504 |
// |
505 |
if ( debug ) printf(" ########## IsInsideAcceptance ######### \n"); |
506 |
// |
507 |
// clear variables |
508 |
// |
509 |
this->Clear(); |
510 |
// |
511 |
section.ToUpper(); |
512 |
// |
513 |
// Count the number of section(s) given as input |
514 |
// |
515 |
Int_t fNumSec = Int_t(section.Contains("XO"))+Int_t(section.Contains("XE"))+Int_t(section.Contains("YO"))+Int_t(section.Contains("YE")); |
516 |
if ( !fNumSec ){ |
517 |
printf(" ERROR: section must be XO or XE or YO or YE while it is %s \n",section.Data()); |
518 |
return false; |
519 |
}; |
520 |
// |
521 |
// If the presampler object exists then use the presampler output instead of the level2 output |
522 |
// |
523 |
CaloLevel2 *cl2 = NULL; |
524 |
if ( cp ){ |
525 |
cl2 = cp->GetCaloLevel2(); |
526 |
} else { |
527 |
cl2 = L2->GetCaloLevel2(); |
528 |
}; |
529 |
// |
530 |
// get the energy for every strip of the calorimeter |
531 |
// |
532 |
for (Int_t ch=0; ch< L2->GetCaloLevel1()->istrip; ch++){ |
533 |
en = L2->GetCaloLevel1()->DecodeEstrip(ch,view,plane,strip); |
534 |
// |
535 |
if ( !usepl18x && view==0 && plane==18 ) en = 0.; |
536 |
// |
537 |
enstrip[view][plane][strip]=en; |
538 |
}; |
539 |
// |
540 |
if ( debug && ((fM1+0.122-0.244*(Float_t)fRad) < 0.) ) printf("Error: (fM1+0.122-0.244*(Float_t)fRad) < 0. fM1 %f fRad %i %f \n",fM1,fRad,(fM1+0.122-0.244*(Float_t)fRad)); |
541 |
// |
542 |
// inclination factor (stolen from Daniele's code) |
543 |
// |
544 |
Float_t ytgx = 0; |
545 |
Float_t ytgy = 0; |
546 |
ytgx = 0.76 * cl2->tanx[0]; |
547 |
ytgy = 0.76 * cl2->tany[0]; |
548 |
X0pl = sqrt( pow(0.76,2.) + pow(ytgx,2.) + pow(ytgy,2.) ); |
549 |
// |
550 |
// sum energy plane by plane for each sections |
551 |
// |
552 |
Float_t fen_xep[11]; |
553 |
Float_t fen_xop[11]; |
554 |
Float_t fen_yep[11]; |
555 |
Float_t fen_yop[11]; |
556 |
memset(fen_xep,0,11*sizeof(Float_t)); |
557 |
memset(fen_xop,0,11*sizeof(Float_t)); |
558 |
memset(fen_yep,0,11*sizeof(Float_t)); |
559 |
memset(fen_yop,0,11*sizeof(Float_t)); |
560 |
// |
561 |
for (Int_t i=0;i<11;i++){ |
562 |
for(strip=0; strip<96; strip++) { |
563 |
fen_xep[i] += enstrip[1][2*i][strip]; |
564 |
fen_yop[i] += enstrip[0][2*i][strip]; |
565 |
fen_xop[i] += enstrip[1][(2*i)+1][strip]; |
566 |
fen_yep[i] += enstrip[0][(2*i)+1][strip]; |
567 |
if ( fRad < 0 ){ |
568 |
// |
569 |
// run over all the strips of the plane |
570 |
// |
571 |
en_xep[i] += enstrip[1][2*i][strip]; |
572 |
en_yop[i] += enstrip[0][2*i][strip]; |
573 |
en_xop[i] += enstrip[1][(2*i)+1][strip]; |
574 |
en_yep[i] += enstrip[0][(2*i)+1][strip]; |
575 |
} else { |
576 |
// |
577 |
// use only the strips inside a cylinder of given radius fRad |
578 |
// |
579 |
if ( cl2->cibar[2*i][1] >= 1 && cl2->cibar[2*i][1] <= 96 && |
580 |
(strip >= (cl2->cibar[2*i][1]-1-fRad)) && (strip <= (cl2->cibar[2*i][1]-1+fRad)) ) en_xep[i] += enstrip[1][2*i][strip]; |
581 |
|
582 |
if ( cl2->cibar[2*i][0] >= 1 && cl2->cibar[2*i][0] <= 96 && |
583 |
(strip >= (cl2->cibar[2*i][0]-1-fRad)) && (strip <= (cl2->cibar[2*i][0]-1+fRad)) ) en_yop[i] += enstrip[0][2*i][strip]; |
584 |
|
585 |
if ( cl2->cibar[(2*i)+1][1] >= 1 && cl2->cibar[(2*i)+1][1] <= 96 && |
586 |
(strip >= (cl2->cibar[(2*i)+1][1]-1-fRad)) && (strip <= (cl2->cibar[(2*i)+1][1]-1+fRad)) ) en_xop[i] += enstrip[1][(2*i)+1][strip]; |
587 |
|
588 |
if ( cl2->cibar[(2*i)+1][0] >= 1 && cl2->cibar[(2*i)+1][0] <= 96 && |
589 |
(strip >= (cl2->cibar[(2*i)+1][0]-1-fRad)) && (strip <= (cl2->cibar[(2*i)+1][0]-1+fRad)) ) en_yep[i] += enstrip[0][(2*i)+1][strip]; |
590 |
}; |
591 |
}; |
592 |
if ( debug ) printf(" ex_xep[%i] %f cibar %i \n",i,en_xep[i],cl2->cibar[2*i][1]); |
593 |
if ( debug ) printf(" ex_xop[%i] %f cibar %i \n",i,en_xop[i],cl2->cibar[(2*i)+1][1]); |
594 |
if ( debug ) printf(" ex_yep[%i] %f cibar %i \n",i,en_yep[i],cl2->cibar[(2*i)+1][0]); |
595 |
if ( debug ) printf(" ex_yop[%i] %f cibar %i \n",i,en_yop[i],cl2->cibar[2*i][0]); |
596 |
energyxe += en_xep[i]; |
597 |
energyyo += en_yop[i]; |
598 |
energyxo += en_xop[i]; |
599 |
energyye += en_yep[i]; |
600 |
}; |
601 |
// |
602 |
// Find the plane of maximum for each section |
603 |
// |
604 |
// |
605 |
Int_t xen = 0; |
606 |
Int_t yon = 0; |
607 |
Int_t xon = 0; |
608 |
Int_t yen = 0; |
609 |
Float_t en = 0.; |
610 |
// |
611 |
if ( section.Contains("XE") ){ |
612 |
yon++; |
613 |
xon++; |
614 |
yen++; |
615 |
for (Int_t ipl =0; ipl < 11; ipl ++) { |
616 |
en = fen_xep[ipl]; |
617 |
if ( !fAllpl ) en = en_xep[ipl]; |
618 |
if(en > xemax_en) { |
619 |
xemax_en = en; |
620 |
fMax_planexe = ipl; |
621 |
}; |
622 |
}; |
623 |
}; |
624 |
// |
625 |
if ( section.Contains("YO") ){ |
626 |
xon++; |
627 |
yen++; |
628 |
for (Int_t ipl =0; ipl < 11; ipl ++) { |
629 |
en = fen_yop[ipl]; |
630 |
if ( !fAllpl ) en = en_yop[ipl]; |
631 |
if(en > yomax_en) { |
632 |
yomax_en = en; |
633 |
fMax_planeyo = ipl; |
634 |
}; |
635 |
}; |
636 |
}; |
637 |
// |
638 |
if ( section.Contains("XO") ){ |
639 |
yen++; |
640 |
for (Int_t ipl =0; ipl < 11; ipl ++) { |
641 |
en = fen_xop[ipl]; |
642 |
if ( !fAllpl ) en = en_xop[ipl]; |
643 |
if(en > xomax_en) { |
644 |
xomax_en = en; |
645 |
fMax_planexo = ipl; |
646 |
}; |
647 |
}; |
648 |
}; |
649 |
// |
650 |
if ( section.Contains("YE") ){ |
651 |
for (Int_t ipl =0; ipl < 11; ipl ++) { |
652 |
en = fen_yep[ipl]; |
653 |
if ( !fAllpl ) en = en_yep[ipl]; |
654 |
if(en > yemax_en) { |
655 |
yemax_en = en; |
656 |
fMax_planeye = ipl; |
657 |
}; |
658 |
}; |
659 |
}; |
660 |
// |
661 |
// the maximum is given externally as X0, convert it to plane and section |
662 |
// |
663 |
if ( x0max > 0. ){ |
664 |
if ( debug ) printf(" CALCULATE MAX PLANE GIVEN X0 ASSUMING PERPENDICULAR TRACK \n"); |
665 |
// Int_t wpl = (Int_t)roundf(x0max/0.76); |
666 |
Int_t wpl = (Int_t)roundf(x0max/X0pl); |
667 |
Bool_t isY = false; |
668 |
// if ( ((x0max/0.76)-(Float_t)wpl) > 0. ) isY = true; |
669 |
if ( ((x0max/X0pl)-(Float_t)wpl) > 0. ) isY = true; |
670 |
xomax_en = 0.; |
671 |
yemax_en = 0.; |
672 |
xemax_en = 0.; |
673 |
yomax_en = 0.; |
674 |
// |
675 |
if ( !(wpl%2) ){ |
676 |
// 0, 2, 4, ... |
677 |
if ( isY ){ |
678 |
if ( section.Contains("YO") ) yomax_en = 1000.; |
679 |
if ( section.Contains("XE") ) xemax_en = 500.; |
680 |
fMax_planeyo=wpl/2; |
681 |
fMax_planexe=wpl/2; |
682 |
if ( section.Contains("XO") ) xomax_en = 10.; |
683 |
if ( section.Contains("YE") ) yemax_en = 5.; |
684 |
} else { |
685 |
if ( section.Contains("YO") ) yomax_en = 500.; |
686 |
if ( section.Contains("XE") ) xemax_en = 1000.; |
687 |
fMax_planeyo=wpl/2; |
688 |
fMax_planexe=wpl/2; |
689 |
if ( section.Contains("XO") ) xomax_en = 5.; |
690 |
if ( section.Contains("YE") ) yemax_en = 10.; |
691 |
}; |
692 |
} else { |
693 |
// 1, 3, 5, ... |
694 |
if ( isY ){ |
695 |
if ( section.Contains("YE") ) yemax_en = 1000.; |
696 |
if ( section.Contains("XO") ) xomax_en = 500.; |
697 |
fMax_planeye=(wpl-1)/2; |
698 |
fMax_planexo=(wpl-1)/2; |
699 |
if ( section.Contains("XE") ) xemax_en = 10.; |
700 |
if ( section.Contains("YO") ) yomax_en = 5.; |
701 |
} else { |
702 |
if ( section.Contains("YE") ) yemax_en = 500.; |
703 |
if ( section.Contains("XO") ) xomax_en = 1000.; |
704 |
fMax_planeye=(wpl-1)/2; |
705 |
fMax_planexo=(wpl-1)/2; |
706 |
if ( section.Contains("XE") ) xemax_en = 5.; |
707 |
if ( section.Contains("YO") ) yomax_en = 10.; |
708 |
}; |
709 |
}; |
710 |
// if ( debug ) printf(" x0max %f x0max/0.76 %f wpl %i isY %i yomax_en %f xemax_en %f yemax_en %f xomax_en %f fMaxplane %i %i %i %i\n",x0max,(x0max/0.76),wpl,isY,yomax_en,xemax_en,yemax_en,xomax_en,fMax_planeyo,fMax_planexe,fMax_planeye,fMax_planexo); |
711 |
if ( debug ) printf(" x0max %f x0max/X0pl %f X0pl %f wpl %i isY %i yomax_en %f xemax_en %f yemax_en %f xomax_en %f fMaxplane %i %i %i %i\n",x0max,(x0max/X0pl),X0pl,wpl,isY,yomax_en,xemax_en,yemax_en,xomax_en,fMax_planeyo,fMax_planexe,fMax_planeye,fMax_planexo); |
712 |
}; |
713 |
// |
714 |
Int_t nPl = fPl; |
715 |
// |
716 |
// Set the maximum in case of coherent mode was selected |
717 |
// |
718 |
if ( !indep ){ |
719 |
nPl = 0; |
720 |
if ( debug ) printf(" A: Check maximum, coherent mode: xoen %f yoen %f xeen %f yeen %f xomax %i yomax %i xemax %i yemax %i\n",xomax_en,yomax_en,xemax_en,yemax_en,fMax_planexo,fMax_planeyo,fMax_planexe,fMax_planeye); |
721 |
Int_t nummod = 0; |
722 |
Int_t numexpl = 0; |
723 |
if ( xomax_en > xemax_en && xomax_en > yemax_en && xomax_en > yomax_en ){ |
724 |
// |
725 |
// Section XO contains the maximum energy release per plane of the whole calorimeter |
726 |
// |
727 |
if ( debug ) printf(" XO is MAX %i %i %i %i\n",xen,yon,xon,yen); |
728 |
// |
729 |
// fMax_plane is the plane of maximum + number of additional dE/dx measurement counting planes from 0 to 43 |
730 |
// |
731 |
fMax_plane = (fNumSec * fMax_planexo) +(Float_t)xon + fPl; |
732 |
// |
733 |
// nummod is the integer part of the number of modules in which the maximum is contained |
734 |
// |
735 |
nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec); |
736 |
// |
737 |
// numexpl is the number of additional planes (0,1,2) inside the module |
738 |
// |
739 |
numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod); |
740 |
// |
741 |
}; |
742 |
if ( xemax_en > xomax_en && xemax_en > yemax_en && xemax_en > yomax_en ){ |
743 |
if ( debug ) printf(" XE is MAX %i %i %i %i\n",xen,yon,xon,yen); |
744 |
fMax_plane = (fNumSec * fMax_planexe) +(Float_t)xen + fPl; |
745 |
nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec); |
746 |
numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod); |
747 |
// |
748 |
}; |
749 |
|
750 |
if ( yemax_en > xomax_en && yemax_en > xemax_en && yemax_en > yomax_en ){ |
751 |
if ( debug ) printf(" YE is MAX %i %i %i %i\n",xen,yon,xon,yen); |
752 |
fMax_plane = (fNumSec * fMax_planeye) +(Float_t)yen + fPl; |
753 |
nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec); |
754 |
numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod); |
755 |
// |
756 |
}; |
757 |
if ( yomax_en > xemax_en && yomax_en > yemax_en && yomax_en > xomax_en ){ |
758 |
if ( debug ) printf(" YO is MAX %i %i %i %i\n",xen,yon,xon,yen); |
759 |
fMax_plane = (fNumSec * fMax_planeyo) +(Float_t)yon + fPl; |
760 |
nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec); |
761 |
numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod); |
762 |
// |
763 |
}; |
764 |
// |
765 |
// find the plane up to which is necessary to integrate the energy for each section |
766 |
// |
767 |
Int_t a = 0; |
768 |
Int_t b = 0; |
769 |
Int_t c = 0; |
770 |
if ( numexpl > xen ) a = 1; |
771 |
if ( numexpl > yon ) b = 1; |
772 |
if ( numexpl > xon ) c = 1; |
773 |
fMax_planexe = nummod; |
774 |
fMax_planeyo = nummod - 1 + a; |
775 |
fMax_planexo = nummod - 1 + b; |
776 |
fMax_planeye = nummod - 1 + c; |
777 |
if ( debug ) printf(" fMax_plane %f nummod %i numexpl %i a %i b %i c %i \n",fMax_plane,nummod,numexpl,a,b,c); |
778 |
if ( debug ) printf(" DONE: Check maximum, coherent mode: xoen %f yoen %f xeen %f yeen %f xomax %i yomax %i xemax %i yemax %i\n",xomax_en,yomax_en,xemax_en,yemax_en,fMax_planexo,fMax_planeyo,fMax_planexe,fMax_planeye); |
779 |
}; |
780 |
// |
781 |
// for each plane of the calorimeter find the position of the track in the direction along the strip (where we do not have a measurement from the selected plane) by looking at the plane above/below of the other view and extrapolating the trajectory to the given plane |
782 |
// |
783 |
// |
784 |
Float_t tgx_cl2; |
785 |
Float_t tgy_cl2; |
786 |
tgx_cl2 = cl2->tanx[0]; |
787 |
tgy_cl2 = cl2->tany[0]; |
788 |
// |
789 |
for (Int_t p=0; p<22; p++){ |
790 |
track_coordy[p][1] = cl2->cbar[p][1]; |
791 |
track_coordx[p][1] = cl2->cbar[p][0] - fabs(trk_z[p][1]-trk_z[p][0])*tgx_cl2; |
792 |
// track_coordx[p][1] = cl2->cbar[p][0] + fabs(trk_z[p][1]-trk_z[p][0])*tgx_cl2; |
793 |
track_coordx[p][0] = cl2->cbar[p][0]; |
794 |
track_coordy[p][0] = cl2->cbar[p][1] - fabs(trk_z[p][1]-trk_z[p][0])*tgy_cl2; |
795 |
// track_coordy[p][0] = cl2->cbar[p][1] + fabs(trk_z[p][1]-trk_z[p][0])*tgy_cl2; |
796 |
if ( debug ) printf(" p %i track_coordy[p][1] %f track_coordx[p][1] %f track_coordx[p][0] %f track_coordy[p][0] %f \n",p,track_coordy[p][1],track_coordx[p][1],track_coordx[p][0],track_coordy[p][0]); |
797 |
}; |
798 |
// |
799 |
if ( debug ) printf(" acceptance fNumSec %i tgx %f tgy %f\n",fNumSec,tgx_cl2,tgy_cl2); |
800 |
// |
801 |
if ( section.Contains("XO") ){ |
802 |
// |
803 |
// find the column hit in the first plane |
804 |
// |
805 |
Int_t ix = -1; |
806 |
Int_t iy = -1; |
807 |
if ( track_coordx[(2*0)+1][1] >= (-12.054+fM) && track_coordx[(2*0)+1][1] <= (-4.246-fM) ) ix = 0; |
808 |
if ( track_coordx[(2*0)+1][1] >= ( -4.004+fM) && track_coordx[(2*0)+1][1] <= ( 3.804-fM) ) ix = 1; |
809 |
if ( track_coordx[(2*0)+1][1] >= ( 4.046+fM) && track_coordx[(2*0)+1][1] <= (11.854-fM) ) ix = 2; |
810 |
if ( cl2->cbar[(2*0)+1][1] >= (xo1 + fM1) && cl2->cbar[(2*0)+1][1] <= (xo2 - fM1) ) iy = 0; |
811 |
if ( cl2->cbar[(2*0)+1][1] >= (xo3 + fM1) && cl2->cbar[(2*0)+1][1] <= (xo4 - fM1) ) iy = 1; |
812 |
if ( cl2->cbar[(2*0)+1][1] >= (xo5 + fM1) && cl2->cbar[(2*0)+1][1] <= (xo6 - fM1) ) iy = 2; |
813 |
if ( ix > -1 && iy > -1 ) fColXO = ix + iy*3; |
814 |
// |
815 |
// check event is inside XO acceptance, if multicol is false (SingleColumn mode) then the track must be contained in a column. |
816 |
// |
817 |
for (Int_t i=0; i<11; i++) { |
818 |
if (( |
819 |
( track_coordx[(2*i)+1][1] >= (-12.054+fM) && track_coordx[(2*i)+1][1] <= (-4.246-fM) && (ix == 0 || multicol) ) || |
820 |
( track_coordx[(2*i)+1][1] >= ( -4.004+fM) && track_coordx[(2*i)+1][1] <= ( 3.804-fM) && (ix == 1 || multicol) ) || |
821 |
( track_coordx[(2*i)+1][1] >= ( 4.046+fM) && track_coordx[(2*i)+1][1] <= (11.854-fM) && (ix == 2 || multicol) ) |
822 |
) && ( |
823 |
( cl2->cbar[(2*i)+1][1] >= (xo1 + fM1) && cl2->cbar[(2*i)+1][1] <= (xo2 - fM1) && (iy == 0 || multicol) ) || |
824 |
( cl2->cbar[(2*i)+1][1] >= (xo3 + fM1) && cl2->cbar[(2*i)+1][1] <= (xo4 - fM1) && (iy == 1 || multicol) ) || |
825 |
( cl2->cbar[(2*i)+1][1] >= (xo5 + fM1) && cl2->cbar[(2*i)+1][1] <= (xo6 - fM1) && (iy == 2 || multicol) ) |
826 |
)){ |
827 |
fXosel = true; |
828 |
fXoout = i; |
829 |
} else { |
830 |
fXosel = false; |
831 |
break; |
832 |
}; |
833 |
}; |
834 |
// |
835 |
// if it goes out of the acceptance BUT the plane up to which we are integrating the energy is contained then the event is "partially" contained |
836 |
// |
837 |
if ( !fXosel && fXoout >= fXomin && fXoout >= (Int_t)(fMax_planexo+nPl) ){ |
838 |
if ( debug ) printf(" XO: this event is only partially contained: fXoout %i fXomin %i fMax_planexo + nPl %i \n",fXoout,fXomin,(Int_t)(fMax_planexo+nPl)); |
839 |
fPartsel = true; |
840 |
fXosel = true; |
841 |
}; |
842 |
// |
843 |
// event is contained (or partially contained) hence we can integrate energy up to the maximum and calculate the energy as measured by this section |
844 |
// |
845 |
if ( fXosel ){ |
846 |
for (Int_t iplm=0; iplm<=TMath::Min(10,(Int_t)(fMax_planexo+nPl)); iplm++){ |
847 |
fXOen_maxplane += en_xop[iplm]; |
848 |
if ( debug ) printf(" XO iplm %i fXOen_maxplane %f en_xop[iplm] %f\n",iplm,fXOen_maxplane,en_xop[iplm]); |
849 |
}; |
850 |
fEnergyxo = fXOen_maxplane/fConv_rxo; |
851 |
// |
852 |
for (Int_t i=0;i<11;i++){ |
853 |
for(strip=0; strip<96; strip++) { |
854 |
// |
855 |
// run over all the strips of the plane |
856 |
// |
857 |
if ( strip >= 0 && strip < 32 ) encol[1][0] += enstrip[1][(2*i)+1][strip]; |
858 |
if ( strip >= 32 && strip < 64 ) encol[1][1] += enstrip[1][(2*i)+1][strip]; |
859 |
if ( strip >= 64 && strip < 96 ) encol[1][2] += enstrip[1][(2*i)+1][strip]; |
860 |
entot[1] += enstrip[1][(2*i)+1][strip]; |
861 |
// |
862 |
}; |
863 |
}; |
864 |
}; |
865 |
}; |
866 |
// |
867 |
if ( section.Contains("XE") ){ |
868 |
// |
869 |
// find the column hit in the first plane |
870 |
// |
871 |
Int_t ix = -1; |
872 |
Int_t iy = -1; |
873 |
if ( track_coordx[(2*0)][1] >= (-11.854+fM) && track_coordx[(2*0)][1] <= (-4.046-fM) ) ix = 0; |
874 |
if ( track_coordx[(2*0)][1] >= ( -3.804+fM) && track_coordx[(2*0)][1] <= ( 4.004-fM) ) ix = 1; |
875 |
if ( track_coordx[(2*0)][1] >= ( 4.246+fM) && track_coordx[(2*0)][1] <= (12.054-fM) ) ix = 2; |
876 |
if ( cl2->cbar[(2*0)][1] >= (xe1 + fM1) && cl2->cbar[(2*0)][1] <= (xe2 - fM1) ) iy = 0; |
877 |
if ( cl2->cbar[(2*0)][1] >= (xe3 + fM1) && cl2->cbar[(2*0)][1] <= (xe4 - fM1) ) iy = 1; |
878 |
if ( cl2->cbar[(2*0)][1] >= (xe5 + fM1) && cl2->cbar[(2*0)][1] <= (xe6 - fM1) ) iy = 2; |
879 |
if ( ix > -1 && iy > -1 ) fColXE = ix + iy*3; |
880 |
// |
881 |
// check event is inside XO acceptance |
882 |
// |
883 |
for (Int_t i=0; i<11; i++) { |
884 |
if (( |
885 |
( track_coordx[(2*i)][1] >= (-11.854+fM) && track_coordx[(2*i)][1] <= (-4.046-fM) && (ix == 0 || multicol) ) || |
886 |
( track_coordx[(2*i)][1] >= ( -3.804+fM) && track_coordx[(2*i)][1] <= ( 4.004-fM) && (ix == 1 || multicol) ) || |
887 |
( track_coordx[(2*i)][1] >= ( 4.246+fM) && track_coordx[(2*i)][1] <= (12.054-fM) && (ix == 2 || multicol) ) |
888 |
) && ( |
889 |
( cl2->cbar[(2*i)][1] >= (xe1 + fM1) && cl2->cbar[(2*i)][1] <= (xe2 - fM1) && (iy == 0 || multicol) ) || |
890 |
( cl2->cbar[(2*i)][1] >= (xe3 + fM1) && cl2->cbar[(2*i)][1] <= (xe4 - fM1) && (iy == 1 || multicol) ) || |
891 |
( cl2->cbar[(2*i)][1] >= (xe5 + fM1) && cl2->cbar[(2*i)][1] <= (xe6 - fM1) && (iy == 2 || multicol) ) |
892 |
)){ |
893 |
fXesel = true; |
894 |
fXeout = i; |
895 |
} else { |
896 |
fXesel = false; |
897 |
break; |
898 |
}; |
899 |
}; |
900 |
// |
901 |
if ( !fXesel && fXeout >= fXemin && fXeout >= (Int_t)(fMax_planexe+nPl) ){ |
902 |
if ( debug ) printf(" XE: this event is only partially contained: fXeout %i fXemin %i fMax_planexe + nPl %i \n",fXeout,fXemin,(Int_t)(fMax_planexe+nPl)); |
903 |
fPartsel = true; |
904 |
fXesel = true; |
905 |
}; |
906 |
if ( fXesel ){ |
907 |
for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planexe+nPl)) ;iplm++){ |
908 |
fXEen_maxplane += en_xep[iplm]; |
909 |
if ( debug ) printf(" XE iplm %i fXOen_maxplane %f en_xop[iplm] %f\n",iplm,fXEen_maxplane,en_xep[iplm]); |
910 |
}; |
911 |
fEnergyxe = fXEen_maxplane/fConv_rxe; |
912 |
// |
913 |
for (Int_t i=0;i<11;i++){ |
914 |
for(strip=0; strip<96; strip++) { |
915 |
// |
916 |
// run over all the strips of the plane |
917 |
// |
918 |
if ( strip >= 0 && strip < 32 ) encol[1][0] += enstrip[1][(2*i)][strip]; |
919 |
if ( strip >= 32 && strip < 64 ) encol[1][1] += enstrip[1][(2*i)][strip]; |
920 |
if ( strip >= 64 && strip < 96 ) encol[1][2] += enstrip[1][(2*i)][strip]; |
921 |
entot[1] += enstrip[1][(2*i)][strip]; |
922 |
// |
923 |
}; |
924 |
}; |
925 |
}; |
926 |
}; |
927 |
// |
928 |
if ( section.Contains("YE") ){ |
929 |
// |
930 |
// find the column hit in the first plane |
931 |
// |
932 |
Int_t ix = -1; |
933 |
Int_t iy = -1; |
934 |
if ( track_coordy[(2*0)+1][0] >= (-12.154+fM) && track_coordy[(2*0)+1][0] <= (-4.346-fM) ) iy = 0; |
935 |
if ( track_coordy[(2*0)+1][0] >= ( -4.104+fM) && track_coordy[(2*0)+1][0] <= ( 3.704-fM) ) iy = 1; |
936 |
if ( track_coordy[(2*0)+1][0] >= ( 3.946+fM) && track_coordy[(2*0)+1][0] <= (11.754-fM) ) iy = 2; |
937 |
if ( cl2->cbar[(2*0)+1][0] >= (ye1 + fM1) && cl2->cbar[(2*0)+1][0] <= (ye2 - fM1) ) ix = 0; |
938 |
if ( cl2->cbar[(2*0)+1][0] >= (ye3 + fM1) && cl2->cbar[(2*0)+1][0] <= (ye4 - fM1) ) ix = 1; |
939 |
if ( cl2->cbar[(2*0)+1][0] >= (ye5 + fM1) && cl2->cbar[(2*0)+1][0] <= (ye6 - fM1) ) ix = 2; |
940 |
if ( ix > -1 && iy > -1 ) fColYE = ix + iy*3; |
941 |
// |
942 |
// check event is inside XO acceptance |
943 |
// |
944 |
for (Int_t i=0; i<11; i++) { |
945 |
if (( |
946 |
( track_coordy[(2*i)+1][0] >= (-12.154+fM) && track_coordy[(2*i)+1][0] <= (-4.346-fM) && (iy == 0 || multicol) ) || |
947 |
( track_coordy[(2*i)+1][0] >= ( -4.104+fM) && track_coordy[(2*i)+1][0] <= ( 3.704-fM) && (iy == 1 || multicol) ) || |
948 |
( track_coordy[(2*i)+1][0] >= ( 3.946+fM) && track_coordy[(2*i)+1][0] <= (11.754-fM) && (iy == 2 || multicol) ) |
949 |
) && ( |
950 |
( cl2->cbar[(2*i)+1][0] >= (ye1 + fM1) && cl2->cbar[(2*i)+1][0] <= (ye2 - fM1) && (ix == 0 || multicol) ) || |
951 |
( cl2->cbar[(2*i)+1][0] >= (ye3 + fM1) && cl2->cbar[(2*i)+1][0] <= (ye4 - fM1) && (ix == 1 || multicol) ) || |
952 |
( cl2->cbar[(2*i)+1][0] >= (ye5 + fM1) && cl2->cbar[(2*i)+1][0] <= (ye6 - fM1) && (ix == 2 || multicol) ) |
953 |
)){ |
954 |
fYesel = true; |
955 |
fYeout = i; |
956 |
} else { |
957 |
fYesel = false; |
958 |
break; |
959 |
}; |
960 |
}; |
961 |
// |
962 |
if ( !fYesel && fYeout >= fYemin && fYeout >= (Int_t)(fMax_planeye+nPl) ){ |
963 |
if ( debug ) printf(" YE: this event is only partially contained: fYeout %i fYemin %i fMax_planeye + nPl %i \n",fYeout,fYemin,(Int_t)(fMax_planeye+nPl)); |
964 |
fPartsel = true; |
965 |
fYesel = true; |
966 |
}; |
967 |
if ( fYesel ){ |
968 |
for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planeye+nPl)) ;iplm++) fYEen_maxplane += en_yep[iplm]; |
969 |
fEnergyye = fYEen_maxplane/fConv_rye; |
970 |
// |
971 |
for (Int_t i=0;i<11;i++){ |
972 |
for(strip=0; strip<96; strip++) { |
973 |
// |
974 |
// run over all the strips of the plane |
975 |
// |
976 |
if ( strip >= 0 && strip < 32 ) encol[0][0] += enstrip[0][(2*i)+1][strip]; |
977 |
if ( strip >= 32 && strip < 64 ) encol[0][1] += enstrip[0][(2*i)+1][strip]; |
978 |
if ( strip >= 64 && strip < 96 ) encol[0][2] += enstrip[0][(2*i)+1][strip]; |
979 |
entot[0] += enstrip[0][(2*i)+1][strip]; |
980 |
// |
981 |
}; |
982 |
}; |
983 |
// |
984 |
}; |
985 |
}; |
986 |
// |
987 |
if ( section.Contains("YO") ){ |
988 |
// |
989 |
// find the column hit in the first plane |
990 |
// |
991 |
Int_t ix = -1; |
992 |
Int_t iy = -1; |
993 |
if ( track_coordy[(2*0)][0] >= (-11.954+fM) && track_coordy[(2*0)][0] <= (-4.146-fM) ) iy = 0; |
994 |
if ( track_coordy[(2*0)][0] >= ( -3.904+fM) && track_coordy[(2*0)][0] <= ( 3.904-fM) ) iy = 1; |
995 |
if ( track_coordy[(2*0)][0] >= ( 4.146+fM) && track_coordy[(2*0)][0] <= (11.954-fM) ) iy = 2; |
996 |
if ( cl2->cbar[(2*0)][0] >= (yo1 + fM1) && cl2->cbar[(2*0)][0] <= (yo2 - fM1) ) ix = 0; |
997 |
if ( cl2->cbar[(2*0)][0] >= (yo3 + fM1) && cl2->cbar[(2*0)][0] <= (yo4 - fM1) ) ix = 1; |
998 |
if ( cl2->cbar[(2*0)][0] >= (yo5 + fM1) && cl2->cbar[(2*0)][0] <= (yo6 - fM1) ) ix = 2; |
999 |
if ( ix > -1 && iy > -1 ) fColYO = ix + iy*3; |
1000 |
// |
1001 |
// check event is inside XO acceptance |
1002 |
// |
1003 |
for (Int_t i=0; i<11; i++) { |
1004 |
if (( |
1005 |
( track_coordy[(2*i)][0] >= (-11.954+fM) && track_coordy[(2*i)][0] <= (-4.146-fM) && (iy == 0 || multicol) ) || |
1006 |
( track_coordy[(2*i)][0] >= ( -3.904+fM) && track_coordy[(2*i)][0] <= ( 3.904-fM) && (iy == 1 || multicol) ) || |
1007 |
( track_coordy[(2*i)][0] >= ( 4.146+fM) && track_coordy[(2*i)][0] <= (11.954-fM) && (iy == 2 || multicol) ) |
1008 |
) && ( |
1009 |
( cl2->cbar[(2*i)][0] >= (yo1 + fM1) && cl2->cbar[(2*i)][0] <= (yo2 - fM1) && (ix == 0 || multicol) ) || |
1010 |
( cl2->cbar[(2*i)][0] >= (yo3 + fM1) && cl2->cbar[(2*i)][0] <= (yo4 - fM1) && (ix == 1 || multicol) ) || |
1011 |
( cl2->cbar[(2*i)][0] >= (yo5 + fM1) && cl2->cbar[(2*i)][0] <= (yo6 - fM1) && (ix == 2 || multicol) ) |
1012 |
)){ |
1013 |
fYosel = true; |
1014 |
fYoout = i; |
1015 |
} else { |
1016 |
fYosel = false; |
1017 |
break; |
1018 |
}; |
1019 |
}; |
1020 |
// |
1021 |
if ( !fYosel && fYoout >= fYomin && fYoout >= (Int_t)(fMax_planeyo+nPl) ){ |
1022 |
if ( debug ) printf(" YO: this event is only partially contained: fYoout %i fYomin %i fMax_planeyo + nPl %i \n",fYoout,fYomin,(Int_t)(fMax_planeyo+nPl)); |
1023 |
fPartsel = true; |
1024 |
fYosel = true; |
1025 |
}; |
1026 |
if ( fYosel ){ |
1027 |
for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planeyo+nPl)) ;iplm++) fYOen_maxplane += en_yop[iplm]; |
1028 |
fEnergyyo = fYOen_maxplane/fConv_ryo; |
1029 |
// |
1030 |
for (Int_t i=0;i<11;i++){ |
1031 |
for(strip=0; strip<96; strip++) { |
1032 |
// |
1033 |
// run over all the strips of the plane |
1034 |
// |
1035 |
if ( strip >= 0 && strip < 32 ) encol[0][0] += enstrip[0][(2*i)][strip]; |
1036 |
if ( strip >= 32 && strip < 64 ) encol[0][1] += enstrip[0][(2*i)][strip]; |
1037 |
if ( strip >= 64 && strip < 96 ) encol[0][2] += enstrip[0][(2*i)][strip]; |
1038 |
entot[0] += enstrip[0][(2*i)][strip]; |
1039 |
// |
1040 |
}; |
1041 |
}; |
1042 |
}; |
1043 |
}; |
1044 |
// |
1045 |
// Count the number of sections in which the event is contained |
1046 |
// |
1047 |
fCount = (Float_t)((Int_t)fXesel+(Int_t)fXosel+(Int_t)fYesel+(Int_t)fYosel); |
1048 |
// |
1049 |
if ( indep ){ |
1050 |
// |
1051 |
// independent mode, average the energy measurement and max plane of the contained sections |
1052 |
// |
1053 |
fSel = ( fXesel || fYesel || fXosel || fYosel ); |
1054 |
fMax_plane = (Float_t)(fMax_planeyo+fMax_planeye+fMax_planexo+fMax_planexe)/fCount; |
1055 |
fEnergy = (fEnergyxe+fEnergyyo+fEnergyye+fEnergyxo)/fCount; |
1056 |
// |
1057 |
} else { |
1058 |
// |
1059 |
// coherent mode, sum the energy [MIP] of the given sections and convert using fConv_rxo. **** NB: it is assumed that the conversion factor is unique and the method SetConvertionFactor(Float_t) has been used**** The event is selected only if it is contained in all the given sections |
1060 |
// |
1061 |
if ( fCount != fNumSec ){ |
1062 |
fSel = false; |
1063 |
} else { |
1064 |
fSel = true; |
1065 |
}; |
1066 |
fEnergy = (fXEen_maxplane+fYOen_maxplane+fYEen_maxplane+fXOen_maxplane)/fConv_rxo; |
1067 |
if ( fSel ){ |
1068 |
if ( fXesel ) fColumn = fColXE; |
1069 |
if ( fXosel ){ |
1070 |
if ( fColXO != fColumn && fColumn > -1 ){ |
1071 |
printf(" ERROR! mismatch in column number between different sections! fColumn %i fColXO %i \n",fColumn,fColXO); |
1072 |
} else { |
1073 |
fColumn = fColXO; |
1074 |
}; |
1075 |
}; |
1076 |
if ( fYesel ){ |
1077 |
if ( fColYE != fColumn && fColumn > -1 ){ |
1078 |
printf(" ERROR! mismatch in column number between different sections! fColumn %i fColYE %i \n",fColumn,fColYE); |
1079 |
} else { |
1080 |
fColumn = fColYE; |
1081 |
}; |
1082 |
}; |
1083 |
if ( fYosel ){ |
1084 |
if ( fColYO != fColumn && fColumn > -1 ){ |
1085 |
printf(" ERROR! mismatch in column number between different sections! fColumn %i fColYO %i \n",fColumn,fColYO); |
1086 |
} else { |
1087 |
fColumn = fColYO; |
1088 |
}; |
1089 |
}; |
1090 |
}; |
1091 |
}; |
1092 |
// |
1093 |
if ( debug ) printf("sel %i indep %i fMax_plane %f conv_r %f en_maxplane %f encalo %f \n",fSel,indep,fMax_plane,fConv_rxo,fXOen_maxplane,fEnergy); |
1094 |
if ( debug ) printf(" IsInside XE %i XO %i YE %i YO %i => SEL %i \n",fXesel,fXosel,fYesel,fYosel,fSel); |
1095 |
// |
1096 |
// finish exit |
1097 |
// |
1098 |
return fSel; |
1099 |
// |
1100 |
} |
1101 |
|
1102 |
void CaloEnergy::Process(){ |
1103 |
TString xo = "XO"; |
1104 |
this->Process(xo); |
1105 |
} |
1106 |
|
1107 |
|
1108 |
void CaloEnergy::Process(TString section){ |
1109 |
// |
1110 |
// process the event |
1111 |
// |
1112 |
TString ntr = section; |
1113 |
if ( !L2 ){ |
1114 |
printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n"); |
1115 |
printf(" ERROR: CaloEnergy variables not filled \n"); |
1116 |
return; |
1117 |
}; |
1118 |
// |
1119 |
Bool_t newentry = false; |
1120 |
// |
1121 |
if ( L2->IsORB() ){ |
1122 |
if ( L2->GetOrbitalInfo()->pkt_num != PKT || L2->GetOrbitalInfo()->OBT != OBT || L2->GetOrbitalInfo()->absTime != atime || strcmp(ntr.Data(),sntr.Data()) ){ |
1123 |
newentry = true; |
1124 |
OBT = L2->GetOrbitalInfo()->OBT; |
1125 |
PKT = L2->GetOrbitalInfo()->pkt_num; |
1126 |
atime = L2->GetOrbitalInfo()->absTime; |
1127 |
sntr = ntr; |
1128 |
}; |
1129 |
} else { |
1130 |
newentry = true; |
1131 |
}; |
1132 |
// |
1133 |
// if we have already called this method for this event and no input changed then return fSel and exit |
1134 |
// |
1135 |
if ( !newentry ) return; |
1136 |
// |
1137 |
// process the event |
1138 |
// |
1139 |
if ( debug ) printf(" Processing event at OBT %u PKT %u time %u section %s\n",OBT,PKT,atime,section.Data()); |
1140 |
// |
1141 |
// check if the cylinder of integration can go out of the sensor given the frame which has been set (if we use all the calorimeter fRad is < 0 and the printout is suppressed) |
1142 |
// |
1143 |
if ( (fM1+0.122-0.244*(Float_t)fRad) < 0.000000001 ) printf("Error: (fM1+0.122-0.244*(Float_t)fRad) < 0. fM1 %f fRad %i %f \n",fM1,fRad,(fM1+0.122-0.244*(Float_t)fRad)); |
1144 |
// |
1145 |
if ( fLong ){ |
1146 |
if ( debug ) printf(" ==================================================================> LONGITUDINAL FIT! \n"); |
1147 |
// |
1148 |
// use long fit to measure energy |
1149 |
// |
1150 |
if ( this->IsInsideAcceptance(section) ){ |
1151 |
// |
1152 |
if ( debug ) printf(" ==================================================================> LONG INSIDE! \n"); |
1153 |
// |
1154 |
Float_t myene[2][22]; |
1155 |
memset(myene,0,(sizeof(Float_t))*2*22); |
1156 |
for (Int_t j=0; j<11; j++){ |
1157 |
if ( section.Contains("XE") ) myene[1][2*j] = en_xep[j]; |
1158 |
if ( section.Contains("YO") ) myene[0][2*j] = en_yop[j]; |
1159 |
if ( section.Contains("XO") ) myene[1][(2*j)+1] = en_xop[j]; |
1160 |
if ( section.Contains("YE") ) myene[0][(2*j)+1] = en_yep[j]; |
1161 |
}; |
1162 |
clong->UnMaskSections(); |
1163 |
if ( !(section.Contains("YE")) ) clong->MaskSection("YE"); |
1164 |
if ( !(section.Contains("YO")) ) clong->MaskSection("YO"); |
1165 |
if ( !(section.Contains("XO")) ) clong->MaskSection("XO"); |
1166 |
if ( !(section.Contains("XE")) ) clong->MaskSection("XE"); |
1167 |
clong->ForceNextFit(); |
1168 |
clong->SetEnergies(myene); |
1169 |
if ( debug ){ |
1170 |
clong->Fit(true); |
1171 |
} else { |
1172 |
clong->Fit(); |
1173 |
}; |
1174 |
if ( clong->GetLowerLimit() != 0. || clong->GetUpperLimit() != 0. ){ |
1175 |
fXOen_maxplane = clong->Get_defE0(); |
1176 |
} else { |
1177 |
fXOen_maxplane = clong->Get_E0(); |
1178 |
}; |
1179 |
fMax_plane = clong->Get_tmax(); |
1180 |
fYOen_maxplane = 0.; |
1181 |
fYEen_maxplane = 0.; |
1182 |
fXEen_maxplane = 0.; |
1183 |
fEnergy=fXOen_maxplane/fConv_rxo; |
1184 |
if ( fEnergy != fEnergy || clong->Get_fitresult() != 0 ) fEnergy = -1.; |
1185 |
// if ( fEnergy != fEnergy ) fEnergy = 1.; |
1186 |
// |
1187 |
} else { |
1188 |
// |
1189 |
// if the event is not in the acceptance, return a negative energy. |
1190 |
// |
1191 |
if ( debug ) printf(" Outside acceptance \n"); |
1192 |
fEnergy *= -1.; |
1193 |
// |
1194 |
}; |
1195 |
// |
1196 |
} else { |
1197 |
// |
1198 |
// use the energy measurement |
1199 |
// |
1200 |
if ( this->IsInsideAcceptance(section) ){ |
1201 |
// |
1202 |
// the event is good |
1203 |
// |
1204 |
if ( debug ) printf(" XE %i XO %i YE %i YO %i \n",fXesel,fXosel,fYesel,fYosel); |
1205 |
// |
1206 |
} else { |
1207 |
// |
1208 |
// if the event is not in the acceptance, return a negative energy. |
1209 |
// |
1210 |
if ( debug ) printf(" Outside acceptance \n"); |
1211 |
fEnergy *= -1.; |
1212 |
// |
1213 |
}; |
1214 |
}; |
1215 |
// |
1216 |
} |