/[PAMELA software]/DarthVader/ToFLevel2/src/ToFLevel2.cpp
ViewVC logotype

Diff of /DarthVader/ToFLevel2/src/ToFLevel2.cpp

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

revision 1.3 by mocchiut, Fri Jun 30 09:22:03 2006 UTC revision 1.34 by mocchiut, Wed Nov 23 21:19:38 2011 UTC
# Line 1  Line 1 
1  #include <TObject.h>  /**
2     * \file ToFLevel2.cpp
3     * \author Gianfranca DeRosa, Wolfgang Menn
4     *
5     * WM dec 2008: Description of "GetdEdx" changed
6     * WM dec 2008: "GetdEdxPaddle" modified: Now includes saturation limit
7     *              PMTs higher than the saturation limit are not used for dEdx
8     * WM apr 2009: bug found by Nicola in method "GetPaddlePlane"
9     */
10    
11  #include <ToFLevel2.h>  #include <ToFLevel2.h>
 #include <iostream>  
12  using namespace std;  using namespace std;
13    ClassImp(ToFPMT);
14    ClassImp(ToFdEdx);
15    ClassImp(ToFGeom);
16  ClassImp(ToFTrkVar);  ClassImp(ToFTrkVar);
17  ClassImp(ToFLevel2);  ClassImp(ToFLevel2);
18    
19  ToFTrkVar::ToFTrkVar() {  ToFPMT::ToFPMT(){
20      pmt_id = 0;
21    trkseqno = 0;    adc = 0.;
22      tdc_tw = 0.;
23    for (Int_t kk=0; kk<13;kk++){    tdc = 0.;
24      beta_a[kk] = 0;  }
   }  
25    
26    for (Int_t kk=0; kk<4;kk++){  ToFPMT::ToFPMT(const ToFPMT &t){
27      for (Int_t hh=0; hh<12;hh++){    pmt_id = t.pmt_id;
28        adc_c[hh][kk] = 0;    adc = t.adc;
29      }    tdc_tw = t.tdc_tw;
30    }    tdc = t.tdc;
31  }  }
32    
33  ToFTrkVar::ToFTrkVar(const ToFTrkVar &t){  void ToFPMT::Clear(Option_t *t){
34      pmt_id = 0;
35      adc = 0.;
36      tdc_tw = 0.;
37      tdc = 0.;
38    }
39    
   trkseqno = t.trkseqno;    
40    
   memcpy(adc_c,t.adc_c,sizeof(adc_c));  
   memcpy(beta_a,t.beta_a,sizeof(beta_a));  
 }  
41    
42  ToFLevel2::ToFLevel2() {      ToFTrkVar::ToFTrkVar() {
43      trkseqno = 0;
44      npmttdc = 0;
45      npmtadc = 0;
46      pmttdc = TArrayI(48);
47      pmtadc = TArrayI(48);
48      tdcflag = TArrayI(48); // gf: 30 Nov 2006
49      adcflag = TArrayI(48); // gf: 30 Nov 2006
50      dedx = TArrayF(48);
51    //    //
   ToFTrk = new TClonesArray("ToFTrkVar",1);  
52    //    //
53    for (Int_t kk=0; kk<3;kk++){    memset(beta,  0, 13*sizeof(Float_t));
54      xtofpos[kk] = 0.;    memset(xtofpos,  0, 3*sizeof(Float_t));
55      ytofpos[kk] = 0.;    memset(ytofpos,  0, 3*sizeof(Float_t));
56    }    memset(xtr_tof,  0, 6*sizeof(Float_t));
57      memset(ytr_tof,  0, 6*sizeof(Float_t));
58    for (Int_t kk=0; kk<6;kk++){    //
59      tof_i_flag[kk] = 0;  };
     tof_j_flag[kk] = 0;  
   }  
     
   for (Int_t kk=0; kk<13;kk++){  
     betatof_a[kk] = 0;  
   }  
60    
61    for (Int_t kk=0; kk<4;kk++){  void ToFTrkVar::Clear(Option_t *t) {
62      for (Int_t hh=0; hh<12;hh++){    trkseqno = 0;
63        adctof_c[hh][kk] = 0;    npmttdc = 0;
64      }    npmtadc = 0;
65    }      pmttdc.Reset();
66      pmtadc.Reset();
67      tdcflag.Reset(); // gf: 30 Nov 2006
68      adcflag.Reset(); // gf: 30 Nov 2006
69      dedx.Reset();
70      //
71      memset(beta,  0, 13*sizeof(Float_t));
72      memset(xtofpos,  0, 3*sizeof(Float_t));
73      memset(ytofpos,  0, 3*sizeof(Float_t));
74      memset(xtr_tof,  0, 6*sizeof(Float_t));
75      memset(ytr_tof,  0, 6*sizeof(Float_t));
76      //
77    };
78    
79    for (Int_t kk=0; kk<4;kk++){  ToFTrkVar::ToFTrkVar(const ToFTrkVar &t){
     for (Int_t hh=0; hh<12;hh++){  
       tdc_c[hh][kk] = 0;  
     }  
   }    
80    
81      trkseqno = t.trkseqno;  
82      //
83      npmttdc = t.npmttdc;
84      npmtadc = t.npmtadc;
85      (t.pmttdc).Copy(pmttdc);
86      (t.pmtadc).Copy(pmtadc);
87      (t.tdcflag).Copy(tdcflag); // gf: 30 Nov 2006
88      (t.adcflag).Copy(adcflag); // gf: 30 Nov 2006
89      (t.dedx).Copy(dedx);
90      //
91      memcpy(beta,t.beta,sizeof(beta));
92      memcpy(xtofpos,t.xtofpos,sizeof(xtofpos));
93      memcpy(ytofpos,t.ytofpos,sizeof(ytofpos));
94      memcpy(xtr_tof,t.xtr_tof,sizeof(xtr_tof));
95      memcpy(ytr_tof,t.ytr_tof,sizeof(ytr_tof));
96      //
97  };  };
98    
99  void ToFLevel2::Clear(){  ToFLevel2::ToFLevel2() {    
100    //    //
101    ToFTrk->Clear();  //  PMT = new TClonesArray("ToFPMT",12); //ELENA
102    //  ToFTrk = new TClonesArray("ToFTrkVar",2); //ELENA
103      PMT = 0; //ELENA
104      ToFTrk = 0; //ELENA
105    //    //
106    for (Int_t kk=0; kk<3;kk++){    this->Clear();
107      xtofpos[kk] = 0.;    //
108      ytofpos[kk] = 0.;  };
   }  
   
   for (Int_t kk=0; kk<6;kk++){  
     tof_i_flag[kk] = 0;  
     tof_j_flag[kk] = 0;  
   }  
     
   for (Int_t kk=0; kk<13;kk++){  
     betatof_a[kk] = 0;  
   }  
109    
110    for (Int_t kk=0; kk<4;kk++){  void ToFLevel2::Set(){//ELENA
111      for (Int_t hh=0; hh<12;hh++){      if(!PMT)PMT = new TClonesArray("ToFPMT",12); //ELENA
112        adctof_c[hh][kk] = 0;      if(!ToFTrk)ToFTrk = new TClonesArray("ToFTrkVar",2); //ELENA
113      }  }//ELENA
   }    
114    
115    for (Int_t kk=0; kk<4;kk++){  void ToFLevel2::Clear(Option_t *t){
116      for (Int_t hh=0; hh<12;hh++){    //
117        tdc_c[hh][kk] = 0;    if(ToFTrk)ToFTrk->Delete(); //ELENA
118      }    if(PMT)PMT->Delete(); //ELENA
119    }      memset(tof_j_flag, 0, 6*sizeof(Int_t));
120      unpackError = 0;
121      //
122  };  };
123    
124    void ToFLevel2::Delete(Option_t *t){ //ELENA
125      //
126      if(ToFTrk){
127          ToFTrk->Delete(); //ELENA
128          delete ToFTrk;  //ELENA
129      }
130      if(PMT){
131          PMT->Delete(); //ELENA
132          delete PMT; //ELENA
133      } //ELENA
134      //
135    }; //ELENA
136    
137  ToFTrkVar *ToFLevel2::GetToFTrkVar(Int_t itrk){  ToFTrkVar *ToFLevel2::GetToFTrkVar(Int_t itrk){
138    //        //    
139    if(itrk >= ntrk()){    if(itrk >= ntrk()){
# Line 99  ToFTrkVar *ToFLevel2::GetToFTrkVar(Int_t Line 142  ToFTrkVar *ToFLevel2::GetToFTrkVar(Int_t
142      return(NULL);      return(NULL);
143    }      }  
144    //    //
145      if(!ToFTrk)return 0; //ELENA
146    TClonesArray &t = *(ToFTrk);    TClonesArray &t = *(ToFTrk);
147    ToFTrkVar *toftrack = (ToFTrkVar*)t[itrk];    ToFTrkVar *toftrack = (ToFTrkVar*)t[itrk];
148    return toftrack;    return toftrack;
149  }  }
150    
151    ToFPMT *ToFLevel2::GetToFPMT(Int_t ihit){
152      //    
153      if(ihit >= npmt()){
154        printf(" ToFLevel2 ERROR: pmt variables set %i does not exists! \n",ihit);
155        printf("                  stored pmt variables = %i \n",npmt());
156        return(NULL);
157      }  
158      //
159      if(!PMT)return 0; //ELENA
160      TClonesArray &t = *(PMT);
161      ToFPMT *tofpmt = (ToFPMT*)t[ihit];
162      return tofpmt;
163    }
164  //--------------------------------------  //--------------------------------------
165  //  //
166  //  //
167  //--------------------------------------  //--------------------------------------
168  /**  /**
169   * Method to get the paddle ID (11 12 21 22 31 32) from the paddle index (0 1 2 3 4 5)   * Method to get the plane ID (11 12 21 22 31 32) from the plane index (0 1 2 3 4 5)
170     * @param Plane index (0,1,2,3,4,5).
171   */   */
172    Int_t  ToFLevel2::GetToFPlaneID(Int_t ip){    Int_t  ToFLevel2::GetToFPlaneID(Int_t ip){
173        if(ip>=0 && ip<6)return 10*((int)(ip/2+1.1))+(ip%2)+1;        if(ip>=0 && ip<6)return 10*((int)(ip/2+1.1))+(ip%2)+1;
174        else return -1;        else return -1;
175    };    };
176  /**  /**
177   * Method to get the paddle index (0 1 2 3 4 5) from the paddle ID (11 12 21 22 31 32)   * Method to get the plane index (0 1 2 3 4 5) from the plane ID (11 12 21 22 31 32)
178     * @param plane Plane ID (11, 12, 21, 22, 31, 32)
179   */   */
180    Int_t  ToFLevel2::GetToFPlaneIndex(Int_t plane_id){    Int_t  ToFLevel2::GetToFPlaneIndex(Int_t plane_id){
181        if(        if(
# Line 129  ToFTrkVar *ToFLevel2::GetToFTrkVar(Int_t Line 189  ToFTrkVar *ToFLevel2::GetToFTrkVar(Int_t
189        else return -1;        else return -1;
190    };    };
191  /**  /**
192   * Method to know if a given ToF paddle was hit, that is there is a TDC signal from both PMTs   * Method to know if a given ToF paddle was hit, that is there is a TDC signal
193     * from both PMTs. The method uses the "tof_j_flag" variable.
194   * @param plane Plane ID (11, 12, 21, 22, 31, 32) or Plane index (0,1,2,3,4,5).   * @param plane Plane ID (11, 12, 21, 22, 31, 32) or Plane index (0,1,2,3,4,5).
195   * @param paddle_id Paddle ID.   * @param paddle_id Paddle ID.
196   * @return 1 if the paddle was hit.   * @return 1 if the paddle was hit.
# Line 159  Int_t ToFLevel2::GetNHitPaddles(Int_t pl Line 220  Int_t ToFLevel2::GetNHitPaddles(Int_t pl
220      for(Int_t i=0; i<8; i++)npad = npad + (int)HitPaddle(plane,i);      for(Int_t i=0; i<8; i++)npad = npad + (int)HitPaddle(plane,i);
221      return npad;      return npad;
222  };  };
223    
224    //wm Nov 08
225    //gf Apr 07
226    /**
227     * Method to get the mean dEdx from a ToF layer - ATTENTION:
228     * It will sum up the dEdx of all the paddles, but since by definition
229     * only the paddle hitted by the track gets a dEdx value and the other
230     * paddles are set to zero, the output is just the dEdx of the hitted
231     * paddle in each layer!
232     * The "adcfl" option is not very useful (an artificial dEdx is per
233     * definition= 1 mip and not a real measurement), anyway left in the code
234     * @param notrack Track Number
235     * @param plane Plane index (0,1,2,3,4,5)
236     * @param adcflag in the plane (100<-> independent of the adcflag; !=0&&!=100 <-> at least one PMT with adcflag!=0; )
237     */
238    Float_t ToFLevel2::GetdEdx(Int_t notrack, Int_t plane, Int_t adcfl){
239    
240      Float_t dedx = 0.;
241      Float_t PadEdx =0.;
242      Int_t SatWarning;
243      Int_t pad=-1;
244      //
245      ToFTrkVar *trk = GetToFTrkVar(notrack);
246      if(!trk) return 0; //ELENA
247      //
248      for (Int_t ii=0; ii<GetNPaddle(plane); ii++){
249        Int_t paddleid=ii;
250        pad = GetPaddleid(plane,paddleid);
251        GetdEdxPaddle(notrack, pad, adcfl, PadEdx, SatWarning);
252        dedx += PadEdx;
253      };
254      //
255      return(dedx);
256    };
257    
258    /**
259     * Method to fill the ADC_C 4x12 matrix with the dEdx values and the TDC 4x12 matrix
260     * with the time-walk corrected TDC values.
261     * @param notrack Track Number
262     * @param adc  ADC_C matrix with dEdx values
263     * @param tdc  TDC matrix
264     */
265    void ToFLevel2::GetMatrix(Int_t notrack, Float_t adc[4][12], Float_t tdc[4][12]){
266      //
267      for (Int_t aa=0; aa<4;aa++){
268        for (Int_t bb=0; bb<12;bb++){
269          adc[aa][bb] = 1000.;
270          tdc[aa][bb] = 4095.;
271        };
272      };
273      //
274      Int_t pmt_id = 0;
275      Int_t hh = 0;
276      Int_t kk = 0;
277      //
278      ToFTrkVar *trk = GetToFTrkVar(notrack);
279      if(!trk)return; //ELENA
280      //
281      for (Int_t i=0; i<trk->npmtadc; i++){
282        //
283        pmt_id = (trk->pmtadc).At(i);
284        //
285        GetPMTIndex(pmt_id,hh,kk);
286        adc[kk][hh] = (trk->dedx).At(i);  
287        //
288      };
289      //
290      for (Int_t i=0; i<npmt(); i++){
291        //
292        ToFPMT *pmt = GetToFPMT(i);
293        if(!pmt)break; //ELENA
294        //
295        GetPMTIndex(pmt->pmt_id,hh,kk);
296        //
297        tdc[kk][hh] = pmt->tdc_tw;  
298        //
299      };
300      //
301      return;
302    };
303    
304    
305    /**
306     * Method to get the plane index (0 - 5) for the PMT_ID as input
307     * @param pmt_id  PMT_ID (0 - 47)
308     */
309    Int_t ToFLevel2::GetPlaneIndex(Int_t pmt_id){
310      TString pmtname = GetPMTName(pmt_id);
311      pmtname.Resize(3);
312      if ( !strcmp(pmtname,"S11") ) return(0);
313      if ( !strcmp(pmtname,"S12") ) return(1);
314      if ( !strcmp(pmtname,"S21") ) return(2);
315      if ( !strcmp(pmtname,"S22") ) return(3);
316      if ( !strcmp(pmtname,"S31") ) return(4);
317      if ( !strcmp(pmtname,"S32") ) return(5);
318      return(-1);
319    };
320    
321    
322    /**
323     * Method to get the PMT_ID if the index (4,12) is given. We have 4 channels on
324     * each of the 12 half-boards, this method decodes which PMT is cables to which
325     * channel.
326     * @param hh Channel
327     * @param kk HalfBoard
328     */
329    Int_t ToFLevel2::GetPMTid(Int_t hh, Int_t kk){
330      //
331      short tof[4][24] = {
332        {4, 4,  4,  4,  1,  1, 2, 2,  3,  3, 3, 3,  3,  3, 1, 1,  1,  1, 2, 3,  3, 3, 3,  4},
333        {1, 3,  5,  7, 10, 12, 2, 4,  2,  4, 6, 8, 10, 12, 1, 5,  3,  9, 7, 9, 11, 1, 5,  9},
334        {2, 2,  2,  2,  1,  1, 1, 1,  4,  4, 4, 4,  4,  4, 2, 1,  2,  1, 2, 2,  2, 3, 3,  4},
335        {6, 8, 12, 10,  8,  6, 4, 2, 12, 10, 8, 6,  4,  2, 9, 7, 11, 11, 5, 3,  1, 3, 7, 11}
336      };
337      //
338      Int_t ind = 0;
339      Int_t k = 0;
340      while (k < 24){
341        Int_t j = 0;
342        while (j < 2){
343          Int_t ch = tof[2*j][k]     - 1;
344          Int_t hb = tof[2*j + 1][k] - 1;      
345          /* tofEvent->tdc[ch][hb] */      
346          if( ch == hh && hb == kk ){      
347            ind = 2*k + j;
348            break;
349          };
350          j++;
351        };
352        k++;
353      };
354      return ind;
355    };
356    
357    
358    /**
359     * Method to get the PMT index if the PMT ID is given. This method is the
360     * "reverse" of method "GetPMTid"
361     * @param ind  PMT_ID (0 - 47)
362     * @param hb   HalfBoard
363     * @param ch   Channel
364     */
365    void ToFLevel2::GetPMTIndex(Int_t ind, Int_t &hb, Int_t &ch){
366      //
367      short tof[4][24] = {
368        {4, 4,  4,  4,  1,  1, 2, 2,  3,  3, 3, 3,  3,  3, 1, 1,  1,  1, 2, 3,  3, 3, 3,  4},
369        {1, 3,  5,  7, 10, 12, 2, 4,  2,  4, 6, 8, 10, 12, 1, 5,  3,  9, 7, 9, 11, 1, 5,  9},
370        {2, 2,  2,  2,  1,  1, 1, 1,  4,  4, 4, 4,  4,  4, 2, 1,  2,  1, 2, 2,  2, 3, 3,  4},
371        {6, 8, 12, 10,  8,  6, 4, 2, 12, 10, 8, 6,  4,  2, 9, 7, 11, 11, 5, 3,  1, 3, 7, 11}
372      };
373      //
374      Int_t k = 0;
375      while (k < 24){
376        Int_t j = 0;
377        while (j < 2){
378          /* tofEvent->tdc[ch][hb] */            
379          if( ind == 2*k + j ){
380            ch = tof[2*j][k]     - 1;
381            hb = tof[2*j + 1][k] - 1;      
382            return;
383          };
384          j++;
385        };
386        k++;
387      };
388      return;
389    };
390    
391    
392    
393    //  wm Nov 08 revision - saturation values included
394    /// gf Apr 07
395    /**
396     * Method to get the dEdx from a given ToF paddle.
397     * If two PMTs are good, the mean dEdx of both PMTs is taken, otherwise
398     * just the dEdx of the "good" PMT. If both PMTs are above saturation => dEdx=1000
399     * @param notrack Track Number
400     * @param Paddle index (0,1,...,23).
401     * @param adcflag in the paddle (100<-> independent of the adcflag; !=0&&!=100 <-> at least one PMT with adcflag!=0; )
402     * @param PadEdx dEdx from a given ToF paddle
403     * @param SatWarning 1 if the PMT ios near saturation region (adcraw ~3000)
404     */
405    void ToFLevel2::GetdEdxPaddle(Int_t notrack, Int_t paddleid, Int_t adcfl, Float_t &PadEdx, Int_t &SatWarning){
406    
407    /*
408    Float_t  PMTsat[48] = {
409    3162.14, 3165.48, 3153.85, 3085.73, 3089.65, 3107.64, 3097.52, 3078.37,
410    3130.05, 3087.07, 3112.22, 3102.92, 3080.58, 3092.55, 3087.94, 3125.03,
411    3094.09, 3143.16, 3125.51, 3181.27, 3092.09, 3124.98, 3069.3, 3095.53,
412    3097.11, 3133.53, 3114.73, 3113.01, 3091.19, 3097.99, 3033.84, 3134.98,
413    3081.37, 3111.04, 3066.77, 3108.17, 3133, 3111.06, 3052.52, 3140.66,
414    3106.33, 3094.85, 3150.85, 3118.8, 3096.24, 3118.47,3111.36, 3117.11 } ;
415    */
416    
417    // new values from Napoli dec 2008
418    Float_t  PMTsat[48] = {
419    3176.35,3178.19,3167.38,3099.73,3117.00,3126.29,3111.44,3092.27,
420    3146.48,3094.41,3132.13,3115.37,3099.32,3110.97,3111.80,3143.14,
421    3106.72,3153.44,3136.00,3188.96,3104.73,3140.45,3073.18,3106.62,
422    3112.48,3146.92,3127.24,3136.52,3109.59,3112.89,3045.15,3147.26,
423    3095.92,3121.05,3083.25,3123.62,3150.92,3125.30,3067.60,3160.18,
424    3119.36,3108.92,3164.77,3133.64,3111.47,3131.98,3128.87,3135.56 };
425    
426    for (Int_t i=0; i<48;i++) PMTsat[i] = PMTsat[i] - 5.;  // safety margin
427    
428    
429      PadEdx = 0.;
430    //  SatWarning = 1000;
431      SatWarning = 0;   // 0=good, increase for each bad PMT
432    
433      Float_t dEdx[48] = {0};
434      Int_t pmt_id = -1;
435      Float_t adcraw[48];
436      //
437      ToFTrkVar *trk = GetToFTrkVar(notrack);
438      if(!trk) return; //ELENA
439      //
440    
441      Int_t pmtleft=-1;
442      Int_t pmtright=-1;
443      GetPaddlePMT(paddleid, pmtleft, pmtright);
444    
445      adcraw[pmtleft] = 4095;
446      adcraw[pmtright] = 4095;
447    
448      
449      for (Int_t jj=0; jj<npmt(); jj++){
450        
451        ToFPMT *pmt = GetToFPMT(jj);
452        if(!pmt)break; //ELENA
453        
454        pmt_id = pmt->pmt_id;
455        if(pmt_id==pmtleft){
456          adcraw[pmtleft] = pmt->adc;
457        }
458        
459        if(pmt_id==pmtright){
460          adcraw[pmtright] = pmt->adc;
461        }
462      }
463    
464      
465      for (Int_t i=0; i<trk->npmtadc; i++){
466    
467        if((trk->adcflag).At(i)==0 || adcfl==100){
468          if((trk->pmtadc).At(i) == pmtleft)dEdx[pmtleft] = (trk->dedx).At(i);
469          if((trk->pmtadc).At(i) == pmtright)dEdx[pmtright] = (trk->dedx).At(i);
470        }else{
471          if((trk->pmtadc).At(i) == pmtleft)dEdx[pmtleft] = 0.;
472          if((trk->pmtadc).At(i) == pmtright)dEdx[pmtright] = 0.;
473        }
474      }
475    
476    
477    //  if( adcraw[pmtleft] >3000 || adcraw[pmtright] >3000)SatWarning=1;  //old version
478    
479    // Increase SatWarning Counter for each PMT>Sat
480      if( adcraw[pmtleft] > PMTsat[pmtleft])SatWarning++;  
481      if( adcraw[pmtright] > PMTsat[pmtright])SatWarning++;
482    
483    // if ADC  > sat set dEdx=1000
484      if( adcraw[pmtleft] > PMTsat[pmtleft]) dEdx[pmtleft] = 1000.;
485      if( adcraw[pmtright] > PMTsat[pmtright]) dEdx[pmtright] = 1000. ;
486    
487    // if two PMT are good, take mean dEdx, otherwise only the good dEdx
488      if(dEdx[pmtleft]<1000 && dEdx[pmtright]<1000) PadEdx = (dEdx[pmtleft]+dEdx[pmtright])*0.5;
489      if(dEdx[pmtleft]==1000 && dEdx[pmtright]<1000) PadEdx = dEdx[pmtright];  
490      if(dEdx[pmtleft]<1000 && dEdx[pmtright]==1000) PadEdx = dEdx[pmtleft];
491      
492    };
493    //
494    
495    
496    // gf Apr 07
497    
498    /**
499     * Method to get the PMT name (like "S11_1A") if the PMT_ID is given.
500     * Indexes of corresponding  plane, paddle and  pmt are also given as output.
501     * @param ind  PMT_ID (0 - 47)
502     * @param iplane plane index (0 - 5)
503     * @param ipaddle paddle index (relative to the plane)
504     * @param ipmt pmt index (0(A), 1(B))
505     */
506    TString ToFLevel2::GetPMTName(Int_t ind, Int_t &iplane, Int_t &ipaddle,Int_t &ipmt){
507      
508      TString pmtname = " ";
509      
510      TString photoS[48] = {
511        "S11_1A", "S11_1B", "S11_2A", "S11_2B", "S11_3A", "S11_3B", "S11_4A",
512        "S11_4B",
513        "S11_5A", "S11_5B", "S11_6A", "S11_6B", "S11_7A", "S11_7B", "S11_8A",
514        "S11_8B",
515        "S12_1A", "S12_1B", "S12_2A", "S12_2B", "S12_3A", "S12_3B", "S12_4A",
516        "S12_4B", "S12_5A",  "S12_5B", "S12_6A", "S12_6B",
517        "S21_1A", "S21_1B", "S21_2A", "S21_2B",
518        "S22_1A", "S22_1B", "S22_2A", "S22_2B",
519        "S31_1A", "S31_1B", "S31_2A", "S31_2B", "S31_3A", "S31_3B",
520        "S32_1A", "S32_1B", "S32_2A", "S32_2B", "S32_3A", "S32_3B"
521      };
522      
523      
524      pmtname = photoS[ind].Data();
525      
526      TString ss = pmtname(1,2);
527      iplane  = (int)(atoi(ss.Data())/10)*2-3+atoi(ss.Data())%10;
528      ss = pmtname(4);
529      ipaddle = atoi(ss.Data())-1 ;
530      if( pmtname.Contains("A") )ipmt=0;
531      if( pmtname.Contains("B") )ipmt=1;
532      
533      return pmtname;
534    };
535    /**
536     * Method to get the PMT name (like "S11_1A") if the PMT_ID is given
537     * @param ind  PMT_ID (0 - 47)
538     */
539    TString ToFLevel2::GetPMTName(Int_t ind){
540    
541      Int_t iplane  = -1;
542      Int_t ipaddle = -1;
543      Int_t ipmt    = -1;
544      return GetPMTName(ind,iplane,ipaddle,ipmt);
545      
546    };
547    
548    // wm jun 08
549    Int_t ToFLevel2::GetPaddleIdOfTrack(Float_t xtr, Float_t ytr, Int_t plane){
550    return GetPaddleIdOfTrack(xtr ,ytr ,plane, 0.4);
551    }
552    
553    // gf Apr 07
554    Int_t ToFLevel2::GetPaddleIdOfTrack(Float_t xtr, Float_t ytr, Int_t plane, Float_t margin){
555      
556      Double_t xt,yt,xl,xh,yl,yh;
557      
558      Float_t tof11_x[8] = {-17.85,-12.75,-7.65,-2.55,2.55,7.65,12.75,17.85};
559      Float_t tof12_y[6] = { -13.75,-8.25,-2.75,2.75,8.25,13.75};
560      Float_t tof21_y[2] = { 3.75,-3.75};
561      Float_t tof22_x[2] = { -4.5,4.5};
562      Float_t tof31_x[3] = { -6.0,0.,6.0};
563      Float_t tof32_y[3] = { -5.0,0.0,5.0};
564      
565      //  S11 8 paddles  33.0 x 5.1 cm
566      //  S12 6 paddles  40.8 x 5.5 cm
567      //  S21 2 paddles  18.0 x 7.5 cm
568      //  S22 2 paddles  15.0 x 9.0 cm
569      //  S31 3 paddles  15.0 x 6.0 cm
570      //  S32 3 paddles  18.0 x 5.0 cm
571      
572      Int_t paddleidoftrack=-1;
573      //
574      
575      //--- S11 ------
576      
577      if(plane==0){
578        xt = xtr;
579        yt = ytr;
580        paddleidoftrack=-1;
581        yl = -33.0/2. ;
582        yh =  33.0/2. ;
583        if ((yt>yl)&&(yt<yh)) {
584          for (Int_t i1=0; i1<8;i1++){
585            xl = tof11_x[i1] - (5.1-margin)/2. ;
586            xh = tof11_x[i1] + (5.1-margin)/2. ;
587            if ((xt>xl)&&(xt<xh))  paddleidoftrack=i1;
588          }
589        }
590      }
591      //      cout<<"S11  "<<paddleidoftrack[0]<<"\n";
592      
593      //--- S12 -------
594      if(plane==1){
595        xt = xtr;
596        yt = ytr;
597        paddleidoftrack=-1;
598        xl = -40.8/2. ;
599        xh =  40.8/2. ;
600        
601        if ((xt>xl)&&(xt<xh)) {
602          for (Int_t i1=0; i1<6;i1++){
603            yl = tof12_y[i1] - (5.5-margin)/2. ;
604            yh = tof12_y[i1] + (5.5-margin)/2. ;
605            if ((yt>yl)&&(yt<yh))  paddleidoftrack=i1;
606          }
607        }
608      }
609      
610      //--- S21 ------
611    
612      if(plane==2){
613        xt = xtr;
614        yt = ytr;
615        paddleidoftrack=-1;
616        xl = -18./2. ;
617        xh =  18./2. ;
618        
619        if ((xt>xl)&&(xt<xh)) {
620          for (Int_t i1=0; i1<2;i1++){
621            yl = tof21_y[i1] - (7.5-margin)/2. ;
622            yh = tof21_y[i1] + (7.5-margin)/2. ;
623            if ((yt>yl)&&(yt<yh))  paddleidoftrack=i1;
624          }
625        }
626      }
627      
628      //--- S22 ------
629      if(plane==3){
630        xt = xtr;
631        yt = ytr;
632        paddleidoftrack=-1;
633        yl = -15./2. ;
634        yh =  15./2. ;
635        
636        if ((yt>yl)&&(yt<yh)) {
637          for (Int_t i1=0; i1<2;i1++){
638            xl = tof22_x[i1] - (9.0-margin)/2. ;
639            xh = tof22_x[i1] + (9.0-margin)/2. ;
640            if ((xt>xl)&&(xt<xh))  paddleidoftrack=i1;
641          }
642        }
643      }  
644      
645      //--- S31 ------
646      if(plane==4){
647        xt = xtr;
648        yt = ytr;
649        paddleidoftrack=-1;
650        yl = -15.0/2. ;
651        yh =  15.0/2. ;
652        
653        if ((yt>yl)&&(yt<yh)) {
654          for (Int_t i1=0; i1<3;i1++){
655            xl = tof31_x[i1] - (6.0-margin)/2. ;
656            xh = tof31_x[i1] + (6.0-margin)/2. ;
657            if ((xt>xl)&&(xt<xh))  paddleidoftrack=i1;
658          }
659        }
660      }  
661      
662      //---  S32 ------
663      if(plane==5){
664        xt = xtr;
665        yt = ytr;
666        paddleidoftrack=-1;
667        xl = -18.0/2. ;
668        xh =  18.0/2. ;
669        
670        if ((xt>xl)&&(xt<xh)) {
671          for (Int_t i1=0; i1<3;i1++){
672            yl = tof32_y[i1] - (5.0-margin)/2. ;
673            yh = tof32_y[i1] + (5.0-margin)/2. ;
674            if ((yt>yl)&&(yt<yh)) paddleidoftrack=i1;
675          }
676        }
677      }
678      
679      return paddleidoftrack;
680    
681    }  
682    
683    //
684    
685    // gf Apr 07
686    
687    void ToFLevel2::GetPMTPaddle(Int_t pmt_id, Int_t &plane, Int_t &paddle){
688      
689      plane = GetPlaneIndex(pmt_id);
690    
691      if(plane == 0){
692        if(pmt_id==0 || pmt_id==1)paddle=0;
693        if(pmt_id==2 || pmt_id==3)paddle=1;
694        if(pmt_id==4 || pmt_id==5)paddle=2;
695        if(pmt_id==6 || pmt_id==7)paddle=3;
696        if(pmt_id==8 || pmt_id==9)paddle=4;
697        if(pmt_id==10 || pmt_id==11)paddle=5;
698        if(pmt_id==12 || pmt_id==13)paddle=6;
699        if(pmt_id==14 || pmt_id==15)paddle=7;
700      }
701      
702      if(plane == 1){
703        if(pmt_id==16 || pmt_id==17)paddle=0;
704        if(pmt_id==18 || pmt_id==19)paddle=1;
705        if(pmt_id==20 || pmt_id==21)paddle=2;
706        if(pmt_id==22 || pmt_id==23)paddle=3;
707        if(pmt_id==24 || pmt_id==25)paddle=4;
708        if(pmt_id==26 || pmt_id==27)paddle=5;
709      }
710      
711      if(plane == 2){
712        if(pmt_id==28 || pmt_id==29)paddle=0;
713        if(pmt_id==30 || pmt_id==31)paddle=1;
714      }
715      
716      if(plane == 3){
717        if(pmt_id==32 || pmt_id==33)paddle=0;
718        if(pmt_id==34 || pmt_id==35)paddle=1;
719      }
720      
721      if(plane == 4){
722        if(pmt_id==36 || pmt_id==37)paddle=0;
723        if(pmt_id==38 || pmt_id==39)paddle=1;
724        if(pmt_id==40 || pmt_id==41)paddle=2;
725      }
726      
727      if(plane == 5){
728        if(pmt_id==42 || pmt_id==43)paddle=0;
729        if(pmt_id==44 || pmt_id==45)paddle=1;
730        if(pmt_id==46 || pmt_id==47)paddle=2;
731      }
732      return;
733    }
734    
735    //
736    
737    // gf Apr 07
738    
739    void ToFLevel2::GetPaddlePMT(Int_t paddle, Int_t &pmtleft, Int_t &pmtright){
740      pmtleft=paddle*2;
741      pmtright= pmtleft+1;  
742      return;
743    }
744    
745    //
746    
747    
748    
749    // // gf Apr 07
750    
751    void ToFLevel2::GetPaddleGeometry(Int_t plane, Int_t paddle, Float_t &xleft, Float_t &xright, Float_t &yleft, Float_t &yright){
752      
753      Int_t i1;
754    
755      Float_t tof11_x[8] = {-17.85,-12.75,-7.65,-2.55,2.55,7.65,12.75,17.85};
756      Float_t tof12_y[6] = { -13.75,-8.25,-2.75,2.75,8.25,13.75};
757      Float_t tof21_y[2] = { 3.75,-3.75};
758      Float_t tof22_x[2] = { -4.5,4.5};
759      Float_t tof31_x[3] = { -6.0,0.,6.0};
760      Float_t tof32_y[3] = { -5.0,0.0,5.0};
761            
762      //  S11 8 paddles  33.0 x 5.1 cm
763      //  S12 6 paddles  40.8 x 5.5 cm
764      //  S21 2 paddles  18.0 x 7.5 cm
765      //  S22 2 paddles  15.0 x 9.0 cm
766      //  S31 3 paddles  15.0 x 6.0 cm
767      //  S32 3 paddles  18.0 x 5.0 cm
768    
769      if(plane==0)
770        {
771          for (i1=0; i1<8;i1++){
772            if(i1 == paddle){
773              xleft = tof11_x[i1] - 5.1/2.;
774              xright = tof11_x[i1] + 5.1/2.;
775              yleft = -33.0/2.;
776              yright = 33.0/2.;
777            }
778          }
779        }
780      
781      if(plane==1)
782        {
783          for (i1=0; i1<6;i1++){
784            if(i1 == paddle){
785              xleft = -40.8/2.;
786              xright = 40.8/2.;
787              yleft = tof12_y[i1] - 5.5/2.;
788              yright = tof12_y[i1] + 5.5/2.;
789            }
790          }
791        }
792    
793      if(plane==2)
794        {
795          for (i1=0; i1<2;i1++){
796            if(i1 == paddle){
797              xleft =  -18./2.;
798              xright = 18./2.;
799              yleft = tof21_y[i1] - 7.5/2.;
800              yright = tof21_y[i1] + 7.5/2.;
801            }
802          }
803        }
804      
805      if(plane==3)
806        {
807          for (i1=0; i1<2;i1++){
808            if(i1 == paddle){
809              xleft = tof22_x[i1] - 9.0/2.;
810              xright = tof22_x[i1] + 9.0/2.;
811              yleft = -15./2.;
812              yright = 15./2.;
813            }
814          }
815        }
816    
817    
818      if(plane==4)
819        {
820          for (i1=0; i1<3;i1++){
821            if(i1 == paddle){
822              xleft = tof31_x[i1] - 6.0/2.;
823              xright = tof31_x[i1] + 6.0/2.;
824              yleft = -15./2.;
825              yright = 15./2.;
826            }
827          }
828        }
829    
830      if(plane==5)
831        {
832          for (i1=0; i1<3;i1++){
833            if(i1 == paddle){
834              xleft = -18.0/2.;
835              xright = 18.0/2.;
836              yleft = tof32_y[i1] - 5.0/2.;
837              yright = tof32_y[i1] + 5.0/2.;
838            }
839          }
840        }
841      return;
842    }
843    
844    // gf Apr 07
845    /**
846     * Method to get the paddle index (0,...23) if the plane ID and the paddle id in the plane is given.
847     * This method is the
848     * "reverse" of method "GetPaddlePlane"
849     * @param plane    (0 - 5)
850     * @param paddle   (plane=0, paddle = 0,...5)
851     * @param padid    (0 - 23)
852     */
853    Int_t ToFLevel2::GetPaddleid(Int_t plane, Int_t paddle)
854    {
855      Int_t padid=-1;
856      Int_t pads[6]={8,6,2,2,3,3};
857    
858      int somma=0;
859      int np=plane;
860      for(Int_t j=0; j<np; j++){
861        somma+=pads[j];
862      }
863      padid=paddle+somma;
864      return padid;
865    
866    }
867    
868    
869    // gf Apr 07
870    /**
871     * Method to get the plane ID and the paddle id in the plane if the paddle index (0,...23) is given.
872     * This method is the
873     * "reverse" of method "GetPaddleid"
874     * @param pad      (0 - 23)
875     * @param plane    (0 - 5)
876     * @param paddle   (plane=0, paddle = 0,...5)
877     */
878    void ToFLevel2::GetPaddlePlane(Int_t pad, Int_t &plane, Int_t &paddle)
879    {
880    
881      Int_t pads11=8;
882      Int_t pads12=6;
883      Int_t pads21=2;
884      Int_t pads22=2;
885      Int_t pads31=3;
886      // Int_t pads32=3;
887    
888      if(pad<8){
889        plane=0;
890        paddle=pad;
891        return;
892      }
893    
894      if((7<pad)&&(pad<14)){
895        plane=1;
896        paddle=pad-pads11;
897        return;
898      }
899      
900      if((13<pad)&&(pad<16)){
901        plane=2;
902        paddle=pad-pads11-pads12;
903        return;
904      }
905    
906      if((15<pad)&&(pad<18)){
907        plane=3;
908        paddle=pad-pads11-pads12-pads21;
909        return;
910      }
911    
912      if((17<pad)&&(pad<21)){
913        plane=4;
914        paddle=pad-pads11-pads12-pads21-pads22;
915        return;
916      }
917    
918      if((20<pad)&&(pad<24)){
919        plane=5;
920        paddle=pad-pads11-pads12-pads21-pads22-pads31;
921        return;
922      }  
923    
924    }
925    
926    
927    Int_t ToFLevel2::GetNPaddle(Int_t plane){
928    
929      Int_t npaddle=-1;
930    
931      Int_t pads11=8;
932      Int_t pads12=6;
933      Int_t pads21=2;
934      Int_t pads22=2;
935      Int_t pads31=3;
936      Int_t pads32=3;
937    
938      if(plane==0)npaddle=pads11;
939      if(plane==1)npaddle=pads12;
940      if(plane==2)npaddle=pads21;
941      if(plane==3)npaddle=pads22;
942      if(plane==4)npaddle=pads31;
943      if(plane==5)npaddle=pads32;
944    
945      return npaddle;
946    
947    }
948    
949    
950    
951    /// wm feb 08
952    
953    /**
954     * Method to calculate Beta from the 12 single measurements
955     * we check the individual weights for artificial TDC values, then calculate
956     * am mean beta for the first time. In a second step we loop again through
957     * the single measurements, checking for the residual from the mean
958     * The cut on the residual reject measurements > "x"-sigma. A chi2 value is
959     * calculated, furthermore a "quality" value by adding the weights which
960     * are finally used. If all measurements are taken, "quality" will be = 22.47.
961     * A chi2 cut around 3-4 and a quality-cut > 20 is needed for clean beta
962     * measurements like antiprotons etc.
963     * The Level2 output is derived in the fortran routines using: 10.,10.,20.
964     * @param notrack Track Number
965     * @param cut on residual: difference between single measurement and mean
966     * @param cut on "quality"
967     * @param cut on chi2
968     */
969    
970    Float_t ToFLevel2::CalcBeta(Int_t notrack, Float_t resmax, Float_t qualitycut, Float_t chi2cut){
971    
972    //  cout<<" in CalcBeta "<<resmax<<" "<<chi2cut<<" "<<qualitycut<<endl;
973    
974      Float_t bxx = 100.;
975      //
976      ToFTrkVar *trk = GetToFTrkVar(notrack);
977      if(!trk) return 0; //ELENA
978    
979    
980      Float_t chi2,xhelp,beta_mean;
981      Float_t w_i[12],quality,sw,sxw,res,betachi,beta_mean_inv;
982      Float_t b[12],tdcfl;
983      Int_t  pmt_id,pmt_plane;
984    
985      for (Int_t i=0; i<12; i++){
986        b[i] = trk->beta[i];
987                                  }
988          
989    
990    //========================================================================
991    //---  Find out ToF layers with artificial TDC values & fill vector    ---
992    //========================================================================
993    
994    Float_t  w_il[6];
995    
996         for (Int_t jj=0; jj<6;jj++) {
997             w_il[jj] = 1000.;
998                                     }
999    
1000    
1001      for (Int_t i=0; i<trk->npmttdc; i++){
1002        //
1003        pmt_id = (trk->pmttdc).At(i);
1004        pmt_plane = GetPlaneIndex(pmt_id);
1005        tdcfl = (trk->tdcflag).At(i);
1006        if (w_il[pmt_plane] != 1.) w_il[pmt_plane] = tdcfl; //tdcflag
1007                                         };
1008      
1009    //========================================================================
1010    //---  Set weights for the 12 measurements using information for top and bottom:
1011    //---  if no measurements: weight = set to very high value=> not used
1012    //---  top or bottom artificial: weight*sqrt(2)
1013    //---  top and bottom artificial: weight*sqrt(2)*sqrt(2)
1014    //========================================================================
1015    
1016    Int_t itop[12] = {0,0,1,1,2,2,3,3,0,0,1,1};
1017    Int_t ibot[12] = {4,5,4,5,4,5,4,5,2,3,2,3};
1018    
1019         xhelp= 1E09;
1020      
1021         for (Int_t jj=0; jj<12;jj++) {
1022         if (jj<4)           xhelp = 0.11;    // S1-S3
1023         if ((jj>3)&&(jj<8)) xhelp = 0.18;    // S2-S3
1024         if (jj>7)           xhelp = 0.28;    // S1-S2
1025         if ((w_il[itop[jj]] == 1000.) && (w_il[ibot[jj]] == 1000.)) xhelp = 1E09;
1026         if ((w_il[itop[jj]] == 1) || (w_il[ibot[jj]] == 1.)) xhelp = xhelp*1.414 ;
1027         if ((w_il[itop[jj]] == 1) && (w_il[ibot[jj]] == 1.)) xhelp = xhelp*2. ;
1028    
1029         w_i[jj] = 1./xhelp;
1030                                      }
1031    
1032    
1033    //========================================================================
1034    //--- Calculate mean beta for the first time -----------------------------
1035    //--- We are using "1/beta" since its error is gaussian ------------------
1036    //========================================================================
1037    
1038          Int_t icount=0;
1039          sw=0.;
1040          sxw=0.;
1041          beta_mean=100.;
1042    
1043              for (Int_t jj=0; jj<12;jj++){
1044            if ((fabs(1./b[jj])>0.1)&&(fabs(1./b[jj])<15.))
1045             {
1046                icount= icount+1;
1047                sxw=sxw + (1./b[jj])*w_i[jj]*w_i[jj] ;
1048                sw =sw + w_i[jj]*w_i[jj] ;
1049    
1050             }
1051             }
1052    
1053          if (icount>0) beta_mean=1./(sxw/sw);
1054          beta_mean_inv = 1./beta_mean;
1055    
1056    //========================================================================
1057    //--- Calculate beta for the second time, use residuals of the single
1058    //--- measurements to get a chi2 value
1059    //========================================================================
1060    
1061          icount=0;
1062          sw=0.;
1063          sxw=0.;
1064          betachi = 100.;
1065          chi2 = 0.;
1066          quality=0.;
1067    
1068    
1069              for (Int_t jj=0; jj<12;jj++){
1070           if ((fabs(1./b[jj])>0.1)&&(fabs(1./b[jj])<15.)&&(w_i[jj]>0.01)) {
1071                res = beta_mean_inv - (1./b[jj]) ;
1072                if (fabs(res*w_i[jj])<resmax)          {;
1073                chi2 = chi2 + pow((res*w_i[jj]),2) ;
1074                icount= icount+1;
1075                sxw=sxw + (1./b[jj])*w_i[jj]*w_i[jj] ;
1076                sw =sw + w_i[jj]*w_i[jj] ;
1077                                                   }
1078                                                                            }
1079                                          }
1080          quality = sqrt(sw) ;
1081    
1082          if (icount==0) chi2 = 1000.;
1083          if (icount>0) chi2 = chi2/(icount) ;
1084          if (icount>0) betachi=1./(sxw/sw);
1085    
1086       bxx = 100.;
1087       if ((chi2 < chi2cut)&&(quality>qualitycut)) bxx = betachi;
1088      //
1089      return(bxx);
1090    };
1091    
1092    
1093    ////////////////////////////////////////////////////
1094    ////////////////////////////////////////////////////
1095    
1096    
1097    /**
1098     * Fills a struct cToFLevel2 with values from a ToFLevel2 object (to put data into a F77 common).
1099     */
1100    void ToFLevel2::GetLevel2Struct(cToFLevel2 *l2) const{
1101    
1102      for(Int_t i=0;i<6;i++)
1103        l2->tof_j_flag[i]=tof_j_flag[i];
1104    
1105      if(ToFTrk){ //ELENA
1106          l2->ntoftrk = ToFTrk->GetEntries();
1107          for(Int_t j=0;j<l2->ntoftrk;j++){
1108              l2->toftrkseqno[j]= ((ToFTrkVar*)ToFTrk->At(j))->trkseqno;
1109              l2->npmttdc[j]= ((ToFTrkVar*)ToFTrk->At(j))->npmttdc;
1110              for(Int_t i=0;i<l2->npmttdc[j];i++){
1111                  l2->pmttdc[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->pmttdc.At(i);
1112                  l2->tdcflag[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->tdcflag.At(i); // gf: 30 Nov 2006
1113              }
1114              for(Int_t i=0;i<13;i++)
1115                  l2->beta[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->beta[i];
1116              
1117              l2->npmtadc[j]= ((ToFTrkVar*)ToFTrk->At(j))->npmtadc;
1118              for(Int_t i=0;i<l2->npmtadc[j];i++){
1119                  l2->pmtadc[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->pmtadc.At(i);
1120                  l2->adcflag[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->adcflag.At(i); // gf: 30 Nov 2006
1121                  l2->dedx[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->dedx.At(i);
1122              }
1123              for(Int_t i=0;i<3;i++){
1124                  l2->xtofpos[i][j]=((ToFTrkVar*)ToFTrk->At(j))->xtofpos[i];
1125                  l2->ytofpos[i][j]=((ToFTrkVar*)ToFTrk->At(j))->ytofpos[i];
1126              }
1127              for(Int_t i=0;i<6;i++){
1128                  l2->xtr_tof[i][j]=((ToFTrkVar*)ToFTrk->At(j))->xtr_tof[i];
1129                  l2->ytr_tof[i][j]=((ToFTrkVar*)ToFTrk->At(j))->ytr_tof[i];
1130              }
1131          }
1132      } //ELENA
1133        
1134      if(PMT){ //ELENA
1135          l2->npmt = PMT->GetEntries();
1136          for(Int_t j=0;j<l2->npmt;j++){
1137              l2->pmt_id[j] = ((ToFPMT*)PMT->At(j))->pmt_id;
1138              l2->adc[j] =((ToFPMT*)PMT->At(j))->adc;
1139              l2->tdc_tw[j] =((ToFPMT*)PMT->At(j))->tdc_tw;
1140          }
1141      } //ELENA
1142    }
1143    
1144    
1145    //
1146    // Reprocessing tool // Emiliano 08/04/07
1147    //
1148    Int_t ToFLevel2::Process(TrkLevel2 *trk, TrigLevel2 *trg, GL_RUN *run, OrbitalInfo *orb, Bool_t force){
1149      //
1150      // Copiare qui qualcosa di simile a calonuclei per evitare di riprocessare sempre tutto
1151      //
1152      printf("\n\n\n ERROR: NOT IMPLEMENTED ANYMORE, write Emiliano if you need this method (Emiliano.Mocchiutti@ts.infn.it) \n\n\n");
1153      return(-1);
1154      //   //
1155      //   // structures to communicate with F77
1156      //   //
1157      //   extern struct ToFInput  tofinput_;
1158    //   extern struct ToFOutput tofoutput_;
1159    //   //
1160    //   // DB connection
1161    //   //
1162    //   TString host;
1163    //   TString user;
1164    //   TString psw;
1165    //   const char *pamdbhost=gSystem->Getenv("PAM_DBHOST");
1166    //   const char *pamdbuser=gSystem->Getenv("PAM_DBUSER");
1167    //   const char *pamdbpsw=gSystem->Getenv("PAM_DBPSW");
1168    //   if ( !pamdbhost ) pamdbhost = "";
1169    //   if ( !pamdbuser ) pamdbuser = "";
1170    //   if ( !pamdbpsw ) pamdbpsw = "";
1171    //   if ( strcmp(pamdbhost,"") ) host = pamdbhost;
1172    //   if ( strcmp(pamdbuser,"") ) user = pamdbuser;
1173    //   if ( strcmp(pamdbpsw,"") ) psw = pamdbpsw;
1174    //   //
1175    //   //
1176    //   TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
1177    //   if ( !dbc->IsConnected() ) return 1;
1178    //   stringstream myquery;
1179    //   myquery.str("");
1180    //   myquery << "SET time_zone='+0:00'";
1181    //   dbc->Query(myquery.str().c_str());
1182    //   GL_PARAM *glparam = new GL_PARAM();
1183    //   glparam->Query_GL_PARAM(1,1,dbc); // parameters stored in DB in GL_PRAM table
1184    //   trk->LoadField(glparam->PATH+glparam->NAME);
1185    //   //
1186    //   Bool_t defcal = true;
1187    //   Int_t error=glparam->Query_GL_PARAM(run->RUNHEADER_TIME,201,dbc); // parameters stored in DB in GL_PRAM table
1188    //   if ( error<0 ) {
1189    //     return(1);
1190    //   };
1191    //   printf(" Reading ToF parameter file: %s \n",(glparam->PATH+glparam->NAME).Data());
1192    //   if ( (UInt_t)glparam->TO_TIME != (UInt_t)4294967295UL ) defcal = false;
1193    //   //
1194    //   Int_t nlen = (Int_t)(glparam->PATH+glparam->NAME).Length();
1195    //   rdtofcal((char *)(glparam->PATH+glparam->NAME).Data(),&nlen);
1196    //   //
1197    //   Int_t adc[4][12];
1198    //   Int_t tdc[4][12];
1199    //   Float_t tdcc[4][12];
1200    //   //
1201    //   // process tof data
1202    //   //
1203    //   for (Int_t hh=0; hh<12;hh++){
1204    //     for (Int_t kk=0; kk<4;kk++){
1205    //            adc[kk][hh] = 4095;
1206    //            tdc[kk][hh] = 4095;
1207    //            tdcc[kk][hh] = 4095.;
1208    //            tofinput_.adc[hh][kk] = 4095;
1209    //            tofinput_.tdc[hh][kk] = 4095;
1210    //     };
1211    //   };
1212    //   Int_t ntrkentry = 0;
1213    //   Int_t npmtentry = 0;
1214    //   Int_t gg = 0;
1215    //   Int_t hh = 0;
1216    //   Int_t adcf[48];
1217    //   memset(adcf, 0, 48*sizeof(Int_t));
1218    //   Int_t tdcf[48];
1219    //   memset(tdcf, 0, 48*sizeof(Int_t));
1220    //   for (Int_t pm=0; pm < this->ntrk() ; pm++){
1221    //      ToFTrkVar *ttf = this->GetToFTrkVar(pm);
1222    //      for ( Int_t nc=0; nc < ttf->npmttdc; nc++){
1223    //             if ( (ttf->tdcflag).At(nc) != 0 ) tdcf[(ttf->pmttdc).At(nc)] = 1;
1224    //      };
1225    //      for ( Int_t nc=0; nc < ttf->npmtadc; nc++){
1226    //             if ( (ttf->adcflag).At(nc) != 0 ) adcf[(ttf->pmtadc).At(nc)] = 1;
1227    //      };
1228    //   };
1229    //   //
1230    //   for (Int_t pm=0; pm < this->npmt() ; pm++){
1231    //      ToFPMT *pmt = this->GetToFPMT(pm);
1232    //      this->GetPMTIndex(pmt->pmt_id, gg, hh);
1233    //      if ( adcf[pmt->pmt_id] == 0 ){
1234    //              tofinput_.adc[gg][hh] = (int)pmt->adc;
1235    //              adc[hh][gg] = (int)pmt->adc;
1236    //      };
1237    //      if ( tdcf[pmt->pmt_id] == 0 ){
1238    //              tofinput_.tdc[gg][hh] = (int)pmt->tdc;
1239    //              tdc[hh][gg] = (int)pmt->tdc;
1240    //      };
1241    //      tdcc[hh][gg] = (float)pmt->tdc_tw;
1242    //      // Int_t pppid = this->GetPMTid(hh,gg);
1243    //      //      printf(" pm %i pmt_id %i pppid %i hh %i gg %i tdcc %f tdc %f adc %f \n",pm,pmt->pmt_id,pppid,hh,gg,pmt->tdc_tw,pmt->tdc,pmt->adc);
1244    //   };
1245    //   //
1246    //   Int_t unpackError = this->unpackError;
1247    //   //
1248    //   for (Int_t hh=0; hh<5;hh++){
1249    //      tofinput_.patterntrig[hh]=trg->patterntrig[hh];
1250    //   };
1251    //   //
1252    //   this->Clear();
1253    //   //
1254    //       Int_t pmt_id = 0;
1255    //       ToFPMT *t_pmt = new ToFPMT();
1256    //       if(!(this->PMT)) this->PMT = new TClonesArray("ToFPMT",12); //ELENA
1257    //       TClonesArray &tpmt = *this->PMT;
1258    //       ToFTrkVar *t_tof = new ToFTrkVar();
1259    //       if(!(this->ToFTrk)) this->ToFTrk = new TClonesArray("ToFTrkVar",2); //ELENA
1260    //       TClonesArray &t = *this->ToFTrk;
1261    //       //
1262    //       //
1263    //       // Here we have calibrated data, ready to be passed to the FORTRAN routine which will extract common and track-related  variables.
1264    //       //
1265    //       npmtentry = 0;
1266    //       //
1267    //       ntrkentry = 0;
1268    //       //
1269    //       // Calculate tracks informations from ToF alone
1270    //       //
1271    //       tofl2com();
1272    //       //
1273    //       memcpy(this->tof_j_flag,tofoutput_.tof_j_flag,6*sizeof(Int_t));
1274    //       //
1275    //       t_tof->trkseqno = -1;
1276    //       //
1277    //       // and now we must copy from the output structure to the level2 class:
1278    //       //
1279    //       t_tof->npmttdc = 0;
1280    //       //
1281    //       for (Int_t hh=0; hh<12;hh++){
1282    //         for (Int_t kk=0; kk<4;kk++){
1283    //           if ( tofoutput_.tofmask[hh][kk] != 0 ){
1284    //             pmt_id = this->GetPMTid(kk,hh);
1285    //             t_tof->pmttdc.AddAt(pmt_id,t_tof->npmttdc);
1286    //             t_tof->tdcflag.AddAt(tofoutput_.tdcflagtof[hh][kk],t_tof->npmttdc); // gf: Jan 09/07
1287    //             t_tof->npmttdc++;
1288    //           };
1289    //         };
1290    //       };
1291    //       for (Int_t kk=0; kk<13;kk++){
1292    //         t_tof->beta[kk] = tofoutput_.betatof_a[kk];
1293    //       }
1294    //       //
1295    //       t_tof->npmtadc = 0;
1296    //       for (Int_t hh=0; hh<12;hh++){
1297    //         for (Int_t kk=0; kk<4;kk++){
1298    //           if ( tofoutput_.adctof_c[hh][kk] < 1000 ){
1299    //             t_tof->dedx.AddAt(tofoutput_.adctof_c[hh][kk],t_tof->npmtadc);
1300    //             pmt_id = this->GetPMTid(kk,hh);
1301    //             t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
1302    //             t_tof->adcflag.AddAt(tofoutput_.adcflagtof[hh][kk],t_tof->npmtadc); // gf: Jan 09/07
1303    //             t_tof->npmtadc++;
1304    //           };
1305    //         };
1306    //       };
1307    //       //
1308    //       memcpy(t_tof->xtofpos,tofoutput_.xtofpos,sizeof(t_tof->xtofpos));
1309    //       memcpy(t_tof->ytofpos,tofoutput_.ytofpos,sizeof(t_tof->ytofpos));
1310    //       memcpy(t_tof->xtr_tof,tofoutput_.xtr_tof,sizeof(t_tof->xtr_tof));
1311    //       memcpy(t_tof->ytr_tof,tofoutput_.ytr_tof,sizeof(t_tof->ytr_tof));
1312    //       //
1313    //       new(t[ntrkentry]) ToFTrkVar(*t_tof);
1314    //       ntrkentry++;
1315    //       t_tof->Clear();
1316    //       //
1317    //       //
1318    //       //
1319    //       t_pmt->Clear();
1320    //       //
1321    //       for (Int_t hh=0; hh<12;hh++){
1322    //         for (Int_t kk=0; kk<4;kk++){
1323    //          // new WM
1324    //           if ( tofoutput_.tdc_c[hh][kk] < 4095 || adc[kk][hh] < 4095  || tdc[kk][hh] < 4095 ){
1325    // //          if ( tdcc[kk][hh] < 4095. || adc[kk][hh] < 4095  || tdc[kk][hh] < 4095 ){
1326    //             //
1327    //             t_pmt->pmt_id = this->GetPMTid(kk,hh);
1328    //             t_pmt->tdc_tw = tofoutput_.tdc_c[hh][kk];
1329    //             t_pmt->adc = (Float_t)adc[kk][hh];
1330    //             t_pmt->tdc = (Float_t)tdc[kk][hh];
1331    //             //
1332    //             new(tpmt[npmtentry]) ToFPMT(*t_pmt);
1333    //             npmtentry++;
1334    //             t_pmt->Clear();
1335    //           };
1336    //         };
1337    //       };
1338    //       //
1339    //       // Calculate track-related variables
1340    //       //
1341    //       if ( trk->ntrk() > 0 ){
1342    //         //
1343    //         // We have at least one track
1344    //         //
1345    //         //
1346    //         // Run over tracks
1347    //         //
1348    //         for(Int_t nt=0; nt < trk->ntrk(); nt++){
1349    //           //
1350    //           TrkTrack *ptt = trk->GetStoredTrack(nt);
1351    //           //
1352    //           // Copy the alpha vector in the input structure
1353    //           //
1354    //           for (Int_t e = 0; e < 5 ; e++){
1355    //             tofinput_.al_pp[e] = ptt->al[e];
1356    //           };
1357    //           //
1358    //           // Get tracker related variables for this track
1359    //           //
1360    //           toftrk();
1361    //           //
1362    //           // Copy values in the class from the structure (we need to use a temporary class to store variables).
1363    //           //
1364    //           t_tof->npmttdc = 0;
1365    //           for (Int_t hh=0; hh<12;hh++){
1366    //             for (Int_t kk=0; kk<4;kk++){
1367    //               if ( tofoutput_.tofmask[hh][kk] != 0 ){
1368    //                 pmt_id = this->GetPMTid(kk,hh);
1369    //                 t_tof->pmttdc.AddAt(pmt_id,t_tof->npmttdc);
1370    //                 t_tof->tdcflag.AddAt(tofoutput_.tdcflag[hh][kk],t_tof->npmttdc); // gf: Jan 09/07
1371    //                 t_tof->npmttdc++;
1372    //               };
1373    //             };
1374    //           };
1375    //           for (Int_t kk=0; kk<13;kk++){
1376    //             t_tof->beta[kk] = tofoutput_.beta_a[kk];
1377    //           };
1378    //           //
1379    //           t_tof->npmtadc = 0;
1380    //           for (Int_t hh=0; hh<12;hh++){
1381    //             for (Int_t kk=0; kk<4;kk++){
1382    //               if ( tofoutput_.adc_c[hh][kk] < 1000 ){
1383    //                 t_tof->dedx.AddAt(tofoutput_.adc_c[hh][kk],t_tof->npmtadc);
1384    //                 pmt_id = this->GetPMTid(kk,hh);
1385    //                 t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
1386    //                 t_tof->adcflag.AddAt(tofoutput_.adcflag[hh][kk],t_tof->npmtadc); // gf: Jan 09/07
1387    //                 t_tof->npmtadc++;
1388    //               };
1389    //             };
1390    //           };
1391    //           //
1392    //           memcpy(t_tof->xtofpos,tofoutput_.xtofpos,sizeof(t_tof->xtofpos));
1393    //           memcpy(t_tof->ytofpos,tofoutput_.ytofpos,sizeof(t_tof->ytofpos));
1394    //           memcpy(t_tof->xtr_tof,tofoutput_.xtr_tof,sizeof(t_tof->xtr_tof));
1395    //           memcpy(t_tof->ytr_tof,tofoutput_.ytr_tof,sizeof(t_tof->ytr_tof));
1396    //           //
1397    //           // Store the tracker track number in order to be sure to have shyncronized data during analysis
1398    //           //
1399    //           t_tof->trkseqno = nt;
1400    //           //
1401    //           // create a new object for this event with track-related variables
1402    //           //
1403    //           new(t[ntrkentry]) ToFTrkVar(*t_tof);
1404    //           ntrkentry++;
1405    //           t_tof->Clear();
1406    //           //
1407    //         }; // loop on all the tracks
1408    //       //
1409    //       this->unpackError = unpackError;
1410    //       if ( defcal ){
1411    //         this->default_calib = 1;
1412    //       } else {
1413    //         this->default_calib = 0;
1414    //       };
1415    //};
1416    //  return(0);
1417    }
1418    
1419    
1420    ToFdEdx::ToFdEdx()
1421    {
1422      memset(conn,0,12*sizeof(Bool_t));
1423      memset(ts,0,12*sizeof(UInt_t));
1424      memset(te,0,12*sizeof(UInt_t));
1425      eDEDXpmt = new TArrayF(48);
1426      Define_PMTsat();
1427      Clear();
1428    }
1429    //------------------------------------------------------------------------
1430    void ToFdEdx::CheckConnectors(UInt_t atime, GL_PARAM *glparam, TSQLServer *dbc)
1431    {
1432      for(int i=0; i<12; i++){
1433        if(atime<=ts[i] || atime>te[i]){
1434          Int_t error=glparam->Query_GL_PARAM(atime,210+i,dbc); // parameters stored in DB in GL_PRAM table
1435          if ( error<0 ) {
1436            conn[i]=false;
1437            ts[i]=0;
1438            te[i]=numeric_limits<UInt_t>::max();
1439          };
1440          if ( !error ){
1441            conn[i]=true;
1442            ts[i]=glparam->FROM_TIME;
1443            te[i]=glparam->TO_TIME;
1444          }
1445          if ( error>0 ){
1446            conn[i]=false;
1447            ts[i]=glparam->TO_TIME;
1448            TSQLResult *pResult;
1449            TSQLRow *row;
1450            TString query= Form("SELECT FROM_TIME FROM GL_PARAM WHERE TYPE=%i AND FROM_TIME>=%i ORDER BY FROM_TIME ASC LIMIT 1;",210+i,atime);
1451            pResult=dbc->Query(query.Data());
1452            if(!pResult->GetRowCount()){
1453              te[i]=numeric_limits<UInt_t>::max();
1454            }else{
1455              row=pResult->Next();
1456              te[i]=(UInt_t)atoll(row->GetField(0));
1457            }
1458          }
1459          //
1460          
1461        }
1462      }
1463    
1464    }
1465    //------------------------------------------------------------------------
1466    void ToFdEdx::Clear(Option_t *option)
1467    {
1468      //
1469      // Set arrays and initialize structure
1470      //  eDEDXpmt.Set(48);    eDEDXpmt.Reset(-1);   // Set array size  and reset structure
1471      eDEDXpmt->Set(48);    eDEDXpmt->Reset(-1);   // Set array size  and reset structure
1472      //
1473    };
1474    
1475    //------------------------------------------------------------------------
1476    void ToFdEdx::Print(Option_t *option)
1477    {
1478      //
1479      printf("========================================================================\n");
1480    
1481    };
1482    
1483    //------------------------------------------------------------------------
1484    void ToFdEdx::Init(pamela::tof::TofEvent *tofl0)
1485    {
1486      //
1487      ToFLevel2 tf;
1488      for (Int_t gg=0; gg<4;gg++){
1489        for (Int_t hh=0; hh<12;hh++){
1490          //          tofinput_.tdc[hh][gg]=tofEvent->tdc[gg][hh];          
1491          int mm = tf.GetPMTid(gg,hh);        
1492          adc[mm]=tofl0->adc[gg][hh];
1493        };      
1494      };
1495      
1496    };
1497    
1498    //------------------------------------------------------------------------
1499    void ToFdEdx::Init(Int_t gg, Int_t hh, Float_t adce)
1500    {
1501      //
1502      ToFLevel2 tf;
1503      //  for (Int_t gg=0; gg<4;gg++){
1504      //    for (Int_t hh=0; hh<12;hh++){
1505      int mm = tf.GetPMTid(gg,hh);    
1506      adc[mm]=adce;
1507      
1508    };
1509    //------------------------------------------------------------------------
1510    void ToFdEdx::Process(UInt_t atime, Float_t betamean, Float_t *xtr_tof, Float_t *ytr_tof, Int_t exitat)
1511    {
1512      // the parameters should be already initialised by InitPar()
1513      //  printf(" in process \n");
1514      Clear();
1515    
1516     // define angle:  
1517      double dx   = xtr_tof[1] - xtr_tof[5];
1518      double dy   = ytr_tof[0] - ytr_tof[4];
1519      double dr   = sqrt(dx*dx+dy*dy);
1520      double theta=atan(dr/76.81);
1521      //
1522      if ( xtr_tof[1] > 99. ||  xtr_tof[5] > 99. || ytr_tof[0] > 99. ||  ytr_tof[4] > 99. ) theta = 0.;
1523      for (Int_t ii=0; ii<6; ii++){
1524        if ( xtr_tof[ii] > 99. ) xtr_tof[ii] = 0.;
1525        if ( ytr_tof[ii] > 99. ) ytr_tof[ii] = 0.;
1526      };
1527      //
1528      
1529    
1530      //--------------------- TABLE OF PERIODS WITH HV PROBLEMS ----------------------------
1531      
1532      int Aconn=conn[0];    // PMT 0,20,22,24
1533      int Bconn=conn[1];    // PMT 6,12,26,34
1534      int Cconn=conn[2];    // PMT 4,14,28,32
1535      int Dconn=conn[3];    // PMT 2,8,10,30
1536      int Econn=conn[4];    // PMT 42,43,44,47
1537      int Fconn=conn[5];    // PMT 7,19,23,27
1538      int Gconn=conn[6];    // PMT 3,11,25,33
1539      int Hconn=conn[7];    // PMT 1,9,13,21
1540      int Iconn=conn[8];    // PMT 5,29,31,35
1541      int Lconn=conn[9];    // PMT 37,40,45,46
1542      int Mconn=conn[10];    // PMT 15,16,17,18
1543      int Nconn=conn[11];    // PMT 36,38,39,41
1544      if( false ) cout << Gconn << Iconn << Lconn <<endl; // to avoid compilation warnings
1545        
1546      //  printf(" size %i \n",eDEDXpmt.GetSize());
1547      for( int ii=0; ii<48; ii++ ) {
1548        //
1549        //    eDEDXpmt.SetAt(-1.,ii);
1550        //    printf(" ii %i beta %f atime %u xtr 1 %f ytr 1 %f adc %f \n",ii,betamean,atime,xtr_tof[0],ytr_tof[0],adc[ii]);
1551    
1552        if( adc[ii] >= 4095. ){
1553          //      eDEDXpmt[ii] = 0.;
1554          eDEDXpmt->AddAt(0.,ii);
1555          continue; // EMILIANO
1556        };
1557    
1558        if( adc[ii] >= (PMTsat[ii]-5.) && adc[ii] < 4095. ){
1559          eDEDXpmt->AddAt(1000.,ii);
1560          continue; // EMILIANO
1561        };
1562    
1563        if( adc[ii] <= 0. ) {
1564          eDEDXpmt->AddAt(1500.,ii);
1565          continue;
1566        };
1567        //
1568        double adcpC   = f_adcPC( adc[ii] );    // - adc conversion in pC
1569        if ( exitat == 0 ){
1570          eDEDXpmt->AddAt((Float_t)adcpC,ii);
1571          continue;
1572        }
1573        //    printf(" e qua? \n");
1574    
1575        double adccorr = adcpC*fabs(cos(theta));    
1576        if(adccorr<=0.)           continue;
1577        if ( exitat == 1 ){
1578          eDEDXpmt->AddAt((Float_t)adccorr,ii);
1579          continue;
1580        }
1581        //    printf(" e quo? \n");
1582    
1583        //    int standard=0;
1584        int S115B_ok=0;
1585        int S115B_break=0;
1586    
1587        if(atime<1158720000)S115B_ok=1;
1588        else S115B_break=1;
1589    
1590    
1591        //------------------------------------------------------------------------
1592        //    printf(" e qui? \n");
1593        //---------------------------------------------------- Z reconstruction
1594    
1595        double adcHe, adcnorm, adclin, dEdx, Zeta;
1596    
1597        adcHe=-2;
1598        adcnorm=-2;
1599        adclin=-2;
1600        dEdx=-2;
1601        Zeta=-2;
1602        Double_t correction = 1.;
1603    
1604        if(Aconn==1 && (ii==0 || ii==20 || ii==22 || ii==24)){
1605          correction = 1.675;
1606        }
1607        else if(Bconn==1 && (ii==6 || ii==12 || ii==26 || ii==34)){
1608          correction = 2.482;
1609        }
1610        else if(Cconn==1 && (ii==4 || ii==14 || ii==28 || ii==32)){
1611          correction = 1.464;
1612        }
1613        else if(Dconn==1 && (ii==2 || ii==8 || ii==10 || ii==30)){
1614          correction = 1.995;
1615        }
1616        else if(Econn==1 && (ii==42 || ii==43 || ii==44 || ii==47)){
1617          correction = 1.273;
1618        }
1619        else if(Fconn==1 && (ii==7 || ii==19 || ii==23 || ii==27)){
1620          correction = 1.565;
1621        }
1622        else if(Mconn==1 && (ii==15 || ii==16 || ii==17 || ii==18)){
1623          correction = 1.565;
1624        }
1625        else if(Nconn==1 && (ii==36 || ii==38 || ii==39 || ii==41)){
1626          correction = 1.018;
1627        }
1628        else if(Hconn==1 && (ii==1 || ii==13 || ii==21 || (ii==9&&S115B_ok==1))){
1629          correction = 1.84;
1630        }
1631        else if(S115B_break==1 && ii==9 && Hconn==1){
1632          correction = 1.64;
1633        }
1634        else correction = 1.;
1635        
1636        if( ii==9 && S115B_break==1 ){
1637          adcHe   = f_att5B( ytr_tof[0] )/correction;
1638        } else {
1639          adcHe   = Get_adc_he(ii, xtr_tof, ytr_tof)/correction;
1640        };
1641        if(adcHe<=0)   continue;
1642        if ( exitat == 2 ){
1643          if(ii==9 && S115B_break==1)  eDEDXpmt->AddAt(36.*(Float_t)adccorr/adcHe,ii);
1644          else  adclin  = 4.*(Float_t)adccorr/adcHe;
1645          continue;
1646        }
1647    
1648        if(ii==9 && S115B_break==1)  adcnorm = f_pos5B(adccorr);
1649        else adcnorm = f_pos( (parPos[ii]), adccorr);
1650        if(adcnorm<=0) continue;
1651        if(ii==9 && S115B_break==1)  adclin  = 36.*adcnorm/adcHe;
1652        else  adclin  = 4.*adcnorm/adcHe;
1653        if(adclin<=0)  continue;
1654        if ( exitat == 3 ){
1655          if(ii==9 && S115B_break==1)  eDEDXpmt->AddAt((Float_t)adclin,ii);
1656          else  eDEDXpmt->AddAt((Float_t)adclin,ii);
1657          continue;
1658        }
1659        //
1660        if ( betamean > 99. ){
1661          //      eDEDXpmt.AddAt((Float_t)adclin,ii);
1662          eDEDXpmt->AddAt((Float_t)adclin,ii);
1663          //      printf(" AAPMT IS %i dedx is %f vector is %f \n",ii,adclin,eDEDXpmt[ii]);
1664          continue;
1665        };
1666        //
1667        double dEdxHe=-2;
1668        if(ii==9 && S115B_break==1){
1669          if( betamean <1. ) dEdxHe = f_BB5B( betamean );
1670          else                       dEdxHe = 33;
1671        } else {
1672          if( betamean <1. ) dEdxHe = f_BB( (parBBneg[ii]), betamean );
1673          else                       dEdxHe = parBBpos[ii];
1674        }
1675        
1676        
1677        if(dEdxHe<=0){
1678          eDEDXpmt->AddAt((Float_t)adclin,ii);
1679          continue;
1680        };
1681    
1682        if(ii==9 && S115B_break==1)  dEdx = f_desatBB5B( adclin );
1683        else  dEdx = f_desatBB((parDesatBB[ii]), adclin );
1684    
1685        if(dEdx<=0){
1686          eDEDXpmt->AddAt((Float_t)adclin,ii);
1687          continue;
1688        };
1689    
1690        eDEDXpmt->AddAt((Float_t)dEdx,ii);
1691        //    eDEDXpmt.AddAt((Float_t)dEdx,ii);
1692    
1693        //    printf(" PMT IS %i dedx is %f vector is %f \n",ii,dEdx,eDEDXpmt[ii]);
1694    
1695      }  //end loop on 48 PMT
1696    
1697    };
1698    
1699    
1700    //------------------------------------------------------------------------
1701    void ToFdEdx::Define_PMTsat()
1702    {
1703      Float_t  sat[48] = {
1704        3176.35,3178.19,3167.38,3099.73,3117.00,3126.29,3111.44,3092.27,
1705        3146.48,3094.41,3132.13,3115.37,3099.32,3110.97,3111.80,3143.14,
1706        3106.72,3153.44,3136.00,3188.96,3104.73,3140.45,3073.18,3106.62,
1707        3112.48,3146.92,3127.24,3136.52,3109.59,3112.89,3045.15,3147.26,
1708        3095.92,3121.05,3083.25,3123.62,3150.92,3125.30,3067.60,3160.18,
1709        3119.36,3108.92,3164.77,3133.64,3111.47,3131.98,3128.87,3135.56 };
1710      PMTsat.Set(48,sat);
1711    }
1712    
1713    //------------------------------------------------------------------------
1714    void ToFdEdx::ReadParBBpos( const char *fname )
1715    {
1716      //  printf("read %s\n",fname);
1717      parBBpos.Set(48);
1718      FILE *fattin = fopen( fname , "r" );
1719      for (int i=0; i<48; i++) {
1720        int   tid=0;
1721        float  tp;
1722        if(fscanf(fattin,"%d %f",
1723                  &tid, &tp )!=2) break;
1724        parBBpos[i]=tp;
1725      }
1726      fclose(fattin);
1727    }
1728    
1729    //------------------------------------------------------------------------
1730    void ToFdEdx::ReadParDesatBB( const char *fname )
1731    {
1732      //  printf("read %s\n",fname);
1733      FILE *fattin = fopen( fname , "r" );
1734      for (int i=0; i<48; i++) {
1735        int   tid=0;
1736        float  tp[3];
1737        if(fscanf(fattin,"%d %f %f %f",
1738                  &tid, &tp[0], &tp[1], &tp[2] )!=4) break;
1739        parDesatBB[i].Set(3,tp);
1740      }
1741      fclose(fattin);
1742    }
1743    
1744    
1745    //------------------------------------------------------------------------
1746    void ToFdEdx::ReadParBBneg( const char *fname )
1747    
1748    {
1749      //  printf("read %s\n",fname);
1750      FILE *fattin = fopen( fname , "r" );
1751      for (int i=0; i<48; i++) {
1752        int   tid=0;
1753        float  tp[3];
1754        if(fscanf(fattin,"%d %f %f %f",
1755                  &tid, &tp[0], &tp[1], &tp[2] )!=4) break;
1756        parBBneg[i].Set(3,tp);
1757      }
1758      fclose(fattin);
1759    }
1760    
1761    //------------------------------------------------------------------------
1762    void ToFdEdx::ReadParPos( const char *fname )
1763    {
1764      //  printf("read %s\n",fname);
1765      FILE *fattin = fopen( fname , "r" );
1766      for (int i=0; i<48; i++) {
1767        int   tid=0;
1768        float  tp[4];
1769        if(fscanf(fattin,"%d %f %f %f %f",
1770                  &tid, &tp[0], &tp[1], &tp[2], &tp[3])!=5) break;
1771        parPos[i].Set(4,tp);
1772      }
1773      fclose(fattin);
1774    }
1775    
1776    //------------------------------------------------------------------------
1777    void ToFdEdx::ReadParAtt( const char *fname )
1778    {
1779      //  printf("read %s\n",fname);
1780      FILE *fattin = fopen( fname , "r" );
1781      for (int i=0; i<48; i++) {
1782        int   tid=0;
1783        float  tp[6];
1784        if(fscanf(fattin,"%d %f %f %f %f %f %f",
1785                  &tid, &tp[0], &tp[1], &tp[2], &tp[3], &tp[4], &tp[5] )!=7) break;
1786        parAtt[i].Set(6,tp);
1787      }
1788      fclose(fattin);
1789    }
1790    
1791    
1792    
1793    
1794    
1795    
1796    double ToFdEdx::f_att( TArrayF &p, float x )
1797    {
1798      return
1799        p[0] +
1800        p[1]*x +
1801        p[2]*x*x +
1802        p[3]*x*x*x +
1803        p[4]*x*x*x*x +
1804        p[5]*x*x*x*x*x;
1805    }
1806    //------------------------------------------------------------------------
1807    double ToFdEdx::f_att5B( float x )
1808    {
1809      return
1810        101.9409 +
1811        6.643781*x +
1812        0.2765518*x*x +
1813        0.004617647*x*x*x +
1814        0.0006195132*x*x*x*x +
1815        0.00002813734*x*x*x*x*x;
1816    }
1817    
1818    
1819    double ToFdEdx::f_pos( TArrayF &p, float x )
1820    {
1821      return
1822        p[0] +
1823        p[1]*x +
1824        p[2]*x*x +
1825        p[3]*x*x*x;
1826    }
1827    
1828    double ToFdEdx::f_pos5B( float x )
1829    {
1830      return
1831        15.45132 +
1832        0.8369721*x +
1833        0.0005*x*x;
1834    }
1835    
1836    
1837    
1838    double ToFdEdx::f_adcPC( float x )
1839    {
1840      return 28.12+0.6312*x-5.647e-05*x*x+3.064e-08*x*x*x;
1841    }
1842    
1843    
1844    float ToFdEdx::Get_adc_he( int id, float pl_x[6], float pl_y[6])
1845    {
1846    
1847      //
1848      // input: id - pmt [0:47}
1849      //             pl_x - coord x of the tof plane
1850      //             pl_y - coord y
1851    
1852      adc_he = 0;
1853      if( eGeom.GetXY(id)==1 )  adc_he = f_att( (parAtt[id]), pl_x[eGeom.GetPlane(id)] );
1854      if( eGeom.GetXY(id)==2 )  adc_he = f_att( (parAtt[id]), pl_y[eGeom.GetPlane(id)] );
1855      return adc_he;
1856    }
1857    
1858    //------------------------------------------------------------------------
1859    double ToFdEdx::f_BB( TArrayF &p, float x )
1860    {
1861      return  p[0]/(x*x)*(log(x*x/(1-x*x)) - p[1]*x*x - p[2]);
1862    }
1863    
1864    //------------------------------------------------------------------------
1865    double ToFdEdx::f_BB5B( float x )
1866    {
1867      return  0.165797/(x*x)*(log(x*x/(1-x*x)) + 140.481*x*x + 52.9258);
1868    }
1869    //------------------------------------------------------------------------
1870    double ToFdEdx::f_desatBB( TArrayF &p, float x )
1871    {
1872      return
1873        p[0] +
1874        p[1]*x +
1875        p[2]*x*x;
1876    }
1877    
1878    //------------------------------------------------------------------------
1879    double ToFdEdx::f_desatBB5B( float x )
1880    {
1881      return
1882        -2.4 +
1883        0.75*x +
1884        0.009*x*x;
1885    }
1886    
1887    
1888    
1889    
1890    
1891    
1892    
1893    
1894    
1895    
1896    
1897    
1898    
1899    
1900    
1901    
Legend:
Removed from v.1.3  
changed lines
  Added in v.1.34
  ViewVC Help
Powered by ViewVC 1.1.23