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

Contents of /DarthVader/ToFLevel2/src/ToFLevel2.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.48 - (show annotations) (download)
Wed May 20 10:57:02 2015 UTC (9 years, 6 months ago) by pam-fi
Branch: MAIN
CVS Tags: HEAD
Changes since 1.47: +10 -9 lines
Error occurred while calculating annotation data.
previous two commits were wrong, sorry... I put everything back

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

  ViewVC Help
Powered by ViewVC 1.1.23