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