/[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.17 - (show annotations) (download)
Fri Nov 27 10:31:08 2009 UTC (15 years ago) by mocchiut
Branch: MAIN
Changes since 1.16: +2 -0 lines
Minor changes

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

  ViewVC Help
Powered by ViewVC 1.1.23