/[PAMELA software]/DarthVader/TrackerLevel2/src/TrkLevel1.cpp
ViewVC logotype

Contents of /DarthVader/TrackerLevel2/src/TrkLevel1.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.27 - (show annotations) (download)
Thu Feb 27 11:24:43 2014 UTC (10 years, 9 months ago) by pam-fi
Branch: MAIN
Changes since 1.26: +10 -4 lines
Added new tracking algorythm

1 /**
2 * \file TrkLevel1.cpp
3 * \author Elena Vannuccini
4 */
5 #include <TrkLevel1.h>
6 #include <iostream>
7 using namespace std;
8 //......................................
9 // F77 routines
10 //......................................
11 extern "C" {
12
13 // int readetaparam_();
14 float cog_(int*,int*);
15 float pfaeta_(int*,float*);
16 float pfaeta2_(int*,float*);
17 float pfaeta3_(int*,float*);
18 float pfaeta4_(int*,float*);
19 float pfaetal_(int*,float*);
20 float digsat_(int*);
21 int npfastrips_(int*,float*);
22
23 }
24 //--------------------------------------
25 //
26 //
27 //--------------------------------------
28 TrkCluster::TrkCluster(){
29
30 // cout << "TrkCluster::TrkCluster()"<<endl;
31 view = -1;
32 maxs = -1;
33 indmax = -1;
34
35 CLlength = 0;
36 clsignal = 0;
37 clsigma = 0;
38 cladc = 0;
39 clbad = 0;
40
41 };
42 //--------------------------------------
43 //
44 //
45 //--------------------------------------
46 TrkCluster::TrkCluster(const TrkCluster& t){
47
48 view = t.view;
49 maxs = t.maxs;
50 indmax = t.indmax;
51
52 CLlength = t.CLlength;
53 if(CLlength){
54 clsignal = new Float_t[CLlength];
55 clsigma = new Float_t[CLlength];
56 cladc = new Int_t[CLlength];
57 clbad = new Bool_t[CLlength];
58 for(Int_t i=0; i<CLlength;i++){
59 clsignal[i] = t.clsignal[i];
60 clsigma[i] = t.clsigma[i];
61 cladc[i] = t.cladc[i];
62 clbad[i] = t.clbad[i];
63 };
64 };
65 };
66 //--------------------------------------
67 //
68 //
69 //--------------------------------------
70 void TrkCluster::Clear(){
71
72 // cout << "void TrkCluster::Clear()"<<endl;
73 if(CLlength){
74 delete [] clsignal;
75 delete [] clsigma;
76 delete [] cladc;
77 delete [] clbad;
78 }
79
80 view = 0;
81 maxs = 0;
82 indmax = 0;
83
84 CLlength = 0;
85 clsignal = 0;
86 clsigma = 0;
87 cladc = 0;
88 clbad = 0;
89
90 };
91 //--------------------------------------
92 //
93 //
94 //--------------------------------------
95 /**
96 * Evaluate the cluster signal including a maximum number of adjacent
97 * strips, around maxs, having a significant signal.
98 * @param nstrip Maximum number of strips.
99 * @param cut Inclusion cut ( s > cut*sigma ).
100 * @param force Falg to force the PFA strip-inclusion pattern (nstrip>0)
101 * If nstrip<=0 only the inclusion cut is used to determine the cluster size.
102 */
103 Float_t TrkCluster::GetSignal(Int_t nstrip, Float_t cut, Bool_t force){
104
105 if(CLlength<=0)return 0;
106
107 Float_t s = 0;
108
109 //-----------------------------------
110 // inlcude strips with s > cut*sigma
111 //-----------------------------------
112
113 if( nstrip<=0 ){
114 // for(Int_t is = 0; is < CLlength; is++){
115 // Float_t scut = cut*clsigma[is];
116 // if(clsignal[is] > scut) s += clsignal[is];
117 // };
118 for(Int_t is = indmax+1; is < CLlength; is++){
119 Float_t scut = cut*clsigma[is];
120 if(clsignal[is] > scut) s += clsignal[is];
121 else break;
122 };
123 for(Int_t is = indmax; is >=0; is--){
124 Float_t scut = cut*clsigma[is];
125 if(clsignal[is] > scut) s += clsignal[is];
126 else break;
127 };
128 return s;
129 };
130
131 //---------------------------------------------------
132 // inlcude strips with s > cut*sigma, up to nstrip.
133 // strips are included in order of decreasing signal
134 //---------------------------------------------------
135 if( !force ){
136
137 Int_t il = indmax;
138 Int_t ir = indmax;
139 Int_t inc = 0;
140
141 if( clsignal[indmax] < cut*clsigma[indmax] ) return 0;
142
143 while ( inc < nstrip ){
144 Float_t sl = -100000;
145 Float_t sr = -100000;
146 if( il >= 0 ) sl = clsignal[il];
147 if( ir < CLlength ) sr = clsignal[ir];
148 if( sl == sr && inc == 0 ){
149 s += clsignal[il]; //cout << inc<<" - "<< clsignal[il]<<" "<<s<<endl;
150 il--;
151 ir++;
152 }else if ( sl >= sr && sl > cut*clsigma[il] && inc !=0 ){
153 s += sl;//cout << inc<<" - "<< clsignal[il]<<" "<<s<<endl;
154 il--;
155 }else if ( sl < sr && sr > cut*clsigma[ir] ){
156 s += sr;//cout << inc<<" - " << clsignal[ir]<<" "<<s<<endl;
157 ir++;
158 }else break;
159
160 inc++;
161 }
162 return s;
163
164 }else{
165 //---------------------------------------------------
166 // evaluate signal using a fixed number of strips,
167 // following the PFA inclusion patters
168 //---------------------------------------------------
169 // --> signal of the central strip
170 Float_t sc = clsignal[indmax];
171 // signal of adjacent strips
172 Float_t sl1 = -9999.;
173 Float_t sl2 = -9999.;
174 Float_t sr1 = -9999.;
175 Float_t sr2 = -9999.;
176 if(indmax-1>=0) sl1 = clsignal[indmax-1];
177 if(indmax-2>=0) sl2 = clsignal[indmax-2];
178 if(indmax+1<CLlength) sr1 = clsignal[indmax+1];
179 if(indmax+2<CLlength) sr2 = clsignal[indmax+2];
180
181 if(nstrip==1){
182 s = sc;
183 }else if(nstrip==2){
184 if( sl1>sr1 && sl1+sc!=0 )s = (sl1+sc);
185 if( sl1<sr1 && sr1+sc!=0 )s = (sc+sr1);
186 if( sl1==sr1 && sl1 != -9999.){
187 if( clsigma[indmax-1] < clsigma[indmax+1] && sl1+sc!=0 )s = (sl1+sc);
188 if( clsigma[indmax-1] > clsigma[indmax+1] && sc+sr1!=0 )s = (sc+sr1);
189 }
190 }else if(nstrip==3){
191 s = (sl1+sc+sr1);
192 }else if(nstrip==4){
193 if( sl2>sr2 && sl2+sl1+sc+sr1!=0 )s = (sl2+sl1+sc+sr1);
194 if( sl2<sr2 && sl1+sc+sr1+sr2!=0 )s = (sl1+sc+sr1+sr2);
195 if( sl2==sr2 && sl2 != -9999.){
196 if( clsigma[indmax-2] < clsigma[indmax+2] && sl2+sl1+sc+sr1!=0 )s = (sl2+sl1+sc+sr1);
197 if( clsigma[indmax-2] > clsigma[indmax+2] && sl1+sc+sr1+sr2!=0 )s = (sl1+sc+sr1+sr2);
198 }
199 }else if(nstrip==5){
200 s = (sl1+sc+sr1);
201 if(sl2 != -9999.)s += sl2;
202 if(sr2 != -9999.)s += sr2;
203 }else{
204 cout << "Float_t TrkCluster::GetSignal("<<nstrip<<","<<cut<<","<<force<<")- not implemented"<<endl;
205 }
206
207 }
208
209 return 0.;
210
211 };
212
213
214 /**
215 * Evaluate the cluster signal-to-noise, as defined by Turchetta, including a
216 * maximum number of adjacent strips, around maxs, having a significant signal.
217 * @param nstrip Maximum number of strips.
218 * @param cut Inclusion cut ( s > cut*sigma ).
219 * If nstrip<=0 only the inclusion cut is used to determine the cluster size.
220 */
221 Float_t TrkCluster::GetSignalToNoise(Int_t nstrip, Float_t cut){
222
223 if(CLlength<=0)return 0;
224
225 Float_t sn = 0;
226
227 if( nstrip<=0 ){
228 for(Int_t is = indmax+1; is < CLlength; is++){
229 Float_t scut = cut*clsigma[is];
230 if(clsignal[is] > scut) sn += clsignal[is]/clsigma[is];
231 else break;
232 };
233 for(Int_t is = indmax; is >=0; is--){
234 Float_t scut = cut*clsigma[is];
235 if(clsignal[is] > scut) sn += clsignal[is]/clsigma[is];
236 else break;
237 };
238 return sn;
239 };
240
241
242 Int_t il = indmax;
243 Int_t ir = indmax;
244 Int_t inc = 0;
245
246 if( clsignal[indmax] < cut*clsigma[indmax] ) return 0;
247
248 while ( inc < nstrip ){
249 Float_t sl = -100000;
250 Float_t sr = -100000;
251 if( il >= 0 ) sl = clsignal[il];
252 if( ir < CLlength ) sr = clsignal[ir];
253 if( sl == sr && inc == 0 ){
254 sn += clsignal[il]/clsigma[il];
255 il--;
256 ir++;
257 }else if ( sl >= sr && sl > cut*clsigma[il] && inc !=0 ){
258 sn += sl/clsigma[il];
259 il--;
260 }else if ( sl < sr && sr > cut*clsigma[ir] ){
261 sn += sr/clsigma[ir];
262 ir++;
263 }else break;
264
265 inc++;
266 }
267 return sn;
268 };
269 /**
270 * Evaluate the cluster multiplicity.
271 * @param cut Inclusion cut.
272 */
273 Int_t TrkCluster::GetMultiplicity(Float_t cut){
274
275 if(CLlength<=0)return 0;
276
277 Int_t m = 0;
278
279 for(Int_t is = indmax+1; is < CLlength; is++){
280 Float_t scut = cut*clsigma[is];
281 if(clsignal[is] > scut) m++;
282 else break;
283 };
284 for(Int_t is = indmax; is >=0; is--){
285 Float_t scut = cut*clsigma[is];
286 if(clsignal[is] > scut) m++;
287 else break;
288 };
289 return m;
290 };
291 /**
292 * True if the cluster contains bad strips.
293 * @param nbad Number of strips around the maximum.
294 */
295 Bool_t TrkCluster::IsBad(Int_t nbad){
296
297 if(CLlength<=0)return 0;
298
299 Int_t il,ir;
300 il = indmax;
301 ir = indmax;
302 for(Int_t i=1; i<nbad; i++){
303 if (ir == CLlength-1 && il == 0)break;
304 else if (ir == CLlength-1 && il != 0)il--;
305 else if (ir != CLlength-1 && il == 0)ir++;
306 else{
307 if(clsignal[il-1] > clsignal[ir+1])il--;
308 else ir++;
309 }
310 }
311 Int_t isbad = 0;
312 for(Int_t i=il; i<=ir; i++)isbad += clbad[i];
313
314 return ( isbad != nbad );
315 };
316 /**
317 * True if the cluster contains saturated strips.
318 * @param nbad Number of strips around the maximum.
319 */
320 Bool_t TrkCluster::IsSaturated(Int_t nbad){
321
322 if(CLlength<=0)return 0;
323
324 Int_t il,ir;
325 il = indmax;
326 ir = indmax;
327 for(Int_t i=1; i<nbad; i++){
328 if (ir == CLlength-1 && il == 0)break;
329 else if (ir == CLlength-1 && il != 0)il--;
330 else if (ir != CLlength-1 && il == 0)ir++;
331 else{
332 if(clsignal[il-1] > clsignal[ir+1])il--;
333 else ir++;
334 }
335 }
336 Int_t isbad = 0;
337 for(Int_t i=il; i<=ir; i++){
338 if( IsX() && cladc[i] > 2980 )isbad++;
339 if( IsY() && cladc[i] < 80 )isbad++;
340 }
341 return ( isbad != 0 );
342
343 }
344 //--------------------------------------
345 //
346 //
347 //--------------------------------------
348 void TrkCluster::Dump(){
349
350 cout << "----- Cluster" << endl;
351 cout << "View "<<view << " - Ladder "<<GetLadder()<<endl;
352 cout << "Position of maximun "<< maxs <<endl;
353 cout << "Multiplicity "<< GetMultiplicity() <<endl;
354 cout << "Tot signal "<< GetSignal() << " (ADC channels)"<<endl ;
355 cout << "Signal/Noise "<< GetSignalToNoise()<<endl;
356 cout << "COG "<< GetCOG(0)<<endl;;
357 cout << "Strip signals ";
358 for(Int_t i =0; i<CLlength; i++)cout << " " <<clsignal[i];
359 cout <<endl<< "Strip sigmas ";
360 for(Int_t i =0; i<CLlength; i++)cout << " " <<clsigma[i];
361 cout <<endl<< "Strip ADC ";
362 for(Int_t i =0; i<CLlength; i++)cout << " " <<cladc[i];
363 cout <<endl<< "Strip BAD ";
364 for(Int_t i =0; i<CLlength; i++){
365 if(i==indmax)cout << " *" <<clbad[i]<<"*";
366 else cout << " " <<clbad[i];
367 }
368 cout << endl;
369
370 }
371 //--------------------------------------
372 //
373 //
374 //--------------------------------------
375 /**
376 * Method to fill a level1 struct with only one cluster (done to use F77 p.f.a. routines on a cluster basis).
377 */
378 void TrkCluster::GetLevel1Struct(cTrkLevel1* l1){
379
380 // cTrkLevel1* l1 = new cTrkLevel1;
381
382 // cTrkLevel1* l1 = &level1event_ ;
383
384 l1->nclstr1 = 1;
385 l1->view[0] = view;
386 l1->ladder[0] = GetLadder();
387 l1->maxs[0] = maxs;
388 l1->mult[0] = GetMultiplicity();
389 l1->dedx[0] = GetSignal();
390 l1->indstart[0] = 1;
391 l1->indmax[0] = indmax+1;
392 l1->totCLlength = CLlength;
393 for(Int_t i=0; i<CLlength; i++){
394 l1->clsignal[i] = clsignal[i];
395 l1->clsigma[i] = clsigma[i];
396 l1->cladc[i] = cladc[i];
397 l1->clbad[i] = clbad[i];
398 };
399
400 // return l1;
401 };
402 //--------------------------------------
403 //
404 //
405 //--------------------------------------
406 /**
407 * Evaluates the Center-Of-Gravity (COG) of the cluster, in strips, relative to the strip with the maximum signal (TrkCluster::maxs).
408 * @param ncog Number of strips to evaluate COG.
409 * If ncog=0, the COG of the cluster is evaluated according to the cluster multiplicity (defined by the inclusion cut).
410 * If ncog>0, the COG is evaluated using ncog strips, even if they have a negative signal (according to G.Landi)
411 *
412 * (NB TrkCluster::GetLevel1Struct() showld be called first, in order to fill the F77 level1 common with this single cluster)
413 */
414 Float_t TrkCluster::GetCOG(Int_t ncog){
415
416 int ic = 1;
417 // GetLevel1Struct(); //Elena: dangerous...
418 return cog_(&ncog,&ic);
419
420 };
421 /**
422 * Evaluates the Center-Of-Gravity (COG) of the cluster, in strips, relative to the strip with the maximum signal (TrkCluster::maxs),
423 * choosing the number of strips according to the angle, as implemented for the eta-algorythm .
424 * @param angle Projected angle in degree.
425 */
426 Float_t TrkCluster::GetCOG(Float_t angle){
427
428 Int_t neta = 0;
429
430 // Float_t eta = GetETA(0,angle);
431 // for(neta=2; neta<10; neta++) if( eta == GetETA(neta,angle) ) break;
432 // if(eta != GetETA(neta,angle) )cout << "Attenzione!! pasticcio "<<endl;
433
434 if( view%2 ){ //Y
435 neta=2;
436 }else{ //X
437 if( fabs(angle) <= 10. ){
438 neta = 2;
439 }else if( fabs(angle) > 10. && fabs(angle) <= 15. ){
440 neta = 3;
441 }else{
442 neta = 4;
443 };
444 };
445
446 return GetCOG(neta);
447
448 };
449 //--------------------------------------
450 //
451 //
452 //--------------------------------------
453 /**
454 * Evaluates the cluster position, in pitch units, relative to the strip
455 * with the maximum signal (TrkCluster::maxs), by applying the non-linear
456 * ETA-algorythm.
457 * @param neta Number of strips to evaluate ETA.
458 * @param angle Projected (effective) angle between particle track and detector plane.
459 * @landi flag to apply Landi correction
460 * Implemented values of neta are 2,3,4. If neta=0, ETA2, ETA3 and ETA4 are applied according to the angle.
461 * (NB TrkCluster::GetLevel1Struct() showld be called first, in order to fill the F77 level1 common with this single cluster)
462 */
463 Float_t TrkCluster::GetETA(Int_t neta, float angle, bool landi){
464
465 // cout << "GetETA(neta,angle) "<< neta << " "<< angle;
466 // LoadPfaParam();
467
468 TrkParams::Load(4);
469 if( !TrkParams::IsLoaded(4) ){
470 cout << "Float_t TrkCluster::GetETA(Int_t neta, float angle, bool landi) --- ERROR --- p.f.a. parameters not loaded"<<endl;
471 return 0;
472 }
473
474 float ax = angle;
475 int ic = 1;
476 //GetLevel1Struct(); //Elena: dangerous...
477 if( neta == 0 && !landi) return pfaeta_(&ic,&ax);
478 else if(neta == 0 && landi ) return pfaetal_(&ic,&ax);
479 else if(neta == 2 ) return pfaeta2_(&ic,&ax);
480 else if(neta == 3 ) return pfaeta3_(&ic,&ax);
481 else if(neta == 4 ) return pfaeta4_(&ic,&ax);
482 else cout << "TrkCluster::GetETA("<<neta<<","<<angle<<","<<landi<<") not implemented\n";
483 return 0;
484
485 };
486
487 /**
488 * Evaluates the cluster position, in pitch units, relative to the strip
489 * with the maximum signal (TrkCluster::maxs), by applying the digital
490 * algorithm for saturated clusters.
491 *
492 * @return The cluster position (0 also if if no saturated strip is found).
493 *
494 * (NB TrkCluster::GetLevel1Struct() showld be called first, in order to fill the F77 level1 common with this single cluster)
495 */
496 Float_t TrkCluster::GetDigSat() {
497
498 // GetLevel1Struct(); //Elena: dangerous...
499 int ic = 1;
500 return digsat_(&ic);
501
502 }
503
504 /**
505 * Evaluates the cluster position, in pitch unit, relative to the strip with
506 * the maximum signal (TrkCluster::maxs), by applying the PFA set as default (see TrkParams).
507 * @param angle Projected (effective) angle between particle track and detector plane.
508 */
509 Float_t TrkCluster::GetPositionPU(float angle){
510
511 if ( TrkParams::GetPFA() == 0 )return GetETA(0,angle,false);
512 else if( TrkParams::GetPFA() == 2 )return GetETA(2,angle,false);
513 else if( TrkParams::GetPFA() == 3 )return GetETA(3,angle,false);
514 else if( TrkParams::GetPFA() == 4 )return GetETA(4,angle,false);
515 else if( TrkParams::GetPFA() == 5 )return GetETA(0,angle,true);
516 else if( TrkParams::GetPFA() == 10 )return GetCOG(0);
517 else if( TrkParams::GetPFA() == 11 )return GetCOG(1);
518 else if( TrkParams::GetPFA() == 12 )return GetCOG(2);
519 else if( TrkParams::GetPFA() == 13 )return GetCOG(3);
520 else if( TrkParams::GetPFA() == 14 )return GetCOG(4);
521 else cout << " TrkCluster::GetPositionPU(float "<<angle<<") -- WARNING -- PFA="<<TrkParams::GetPFA()<<" not implemented"<<endl;
522
523 return 0.;
524
525 }
526
527 /**
528 * Give the number of strip used to evaluate the cluster coordinate
529 * according to the p.f.a.
530 * It returns 0 when the COG is used (in this case the number of strip used
531 * equals the multiplicity).
532 * (NB TrkCluster::GetLevel1Struct() showld be called first, in order to fill the F77 level1 common with this single cluster)
533 */
534 Int_t TrkCluster::GetPFAstrips(float angle){
535
536 float ax = angle;
537 int ic = 1;
538 // GetLevel1Struct(); //Elena: dangerous...
539 return npfastrips_(&ic,&ax);
540
541 }
542
543 //--------------------------------------
544 //
545 //
546 //--------------------------------------
547 TrkLevel1::TrkLevel1(){
548
549 // cout << "TrkLevel1::TrkLevel1()"<<endl;
550 // Cluster = new TClonesArray("TrkCluster");
551 Cluster = 0;
552 for(Int_t i=0; i<12 ; i++){
553 good[i] = -1;
554 for(Int_t j=0; j<24 ; j++){
555 cn[j][i]=0;
556 cnn[j][i]=0;
557 };
558 };
559 // TrkParams::SetTrackingMode();
560 // TrkParams::SetPrecisionFactor();
561 // TrkParams::SetStepMin();
562 TrkParams::SetMiniDefault();
563 TrkParams::SetPFA();
564 }
565 //--------------------------------------
566 //
567 //
568 //--------------------------------------
569 void TrkLevel1::Set(){
570 if(!Cluster)Cluster = new TClonesArray("TrkCluster");
571 }
572 //--------------------------------------
573 //
574 //
575 //--------------------------------------
576 void TrkLevel1::Dump(){
577
578 cout<<"DSP status: ";
579 for(Int_t i=0; i<12 ; i++)cout<<good[i]<<" ";
580 cout<<endl;
581 cout<<"VA1 mask : "<<endl;
582 for(Int_t i=0; i<12 ; i++){
583 for(Int_t ii=0; ii<24 ; ii++){
584 Int_t mask = cnn[ii][i];
585 if(mask>0)mask=1;
586 cout<<mask<<" ";
587 }
588 cout <<endl;
589 }
590
591 if(!Cluster)return;
592 TClonesArray &t = *Cluster;
593 for(int i=0; i<this->nclstr(); i++) ((TrkCluster *)t[i])->Dump();
594
595 }
596 /**
597 * \brief Dump processing status
598 */
599 void TrkLevel1::StatusDump(int view){
600 cout << "DSP n. "<<view+1<<" (level1-)status: "<<hex<<showbase<<good[view]<<dec<<endl;
601 };
602 /**
603 * \brief Check event status
604 *
605 * Check the event status, according to a flag-mask given as input.
606 * Return true if the view passes the check.
607 *
608 * @param view View number (0-11)
609 * @param flagmask Mask of flags to check (eg. flagmask=0x111 no missing packet,
610 * no crc error, no software alarm)
611 *
612 * @see TrkLevel2 class definition to know how the status flag is defined
613 *
614 */
615 Bool_t TrkLevel1::StatusCheck(int view, int flagmask){
616
617 if( view<0 || view >= 12)return false;
618 return !(good[view]&flagmask);
619
620 };
621
622
623 //--------------------------------------
624 //
625 //
626 //--------------------------------------
627 /**
628 * Fills a TrkLevel1 object with values from a struct cTrkLevel1 (to get data from F77 common).
629 */
630 void TrkLevel1::SetFromLevel1Struct(cTrkLevel1 *l1, Bool_t full){
631
632 // cout << "void TrkLevel1::SetFromLevel1Struct(cTrkLevel1 *l1, Bool_t full)"<<endl;
633
634 Clear();
635 // ---------------
636 // *** CLUSTER ***
637 // ---------------
638 TrkCluster* t_cl = new TrkCluster();
639 if(!Cluster)Cluster = new TClonesArray("TrkCluster");
640 TClonesArray &t = *Cluster;
641 for(int i=0; i<l1->nclstr1; i++){
642
643 t_cl->Clear();
644 // if( full || (!full && l1->whichtrack[i]) ){
645
646 t_cl->view = l1->view[i];
647 t_cl->maxs = l1->maxs[i];
648
649 if( full || (!full && l1->whichtrack[i]) ){
650 t_cl->indmax = l1->indmax[i] - l1->indstart[i];
651 Int_t from = l1->indstart[i] -1;
652 Int_t to = l1->totCLlength ;
653 if(i != l1->nclstr1-1)to = l1->indstart[i+1] -1 ;
654 t_cl->CLlength = to - from ;
655
656 t_cl->clsignal = new Float_t[t_cl->CLlength];
657 t_cl->clsigma = new Float_t[t_cl->CLlength];
658 t_cl->cladc = new Int_t[t_cl->CLlength];
659 t_cl->clbad = new Bool_t[t_cl->CLlength];
660
661 Int_t index = 0;
662 for(Int_t is = from; is < to; is++ ){
663 t_cl->clsignal[index] = (Float_t) l1->clsignal[is];
664 t_cl->clsigma[index] = (Float_t) l1->clsigma[is];
665 t_cl->cladc[index] = (Int_t) l1->cladc[is];
666 t_cl->clbad[index] = (Bool_t) l1->clbad[is];
667 index++;
668 };
669 }
670 new(t[i]) TrkCluster(*t_cl); // <<< store cluster
671 };
672
673 delete t_cl;
674
675 // -------------------------
676 // ****general variables****
677 // -------------------------
678 for(Int_t i=0; i<12 ; i++){
679 good[i] = l1->good[i];
680 for(Int_t j=0; j<24 ; j++){
681 cn[j][i] = l1->cnev[j][i];
682 // cnrms[j][i] = l1->cnrmsev[j][i];
683 cnn[j][i] = l1->cnnev[j][i];
684 };
685 };
686
687 }
688 /**
689 * Fills a struct cTrkLevel1 with values from a TrkLevel1 object (to put data into a F77 common).
690 */
691
692 void TrkLevel1::GetLevel1Struct(cTrkLevel1* l1) {
693
694 // cTrkLevel1* l1 = &level1event_ ;
695
696 for(Int_t i=0; i<12 ; i++){
697 l1->good[i] = good[i];
698 for(Int_t j=0; j<24 ; j++){
699 l1->cnev[j][i] = cn[j][i] ;
700 l1->cnnev[j][i] = cnn[j][i] ;
701 l1->cnrmsev[j][i] = 0. ;
702 };
703 l1->fshower[i] = 0;
704 };
705
706 l1->nclstr1=0;
707 l1->totCLlength=0;
708 Int_t index=0;
709 if(Cluster){
710 Int_t i=0;
711 for(Int_t ii=0;ii<Cluster->GetEntries();ii++){
712 TrkCluster *clu = GetCluster(ii);
713 // ----------------------------------------
714 // attenzione!!
715 // se il cluster non e` salvato (view = 0)
716 // DEVE essere escluso dal common F77
717 // ----------------------------------------
718 if(clu->view != 0 ){
719 l1->view[i] = clu->view;
720 l1->ladder[i] = clu->GetLadder();
721 l1->maxs[i] = clu->maxs;
722 l1->mult[i] = clu->GetMultiplicity();
723 l1->dedx[i] = clu->GetSignal();
724 l1->indstart[i] = index+1;
725 l1->indmax[i] = l1->indstart[i] + clu->indmax;
726 l1->totCLlength += clu->CLlength;
727 for(Int_t iw=0; iw < clu->CLlength; iw++){
728 l1->clsignal[index] = clu->clsignal[iw];
729 l1->clsigma[index] = clu->clsigma[iw];
730 l1->cladc[index] = clu->cladc[iw];
731 l1->clbad[index] = clu->clbad[iw];
732 index++;
733 }
734 i++;
735 }
736 }
737 l1->nclstr1 = i;
738 }
739
740 // return l1;
741 }
742 //--------------------------------------
743 //
744 //
745 //--------------------------------------
746 void TrkLevel1::Clear(){
747
748 for(Int_t i=0; i<12 ; i++){
749 good[i] = -1;
750 for(Int_t j=0; j<24 ; j++){
751 cn[j][i] = 0;
752 cnn[j][i] = 0;
753 };
754 };
755 // if(Cluster)Cluster->Clear("C");
756 if(Cluster)Cluster->Delete();
757
758 }
759 //--------------------------------------
760 //
761 //
762 //--------------------------------------
763 void TrkLevel1::Delete(){
764
765 // Clear();
766 if(Cluster)Cluster->Delete();
767 if(Cluster)delete Cluster;
768
769 }
770 //--------------------------------------
771 //
772 //
773 //--------------------------------------
774 TrkCluster *TrkLevel1::GetCluster(int is){
775
776 if(!Cluster)return 0;
777 if(is >= nclstr()){
778 cout << "** TrkLevel1::GetCluster(int) ** Cluster "<< is << " does not exits! " << endl;
779 cout << "( Stored clusters nclstr() = "<< this->nclstr()<<" )" << endl;
780 return 0;
781 }
782
783 TClonesArray &t = *(Cluster);
784 TrkCluster *cluster = (TrkCluster*)t[is];
785 return cluster;
786 }
787
788
789 // int TrkLevel1::GetPfaNbinsAngle(){
790 // TrkParams::Load(4);
791 // if( !TrkParams::IsLoaded(4) ){
792 // cout << "int TrkLevel1::GetPfaNbinsAngle() --- ERROR --- p.f.a. parameters not loaded"<<endl;
793 // return 0;
794 // }
795 // return pfa_.nangbin;
796 // };
797
798 // int TrkLevel1::GetPfaNbinsETA(){
799 // TrkParams::Load(4);
800 // if( !TrkParams::IsLoaded(4) ){
801 // cout << "int TrkLevel1::GetPfaNbinsETA() --- ERROR --- p.f.a. parameters not loaded"<<endl;
802 // return 0;
803 // }
804 // return pfa_.netaval;
805 // };
806
807 // /**
808 // *
809 // *
810 // */
811 // float* TrkLevel1::GetPfaCoord(TString pfa, int nview, int nladder, int nang){
812
813 // TrkParams::Load(4);
814 // if( !TrkParams::IsLoaded(4) ){
815 // cout << "float* TrkLevel1::GetPfaCoord(TString pfa, int nview, int nladder, int nang) --- ERROR --- p.f.a. parameters not loaded"<<endl;
816 // return 0;
817 // }
818
819 // int nbins = GetPfaNbinsETA();
820 // if(!nbins)return 0;
821
822 // float *fcorr = new float [nbins];
823
824 // if(!pfa.CompareTo("ETA2",TString::kIgnoreCase)){
825 // for(int ib=0; ib<nbins; ib++){
826 // fcorr[ib] = pfa_.feta2[nang][nladder][nview][ib];
827 // cout << pfa_.eta2[nang][ib] << " - " << pfa_.feta2[nang][nladder][nview][ib]<<endl;;
828 // }
829 // }else if (!pfa.CompareTo("ETA3",TString::kIgnoreCase)){
830 // for(int ib=0; ib<nbins; ib++)fcorr[ib] = pfa_.feta3[nang][nladder][nview][ib];
831 // }else if (!pfa.CompareTo("ETA4",TString::kIgnoreCase)){
832 // for(int ib=0; ib<nbins; ib++)fcorr[ib] = pfa_.feta4[nang][nladder][nview][ib];
833 // }else{
834 // cout << pfa<<" pfa parameters not implemented "<<endl;
835 // return 0;
836 // }
837
838 // return fcorr;
839
840 // };
841
842 // float* TrkLevel1::GetPfaAbs(TString pfa, int nang){
843
844 // TrkParams::Load(4);
845 // if( !TrkParams::IsLoaded(4) ){
846 // cout << "float* TrkLevel1::GetPfaAbs(TString pfa, int nang) --- ERROR --- p.f.a. parameters not loaded"<<endl;
847 // return 0;
848 // }
849
850 // int nbins = GetPfaNbinsETA();
851 // if(!nbins)return 0;
852
853 // float *fcorr = new float [nbins];
854
855 // if(!pfa.CompareTo("ETA2",TString::kIgnoreCase)){
856 // for(int ib=0; ib<nbins; ib++)fcorr[ib] = pfa_.eta2[nang][ib];
857 // }else if (!pfa.CompareTo("ETA3",TString::kIgnoreCase)){
858 // for(int ib=0; ib<nbins; ib++)fcorr[ib] = pfa_.eta3[nang][ib];
859 // }else if (!pfa.CompareTo("ETA4",TString::kIgnoreCase)){
860 // for(int ib=0; ib<nbins; ib++)fcorr[ib] = pfa_.eta4[nang][ib];
861 // }else{
862 // cout << pfa<<" pfa parameters not implemented "<<endl;
863 // return 0;
864 // }
865
866 // return fcorr;
867
868 // };
869
870 /**
871 * Method to call the F77 routine that performs level1->level2 processing.
872 * The level2 output is stored in a common block, which can be retrieved
873 * by mean of the method TrkLevel2::SetFromLevel2Struct().
874 * NB If the TrkLevel1 object is readout from a tree, and the
875 * TrkLevel1::ProcessEvent(int pfa) is used to reprocess the event, attention
876 * should be payed to the fact that single clusters (clusters not associated
877 * with any track) might not be stored. Full reprocessing should be done starting
878 * from level0 data.
879 */
880 //int TrkLevel1::ProcessEvent(int pfa){
881 int TrkLevel1::ProcessEvent(){
882
883 // cout << "int TrkLevel1::ProcessEvent()" << endl;
884 TrkParams::Load( );
885 if( !TrkParams::IsLoaded() )return 0;
886
887 GetLevel1Struct();
888
889 // analysisflight_(&pfa);
890 // TrkParams::SetPFA(pfa);
891 analysisflight_();
892
893 return 1;
894
895 }
896
897
898 ClassImp(TrkLevel1);
899 ClassImp(TrkCluster);

  ViewVC Help
Powered by ViewVC 1.1.23