HLT/exa/fill_pp.C

   1 // $Id$
   2
   3
   4 /**
   5    fill_pp -> store track info in ntuple root file
   6    fill_good_tracks_pp -> store good track info in ntuple root file
   7    print_pp -> print selected events from ntuple root file
   8    plot_pp -> plot parts of selected events from ntuple root file
   9    eval_pp -> simple evaluation of impact parameter cut effects
  10    Author: C. Loizides
  11 */
  12
  13
  14 /* ---fill_pp---
  15    Fill ntuple with relevant pp event data
  16    from the tracks files:
  17
  18    Ev is the event number to use
  19    Path gives the path where you find the digit files (and possibly also l3transform.config)
  20    Rfile gives the root file where the information will be added
  21    Saveev will be Ev unless you specify something else (useful for storing the trigger event of pilup events)
  22 */
  23
  24 void fill_pp(Int_t ev=0,Char_t *path="./",Char_t *rfile="pptracks.root",Int_t saveev=-1)
  25 {
  26   AliL3Logger l;
  27   l.Set(AliL3Logger::kAll);
  28   l.UseStdout();
  29   //l.UseStream();
  30
  31   if(getenv("TRANSFORMCONFIGPATH")){
  32     AliL3Transform::Init(getenv("TRANSFORMCONFIGPATH"));
  33   } else AliL3Transform::Init(path);
  34
  35   TNtuple *ntuppel = new TNtuple("ntuppel","pptracks","pt:phi:eta:xvert:yvert:zvert:imp:nhits:px:py:pz:event:mc");
  36   Float_t meas[13];
  37
  38   TFile file(rfile,"UPDATE");
  39   Char_t name[1024];
  40
  41   AliL3Evaluate *eval=new AliL3Evaluate(path,63,63);
  42   eval->LoadData(ev,-1);
  43   eval->AssignIDs();
  44
  45   AliL3TrackArray *a=eval->GetTracks();
  46   a->Compress();
  47
  48   AliL3Track t;
  49   Float_t phi,pt,eta;
  50   Int_t id;
  51   Int_t nhits;
  52
  53   AliL3Vertex vertex;
  54   Double_t xc,yc,zc,impact;
  55
  56   Int_t ntracks=a.GetNTracks();
  57   for(Int_t i=0;i<ntracks;i++){
  58     t.Set(a.GetTrack(i));
  59
  60     t.CalculateHelix();
  61     t.GetClosestPoint(&vertex,xc,yc,zc);
  62
  63     nhits=t.GetNHits();
  64     if(nhits<63) continue;
  65     pt=t.GetPt();
  66     eta=t.GetPseudoRapidity();
  67     phi=t.GetPsi();
  68     id=t.GetMCid();
  69
  70     //impact=TMath::fabs(zc);
  71     impact=TMath::Sqrt(xc*xc+yc*yc+zc*zc);
  72
  73     meas[0]=pt;
  74     meas[1]=phi;
  75     meas[2]=eta;
  76     meas[3]=xc;
  77     meas[4]=yc;
  78     meas[5]=zc;
  79     meas[6]=impact;
  80     meas[7]=nhits;
  81     meas[8]=t.GetPx();
  82     meas[9]=t.GetPy();
  83     meas[10]=t.GetPz();
  84     meas[11]=ev;
  85     meas[12]=id;
  86
  87     ntuppel->Fill(meas);
  88
  89     cout << i << ": pt " << pt << " " << phi << " " << eta << " " << id << " nhits " << nhits << " impact " << impact << endl;
  90
  91     //if(id!=-1) cout << "mc " << id << " " << t.GetCenterX() << " " << t.GetCenterY() << " " << t.GetRadius() << " fst point " << t.GetFirstPointX() << " " << t.GetFirstPointY() << " " << t.GetFirstPointZ() << " tgl " << t.GetTgl() << endl;
  92   }
  93
  94   if(saveev==-1) saveev=ev;
  95   sprintf(name,"pptracks-%d",saveev);
  96   ntuppel->Write(name);
  97   file.Close();
  98
  99   delete ntuppel;
 100   delete eval;
 101 }
 102
 103 /*--fill_good_tracks_pp--
 104
 105   Fill good tracks taken from offline into ntuple root file for later evaluation.
 106
 107   Evs is the event to start with
 108   Eve is the last event to process
 109   Path gives the path where the offline good_particles_tpc files are
 110   Rfile specifies the recreated root file where the ntuples will be stored
 111 */
 112
 113 void fill_good_tracks_pp(Int_t evs=0,Int_t eve=100,Char_t *path="./",Char_t *rfile="pp-good-tracks.root")
 114 {
 115   TFile file(rfile,"RECREATE");
 116   if(!file.IsOpen()) return;
 117
 118   Float_t meas[10];
 119   Char_t name[1024];
 120
 121   for(Int_t i=evs;i<eve;i++){
 122     TNtuple *ntuppel = new TNtuple("ntuppel-goodtracks","good-pptracks","mc:pdg:px:py:pz:xvert:yvert:zvert:nhits:sector");
 123
 124     sprintf(name,"%s/good_tracks_tpc_%d",path,i);
 125     FILE *f=fopen(name,"r");
 126     cout << "Event " << i << " " << name << endl;
 127
 128     while(!feof(f)){
 129       Int_t mc,pdg,nhits,sector;
 130       Float_t px,py,pz,xv,yv,zv;
 131       fscanf(f,"%d %d %f %f %f %f %f %f %d %d\n",&mc,&pdg,&px,&py,&pz,&xv,&yv,&zv,&nhits,&sector);
 132
 133       if(nhits<63) continue;
 134       if((px==0)&&(py==0)&&(pz==0)) continue;
 135       if((xv==0)&&(yv==0)&&(zv==0)) continue;
 136
 137       //cout << mc << " " << pdg << " " << px << " " << py << " " << pz << " " << xv << " " << yv << " " << zv << " " << nhits << " " << sector << endl;
 138
 139       meas[0]=mc;
 140       meas[1]=pdg;
 141       meas[2]=px;
 142       meas[3]=py;
 143       meas[4]=pz;
 144       meas[5]=xv;
 145       meas[6]=yv;
 146       meas[7]=zv;
 147       meas[8]=nhits;
 148       meas[9]=sector;
 149
 150       ntuppel->Fill(meas);
 151     }
 152     fclose(f);
 153     sprintf(name,"good-pptracks-%d",i);
 154     ntuppel->Write(name);
 155     delete ntuppel;
 156   }
 157
 158   file.Close();
 159 }
 160
 161 /* ---print_pp---
 162    Looping over event reading root file containing ntuples
 163    and printing wanted information:
 164
 165    Pfile specifies the root file created with fill_pp
 166    Evi gives the first event to print
 167    Eve the last event to print
 168 */
 169
 170 void print_pp(Char_t *pfile, Int_t evi=0, Int_t eve=-1)
 171 {
 172   TFile file1=TFile(pfile,"READ");
 173
 174   if(!file1.IsOpen()) return;
 175
 176   //"pt:phi:eta:xvert:yvert:zvert:imp:nhits:px:py:pz:event:mc")
 177   Float_t pt1,phi1,eta1,xvert1,yvert1,zvert1,imp1,nhits1,px1,py1,pz1,event1,mc1; //pileup
 178   Float_t pt2,phi2,eta2,xvert2,yvert2,zvert2,imp2,nhits2,px2,py2,pz2,event2,mc2; //pp
 179
 180   Int_t nent1;
 181   Char_t dummy[1024];
 182   if(eve==-1) eve=evi;
 183
 184   for(Int_t i=evi;i<=eve;i++){ //loop over pileup-events
 185     sprintf(dummy,"pptracks-%d",i);
 186     TNtuple *ntuppel1 = (TNtuple*)file1.Get(dummy);
 187     ntuppel1->SetBranchAddress("pt",&pt1);
 188     ntuppel1->SetBranchAddress("phi",&phi1);
 189     ntuppel1->SetBranchAddress("eta",&eta1);
 190     ntuppel1->SetBranchAddress("xvert",&xvert1);
 191     ntuppel1->SetBranchAddress("yvert",&yvert1);
 192     ntuppel1->SetBranchAddress("zvert",&zvert1);
 193     ntuppel1->SetBranchAddress("imp",&imp1);
 194     ntuppel1->SetBranchAddress("px",&px1);
 195     ntuppel1->SetBranchAddress("py",&py1);
 196     ntuppel1->SetBranchAddress("pz",&pz1);
 197     ntuppel1->SetBranchAddress("event",&event1);
 198     ntuppel1->SetBranchAddress("mc",&mc1);
 199
 200     nent1=ntuppel1->GetEntries();
 201     cout << "Event: " << i << " " << nent1 << " tracks" << endl;
 202     for(Int_t j=0;j<nent1;j++){
 203       ntuppel1->GetEvent(j);
 204
 205       Float_t bt=TMath::Sqrt(xvert1*xvert1+yvert1*yvert1);
 206       Float_t bz=TMath::fabs(zvert1);
 207
 208       if(mc1!=-1)
 209       cout << j << ": " << event1 << " " << pt1 << " " << phi1 << " " << eta1 << " mc " << mc1 << " bt " << bt << " bz " << bz << endl;
 210
 211     }//pp track loop
 212   }//pileup loop
 213
 214   file1.Close();
 215 }
 216
 217 /* ---plot_pp---
 218    Plot some variables from the ntuples (currently impact parameter
 219    but easily changeable.
 220
 221    Pfile is the root file containing the ntuples.
 222    Evi and Eve are start and end event number to process
 223 */
 224
 225 void plot_pp(Char_t *pfile, Int_t evi=0, Int_t eve=100)
 226 {
 227   TH1F *hist1=new TH1F("h1","test",100,0,1000);
 228
 229   TFile file1=TFile(pfile,"READ");
 230
 231   if(!file1.IsOpen()) return;
 232
 233   //"pt:phi:eta:xvert:yvert:zvert:imp:nhits:px:py:pz:event:mc")
 234   Float_t pt1,phi1,eta1,xvert1,yvert1,zvert1,imp1,nhits1,px1,py1,pz1,event1,mc1; //pileup
 235   Float_t pt2,phi2,eta2,xvert2,yvert2,zvert2,imp2,nhits2,px2,py2,pz2,event2,mc2; //pp
 236
 237   Int_t nent1;
 238   Char_t dummy[1024];
 239
 240   for(Int_t i=evi;i<eve;i++){ //loop over pileup-events
 241     sprintf(dummy,"pptracks-%d",i);
 242     TNtuple *ntuppel1 = (TNtuple*)file1.Get(dummy);
 243     if(!ntuppel1){
 244       cerr << "Error getting ntuple from " << dummy<<endl;
 245       continue;
 246     }
 247     ntuppel1->SetBranchAddress("pt",&pt1);
 248     ntuppel1->SetBranchAddress("phi",&phi1);
 249     ntuppel1->SetBranchAddress("eta",&eta1);
 250     ntuppel1->SetBranchAddress("xvert",&xvert1);
 251     ntuppel1->SetBranchAddress("yvert",&yvert1);
 252     ntuppel1->SetBranchAddress("zvert",&zvert1);
 253     ntuppel1->SetBranchAddress("imp",&imp1);
 254     ntuppel1->SetBranchAddress("px",&px1);
 255     ntuppel1->SetBranchAddress("py",&py1);
 256     ntuppel1->SetBranchAddress("pz",&pz1);
 257     ntuppel1->SetBranchAddress("event",&event1);
 258     ntuppel1->SetBranchAddress("mc",&mc1);
 259
 260     nent1=ntuppel1->GetEntries();
 261     cout << "Event: " << i << " " << nent1 << " tracks" << endl;
 262     for(Int_t j=0;j<nent1;j++){
 263       ntuppel1->GetEvent(j);
 264
 265       if(TMath::fabs(eta1)>0.9) continue;
 266       if(pt1<0.1) continue;
 267
 268       Float_t bt=TMath::Sqrt(xvert1*xvert1+yvert1*yvert1);
 269       Float_t bz=TMath::fabs(zvert1);
 270       hist1->Fill(bz,1);
 271     }//pp track loop
 272   }//pileup loop
 273
 274   file1.Close();
 275
 276   hist1->Sumw2();
 277   hist1->Draw("e1p");
 278 }
 279
 280 /* ---eval_pp---
 281    Little macro to play with pp pileup evaluation and do some
 282    counting. See real trigger macro trigger.C instead.
 283
 284    Ev is the event number
 285    Path is the path to the track arrays.
 286 */
 287
 288 void eval_pp(Int_t ev=0,Char_t *path="./")
 289 {
 290   AliL3Logger l;
 291   l.Set(AliL3Logger::kAll);
 292   l.UseStdout();
 293   //l.UseStream();
 294
 295   AliL3Evaluate *eval=new AliL3Evaluate(path,63,63);
 296   eval->LoadData(ev,-1);
 297   eval->AssignIDs();
 298
 299   AliL3TrackArray *a=eval->GetTracks();
 300
 301   AliL3Track t;
 302   Float_t phi,pt,tgl,eta;
 303   Int_t id;
 304   Int_t nhits;
 305   Int_t abspos=0,pos=0,absneg=0,neg=0;
 306   //assume vertex around zero
 307   AliL3Vertex vertex;
 308   Double_t xc,yc,zc,impact;
 309
 310   Int_t nents=0;
 311   Int_t ntracks=a.GetNTracks();
 312   for(Int_t i=0;i<ntracks;i++){
 313     t.Set(a.GetTrack(i));
 314
 315     t.CalculateHelix();
 316     t.GetClosestPoint(&vertex,xc,yc,zc);
 317
 318     nhits=t.GetNHits();
 319     if(nhits<63) continue;
 320     pt=t.GetPt();
 321     tgl=t.GetTgl();
 322     phi=100;
 323     id=t.GetMCid();
 324     cout << id << endl;
 325     if(t.GetPy()!=0) phi=TMath::ATan(t.GetPy()/t.GetPx());
 326
 327     //impact=TMath::fabs(zc);
 328     impact=TMath::Sqrt(xc*xc+yc*yc+zc*zc);
 329
 330     cout << i << ": pt " << pt << " " << phi << " " << tgl << " " << id << " nhits " << nhits << " impact " << impact << endl;
 331
 332     if(id<0) absneg++;
 333     else abspos++;
 334     nents++;
 335     if(impact>10) continue;
 336     if(id<0) neg++;
 337     else pos++;
 338   }
 339
 340   cout << "Result: " << nents << " " << abspos << " " << absneg << endl;
 341   cout << "After cut: " << pos << " " << neg << endl;
 342 }
 343
 344 void simple_eff_pp(Char_t *pfile, Char_t *reffile, Int_t evi=0, Int_t eve=100,Int_t cut=-1)
 345 {
 346   //"pt:phi:eta:xvert:yvert:zvert:imp:nhits:px:py:pz:event:mc")
 347   Float_t pt1,phi1,eta1,xvert1,yvert1,zvert1,imp1,nhits1,px1,py1,pz1,event1,mc1; //pileup
 348   Float_t pt2,phi2,eta2,xvert2,yvert2,zvert2,imp2,nhits2,px2,py2,pz2,event2,mc2; //pp
 349
 350   Int_t countgoodevents=0;
 351   Int_t nent1,nent2;
 352   Char_t dummy[1024];
 353
 354   Int_t tnorm;
 355   Int_t fnorm;
 356   Int_t pnorm;
 357   Int_t tcounts=0;
 358   Int_t fcounts=0;
 359   Int_t pcounts=0;
 360   Float_t fmean=0;
 361   Float_t pmean=0;
 362
 363   TFile file1=TFile(pfile,"READ");
 364   TFile file2=TFile(reffile,"READ");
 365   if(!file1.IsOpen() || !file2.IsOpen()) return;
 366
 367   for(Int_t i=evi;i<eve;i++){ //loop over pileup-events
 368     sprintf(dummy,"pptracks-%d",i);
 369     TNtuple *ntuppel1 = (TNtuple*)file1.Get(dummy);
 370     ntuppel1->SetBranchAddress("pt",&pt1);
 371     ntuppel1->SetBranchAddress("phi",&phi1);
 372     ntuppel1->SetBranchAddress("eta",&eta1);
 373     ntuppel1->SetBranchAddress("xvert",&xvert1);
 374     ntuppel1->SetBranchAddress("yvert",&yvert1);
 375     ntuppel1->SetBranchAddress("zvert",&zvert1);
 376     ntuppel1->SetBranchAddress("imp",&imp1);
 377     ntuppel1->SetBranchAddress("px",&px1);
 378     ntuppel1->SetBranchAddress("py",&py1);
 379     ntuppel1->SetBranchAddress("pz",&pz1);
 380     ntuppel1->SetBranchAddress("event",&event1);
 381     ntuppel1->SetBranchAddress("mc",&mc1);
 382
 383     tnorm=0;
 384     fnorm=0;
 385     pnorm=0;
 386     nent1=ntuppel1->GetEntries();
 387     for(Int_t j=0;j<nent1;j++){
 388       ntuppel1->GetEvent(j);
 389
 390       if(j==0){ //get info of triggered event from file2
 391
 392         sprintf(dummy,"good-pptracks-%d",event1);
 393         TNtuple *ntuppel2 = (TNtuple*)file2.Get(dummy);
 394         Float_t nhits2,pdg2,sector2;
 395         if(ntuppel2){
 396           ntuppel2->SetBranchAddress("mc",&mc2);
 397           ntuppel2->SetBranchAddress("pdg",&pdg2);
 398           ntuppel2->SetBranchAddress("px",&px2);
 399           ntuppel2->SetBranchAddress("py",&py2);
 400           ntuppel2->SetBranchAddress("pz",&pz2);
 401           ntuppel2->SetBranchAddress("xvert",&xvert2);
 402           ntuppel2->SetBranchAddress("yvert",&yvert2);
 403           ntuppel2->SetBranchAddress("zvert",&zvert2);
 404           ntuppel2->SetBranchAddress("nhits",&nhits2);
 405           ntuppel2->SetBranchAddress("sector",&sector2);
 406           nent2=ntuppel2->GetEntries();
 407           for(Int_t k=0;k<nent2;k++){
 408             ntuppel2->GetEvent(k);
 409             if(mc2<0)continue; //dont count fake tracks
 410             //if(TMath::fabs(zvert2)>30) continue;
 411             pt2=TMath::Sqrt(px2*px2+py2*py2);
 412             if(nhits2<63) continue;
 413             if(pt2<0.1) continue;
 414             if((px2==0)&&(py2==0)&&(pz2==0)) continue;
 415             Float_t p=TMath::Sqrt(px2*px2+py2*py2+pz2*pz2);
 416             eta2=0.5*TMath::log((p+pz2)/(p-pz2));
 417             if(TMath::fabs(eta2)>0.9) continue;
 418             tnorm++;
 419           }
 420         } else {
 421           sprintf(dummy,"pptracks-%d",event1);
 422           TNtuple *ntuppel2 = (TNtuple*)file2.Get(dummy);
 423           if(ntuppel2){
 424             ntuppel2->SetBranchAddress("pt",&pt2);
 425             ntuppel2->SetBranchAddress("phi",&phi2);
 426             ntuppel2->SetBranchAddress("eta",&eta2);
 427             ntuppel2->SetBranchAddress("xvert",&xvert2);
 428             ntuppel2->SetBranchAddress("yvert",&yvert2);
 429             ntuppel2->SetBranchAddress("zvert",&zvert2);
 430             ntuppel2->SetBranchAddress("imp",&imp2);
 431             ntuppel2->SetBranchAddress("px",&px2);
 432             ntuppel2->SetBranchAddress("py",&py2);
 433             ntuppel2->SetBranchAddress("pz",&pz2);
 434             ntuppel2->SetBranchAddress("event",&event2);
 435             ntuppel2->SetBranchAddress("mc",&mc2);
 436
 437             nent2=ntuppel2->GetEntries();
 438             for(Int_t k=0;k<nent2;k++){
 439               ntuppel2->GetEvent(k);
 440               if(mc2<0)continue; //dont count fake tracks
 441               if(pt2<0.1) continue;
 442               //if(TMath::fabs(zvert2)>30) continue;
 443               if(TMath::fabs(eta2)>0.9) continue;
 444               tnorm++;
 445             }
 446           }
 447         }
 448       }
 449       if(tnorm==0){
 450         tnorm=0;
 451         break; //triggered event has to have one track at least
 452       }
 453
 454       if(pt1<0.1) continue;
 455       if(TMath::fabs(eta1>0.9)) continue;
 456
 457       if(mc1<-1) continue; //dont count fake tracks
 458       else if(mc1==-1) fnorm++; //count tracks from pileup
 459       else {
 460         if((cut>0)&&(zvert1>cut)) continue;
 461         pnorm++; //count real tracks
 462       }
 463     }//pp track loop
 464
 465     if(tnorm==0) continue; //break from above
 466     if(pnorm>tnorm) pnorm=tnorm; //prob. multiple found tracks
 467     countgoodevents++; //for normalization
 468     fcounts+=fnorm;
 469     pcounts+=pnorm;
 470     tcounts+=tnorm;
 471     pmean+=(Float_t)pnorm/(Float_t)tnorm;
 472     fmean+=(Float_t)fnorm/(Float_t)tnorm;
 473   }//pileup loop
 474
 475   file1.Close();
 476   file2.Close();
 477
 478   Float_t peff,feff;
 479   peff=(Float_t)pcounts/tcounts;
 480   feff=(Float_t)fcounts/tcounts;
 481   pmean/=countgoodevents;
 482   fmean/=countgoodevents;
 483   cout << "Events used: " << countgoodevents << endl;
 484   cout << "Good Eff " << peff << " " << pcounts << " " << tcounts << endl;
 485   cout << "Fake Eff " << feff << " " << fcounts << " " << tcounts << endl;
 486   cout << "Means are " << pmean << " " << fmean << endl;
 487 }
 488
 489 //////////////////////////
 490 // dont look below here //
 491 //////////////////////////
 492
 493
 494 /**
 495     Store certain amount of track
 496     arrays in a root file for later
 497     access - okay needs to many changes
 498     in AliL3TrackArray and AliL3Track
 499     classes so stopped working on it
 500     (see fill_pp instead).
 501 */
 502
 503 void store_pp(Int_t evs=0,Int_t eve=100,Char_t *path="./",Char_t *rfile="tracks.root")
 504 {
 505   AliL3Logger l;
 506   l.Set(AliL3Logger::kAll);
 507   l.UseStdout();
 508   //l.UseStream();
 509
 510   if(getenv("TRANSFORMCONFIGPATH")){
 511     AliL3Transform::Init(getenv("TRANSFORMCONFIGPATH"));
 512   }
 513
 514   TFile file(rfile,"UPDATE");
 515   Char_t name[1024];
 516
 517   for(Int_t i=evs;i<eve;i++){
 518     AliL3Evaluate *eval=new AliL3Evaluate(path,63,63);
 519     eval->LoadData(i,-1);
 520     eval->AssignIDs();
 521
 522     AliL3TrackArray *a=eval->GetTracks();
 523     sprintf(name,"tracks-%d",i);
 524     a->Compress();
 525     a->Write(name);
 526
 527     delete eval;
 528   }
 529
 530   file.Close();
 531 }
 532
 533 /**
 534     Loop over stored track arrays
 535     and print content. Does not properly
 536     work (see store_pp).
 537 */
 538
 539 void load_pp(Int_t evs=0,Int_t eve=100,Char_t *path="./")
 540 {
 541   AliL3Logger l;
 542   l.Set(AliL3Logger::kAll);
 543   l.UseStdout();
 544   //l.UseStream();
 545
 546   if(getenv("TRANSFORMCONFIGPATH")){
 547     AliL3Transform::Init(getenv("TRANSFORMCONFIGPATH"));
 548   } else {
 549     AliL3Transform::Init(path);
 550   }
 551
 552   for(Int_t j=evs;j<eve;j++){
 553
 554     AliL3Evaluate *eval=new AliL3Evaluate(path,63,63);
 555     eval->LoadData(j,-1); //load whole patch
 556     eval->AssignIDs();
 557
 558     AliL3TrackArray *a=eval->GetTracks();
 559
 560     AliL3Track t;
 561     Float_t phi,pt,tgl,eta;
 562     Int_t id;
 563     Int_t nhits;
 564
 565     AliL3Vertex vertex;
 566     Double_t xc,yc,zc,impact;
 567
 568     Int_t nents=0;
 569     Int_t ntracks=a.GetNTracks();
 570     for(Int_t i=0;i<ntracks;i++){
 571       t.Set(a.GetTrack(i));
 572
 573       t.CalculateHelix();
 574       t.GetClosestPoint(&vertex,xc,yc,zc);
 575
 576       nhits=t.GetNHits();
 577       if(nhits<63) continue;
 578       pt=t.GetPt();
 579       tgl=t.GetTgl();
 580       phi=100;
 581       id=t.GetMCid();
 582       if(t.GetPy()!=0) phi=TMath::ATan(t.GetPy()/t.GetPx());
 583
 584       impact=TMath::fabs(zc);
 585       //impact=TMath::Sqrt(xc*xc+yc*yc+zc*zc);
 586
 587       cout << i << ": pt " << pt << " " << phi << " " << tgl << " " << id << " nhits " << nhits << " impact " << impact << endl;
 588
 589       nents++;
 590     }
 591
 592     delete a;
 593   }
 594 }