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

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

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.1.1.1 by mocchiut, Mon Jun 29 09:42:10 2009 UTC revision 1.12 by mocchiut, Tue Aug 18 09:24:51 2009 UTC
# Line 17  CaloEnergy::CaloEnergy(PamLevel2 *l2p){ Line 17  CaloEnergy::CaloEnergy(PamLevel2 *l2p){
17    //    //
18    if ( !L2->IsORB() ) printf(" WARNING: OrbitalInfo Tree is needed, the plugin could not work properly without it \n");    if ( !L2->IsORB() ) printf(" WARNING: OrbitalInfo Tree is needed, the plugin could not work properly without it \n");
19    //    //
   OBT = 0;  
   PKT = 0;  
   atime = 0;  
   sntr = "start";  
   AOBT = 0;  
   APKT = 0;  
   aatime = 0;  
   asntr = "start";  
   //  
   debug = false;  
   //  
   fLong = false;  
20    fSimu = false;    fSimu = false;
21    fPl = 1;    this->Set();
   //  
   this->DefineGeometry();  
   fXosel = true;  
   fXesel = true;  
   fYosel = true;  
   fYesel = true;  
   fConv_r = 44.4;  
22    //    //
23  }  }
24    
# Line 49  CaloEnergy::CaloEnergy(PamLevel2 *l2p, B Line 30  CaloEnergy::CaloEnergy(PamLevel2 *l2p, B
30    //    //
31    if ( !L2->IsORB() ) printf(" WARNING: OrbitalInfo Tree is needed, the plugin could not work properly without it \n");    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::UseCaloPreSampler(){
39      //
40      // 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
41      //
42      if ( !cp ) cp = new CaloPreSampler(L2);
43      cp->SplitInto(0,22);
44      cp->SetForceFitMode(1000);
45    //   cp->UseTracker(false);
46    //   cp->ForceCaloFit();
47    //   cp->SetDebug(true);
48    //   cp->Process();
49      if ( clong ) clong->SetCaloLevel2Pointer(cp->GetLevel2Pointer());
50    }
51    
52    
53    void CaloEnergy::UseLongFit(){
54      fPl = 0;
55      fLong = true;
56      if ( !clong ){
57        clong = new CaloLong(L2);
58        if ( cp ) clong->SetCaloLevel2Pointer(cp->GetLevel2Pointer());
59        clong->SplitInto(0,22);
60      };
61      //
62    }
63    
64    void CaloEnergy::Set(){
65      //
66      // set default values, NB default conversion factor for energy is just very approximated!
67      //
68    OBT = 0;    OBT = 0;
69    PKT = 0;    PKT = 0;
70    atime = 0;    atime = 0;
# Line 61  CaloEnergy::CaloEnergy(PamLevel2 *l2p, B Line 77  CaloEnergy::CaloEnergy(PamLevel2 *l2p, B
77    //    //
78    debug = false;    debug = false;
79    //    //
80      indep = false;
81      //
82    fLong = false;    fLong = false;
   fSimu = simulation;  
83    fPl = 1;    fPl = 1;
84      fRad = -1;
85      cp = NULL;
86      clong = NULL;
87      x0max = -1.;
88    //    //
89    this->DefineGeometry();    this->DefineGeometry();
90    fXosel =true;    fXosel =true;
91    fXesel = true;    fXesel = true;
92    fYosel = true;    fYosel = true;
93    fYesel = true;    fYesel = true;
94    fConv_r = 44.4;    fConv_rxe = 44.4;
95      fConv_rxo = 44.4;
96      fConv_ryo = 44.4;
97      fConv_rye = 44.4;
98      fXomin = 1000;
99      fXemin = 1000;
100      fYomin = 1000;
101      fYemin = 1000;
102    //    //
103  }  }
104    
105    void CaloEnergy::SetMinimumContainment(TString section, Int_t plane){
106      section.ToUpper();
107      if ( section.Contains("XO") ) fXomin = plane;
108      if ( section.Contains("XE") ) fXemin = plane;
109      if ( section.Contains("YO") ) fYomin = plane;
110      if ( section.Contains("YE") ) fYemin = plane;
111    }
112    
113    void CaloEnergy::SetMinimumContainment(Int_t plane){
114      this->SetMinimumContainment("XEXOYEYO",plane);
115    }
116    
117    Int_t CaloEnergy::GetMinimumContainment(TString section){
118      section.ToUpper();
119      if ( section.Contains("XO") ) return(fXomin);
120      if ( section.Contains("XE") ) return(fXemin);
121      if ( section.Contains("YE") ) return(fYemin);
122      if ( section.Contains("YO") ) return(fYomin);
123      printf(" ERROR: section not recognized \n");
124      return(-1000);
125    }
126    
127    void CaloEnergy::SetConversionFactor(TString section, Float_t conv){
128      section.ToUpper();
129      if ( section.Contains("XO") ) fConv_rxo = conv;
130      if ( section.Contains("XE") ) fConv_rxe = conv;
131      if ( section.Contains("YO") ) fConv_ryo = conv;
132      if ( section.Contains("YE") ) fConv_rye = conv;
133    }
134    
135    void CaloEnergy::SetConversionFactor(Float_t conv){
136      this->SetConversionFactor("XEXOYEYO",conv);
137    }
138    
139    Float_t CaloEnergy::GetConversionFactor(TString section){
140      section.ToUpper();
141      if ( section.Contains("XO") ) return(fConv_rxo);
142      if ( section.Contains("XE") ) return(fConv_rxe);
143      if ( section.Contains("YO") ) return(fConv_ryo);
144      if ( section.Contains("YE") ) return(fConv_rye);
145      printf(" ERROR: section not recognized \n");
146      return(-1000.);
147    }
148    
149    Int_t CaloEnergy::GetMaxplane(TString section){
150      section.ToUpper();
151      if ( section.Contains("XO") ) return fMax_planexo;
152      if ( section.Contains("XE") ) return fMax_planexe;
153      if ( section.Contains("YO") ) return fMax_planeyo;
154      if ( section.Contains("YE") ) return fMax_planeye;
155      return(-1);
156    }
157    
158    Float_t CaloEnergy::GetEnergyAtMaxplane(TString section){
159      section.ToUpper();
160      if ( section.Contains("XO") ) return xomax_en;
161      if ( section.Contains("XE") ) return xemax_en;
162      if ( section.Contains("YO") ) return yomax_en;
163      if ( section.Contains("YE") ) return yemax_en;
164      return(-1);
165    }
166    
167    Float_t CaloEnergy::GetMaxEnergy(TString section){
168      section.ToUpper();
169      if ( fLong ){
170        this->Process(section);    
171        return fXOen_maxplane;
172      } else {
173        if ( section.Contains("XO") ) return fXOen_maxplane;
174        if ( section.Contains("XE") ) return fXEen_maxplane;
175        if ( section.Contains("YO") ) return fYOen_maxplane;
176        if ( section.Contains("YE") ) return fYEen_maxplane;
177      };
178      return(-1);
179    }
180    
181    Float_t CaloEnergy::GetMaxEnergy(){  
182      if ( fLong ){
183        if ( debug ) printf(" oh! call process! with asntr %s and sntr %s \n",asntr.Data(),sntr.Data());
184        this->Process(asntr);
185      };
186      return((fXEen_maxplane+fYOen_maxplane+fYEen_maxplane+fXOen_maxplane));
187    }
188    
189  void CaloEnergy::DefineGeometry(){  void CaloEnergy::DefineGeometry(){
190    //    //
191    fM = 2.;    // Use CaloStrip to determine once the position of border strips for each section
192    fM1 = fM - 0.122 -0.096;    //
193      //
194      //  fM = 2. + 0.096; // real position from cbar BUG the 0.096 is already taken into account in the border calculation made by Giovanna
195      fM = 2. ; // real position from cbar
196      //  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
197      fM1 = 2. - 0.122 - 0.096 + 0.096; // due to calculation of xe1 etc.
198    if ( fM1 < 0. ) fM1 = 0.;    if ( fM1 < 0. ) fM1 = 0.;
199    //    //
200    CaloStrip *cs = new CaloStrip(fSimu);    CaloStrip *cs = new CaloStrip(fSimu);
201      //
202    // vista y piano 0 striscia 0    // view y plane 0 strip 0
203    cs->Set(1,0,0);    cs->Set(1,0,0);
204    xe1= cs->GetY();    xe1= cs->GetY();
205    z1= cs->GetZ();    //  z1= cs->GetZ();
206    // vista y piano 0 striscia 31    // view y plane 0 strip 31
207    cs->Set(1,0,31);    cs->Set(1,0,31);
208    xe2= cs->GetY();    xe2= cs->GetY();
209    // vista y piano 0 striscia 32    // view y plane 0 strip 32
210    cs->Set(1,0,32);    cs->Set(1,0,32);
211    xe3= cs->GetY();    xe3= cs->GetY();
212    // vista y piano 0 striscia 63    // view y plane 0 strip 63
213    cs->Set(1,0,63);    cs->Set(1,0,63);
214    xe4= cs->GetY();    xe4= cs->GetY();
215    // vista y piano 0 striscia 64    // view y plane 0 strip 64
216    cs->Set(1,0,64);    cs->Set(1,0,64);
217    xe5= cs->GetY();    xe5= cs->GetY();
218    // vista y piano 0 striscia 95    // view y plane 0 strip 95
219    cs->Set(1,0,95);    cs->Set(1,0,95);
220    xe6= cs->GetY();    xe6= cs->GetY();
221    
222    // vista x piano 0 striscia 0    //  printf(" xe1 %f xe2 %f xe3 %f xe4 %f xe5 %f xe6 %f \n",xe1,xe2,xe3,xe4,xe5,xe6);
223    
224      // view x plane 0 strip 0
225    cs->Set(0,0,0);    cs->Set(0,0,0);
226    yo1= cs->GetX();    yo1= cs->GetX();
227    z2= cs->GetZ();    //  z2= cs->GetZ();
228    // vista x piano 0 striscia 31    // view x plane 0 strip 31
229    cs->Set(0,0,31);    cs->Set(0,0,31);
230    yo2= cs->GetX();    yo2= cs->GetX();
231    // vista x piano 0 striscia 32    // view x plane 0 strip 32
232    cs->Set(0,0,32);    cs->Set(0,0,32);
233    yo3= cs->GetX();    yo3= cs->GetX();
234    // vista x piano 0 striscia 63    // view x plane 0 strip 63
235    cs->Set(0,0,63);    cs->Set(0,0,63);
236    yo4= cs->GetX();    yo4= cs->GetX();
237    // vista x piano 0 striscia 64    // view x plane 0 strip 64
238    cs->Set(0,0,64);    cs->Set(0,0,64);
239    yo5= cs->GetX();    yo5= cs->GetX();
240    // vista x piano 0 striscia 95    // view x plane 0 strip 95
241    cs->Set(0,0,95);    cs->Set(0,0,95);
242    yo6= cs->GetX();    yo6= cs->GetX();
243    
244        // view y plane 1 strip 0
   // vista y piano 1 striscia 0  
245    cs->Set(1,1,0);    cs->Set(1,1,0);
246    xo1= cs->GetY();    xo1= cs->GetY();
247    z3= cs->GetZ();    //  z3= cs->GetZ();
248    // vista y piano 1 striscia 31    // view y plane 1 strip 31
249    cs->Set(1,1,31);    cs->Set(1,1,31);
250    xo2= cs->GetY();    xo2= cs->GetY();
251    // vista y piano 1 striscia 32    // view y plane 1 strip 32
252    cs->Set(1,1,32);    cs->Set(1,1,32);
253    xo3= cs->GetY();    xo3= cs->GetY();
254    // vista y piano 1 striscia 63    // view y plane 1 strip 63
255    cs->Set(1,1,63);    cs->Set(1,1,63);
256    xo4= cs->GetY();    xo4= cs->GetY();
257    // vista y piano 1 striscia 64    // view y plane 1 strip 64
258    cs->Set(1,1,64);    cs->Set(1,1,64);
259    xo5= cs->GetY();    xo5= cs->GetY();
260    // vista y piano 1 striscia 95    // view y plane 1 strip 95
261    cs->Set(1,1,95);    cs->Set(1,1,95);
262    xo6= cs->GetY();    xo6= cs->GetY();
263        
264    // vista x piano 1 striscia 0    // view x plane 1 strip 0
265    cs->Set(0,1,0);    cs->Set(0,1,0);
266    ye1= cs->GetX();    ye1= cs->GetX();
267    z4= cs->GetZ();    //  z4= cs->GetZ();
268    // vista x piano 1 striscia 31    // view x plane 1 strip 31
269    cs->Set(0,1,31);    cs->Set(0,1,31);
270    ye2= cs->GetX();    ye2= cs->GetX();
271    // vista x piano 1 striscia 32    // view x plane 1 strip 32
272    cs->Set(0,1,32);    cs->Set(0,1,32);
273    ye3= cs->GetX();    ye3= cs->GetX();
274    // vista x piano 1 striscia 63    // view x plane 1 strip 63
275    cs->Set(0,1,63);    cs->Set(0,1,63);
276    ye4= cs->GetX();    ye4= cs->GetX();
277    // vista x piano 1 striscia 64    // view x plane 1 strip 64
278    cs->Set(0,1,64);    cs->Set(0,1,64);
279    ye5= cs->GetX();    ye5= cs->GetX();
280    // vista x piano 1 striscia 95    // view x plane 1 strip 95
281    cs->Set(0,1,95);    cs->Set(0,1,95);
282    ye6= cs->GetX();    ye6= cs->GetX();
283    
284      //printf(" ye1 %f ye2 %f ye3 %f ye4 %f ye5 %f ye6 %f \n",ye1,ye2,ye3,ye4,ye5,ye6);
285      
286    //    //
287    for (Int_t p = 0; p<22; p ++){    for (Int_t p = 0; p<22; p ++){
288      for (Int_t v = 0; v<2; v++ ){        for (Int_t v = 0; v<2; v++ ){  
289        cs->Set(v,p,0);        cs->Set(v,p,0);
290        trk_z[p][v]= cs->GetZ();                  // << cooordinata Z piano        trk_z[p][v]= cs->GetZ();                  //  Z coord for each plane
291      };      };
292    };    };
293    //    //
# Line 177  void CaloEnergy::DefineGeometry(){ Line 297  void CaloEnergy::DefineGeometry(){
297    
298  void CaloEnergy::Clear(){  void CaloEnergy::Clear(){
299    //    //
300      // clear variables
301      //
302      fPartsel = false;
303      fSel = false;
304      fXosel = false;
305      fXesel = false;
306      fYosel = false;
307      fYesel = false;
308      fCount = 0.;
309    fEnergy = 0.;    fEnergy = 0.;
310      fEnergyxe = 0.;
311      fEnergyxo = 0.;
312      fEnergyye = 0.;
313      fEnergyyo = 0.;
314      fMax_plane = 0;
315      fMax_planexo = 0;
316      fMax_planexe = 0;
317      fMax_planeyo = 0;
318      fMax_planeye = 0;
319      xomax_en= 0.;
320      xemax_en= 0.;
321      yomax_en= 0.;
322      yemax_en= 0.;
323      //
324    memset(enstrip,0,2*22*96*(sizeof(Float_t)));      memset(enstrip,0,2*22*96*(sizeof(Float_t)));  
325    en = 0.;    en = 0.;
326    view = 0;    view = 0;
# Line 187  void CaloEnergy::Clear(){ Line 330  void CaloEnergy::Clear(){
330    energyyo = 0.;    energyyo = 0.;
331    energyxo = 0.;    energyxo = 0.;
332    energyye = 0.;    energyye = 0.;
333      fYoout = 0;
334      fYeout = 0;
335      fXoout = 0;
336      fXeout = 0;
337      fXEen_maxplane = 0.;  
338      fXOen_maxplane = 0.;  
339      fYEen_maxplane = 0.;  
340      fYOen_maxplane = 0.;  
341    memset(en_xep,0,11*sizeof(Float_t));    memset(en_xep,0,11*sizeof(Float_t));
342    memset(en_yep,0,11*sizeof(Float_t));    memset(en_yep,0,11*sizeof(Float_t));
343    memset(en_xop,0,11*sizeof(Float_t));    memset(en_xop,0,11*sizeof(Float_t));
344    memset(en_yop,0,11*sizeof(Float_t));    memset(en_yop,0,11*sizeof(Float_t));
345    //    //
346      fColumn = -1;
347      memset(encol,0,2*3*sizeof(Float_t));
348      //
349  }  }
350    
351  void CaloEnergy::Print(){  void CaloEnergy::Print(){
# Line 199  void CaloEnergy::Print(){ Line 353  void CaloEnergy::Print(){
353    printf("========================================================================\n");    printf("========================================================================\n");
354    printf(" OBT: %u PKT: %u ATIME: %u \n",OBT,PKT,atime);    printf(" OBT: %u PKT: %u ATIME: %u \n",OBT,PKT,atime);
355    printf(" fEnergy :.............. %f \n",fEnergy);      printf(" fEnergy :.............. %f \n",fEnergy);  
356    printf(" fMax_plane :........... %i \n",fMax_plane);    printf(" fMax_plane :........... %f \n",fMax_plane);
357      printf(" fMax_planexo :......... %i \n",fMax_planexo);
358      printf(" fMax_planexe :......... %i \n",fMax_planexe);
359      printf(" fMax_planeyo :......... %i \n",fMax_planeyo);
360      printf(" fMax_planeye :......... %i \n",fMax_planeye);
361    printf(" fCount  :.............. %f \n",fCount);    printf(" fCount  :.............. %f \n",fCount);
362      printf(" fSel    :.............. %i \n",fSel);
363      printf(" fPartsel:.............. %i \n",fPartsel);
364    printf(" fXesel  :.............. %i \n",fXesel);    printf(" fXesel  :.............. %i \n",fXesel);
365    printf(" fXosel  :.............. %i \n",fXosel);    printf(" fXosel  :.............. %i \n",fXosel);
366    printf(" fYesel  :.............. %i \n",fYesel);    printf(" fYesel  :.............. %i \n",fYesel);
367    printf(" fYosel  :.............. %i \n",fYosel);    printf(" fYosel  :.............. %i \n",fYosel);
368      printf(" fXemin  :.............. %i \n",fXemin);
369      printf(" fXomin  :.............. %i \n",fXomin);
370      printf(" fYemin  :.............. %i \n",fYemin);
371      printf(" fYomin  :.............. %i \n",fYomin);
372      printf(" fXeout  :.............. %i \n",fXeout);
373      printf(" fXoout  :.............. %i \n",fXoout);
374      printf(" fYeout  :.............. %i \n",fYeout);
375      printf(" fYoout  :.............. %i \n",fYoout);
376    printf(" fSimu   :.............. %i \n",fSimu);    printf(" fSimu   :.............. %i \n",fSimu);
377    printf(" fM      :.............. %f \n",fM);    printf(" fM      :.............. %f \n",fM);
378    printf(" fM1     :.............. %f \n",fM1);    printf(" fM1     :.............. %f \n",fM1);
379    printf(" fPl     :.............. %f \n",fPl);    printf(" fRad    :.............. %i \n",fRad);
380    printf(" fConv_r :.............. %f \n",fConv_r);    printf(" fPl     :.............. %i \n",fPl);
381      printf(" fColumn :.............. %i \n",fColumn);
382      printf(" fConv_rxe ............. %f \n",fConv_rxe);
383      printf(" fConv_rxo ............. %f \n",fConv_rxo);
384      printf(" fConv_ryo ............. %f \n",fConv_ryo);
385      printf(" fConv_rye ............. %f \n",fConv_rye);
386    printf(" fLong   :.............. %i \n",fLong);    printf(" fLong   :.............. %i \n",fLong);
387    printf(" energyxe:.............. %f \n",energyxe);    printf(" energyxe:.............. %f \n",energyxe);
388    printf(" energyxo:.............. %f \n",energyxo);    printf(" energyxo:.............. %f \n",energyxo);
389    printf(" energyye:.............. %f \n",energyye);    printf(" energyye:.............. %f \n",energyye);
390    printf(" energyyo:.............. %f \n",energyyo);    printf(" energyyo:.............. %f \n",energyyo);
391      printf(" fXEen_maxplane:........ %f \n",fXEen_maxplane);
392      printf(" fXOen_maxplane:........ %f \n",fXOen_maxplane);
393      printf(" fYEen_maxplane:........ %f \n",fYEen_maxplane);
394      printf(" fYOen_maxplane:........ %f \n",fYOen_maxplane);
395      printf(" x0max   :.............. %f \n",x0max);
396    printf(" debug   :.............. %i \n",debug);    printf(" debug   :.............. %i \n",debug);
397    
398    printf("========================================================================\n");    printf("========================================================================\n");
# Line 223  void CaloEnergy::Print(){ Line 401  void CaloEnergy::Print(){
401    
402  void CaloEnergy::Delete(){  void CaloEnergy::Delete(){
403    Clear();    Clear();
404      delete this;
405  }  }
406    
407    
408  Bool_t CaloEnergy::IsInsideAcceptance(TString section){  Bool_t CaloEnergy::IsInsideAcceptance(TString section){
   return(this->IsInsideAcceptance(section,true));    
 }  
   
 Bool_t CaloEnergy::IsInsideAcceptance(TString section, Bool_t fast){  
409    //    //
410      // check if the event is inside the acceptance of the given section(s)
411    //      //  
412    TString ntr = section;    TString ntr = section;
413    if ( !L2 ){    if ( !L2 ){
# Line 254  Bool_t CaloEnergy::IsInsideAcceptance(TS Line 430  Bool_t CaloEnergy::IsInsideAcceptance(TS
430      newentry = true;      newentry = true;
431    };    };
432    //    //
433      // if we have already called this method for this event and no input changed then return fSel and exit
434      //
435    if ( !newentry ) return fSel;    if ( !newentry ) return fSel;
436    //    //
437      // process the event
438      //
439    if ( debug ) printf(" ########## IsInsideAcceptance ######### \n");    if ( debug ) printf(" ########## IsInsideAcceptance ######### \n");
440    //    //
441    fSel = false;    // clear variables
442      //
443      this->Clear();
444      //
445    section.ToUpper();    section.ToUpper();
   fXosel = false;  
   fXesel = false;  
   fYosel = false;  
   fYesel = false;  
   fCount = 0.;  
446    //    //
447    Int_t ss = Int_t(section.Contains("XO"))+Int_t(section.Contains("XE"))+Int_t(section.Contains("YO"))+Int_t(section.Contains("YE"));    // Count the number of section(s) given as input
448    if ( ss > 1 ) fast = false;    //
449    if ( !ss ){    Int_t fNumSec = Int_t(section.Contains("XO"))+Int_t(section.Contains("XE"))+Int_t(section.Contains("YO"))+Int_t(section.Contains("YE"));
450      if ( !fNumSec ){
451      printf(" ERROR: section must be XO or XE or YO or YE while it is %s \n",section.Data());      printf(" ERROR: section must be XO or XE or YO or YE while it is %s \n",section.Data());
452        return false;
453      };
454      //
455      // If the presampler object exists then use the presampler output instead of the level2 output
456      //
457      CaloLevel2 *cl2 = NULL;
458      if ( cp ){
459        cl2 = cp->GetCaloLevel2();
460      } else {
461        cl2 = L2->GetCaloLevel2();
462      };
463      //
464      // get the energy for every strip of the calorimeter
465      //
466      for (Int_t ch=0; ch< L2->GetCaloLevel1()->istrip; ch++){
467        en = L2->GetCaloLevel1()->DecodeEstrip(ch,view,plane,strip);
468        enstrip[view][plane][strip]=en;
469      };
470      //
471      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));
472      //
473      // sum energy plane by plane for each sections
474      //
475      for (Int_t i=0;i<11;i++){
476        for(strip=0; strip<96; strip++) {
477          if ( fRad < 0 ){
478            //
479            // run over all the strips of the plane
480            //
481            en_xep[i] += enstrip[1][2*i][strip];
482            en_yop[i] += enstrip[0][2*i][strip];
483            en_xop[i] += enstrip[1][(2*i)+1][strip];
484            en_yep[i] += enstrip[0][(2*i)+1][strip];
485          } else {
486            //
487            // use only the strips inside a cylinder of given radius fRad
488            //
489            if ( cl2->cibar[2*i][1] >= 1 && cl2->cibar[2*i][1] <= 96 && (strip >= (cl2->cibar[2*i][1]-1-fRad)) &&  (strip <= (cl2->cibar[2*i][1]-1+fRad)) ) en_xep[i] += enstrip[1][2*i][strip];
490            if ( cl2->cibar[2*i][0] >= 1 && cl2->cibar[2*i][0] <= 96 && (strip >= (cl2->cibar[2*i][0]-1-fRad)) &&  (strip <= (cl2->cibar[2*i][0]-1+fRad)) ) en_yop[i] += enstrip[0][2*i][strip];
491            if ( cl2->cibar[(2*i)+1][1] >= 1 && cl2->cibar[(2*i)+1][1] <= 96 && (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];
492            if ( cl2->cibar[(2*i)+1][0] >= 1 && cl2->cibar[(2*i)+1][0] <= 96 && (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];
493          };
494        };
495        energyxe += en_xep[i];
496        energyyo += en_yop[i];
497        energyxo += en_xop[i];
498        energyye += en_yep[i];
499      };
500      //
501      // Find the plane of maximum for each section
502      //
503      //
504      Int_t xen = 0;
505      Int_t yon = 0;
506      Int_t xon = 0;
507      Int_t yen = 0;
508      //
509      if ( section.Contains("XE") ){
510        yon++;
511        xon++;
512        yen++;
513        for (Int_t ipl =0; ipl < 11; ipl ++) {
514          if(en_xep[ipl] > xemax_en) {
515            xemax_en = en_xep[ipl];
516            fMax_planexe = ipl;
517          };
518        };
519      };
520      //
521      if ( section.Contains("YO") ){
522        xon++;
523        yen++;
524        for (Int_t ipl =0; ipl < 11; ipl ++) {
525          if(en_yop[ipl] > yomax_en) {
526            yomax_en = en_yop[ipl];
527            fMax_planeyo = ipl;
528          };
529        };
530      };
531      //
532      if ( section.Contains("XO") ){
533        yen++;
534        for (Int_t ipl =0; ipl < 11; ipl ++) {
535          if(en_xop[ipl] > xomax_en) {
536            xomax_en = en_xop[ipl];
537            fMax_planexo = ipl;
538          };
539        };
540      };
541      //
542      if ( section.Contains("YE") ){
543        for (Int_t ipl =0; ipl < 11; ipl ++) {
544          if(en_yep[ipl] > yemax_en) {
545            yemax_en = en_yep[ipl];
546            fMax_planeye = ipl;
547          };
548        };
549      };
550      //
551      // the maximum is given externally as X0, convert it to plane and section
552      //
553      if ( x0max > 0. ){
554        if ( debug ) printf(" CALCULATE MAX PLANE GIVEN X0 ASSUMING PERPENDICULAR TRACK \n");
555        Int_t wpl = (Int_t)roundf(x0max/0.76);
556        Bool_t isY = false;
557        if ( ((x0max/0.76)-(Float_t)wpl) > 0. ) isY = true;
558        xomax_en = 0.;
559        yemax_en = 0.;
560        xemax_en = 0.;
561        yomax_en = 0.;
562        //
563        if ( !(wpl%2) ){
564          // 0, 2, 4, ...
565          if ( isY ){
566            if ( section.Contains("YO") ) yomax_en = 1000.;
567            if ( section.Contains("XE") ) xemax_en = 500.;
568            fMax_planeyo=wpl/2;
569            fMax_planexe=wpl/2;
570            if ( section.Contains("XO") ) xomax_en = 10.;
571            if ( section.Contains("YE") ) yemax_en = 5.;
572          } else {
573            if ( section.Contains("YO") ) yomax_en = 500.;
574            if ( section.Contains("XE") ) xemax_en = 1000.;
575            fMax_planeyo=wpl/2;
576            fMax_planexe=wpl/2;
577            if ( section.Contains("XO") ) xomax_en = 5.;
578            if ( section.Contains("YE") ) yemax_en = 10.;
579          };
580        } else {
581          // 1, 3, 5, ...
582          if ( isY ){
583            if ( section.Contains("YE") ) yemax_en = 1000.;
584            if ( section.Contains("XO") ) xomax_en = 500.;
585            fMax_planeye=(wpl-1)/2;
586            fMax_planexo=(wpl-1)/2;
587            if ( section.Contains("XE") ) xemax_en = 10.;
588            if ( section.Contains("YO") ) yomax_en = 5.;
589          } else {
590            if ( section.Contains("YE") ) yemax_en = 500.;
591            if ( section.Contains("XO") ) xomax_en = 1000.;
592            fMax_planeye=(wpl-1)/2;
593            fMax_planexo=(wpl-1)/2;
594            if ( section.Contains("XE") ) xemax_en = 5.;
595            if ( section.Contains("YO") ) yomax_en = 10.;
596          };
597        };
598        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);
599    };    };
600    //    //
601    Float_t  track_coordx[22][2];    Int_t nPl = fPl;
602    Float_t  track_coordy[22][2];    //
603      // Set the maximum in case of coherent mode was selected
604      //
605      if ( !indep ){
606        nPl = 0;
607        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);
608        Int_t nummod = 0;
609        Int_t numexpl = 0;
610        if ( xomax_en > xemax_en && xomax_en > yemax_en && xomax_en > yomax_en ){
611          //
612          // Section XO contains the maximum energy release per plane of the whole calorimeter
613          //
614          if ( debug ) printf(" XO is MAX %i %i %i %i\n",xen,yon,xon,yen);
615          //
616          // fMax_plane is the plane of maximum + number of additional dE/dx measurement counting planes from 0 to 43
617          //
618          fMax_plane = (fNumSec * fMax_planexo) +(Float_t)xon + fPl;
619          //
620          // nummod is the integer part of the number of modules in which the maximum is contained
621          //
622          nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec);
623          //
624          // numexpl is the number of additional planes (0,1,2) inside the module
625          //
626          numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod);
627          //
628        };
629        if ( xemax_en > xomax_en && xemax_en > yemax_en && xemax_en > yomax_en ){
630          if ( debug ) printf(" XE is MAX %i %i %i %i\n",xen,yon,xon,yen);
631          fMax_plane = (fNumSec * fMax_planexe) +(Float_t)xen + fPl;
632          nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec);
633          numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod);
634          //
635       };
636    
637        if ( yemax_en > xomax_en && yemax_en > xemax_en && yemax_en > yomax_en ){
638          if ( debug ) printf(" YE is MAX %i %i %i %i\n",xen,yon,xon,yen);
639          fMax_plane = (fNumSec * fMax_planeye) +(Float_t)yen + fPl;
640          nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec);
641          numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod);
642          //
643        };
644        if ( yomax_en > xemax_en && yomax_en > yemax_en && yomax_en > xomax_en ){
645          if ( debug ) printf(" YO is MAX %i %i %i %i\n",xen,yon,xon,yen);
646          fMax_plane = (fNumSec * fMax_planeyo) +(Float_t)yon + fPl;
647          nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec);
648          numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod);
649          //
650        };
651        //
652        // find the plane up to which is necessary to integrate the energy for each section
653        //
654        Int_t a = 0;
655        Int_t b = 0;
656        Int_t c = 0;
657        if ( numexpl > xen ) a = 1;
658        if ( numexpl > yon ) b = 1;
659        if ( numexpl > xon ) c = 1;
660        fMax_planexe = nummod;
661        fMax_planeyo = nummod - 1 + a;
662        fMax_planexo = nummod - 1 + b;
663        fMax_planeye = nummod - 1 + c;      
664        if ( debug ) printf(" fMax_plane %f nummod %i numexpl %i a %i b %i c %i \n",fMax_plane,nummod,numexpl,a,b,c);
665        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);
666      };
667      //
668      // 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
669      //
670    //    //
671    Float_t tgx_cl2;    Float_t tgx_cl2;
672    Float_t tgy_cl2;      Float_t tgy_cl2;  
673    tgx_cl2 = L2->GetCaloLevel2()->tanx[0];    tgx_cl2 = cl2->tanx[0];
674    tgy_cl2 = L2->GetCaloLevel2()->tany[0];    tgy_cl2 = cl2->tany[0];
675    //    //
676    for (Int_t p=0; p<22; p++){        for (Int_t p=0; p<22; p++){    
677      track_coordy[p][1] = L2->GetCaloLevel2()->cbar[p][1];      track_coordy[p][1] = cl2->cbar[p][1];
678      track_coordx[p][1] = L2->GetCaloLevel2()->cbar[p][0] - fabs(trk_z[p][1]-trk_z[p][0])*tgx_cl2;      track_coordx[p][1] = cl2->cbar[p][0] - fabs(trk_z[p][1]-trk_z[p][0])*tgx_cl2;
679      track_coordx[p][0] = L2->GetCaloLevel2()->cbar[p][0];      //    track_coordx[p][1] = cl2->cbar[p][0] + fabs(trk_z[p][1]-trk_z[p][0])*tgx_cl2;
680      track_coordy[p][0] = L2->GetCaloLevel2()->cbar[p][1] - fabs(trk_z[p][1]-trk_z[p][0])*tgy_cl2;      track_coordx[p][0] = cl2->cbar[p][0];
681        track_coordy[p][0] = cl2->cbar[p][1] - fabs(trk_z[p][1]-trk_z[p][0])*tgy_cl2;
682        //    track_coordy[p][0] = cl2->cbar[p][1] + fabs(trk_z[p][1]-trk_z[p][0])*tgy_cl2;
683      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]);      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]);
684    };    };
685    //    //
686    if ( debug ) printf(" acceptance fast %i ss %i tgx %f tgy %f\n",fast,ss,tgx_cl2,tgy_cl2);    if ( debug ) printf(" acceptance fNumSec %i tgx %f tgy %f\n",fNumSec,tgx_cl2,tgy_cl2);
687    //    //
688    if ( section.Contains("XO") ){    if ( section.Contains("XO") ){
689        //
690        // check event is inside XO acceptance
691        //
692      for (Int_t i=0; i<11; i++) {      for (Int_t i=0; i<11; i++) {
693        if        if
694          (          (
695           (((track_coordx[2*i+1][1]>=(-12.054+fM))&&           (((track_coordx[(2*i)+1][1]>=(-12.054+fM))&&
696             (track_coordx[2*i+1][1]<=(-4.246-fM)))&&             (track_coordx[(2*i)+1][1]<=(-4.246-fM)))&&
697            (((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo1 + fM1)&&            (((cl2->cbar[(2*i)+1][1]>=xo1 + fM1)&&
698              (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo2 - fM1 ))||              (cl2->cbar[(2*i)+1][1]<=xo2 - fM1 ))||
699             ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo3 + fM1)&&             ((cl2->cbar[(2*i)+1][1]>=xo3 + fM1)&&
700              (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo4 - fM1 ))||              (cl2->cbar[(2*i)+1][1]<=xo4 - fM1 ))||
701             ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo5 + fM1)&&             ((cl2->cbar[(2*i)+1][1]>=xo5 + fM1)&&
702              (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo6 - fM1 ))))||              (cl2->cbar[(2*i)+1][1]<=xo6 - fM1 ))))||
703                    
704           (((track_coordx[2*i+1][1]>=(-4.004+fM))&&           (((track_coordx[(2*i)+1][1]>=(-4.004+fM))&&
705             (track_coordx[2*i+1][1]<=(3.804-fM)))&&             (track_coordx[(2*i)+1][1]<=(3.804-fM)))&&
706            (((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo1 + fM1)&&            (((cl2->cbar[(2*i)+1][1]>=xo1 + fM1)&&
707              (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo2 - fM1 ))||              (cl2->cbar[(2*i)+1][1]<=xo2 - fM1 ))||
708             ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo3 + fM1)&&             ((cl2->cbar[(2*i)+1][1]>=xo3 + fM1)&&
709              (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo4 - fM1))||              (cl2->cbar[(2*i)+1][1]<=xo4 - fM1))||
710             ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo5 + fM1)&&             ((cl2->cbar[(2*i)+1][1]>=xo5 + fM1)&&
711              (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo6 - fM1 ))))||              (cl2->cbar[(2*i)+1][1]<=xo6 - fM1 ))))||
712                    
713           (((track_coordx[2*i+1][1]>=(4.046+fM))&&           (((track_coordx[(2*i)+1][1]>=(4.046+fM))&&
714             (track_coordx[2*i+1][1]<=(11.854-fM)))&&             (track_coordx[(2*i)+1][1]<=(11.854-fM)))&&
715            (((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo1 + fM1)&&            (((cl2->cbar[(2*i)+1][1]>=xo1 + fM1)&&
716              (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo2 - fM1))||              (cl2->cbar[(2*i)+1][1]<=xo2 - fM1))||
717             ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo3 + fM1)&&             ((cl2->cbar[(2*i)+1][1]>=xo3 + fM1)&&
718              (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo4 - fM1 ))||              (cl2->cbar[(2*i)+1][1]<=xo4 - fM1 ))||
719             ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo5 + fM1)&&             ((cl2->cbar[(2*i)+1][1]>=xo5 + fM1)&&
720              (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo6 - fM1 ))))              (cl2->cbar[(2*i)+1][1]<=xo6 - fM1 ))))
721           ){           ){
722          fXosel = true;          fXosel = true;
723            fXoout = i;
724        } else {        } else {
         if ( fast ) return false;          
725          fXosel = false;          fXosel = false;
726          break;          break;
727        };        };
728      };      };
729        //
730        // 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
731        //
732        if ( !fXosel && fXoout >= fXomin && fXoout >= (Int_t)(fMax_planexo+nPl) ){
733          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));
734          fPartsel = true;
735          fXosel = true;
736        };
737        //
738        // event is contained (or partially contained) hence we can integrate energy up to the maximum and calculate the energy as measured by this section
739        //
740        if ( fXosel ){
741          for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planexo+nPl)) ;iplm++){        
742            fXOen_maxplane += en_xop[iplm];
743            if ( debug ) printf(" XO iplm %i fXOen_maxplane %f  en_xop[iplm] %f\n",iplm,fXOen_maxplane,en_xop[iplm]);
744          };
745          fEnergyxo = fXOen_maxplane/fConv_rxo;
746        };
747    };    };
748    //    //
749    if ( section.Contains("XE") ){    if ( section.Contains("XE") ){
# Line 336  Bool_t CaloEnergy::IsInsideAcceptance(TS Line 752  Bool_t CaloEnergy::IsInsideAcceptance(TS
752          (          (
753           (((track_coordx[2*i][1]>=(-11.854+fM))&&           (((track_coordx[2*i][1]>=(-11.854+fM))&&
754             (track_coordx[2*i][1]<=(-4.046-fM)))&&             (track_coordx[2*i][1]<=(-4.046-fM)))&&
755            (((L2->GetCaloLevel2()->cbar[2*i][1]>=xe1 + fM1)&&            (((cl2->cbar[2*i][1]>=xe1 + fM1)&&
756              (L2->GetCaloLevel2()->cbar[2*i][1]<=xe2 - fM1 ))||              (cl2->cbar[2*i][1]<=xe2 - fM1 ))||
757             ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe3 + fM1)&&             ((cl2->cbar[2*i][1]>=xe3 + fM1)&&
758              (L2->GetCaloLevel2()->cbar[2*i][1]<=xe4 - fM1 ))||              (cl2->cbar[2*i][1]<=xe4 - fM1 ))||
759             ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe5 + fM1)&&             ((cl2->cbar[2*i][1]>=xe5 + fM1)&&
760              (L2->GetCaloLevel2()->cbar[2*i][1]<=xe6 - fM1 ))))||              (cl2->cbar[2*i][1]<=xe6 - fM1 ))))||
761                    
762           (((track_coordx[2*i][1]>=(-3.804+fM))&&           (((track_coordx[2*i][1]>=(-3.804+fM))&&
763             (track_coordx[2*i][1]<=(4.004-fM)))&&             (track_coordx[2*i][1]<=(4.004-fM)))&&
764            (((L2->GetCaloLevel2()->cbar[2*i][1]>=xe1 + fM1)&&            (((cl2->cbar[2*i][1]>=xe1 + fM1)&&
765              (L2->GetCaloLevel2()->cbar[2*i][1]<=xe2 - fM1 ))||              (cl2->cbar[2*i][1]<=xe2 - fM1 ))||
766             ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe3 + fM1)&&             ((cl2->cbar[2*i][1]>=xe3 + fM1)&&
767              (L2->GetCaloLevel2()->cbar[2*i][1]<=xe4 - fM1))||              (cl2->cbar[2*i][1]<=xe4 - fM1))||
768             ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe5 + fM1)&&             ((cl2->cbar[2*i][1]>=xe5 + fM1)&&
769              (L2->GetCaloLevel2()->cbar[2*i][1]<=xe6 - fM1 ))))||              (cl2->cbar[2*i][1]<=xe6 - fM1 ))))||
770                    
771           (((track_coordx[2*i][1]>=(4.246+fM))&&           (((track_coordx[2*i][1]>=(4.246+fM))&&
772             (track_coordx[2*i][1]<=(12.054-fM)))&&             (track_coordx[2*i][1]<=(12.054-fM)))&&
773            (((L2->GetCaloLevel2()->cbar[2*i][1]>=xe1 + fM1)&&            (((cl2->cbar[2*i][1]>=xe1 + fM1)&&
774              (L2->GetCaloLevel2()->cbar[2*i][1]<=xe2 - fM1))||              (cl2->cbar[2*i][1]<=xe2 - fM1))||
775             ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe3 + fM1)&&             ((cl2->cbar[2*i][1]>=xe3 + fM1)&&
776              (L2->GetCaloLevel2()->cbar[2*i][1]<=xe4 - fM1 ))||              (cl2->cbar[2*i][1]<=xe4 - fM1 ))||
777             ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe5 + fM1)&&             ((cl2->cbar[2*i][1]>=xe5 + fM1)&&
778              (L2->GetCaloLevel2()->cbar[2*i][1]<=xe6 - fM1 ))))              (cl2->cbar[2*i][1]<=xe6 - fM1 ))))
779           ){           ){
780          fXesel = true;          fXesel = true;
781            fXeout = i;
782        } else {        } else {
         if ( fast ) return false;          
783          fXesel = false;          fXesel = false;
784          break;          break;
785        };        };
786      };      };
787        if ( !fXesel && fXeout >= fXemin && fXeout >= (Int_t)(fMax_planexe+nPl) ){
788          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));
789          fPartsel = true;
790          fXesel = true;
791        };
792        if ( fXesel ){
793          for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planexe+nPl)) ;iplm++){
794            fXEen_maxplane += en_xep[iplm];
795            if ( debug ) printf(" XE iplm %i fXOen_maxplane %f  en_xop[iplm] %f\n",iplm,fXEen_maxplane,en_xep[iplm]);
796          };
797          fEnergyxe = fXEen_maxplane/fConv_rxe;
798        };
799    };      };  
800    //    //
801    if ( section.Contains("YE") ){    if ( section.Contains("YE") ){
802      for (Int_t i=0; i<11; i++) {      for (Int_t i=0; i<11; i++) {
803        if        if
804          (          (
805           (((track_coordy[2*i+1][0]>=(-12.154+fM))&&           (((track_coordy[(2*i)+1][0]>=(-12.154+fM))&&
806             (track_coordy[2*i+1][0]<=(-4.346-fM)))&&             (track_coordy[(2*i)+1][0]<=(-4.346-fM)))&&
807            (((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye1 + fM1)&&            (((cl2->cbar[(2*i)+1][0]>=ye1 + fM1)&&
808              (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye2 - fM1 ))||              (cl2->cbar[(2*i)+1][0]<=ye2 - fM1 ))||
809             ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye3 + fM1)&&             ((cl2->cbar[(2*i)+1][0]>=ye3 + fM1)&&
810              (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye4 - fM1 ))||              (cl2->cbar[(2*i)+1][0]<=ye4 - fM1 ))||
811             ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye5 + fM1)&&             ((cl2->cbar[(2*i)+1][0]>=ye5 + fM1)&&
812              (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye6 - fM1 ))))||              (cl2->cbar[(2*i)+1][0]<=ye6 - fM1 ))))||
813                    
814           (((track_coordy[2*i+1][0]>=(-4.104+fM))&&           (((track_coordy[(2*i)+1][0]>=(-4.104+fM))&&
815             (track_coordy[2*i+1][0]<=(3.704-fM)))&&             (track_coordy[(2*i)+1][0]<=(3.704-fM)))&&
816            (((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye1 + fM1)&&            (((cl2->cbar[(2*i)+1][0]>=ye1 + fM1)&&
817              (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye2 - fM1 ))||              (cl2->cbar[(2*i)+1][0]<=ye2 - fM1 ))||
818             ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye3 + fM1)&&             ((cl2->cbar[(2*i)+1][0]>=ye3 + fM1)&&
819              (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye4 - fM1))||              (cl2->cbar[(2*i)+1][0]<=ye4 - fM1))||
820             ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye5 + fM1)&&             ((cl2->cbar[(2*i)+1][0]>=ye5 + fM1)&&
821              (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye6 - fM1 ))))||              (cl2->cbar[(2*i)+1][0]<=ye6 - fM1 ))))||
822                    
823           (((track_coordy[2*i+1][0]>=(3.946+fM))&&           (((track_coordy[(2*i)+1][0]>=(3.946+fM))&&
824             (track_coordy[2*i+1][0]<=(11.754-fM)))&&             (track_coordy[(2*i)+1][0]<=(11.754-fM)))&&
825            (((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye1 + fM1)&&            (((cl2->cbar[(2*i)+1][0]>=ye1 + fM1)&&
826              (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye2 - fM1))||              (cl2->cbar[(2*i)+1][0]<=ye2 - fM1))||
827             ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye3 + fM1)&&             ((cl2->cbar[(2*i)+1][0]>=ye3 + fM1)&&
828              (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye4 - fM1 ))||              (cl2->cbar[(2*i)+1][0]<=ye4 - fM1 ))||
829             ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye5 + fM1)&&             ((cl2->cbar[(2*i)+1][0]>=ye5 + fM1)&&
830              (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye6 - fM1 ))))              (cl2->cbar[(2*i)+1][0]<=ye6 - fM1 ))))
831           ){           ){
832          fYesel = true;          fYesel = true;
833            fYeout = i;
834        } else {        } else {
         if ( fast ) return false;          
835          fYesel = false;          fYesel = false;
836          break;          break;
837        };        };
838      };      };
839        if ( !fYesel && fYeout >= fYemin && fYeout >= (Int_t)(fMax_planeye+nPl) ){
840          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));
841          fPartsel = true;
842          fYesel = true;
843        };
844        if ( fYesel ){
845          for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planeye+nPl)) ;iplm++) fYEen_maxplane += en_yep[iplm];
846          fEnergyye = fYEen_maxplane/fConv_rye;
847        };
848    };      };  
849    //    //
850    if ( section.Contains("YO") ){    if ( section.Contains("YO") ){
851      for (Int_t i=0; i<11; i++) {      for (Int_t i=0; i<11; i++) {
852        if ( debug ) printf(" i %i track_coordy[2*i][0] %f L2->GetCaloLevel2()->cbar[2*i][0] %f \n",i,track_coordy[2*i][0],L2->GetCaloLevel2()->cbar[2*i][0]);        if ( debug ) printf(" i %i track_coordy[2*i][0] %f cl2->cbar[2*i][0] %f \n",i,track_coordy[2*i][0],cl2->cbar[2*i][0]);
853        if ( debug ) printf(" i %i fm %f fm1 %f yo1 %g yo2 %f yo3 %f yo4 %f yo5 %f yo6 %f \n",i,fM,fM1,yo1,yo2,yo3,yo4,yo5,yo6);        if ( debug ) printf(" i %i fm %f fm1 %f yo1 %g yo2 %f yo3 %f yo4 %f yo5 %f yo6 %f \n",i,fM,fM1,yo1,yo2,yo3,yo4,yo5,yo6);
854        if        if
855          (          (
856           (((track_coordy[2*i][0]>=(-11.954+fM))&&           (((track_coordy[2*i][0]>=(-11.954+fM))&&
857             (track_coordy[2*i][0]<=(-4.146-fM)))&&             (track_coordy[2*i][0]<=(-4.146-fM)))&&
858            (((L2->GetCaloLevel2()->cbar[2*i][0]>=yo1 + fM1)&&            (((cl2->cbar[2*i][0]>=yo1 + fM1)&&
859              (L2->GetCaloLevel2()->cbar[2*i][0]<=yo2 - fM1 ))||              (cl2->cbar[2*i][0]<=yo2 - fM1 ))||
860             ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo3 + fM1)&&             ((cl2->cbar[2*i][0]>=yo3 + fM1)&&
861              (L2->GetCaloLevel2()->cbar[2*i][0]<=yo4 - fM1 ))||              (cl2->cbar[2*i][0]<=yo4 - fM1 ))||
862             ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo5 + fM1)&&             ((cl2->cbar[2*i][0]>=yo5 + fM1)&&
863              (L2->GetCaloLevel2()->cbar[2*i][0]<=yo6 - fM1 ))))||              (cl2->cbar[2*i][0]<=yo6 - fM1 ))))||
864                    
865           (((track_coordy[2*i][0]>=(-3.904+fM))&&           (((track_coordy[2*i][0]>=(-3.904+fM))&&
866             (track_coordy[2*i][0]<=(+3.904-fM)))&&             (track_coordy[2*i][0]<=(+3.904-fM)))&&
867            (((L2->GetCaloLevel2()->cbar[2*i][0]>=yo1 + fM1)&&            (((cl2->cbar[2*i][0]>=yo1 + fM1)&&
868              (L2->GetCaloLevel2()->cbar[2*i][0]<=yo2 - fM1 ))||              (cl2->cbar[2*i][0]<=yo2 - fM1 ))||
869             ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo3 + fM1)&&             ((cl2->cbar[2*i][0]>=yo3 + fM1)&&
870              (L2->GetCaloLevel2()->cbar[2*i][0]<=yo4 - fM1))||              (cl2->cbar[2*i][0]<=yo4 - fM1))||
871             ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo5 + fM1)&&             ((cl2->cbar[2*i][0]>=yo5 + fM1)&&
872              (L2->GetCaloLevel2()->cbar[2*i][0]<=yo6 - fM1 ))))||              (cl2->cbar[2*i][0]<=yo6 - fM1 ))))||
873                    
874           (((track_coordy[2*i][0]>=(4.146+fM))&&           (((track_coordy[2*i][0]>=(4.146+fM))&&
875             (track_coordy[2*i][0]<=(11.954-fM)))&&             (track_coordy[2*i][0]<=(11.954-fM)))&&
876            (((L2->GetCaloLevel2()->cbar[2*i][0]>=yo1 + fM1)&&            (((cl2->cbar[2*i][0]>=yo1 + fM1)&&
877              (L2->GetCaloLevel2()->cbar[2*i][0]<=yo2 - fM1))||              (cl2->cbar[2*i][0]<=yo2 - fM1))||
878             ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo3 + fM1)&&             ((cl2->cbar[2*i][0]>=yo3 + fM1)&&
879              (L2->GetCaloLevel2()->cbar[2*i][0]<=yo4 - fM1 ))||              (cl2->cbar[2*i][0]<=yo4 - fM1 ))||
880             ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo5 + fM1)&&             ((cl2->cbar[2*i][0]>=yo5 + fM1)&&
881              (L2->GetCaloLevel2()->cbar[2*i][0]<=yo6 - fM1 ))))              (cl2->cbar[2*i][0]<=yo6 - fM1 ))))
882           ){           ){
883          fYosel = true;          fYosel = true;
884            fYoout = i;
885        } else {        } else {
         if ( fast ) return false;          
886          fYosel = false;          fYosel = false;
887          break;          break;
888        };        };
889      };      };
890        if ( !fYosel && fYoout >= fYomin && fYoout >= (Int_t)(fMax_planeyo+nPl) ){
891          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));
892          fPartsel = true;
893          fYosel = true;
894        };
895        if ( fYosel ){
896          for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planeyo+nPl)) ;iplm++) fYOen_maxplane += en_yop[iplm];
897          fEnergyyo = fYOen_maxplane/fConv_ryo;
898        };
899    };      };  
900    //    //
901    fSel = ( fXesel || fYesel || fXosel || fYosel );    // Count the number of sections in which the event is contained
902      //
903    fCount = (Float_t)((Int_t)fXesel+(Int_t)fXosel+(Int_t)fYesel+(Int_t)fYosel);    fCount = (Float_t)((Int_t)fXesel+(Int_t)fXosel+(Int_t)fYesel+(Int_t)fYosel);
904      //
905      if ( indep ){
906        //
907        // independent mode, average the energy measurement and max plane of the contained sections
908        //
909        fSel = ( fXesel || fYesel || fXosel || fYosel );
910        fMax_plane = (Float_t)(fMax_planeyo+fMax_planeye+fMax_planexo+fMax_planexe)/fCount;
911        fEnergy = (fEnergyxe+fEnergyyo+fEnergyye+fEnergyxo)/fCount;
912        //
913      } else {
914        //
915        // 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
916        //
917        if ( fCount != fNumSec ){
918          fSel = false;
919        } else {
920          fSel = true;
921        };
922        fEnergy = (fXEen_maxplane+fYOen_maxplane+fYEen_maxplane+fXOen_maxplane)/fConv_rxo;
923      };
924      //
925      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);
926    if ( debug ) printf(" IsInside XE %i XO %i YE %i YO %i => SEL %i \n",fXesel,fXosel,fYesel,fYosel,fSel);    if ( debug ) printf(" IsInside XE %i XO %i YE %i YO %i => SEL %i \n",fXesel,fXosel,fYesel,fYosel,fSel);
927    //    //
928      // finish exit
929      //
930    return fSel;    return fSel;
931    //    //
932  }  }
# Line 468  void CaloEnergy::Process(){ Line 939  void CaloEnergy::Process(){
939    
940  void CaloEnergy::Process(TString section){  void CaloEnergy::Process(TString section){
941    //      //  
942      // process the event
943      //
944    TString ntr = section;    TString ntr = section;
945    if ( !L2 ){    if ( !L2 ){
946      printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n");      printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n");
# Line 489  void CaloEnergy::Process(TString section Line 962  void CaloEnergy::Process(TString section
962      newentry = true;      newentry = true;
963    };    };
964    //    //
965      // if we have already called this method for this event and no input changed then return fSel and exit
966      //
967    if ( !newentry ) return;    if ( !newentry ) return;
968    //    //
969      // process the event
970      //
971    if ( debug ) printf(" Processing event at OBT %u PKT %u time %u section %s\n",OBT,PKT,atime,section.Data());    if ( debug ) printf(" Processing event at OBT %u PKT %u time %u section %s\n",OBT,PKT,atime,section.Data());
972    //    //
973    Clear();    // 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)
974      //
975      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));
976    //    //
977    if ( fLong ){    if ( fLong ){
978      // use long fit to measure energy      if ( debug ) printf(" ==================================================================> LONGITUDINAL FIT! \n");
979    } else {      //
980        // use long fit to measure energy
981      //      //
982      if ( this->IsInsideAcceptance(section) ){      if ( this->IsInsideAcceptance(section) ){
983        //        //
984        for (Int_t ch=0; ch< L2->GetCaloLevel1()->istrip; ch++){        if ( debug ) printf(" ==================================================================> LONG INSIDE! \n");
985          en = L2->GetCaloLevel1()->DecodeEstrip(ch,view,plane,strip);        //
986          enstrip[view][plane][strip]=en;        Float_t myene[2][22];
987        };        memset(myene,0,(sizeof(Float_t))*2*22);
988        for (Int_t i=0;i<11;i++){        for (Int_t j=0; j<11; j++){
989          for(strip=0; strip<96; strip++) {          if ( section.Contains("XE") ) myene[1][2*j] = en_xep[j];
990            en_xep[i] += enstrip[1][2*i][strip];          if ( section.Contains("YO") ) myene[0][2*j] = en_yop[j];
991            en_yop[i] += enstrip[0][2*i][strip];          if ( section.Contains("XO") ) myene[1][(2*j)+1] = en_xop[j];
992            en_xop[i] += enstrip[1][2*i+1][strip];          if ( section.Contains("YE") ) myene[0][(2*j)+1] = en_yep[j];    
           en_yep[i] += enstrip[0][2*i+1][strip];  
         };  
         energyxe += en_xep[i];  
         energyyo += en_yop[i];  
         energyxo += en_xop[i];  
         energyye += en_yep[i];  
       };  
       fMax_plane = 0;  
       Float_t en_maxplane = 0.;    
       Float_t max_en= 0.;  
       if ( fXosel ){  
         max_en = 0.;  
         for (Int_t ipl =0; ipl < 11; ipl ++) {  
           if(en_xop[ipl] > max_en) {  
             max_en = en_xop[ipl];  
             fMax_plane = ipl;  
           };  
         };  
         //  
         for (Int_t iplm=0;iplm<=(fMax_plane+fPl) ;iplm++) en_maxplane += en_xop[iplm];  
       };  
       if ( fXesel ){  
         max_en = 0.;  
         for (Int_t ipl =0; ipl < 11; ipl ++) {  
           if(en_xep[ipl] > max_en) {  
             max_en = en_xep[ipl];  
             fMax_plane = ipl;  
           };  
         };  
         //  
         for (Int_t iplm=0;iplm<=(fMax_plane+fPl) ;iplm++) en_maxplane += en_xep[iplm];  
       };  
       if ( fYosel ){  
         max_en = 0.;  
         for (Int_t ipl =0; ipl < 11; ipl ++) {  
           if(en_yop[ipl] > max_en) {  
             max_en = en_yop[ipl];  
             fMax_plane = ipl;  
           };  
         };  
         //  
         for (Int_t iplm=0;iplm<=(fMax_plane+fPl) ;iplm++) en_maxplane += en_yop[iplm];  
993        };        };
994        if ( fYesel ){        clong->UnMaskSections();
995          max_en = 0.;        if ( !(section.Contains("YE")) ) clong->MaskSection("YE");
996          for (Int_t ipl =0; ipl < 11; ipl ++) {        if ( !(section.Contains("YO")) ) clong->MaskSection("YO");
997            if(en_yep[ipl] > max_en) {        if ( !(section.Contains("XO")) ) clong->MaskSection("XO");
998              max_en = en_yep[ipl];        if ( !(section.Contains("XE")) ) clong->MaskSection("XE");
999              fMax_plane = ipl;        clong->ForceNextFit();
1000            };        clong->SetEnergies(myene);
1001          };        if ( debug ){
1002          //          clong->Fit(true);
1003          for (Int_t iplm=0;iplm<=(fMax_plane+fPl) ;iplm++) en_maxplane += en_yep[iplm];        } else {
1004            clong->Fit();
1005          };      
1006          if ( clong->GetLowerLimit() != 0. || clong->GetUpperLimit() != 0. ){
1007            fXOen_maxplane = clong->Get_defE0();
1008          } else {
1009            fXOen_maxplane = clong->Get_E0();
1010        };        };
1011          fMax_plane = clong->Get_tmax();
1012          fYOen_maxplane  = 0.;
1013          fYEen_maxplane = 0.;
1014          fXEen_maxplane = 0.;
1015          fEnergy=fXOen_maxplane/fConv_rxo;
1016          if ( fEnergy != fEnergy || clong->Get_fitresult() != 0 ) fEnergy = -1.;
1017          //      if ( fEnergy != fEnergy ) fEnergy = 1.;
1018          //
1019        } else {
1020          //
1021          // if the event is not in the acceptance, return a negative energy.
1022          //
1023          if ( debug ) printf(" Outside acceptance \n");
1024          fEnergy *= -1.;
1025          //
1026        };
1027        //
1028      } else {
1029        //
1030        // use the energy measurement
1031        //
1032        if ( this->IsInsideAcceptance(section) ){
1033        //        //
1034        fEnergy = en_maxplane/fCount/fConv_r;        // the event is good
1035        //        //
       if ( debug ) printf(" fMax_plane %i conv_r %f en_maxplane %f encalo %f  \n",fMax_plane,fConv_r,en_maxplane,fEnergy);  
1036        if ( debug ) printf(" XE %i XO %i YE %i YO %i \n",fXesel,fXosel,fYesel,fYosel);        if ( debug ) printf(" XE %i XO %i YE %i YO %i \n",fXesel,fXosel,fYesel,fYosel);
1037        //        //
1038      } else {      } else {
1039          //
1040          // if the event is not in the acceptance, return a negative energy.
1041          //
1042        if ( debug ) printf(" Outside acceptance \n");        if ( debug ) printf(" Outside acceptance \n");
1043        fEnergy = -1.;        fEnergy *= -1.;
1044          //
1045      };      };
1046    };    };
1047    //    //
Legend:
Removed from v.1.1.1.1  
changed lines
  Added in v.1.12
  ViewVC Help
Powered by ViewVC 1.1.23