/[PAMELA software]/calo/flight/CaloEnergy/src/CaloEnergy.cpp
ViewVC logotype

Contents of /calo/flight/CaloEnergy/src/CaloEnergy.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.15 - (show annotations) (download)
Thu Sep 10 12:53:30 2009 UTC (15 years, 4 months ago) by mocchiut
Branch: MAIN
Changes since 1.14: +17 -3 lines
Bugs fixed and new methods implemented

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

  ViewVC Help
Powered by ViewVC 1.1.23