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

Annotation of /DarthVader/ToFLevel2/src/ToFLevel2.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.45 - (hide annotations) (download)
Tue Mar 24 11:09:54 2015 UTC (9 years, 8 months ago) by pam-fi
Branch: MAIN
Changes since 1.44: +23 -9 lines
code modified to implement extended tracks

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

  ViewVC Help
Powered by ViewVC 1.1.23