HLT/TPCLib/tracking-ca/AliHLTTPCCAPerformance.cxx

   1 // $Id$
   2 //***************************************************************************
   3 // This file is property of and copyright by the ALICE HLT Project          *
   4 // ALICE Experiment at CERN, All rights reserved.                           *
   5 //                                                                          *
   6 // Primary Authors: Sergey Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de> *
   7 //                  Ivan Kisel <kisel@kip.uni-heidelberg.de>                *
   8 //                  for The ALICE HLT Project.                              *
   9 //                                                                          *
  10 // Permission to use, copy, modify and distribute this software and its     *
  11 // documentation strictly for non-commercial purposes is hereby granted     *
  12 // without fee, provided that the above copyright notice appears in all     *
  13 // copies and that both the copyright notice and this permission notice     *
  14 // appear in the supporting documentation. The authors make no claims       *
  15 // about the suitability of this software for any purpose. It is            *
  16 // provided "as is" without express or implied warranty.                    *
  17 //***************************************************************************
  18
  19 #include "AliHLTTPCCAPerformance.h"
  20 #include "AliHLTTPCCAGBHit.h"
  21 #include "AliHLTTPCCAMCTrack.h"
  22 #include "AliHLTTPCCAMCPoint.h"
  23 #include "AliHLTTPCCAOutTrack.h"
  24 #include "AliHLTTPCCAGBTrack.h"
  25 #include "AliHLTTPCCAGBTracker.h"
  26 #include "AliHLTTPCCATracker.h"
  27
  28
  29 #include "TMath.h"
  30 #include "TROOT.h"
  31 #include "Riostream.h"
  32 #include "TFile.h"
  33 #include "TH1.h"
  34 #include "TH2.h"
  35 #include "TProfile.h"
  36
  37
  38 ClassImp(AliHLTTPCCAPerformance)
  39
  40
  41
  42 AliHLTTPCCAPerformance::AliHLTTPCCAPerformance()
  43   : TObject(),
  44     fTracker(0),
  45     fHitLabels(0),
  46     fNHits(0),
  47     fMCTracks(0),
  48     fNMCTracks(0),
  49     fMCPoints(0),
  50     fNMCPoints(0),
  51     fDoClusterPulls(1),
  52     fStatNEvents(0),
  53     fStatNRecTot(0),
  54     fStatNRecOut(0),
  55     fStatNGhost(0),
  56     fStatNMCAll(0),
  57     fStatNRecAll(0),
  58     fStatNClonesAll(0),
  59     fStatNMCRef(0),
  60     fStatNRecRef(0),
  61     fStatNClonesRef(0),
  62     fStatGBNRecTot(0),
  63     fStatGBNRecOut(0),
  64     fStatGBNGhost(0),
  65     fStatGBNMCAll(0),
  66     fStatGBNRecAll(0),
  67     fStatGBNClonesAll(0),
  68     fStatGBNMCRef(0),
  69     fStatGBNRecRef(0),
  70     fStatGBNClonesRef(0),
  71     fHistoDir(0),
  72     fhResY(0),
  73     fhResZ(0),
  74     fhResSinPhi(0),
  75     fhResDzDs(0),
  76     fhResPt(0),
  77     fhPullY(0),
  78     fhPullZ(0),
  79     fhPullSinPhi(0),
  80     fhPullDzDs(0),
  81     fhPullQPt(0),
  82     fhHitErrY(0),
  83     fhHitErrZ(0),
  84     fhHitResY(0),
  85     fhHitResZ(0),
  86     fhHitPullY(0),
  87     fhHitPullZ(0),
  88     fhHitResY1(0),
  89     fhHitResZ1(0),
  90     fhHitPullY1(0),
  91     fhHitPullZ1(0),
  92     fhCellPurity(0),
  93     fhCellNHits(0),
  94     fhCellPurityVsN(0),
  95     fhCellPurityVsPt(0),
  96     fhEffVsP(0),
  97     fhGBEffVsP(0)
  98 {
  99   //* constructor
 100 }
 101
 102
 103 AliHLTTPCCAPerformance::AliHLTTPCCAPerformance(const AliHLTTPCCAPerformance&)
 104   : TObject(),
 105     fTracker(0),
 106     fHitLabels(0),
 107     fNHits(0),
 108     fMCTracks(0),
 109     fNMCTracks(0),
 110     fMCPoints(0),
 111     fNMCPoints(0),
 112     fDoClusterPulls(1),
 113     fStatNEvents(0),
 114     fStatNRecTot(0),
 115     fStatNRecOut(0),
 116     fStatNGhost(0),
 117     fStatNMCAll(0),
 118     fStatNRecAll(0),
 119     fStatNClonesAll(0),
 120     fStatNMCRef(0),
 121     fStatNRecRef(0),
 122     fStatNClonesRef(0),
 123     fStatGBNRecTot(0),
 124     fStatGBNRecOut(0),
 125     fStatGBNGhost(0),
 126     fStatGBNMCAll(0),
 127     fStatGBNRecAll(0),
 128     fStatGBNClonesAll(0),
 129     fStatGBNMCRef(0),
 130     fStatGBNRecRef(0),
 131     fStatGBNClonesRef(0),
 132     fHistoDir(0),
 133     fhResY(0),
 134     fhResZ(0),
 135     fhResSinPhi(0),
 136     fhResDzDs(0),
 137     fhResPt(0),
 138     fhPullY(0),
 139     fhPullZ(0),
 140     fhPullSinPhi(0),
 141     fhPullDzDs(0),
 142     fhPullQPt(0),
 143     fhHitErrY(0),
 144     fhHitErrZ(0),
 145     fhHitResY(0),
 146     fhHitResZ(0),
 147     fhHitPullY(0),
 148     fhHitPullZ(0),
 149     fhHitResY1(0),
 150     fhHitResZ1(0),
 151     fhHitPullY1(0),
 152     fhHitPullZ1(0),
 153     fhCellPurity(0),
 154     fhCellNHits(0),
 155     fhCellPurityVsN(0),
 156     fhCellPurityVsPt(0),
 157     fhEffVsP(0),
 158     fhGBEffVsP(0)
 159 {
 160   //* dummy
 161 }
 162
 163 AliHLTTPCCAPerformance &AliHLTTPCCAPerformance::operator=(const AliHLTTPCCAPerformance&)
 164 {
 165   //* dummy
 166   return *this;
 167 }
 168
 169 AliHLTTPCCAPerformance::~AliHLTTPCCAPerformance()
 170 {
 171   //* destructor
 172   StartEvent();
 173 }
 174
 175 void AliHLTTPCCAPerformance::SetTracker( AliHLTTPCCAGBTracker *Tracker )
 176 {
 177   //* set pointer to HLT CA Global tracker
 178   fTracker = Tracker;
 179 }
 180
 181 void AliHLTTPCCAPerformance::StartEvent()
 182 {
 183   //* clean up arrays
 184   if( !fHistoDir )  CreateHistos();
 185   if( fHitLabels ) delete[] fHitLabels;
 186   fHitLabels = 0;
 187   fNHits = 0;
 188   if( fMCTracks ) delete[] fMCTracks;
 189   fMCTracks = 0;
 190   fNMCTracks = 0;
 191   if( fMCPoints ) delete[] fMCPoints;
 192   fMCPoints = 0;
 193   fNMCPoints = 0;
 194 }
 195
 196 void AliHLTTPCCAPerformance::SetNHits( Int_t NHits )
 197 {
 198   //* set number of hits
 199   if( fHitLabels ) delete[] fHitLabels;
 200   fHitLabels = 0;
 201   fHitLabels = new AliHLTTPCCAHitLabel[ NHits ];
 202   fNHits = NHits;
 203 }
 204
 205 void AliHLTTPCCAPerformance::SetNMCTracks( Int_t NMCTracks )
 206 {
 207   //* set number of MC tracks
 208   if( fMCTracks ) delete[] fMCTracks;
 209   fMCTracks = 0;
 210   fMCTracks = new AliHLTTPCCAMCTrack[ NMCTracks ];
 211   fNMCTracks = NMCTracks;
 212 }
 213
 214 void AliHLTTPCCAPerformance::SetNMCPoints( Int_t NMCPoints )
 215 {
 216   //* set number of MC points
 217   if( fMCPoints ) delete[] fMCPoints;
 218   fMCPoints = 0;
 219   fMCPoints = new AliHLTTPCCAMCPoint[ NMCPoints ];
 220   fNMCPoints = 0;
 221 }
 222
 223 void AliHLTTPCCAPerformance::ReadHitLabel( Int_t HitID,
 224                                            Int_t lab0, Int_t lab1, Int_t lab2 )
 225 {
 226   //* read the hit labels
 227   AliHLTTPCCAHitLabel hit;
 228   hit.fLab[0] = lab0;
 229   hit.fLab[1] = lab1;
 230   hit.fLab[2] = lab2;
 231   fHitLabels[HitID] = hit;
 232 }
 233
 234 void AliHLTTPCCAPerformance::ReadMCTrack( Int_t index, const TParticle *part )
 235 {
 236   //* read mc track to the local array
 237   fMCTracks[index] = AliHLTTPCCAMCTrack(part);
 238 }
 239
 240 void AliHLTTPCCAPerformance::ReadMCTPCTrack( Int_t index, Float_t X, Float_t Y, Float_t Z,
 241                                              Float_t Px, Float_t Py, Float_t Pz )
 242 {
 243   //* read mc track parameters at TPC
 244   fMCTracks[index].SetTPCPar(X,Y,Z,Px,Py,Pz);
 245 }
 246
 247 void AliHLTTPCCAPerformance::ReadMCPoint( Int_t TrackID, Float_t X, Float_t Y, Float_t Z, Float_t Time, Int_t iSlice )
 248 {
 249   //* read mc point to the local array
 250   AliHLTTPCCAMCPoint &p = fMCPoints[fNMCPoints];
 251   p.TrackID() = TrackID;
 252   p.X() = X;
 253   p.Y() = Y;
 254   p.Z() = Z;
 255   p.Time() = Time;
 256   p.ISlice() = iSlice;
 257   fTracker->Slices()[iSlice].Param().Global2Slice( X, Y, Z,
 258                                                    &p.Sx(), &p.Sy(), &p.Sz() );
 259   if( X*X + Y*Y>10.) fNMCPoints++;
 260 }
 261
 262 void AliHLTTPCCAPerformance::CreateHistos()
 263 {
 264   //* create performance histogramms
 265   TDirectory *curdir = gDirectory;
 266   fHistoDir = gROOT->mkdir("HLTTPCCATrackerPerformance");
 267   fHistoDir->cd();
 268   gDirectory->mkdir("TrackFit");
 269   gDirectory->cd("TrackFit");
 270
 271   fhResY = new TH1D("resY", "track Y resoltion [cm]", 30, -.5, .5);
 272   fhResZ = new TH1D("resZ", "track Z resoltion [cm]", 30, -.5, .5);
 273   fhResSinPhi = new TH1D("resSinPhi", "track SinPhi resoltion ", 30, -.03, .03);
 274   fhResDzDs = new TH1D("resDzDs", "track DzDs resoltion ", 30, -.01, .01);
 275   fhResPt = new TH1D("resPt", "track telative Pt resoltion", 30, -.2, .2);
 276   fhPullY = new TH1D("pullY", "track Y pull", 30, -10., 10.);
 277   fhPullZ = new TH1D("pullZ", "track Z pull", 30, -10., 10.);
 278   fhPullSinPhi = new TH1D("pullSinPhi", "track SinPhi pull", 30, -10., 10.);
 279   fhPullDzDs = new TH1D("pullDzDs", "track DzDs pull", 30, -10., 10.);
 280   fhPullQPt = new TH1D("pullQPt", "track Q/Pt pull", 30, -10., 10.);
 281
 282   gDirectory->cd("..");
 283
 284   fhEffVsP = new TProfile("EffVsP", "Eff vs P", 100, 0., 5.);
 285   fhGBEffVsP = new TProfile("GBEffVsP", "Global tracker: Eff vs P", 100, 0., 5.);
 286
 287   gDirectory->mkdir("Clusters");
 288   gDirectory->cd("Clusters");
 289   fhHitResY = new TH1D("resHitY", "Y cluster resoltion [cm]", 100, -2., 2.);
 290   fhHitResZ = new TH1D("resHitZ", "Z cluster resoltion [cm]", 100, -2., 2.);
 291   fhHitPullY = new TH1D("pullHitY", "Y cluster pull", 50, -10., 10.);
 292   fhHitPullZ = new TH1D("pullHitZ", "Z cluster pull", 50, -10., 10.);
 293
 294   fhHitResY1 = new TH1D("resHitY1", "Y cluster resoltion [cm]", 100, -2., 2.);
 295   fhHitResZ1 = new TH1D("resHitZ1", "Z cluster resoltion [cm]", 100, -2., 2.);
 296   fhHitPullY1 = new TH1D("pullHitY1", "Y cluster pull", 50, -10., 10.);
 297   fhHitPullZ1 = new TH1D("pullHitZ1", "Z cluster pull", 50, -10., 10.);
 298
 299   fhHitErrY = new TH1D("HitErrY", "Y cluster error [cm]", 100, 0., 1.);
 300   fhHitErrZ = new TH1D("HitErrZ", "Z cluster error [cm]", 100, 0., 1.);
 301
 302   gDirectory->cd("..");
 303
 304   gDirectory->mkdir("Cells");
 305   gDirectory->cd("Cells");
 306   fhCellPurity = new TH1D("CellPurity", "Cell Purity", 100, -0.1, 1.1);
 307   fhCellNHits = new TH1D("CellNHits", "Cell NHits", 40, 0., 40.);
 308   fhCellPurityVsN = new TProfile("CellPurityVsN", "Cell purity Vs N hits", 40, 2., 42.);
 309   fhCellPurityVsPt = new TProfile("CellPurityVsPt", "Cell purity Vs Pt", 100, 0., 5.);
 310   gDirectory->cd("..");
 311
 312   curdir->cd();
 313 }
 314
 315 void AliHLTTPCCAPerformance::WriteDir2Current( TObject *obj )
 316 {
 317   //* recursive function to copy the directory 'obj' to the current one
 318   if( !obj->IsFolder() ) obj->Write();
 319   else{
 320     TDirectory *cur = gDirectory;
 321     TDirectory *sub = cur->mkdir(obj->GetName());
 322     sub->cd();
 323     TList *listSub = ((TDirectory*)obj)->GetList();
 324     TIter it(listSub);
 325     while( TObject *obj1=it() ) WriteDir2Current(obj1);
 326     cur->cd();
 327   }
 328 }
 329
 330 void AliHLTTPCCAPerformance::WriteHistos()
 331 {
 332   //* write histograms to the file
 333   TDirectory *curr = gDirectory;
 334   // Open output file and write histograms
 335   TFile* outfile = new TFile("HLTTPCCATrackerPerformance.root","RECREATE");
 336   outfile->cd();
 337   WriteDir2Current(fHistoDir);
 338   outfile->Close();
 339   curr->cd();
 340 }
 341
 342
 343 void AliHLTTPCCAPerformance::SlicePerformance( Int_t iSlice, Bool_t PrintFlag )
 344 {
 345   //* calculate slice tracker performance
 346   if( !fTracker ) return;
 347
 348   Int_t nRecTot = 0, nGhost=0, nRecOut=0;
 349   Int_t nMCAll = 0, nRecAll=0, nClonesAll=0;
 350   Int_t nMCRef = 0, nRecRef=0, nClonesRef=0;
 351   AliHLTTPCCATracker &slice = fTracker->Slices()[iSlice];
 352
 353   Int_t firstSliceHit = 0;
 354   for( ; firstSliceHit<fTracker->NHits(); firstSliceHit++){
 355     if( fTracker->Hits()[firstSliceHit].ISlice()==iSlice ) break;
 356   }
 357   Int_t endSliceHit = firstSliceHit;
 358
 359   for( ; endSliceHit<fTracker->NHits(); endSliceHit++){
 360     if( fTracker->Hits()[endSliceHit].ISlice()!=iSlice ) break;
 361   }
 362
 363   { // Cell construction performance
 364
 365     for( Int_t iRow=0; iRow<slice.Param().NRows(); iRow++ ){
 366       AliHLTTPCCARow &row = slice.Rows()[iRow];
 367       for (Int_t ic = 0; ic<row.NCells(); ic++){
 368         AliHLTTPCCACell &c  = row.Cells()[ic];
 369         Int_t *lb = new Int_t[c.NHits()*3];
 370         Int_t nla = 0;
 371         for( Int_t j=0; j<c.NHits(); j++){
 372           AliHLTTPCCAHit &h = row.GetCellHit(c,j);
 373           //cout<<"hit ID="<<h.ID()<<" of"<<fTracker->NHits()<<endl;
 374           //cout<<"gb hit ID="<<fTracker->Hits()[h.ID()].ID()<<endl;
 375           AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[h.ID()].ID()];
 376           if( l.fLab[0]>=0 ) lb[nla++]= l.fLab[0];
 377           if( l.fLab[1]>=0 ) lb[nla++]= l.fLab[1];
 378           if( l.fLab[2]>=0 ) lb[nla++]= l.fLab[2];
 379         }
 380         sort( lb, lb+nla );
 381         Int_t labmax = -1, labcur=-1, lmax = 0, lcurr=0;
 382         for( Int_t i=0; i<nla; i++ ){
 383           if( lb[i]!=labcur ){
 384             if( labcur>=0 && lmax<lcurr ){
 385               lmax = lcurr;
 386               labmax = labcur;
 387             }
 388             labcur = lb[i];
 389             lcurr = 0;
 390           }
 391           lcurr++;
 392         }
 393         if( labcur>=0 && lmax<lcurr ){
 394           lmax = lcurr;
 395           labmax = labcur;
 396         }
 397
 398         Int_t label = labmax;
 399         lmax = 0;
 400         for( Int_t j=0; j<c.NHits(); j++){
 401           AliHLTTPCCAHit &h = row.GetCellHit(c,j);
 402           AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[h.ID()].ID()];
 403           if( l.fLab[0]==label || l.fLab[1]==label || l.fLab[2]==label ) lmax++;
 404         }
 405
 406         nla = c.NHits();
 407         if( nla>0 && label>=0 ){
 408           double purity = double(lmax)/double(nla);
 409           fhCellPurity->Fill(purity);
 410           fhCellPurityVsN->Fill(c.NHits(),purity);
 411           fhCellPurityVsPt->Fill(fMCTracks[label].Pt(),purity);
 412         }
 413         fhCellNHits->Fill(c.NHits());
 414         if(lb) delete[] lb;
 415       }
 416     }
 417   }
 418
 419   // Select reconstructable MC tracks
 420
 421   {
 422     for (Int_t imc=0; imc<fNMCTracks; imc++) fMCTracks[imc].NHits() = 0;
 423
 424     for( Int_t ih=firstSliceHit; ih<endSliceHit; ih++){
 425       AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[ih].ID()];
 426       if( l.fLab[0]>=0 ) fMCTracks[l.fLab[0]].NHits()++;
 427       if( l.fLab[1]>=0 ) fMCTracks[l.fLab[1]].NHits()++;
 428       if( l.fLab[2]>=0 ) fMCTracks[l.fLab[2]].NHits()++;
 429     }
 430
 431     for (Int_t imc=0; imc<fNMCTracks; imc++) {
 432       AliHLTTPCCAMCTrack &mc = fMCTracks[imc];
 433       mc.Set() = 0;
 434       mc.NReconstructed() = 0;
 435       mc.NTurns() = 1;
 436       if( mc.NHits() >=  10 && mc.P()>=.05 ){
 437         mc.Set() = 1;
 438         nMCAll++;
 439         if( mc.P()>=1. ){
 440           mc.Set() = 2;
 441           nMCRef++;
 442         }
 443       }
 444     }
 445   }
 446
 447   Int_t traN = slice.NOutTracks();
 448   Int_t *traLabels = 0;
 449   Double_t *traPurity = 0;
 450   traLabels = new Int_t[traN];
 451   traPurity = new Double_t[traN];
 452   {
 453     for (Int_t itr=0; itr<traN; itr++) {
 454       traLabels[itr]=-1;
 455       traPurity[itr]= 0;
 456       AliHLTTPCCAOutTrack &tCA = slice.OutTracks()[itr];
 457       Int_t nhits = tCA.NHits();
 458       Int_t *lb = new Int_t[nhits*3];
 459       Int_t nla=0;
 460       //cout<<"\nHit labels:"<<endl;
 461       for( Int_t ihit=0; ihit<nhits; ihit++){
 462         Int_t index = slice.OutTrackHits()[tCA.FirstHitRef()+ihit];
 463         AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[index].ID()];
 464         //cout<<l.fLab[0]<<" "<<l.fLab[1]<<" "<<l.fLab[2]<<endl;
 465         if(l.fLab[0]>=0 ) lb[nla++]= l.fLab[0];
 466         if(l.fLab[1]>=0 ) lb[nla++]= l.fLab[1];
 467         if(l.fLab[2]>=0 ) lb[nla++]= l.fLab[2];
 468       }
 469       sort( lb, lb+nla );
 470       Int_t labmax = -1, labcur=-1, lmax = 0, lcurr=0;
 471       for( Int_t i=0; i<nla; i++ ){
 472         if( lb[i]!=labcur ){
 473           if( labcur>=0 && lmax<lcurr ){
 474             lmax = lcurr;
 475             labmax = labcur;
 476           }
 477           labcur = lb[i];
 478           lcurr = 0;
 479         }
 480         lcurr++;
 481       }
 482       if( labcur>=0 && lmax<lcurr ){
 483         lmax = lcurr;
 484         labmax = labcur;
 485       }
 486       lmax = 0;
 487       for( Int_t ihit=0; ihit<nhits; ihit++){
 488         Int_t index = slice.OutTrackHits()[tCA.FirstHitRef()+ihit];
 489         AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[index].ID()];
 490         if( l.fLab[0] == labmax || l.fLab[1] == labmax || l.fLab[2] == labmax
 491             ) lmax++;
 492       }
 493       traLabels[itr] = labmax;
 494       traPurity[itr] = ( (nhits>0) ?double(lmax)/double(nhits) :0 );
 495       //cout<<"perf track "<<itr<<": "<<nhits<<" "<<labmax<<" "<<traPurity[itr]<<endl;
 496       if( lb ) delete[] lb;
 497     }
 498   }
 499
 500   nRecTot+= traN;
 501
 502   for(Int_t itr=0; itr<traN; itr++){
 503     if( traPurity[itr]<.9 || traLabels[itr]<0 || traLabels[itr]>=fNMCTracks){
 504       nGhost++;
 505       continue;
 506     }
 507
 508     AliHLTTPCCAMCTrack &mc = fMCTracks[traLabels[itr]];
 509     mc.NReconstructed()++;
 510     if( mc.Set()== 0 ) nRecOut++;
 511     else{
 512       if( mc.NReconstructed()==1 ) nRecAll++;
 513       else if(mc.NReconstructed() > mc.NTurns() ) nClonesAll++;
 514       if( mc.Set()==2 ){
 515         if( mc.NReconstructed()==1 ) nRecRef++;
 516         else if(mc.NReconstructed() > mc.NTurns() ) nClonesRef++;
 517       }
 518     }
 519   }
 520
 521   for (Int_t ipart=0; ipart<fNMCTracks; ipart++) {
 522     AliHLTTPCCAMCTrack &mc = fMCTracks[ipart];
 523     if( mc.Set()>0 ) fhEffVsP->Fill(mc.P(), ( mc.NReconstructed()>0 ?1 :0));
 524   }
 525
 526   if( traLabels ) delete[] traLabels;
 527   if( traPurity ) delete[] traPurity;
 528
 529   fStatNRecTot += nRecTot;
 530   fStatNRecOut += nRecOut;
 531   fStatNGhost  += nGhost;
 532   fStatNMCAll  += nMCAll;
 533   fStatNRecAll  += nRecAll;
 534   fStatNClonesAll  += nClonesAll;
 535   fStatNMCRef  += nMCRef;
 536   fStatNRecRef  += nRecRef;
 537   fStatNClonesRef  += nClonesRef;
 538
 539   if( nMCAll ==0 ) return;
 540
 541   if( PrintFlag ){
 542     cout<<"Performance for slice "<<iSlice<<" : "<<endl;
 543     cout<<" N tracks : "
 544         <<nMCAll<<" mc all, "
 545         <<nMCRef<<" mc ref, "
 546         <<nRecTot<<" rec total, "
 547         <<nRecAll<<" rec all, "
 548         <<nClonesAll<<" clones all, "
 549         <<nRecRef<<" rec ref, "
 550         <<nClonesRef<<" clones ref, "
 551         <<nRecOut<<" out, "
 552         <<nGhost<<" ghost"<<endl;
 553
 554     Int_t nRecExtr = nRecAll - nRecRef;
 555     Int_t nMCExtr = nMCAll - nMCRef;
 556     Int_t nClonesExtr = nClonesAll - nClonesRef;
 557
 558     Double_t dRecTot = (nRecTot>0 ) ? nRecTot :1;
 559     Double_t dMCAll = (nMCAll>0 ) ? nMCAll :1;
 560     Double_t dMCRef = (nMCRef>0 ) ? nMCRef :1;
 561     Double_t dMCExtr = (nMCExtr>0 ) ? nMCExtr :1;
 562     Double_t dRecAll = (nRecAll+nClonesAll>0 ) ? nRecAll+nClonesAll :1;
 563     Double_t dRecRef = (nRecRef+nClonesRef>0 ) ? nRecRef+nClonesRef :1;
 564     Double_t dRecExtr = (nRecExtr+nClonesExtr>0 ) ? nRecExtr+nClonesExtr :1;
 565
 566     cout<<" EffRef = ";
 567     if( nMCRef>0 ) cout<<nRecRef/dMCRef; else cout<<"_";
 568     cout<<", CloneRef = ";
 569     if( nRecRef >0 ) cout << nClonesRef/dRecRef; else cout<<"_";
 570     cout<<endl;
 571     cout<<" EffExtra = ";
 572     if( nMCExtr>0 ) cout << nRecExtr/dMCExtr; else cout<<"_";
 573     cout <<", CloneExtra = ";
 574     if( nRecExtr>0 ) cout << nClonesExtr/dRecExtr; else cout<<"_";
 575     cout<<endl;
 576     cout<<" EffAll = ";
 577     if( nMCAll>0 ) cout<<nRecAll/dMCAll; else cout<<"_";
 578     cout <<", CloneAll = ";
 579     if( nRecAll>0 ) cout << nClonesAll/dRecAll; else cout<<"_";
 580     cout <<endl;
 581     cout<<" Out = ";
 582     if( nRecTot>0 ) cout <<nRecOut/dRecTot; else cout<<"_";
 583     cout <<", Ghost = ";
 584     if( nRecTot>0 ) cout<<nGhost/dRecTot; else cout<<"_";
 585     cout<<endl;
 586   }
 587 }
 588
 589
 590 void AliHLTTPCCAPerformance::Performance()
 591 {
 592   // main routine for performance calculation
 593   //SG!!!
 594   /*
 595   fStatNEvents=0;
 596     fStatNRecTot=0;
 597     fStatNRecOut=0;
 598     fStatNGhost=0;
 599     fStatNMCAll=0;
 600     fStatNRecAll=0;
 601     fStatNClonesAll=0;
 602     fStatNMCRef=0;
 603     fStatNRecRef=0;
 604     fStatNClonesRef=0;
 605   */
 606   fStatNEvents++;
 607   for( Int_t islice=0; islice<fTracker->NSlices(); islice++){
 608     SlicePerformance(islice,0);
 609   }
 610
 611   // global tracker performance
 612   {
 613       if( !fTracker ) return;
 614
 615       Int_t nRecTot = 0, nGhost=0, nRecOut=0;
 616       Int_t nMCAll = 0, nRecAll=0, nClonesAll=0;
 617       Int_t nMCRef = 0, nRecRef=0, nClonesRef=0;
 618
 619       // Select reconstructable MC tracks
 620
 621       {
 622         for (Int_t imc=0; imc<fNMCTracks; imc++) fMCTracks[imc].NHits() = 0;
 623
 624         for( Int_t ih=0; ih<fNHits; ih++){
 625           AliHLTTPCCAHitLabel &l = fHitLabels[ih];
 626           if( l.fLab[0]>=0 ) fMCTracks[l.fLab[0]].NHits()++;
 627           if( l.fLab[1]>=0 ) fMCTracks[l.fLab[1]].NHits()++;
 628           if( l.fLab[2]>=0 ) fMCTracks[l.fLab[2]].NHits()++;
 629         }
 630
 631         for (Int_t imc=0; imc<fNMCTracks; imc++) {
 632           AliHLTTPCCAMCTrack &mc = fMCTracks[imc];
 633           mc.Set() = 0;
 634           mc.NReconstructed() = 0;
 635           mc.NTurns() = 1;
 636           if( mc.NHits() >=  50 && mc.P()>=.05 ){
 637             mc.Set() = 1;
 638             nMCAll++;
 639             if( mc.P()>=1. ){
 640               mc.Set() = 2;
 641               nMCRef++;
 642             }
 643           }
 644         }
 645       }
 646
 647       Int_t traN = fTracker->NTracks();
 648       Int_t *traLabels = 0;
 649       Double_t *traPurity = 0;
 650       traLabels = new Int_t[traN];
 651       traPurity = new Double_t[traN];
 652       {
 653         for (Int_t itr=0; itr<traN; itr++) {
 654           traLabels[itr]=-1;
 655           traPurity[itr]= 0;
 656           AliHLTTPCCAGBTrack &tCA = fTracker->Tracks()[itr];
 657           Int_t nhits = tCA.NHits();
 658           Int_t *lb = new Int_t[nhits*3];
 659           Int_t nla=0;
 660           for( Int_t ihit=0; ihit<nhits; ihit++){
 661             Int_t index = fTracker->TrackHits()[tCA.FirstHitRef()+ihit];
 662             AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[index].ID()];
 663             if(l.fLab[0]>=0 ) lb[nla++]= l.fLab[0];
 664             if(l.fLab[1]>=0 ) lb[nla++]= l.fLab[1];
 665             if(l.fLab[2]>=0 ) lb[nla++]= l.fLab[2];
 666           }
 667           sort( lb, lb+nla );
 668           Int_t labmax = -1, labcur=-1, lmax = 0, lcurr=0;
 669           for( Int_t i=0; i<nla; i++ ){
 670             if( lb[i]!=labcur ){
 671               if( labcur>=0 && lmax<lcurr ){
 672                 lmax = lcurr;
 673                 labmax = labcur;
 674               }
 675               labcur = lb[i];
 676               lcurr = 0;
 677             }
 678             lcurr++;
 679           }
 680           if( labcur>=0 && lmax<lcurr ){
 681             lmax = lcurr;
 682             labmax = labcur;
 683           }
 684           lmax = 0;
 685           for( Int_t ihit=0; ihit<nhits; ihit++){
 686             Int_t index = fTracker->TrackHits()[tCA.FirstHitRef()+ihit];
 687             AliHLTTPCCAHitLabel &l = fHitLabels[fTracker->Hits()[index].ID()];
 688             if( l.fLab[0] == labmax || l.fLab[1] == labmax || l.fLab[2] == labmax
 689                 ) lmax++;
 690           }
 691           traLabels[itr] = labmax;
 692           traPurity[itr] = ( (nhits>0) ?double(lmax)/double(nhits) :0 );
 693           if( lb ) delete[] lb;
 694         }
 695       }
 696
 697       nRecTot+= traN;
 698       for(Int_t itr=0; itr<traN; itr++){
 699         if( traPurity[itr]<.9 || traLabels[itr]<0 || traLabels[itr]>=fNMCTracks){
 700           nGhost++;
 701           continue;
 702         }
 703
 704         AliHLTTPCCAMCTrack &mc = fMCTracks[traLabels[itr]];
 705
 706         mc.NReconstructed()++;
 707         if( mc.Set()== 0 ) nRecOut++;
 708         else{
 709           if( mc.NReconstructed()==1 ) nRecAll++;
 710           else if(mc.NReconstructed() > mc.NTurns() ) nClonesAll++;
 711           if( mc.Set()==2 ){
 712             if( mc.NReconstructed()==1 ) nRecRef++;
 713             else if(mc.NReconstructed() > mc.NTurns() ) nClonesRef++;
 714           }
 715         }
 716
 717         // track resolutions
 718         while( TMath::Abs(mc.TPCPar()[0]) + TMath::Abs(mc.TPCPar()[1])>1 ){
 719           if( traPurity[itr]<.90 ) break;
 720           AliHLTTPCCAGBTrack &t = fTracker->Tracks()[itr];
 721           AliHLTTPCCATrackParam p = t.Param();
 722           Double_t cosA = TMath::Cos( t.Alpha() );
 723           Double_t sinA = TMath::Sin( t.Alpha() );
 724           Double_t mcX =  mc.TPCPar()[0]*cosA + mc.TPCPar()[1]*sinA;
 725           Double_t mcY = -mc.TPCPar()[0]*sinA + mc.TPCPar()[1]*cosA;
 726           Double_t mcZ =  mc.TPCPar()[2];
 727           Double_t mcEx =  mc.TPCPar()[3]*cosA + mc.TPCPar()[4]*sinA;
 728           Double_t mcEy = -mc.TPCPar()[3]*sinA + mc.TPCPar()[4]*cosA;
 729           Double_t mcEz =  mc.TPCPar()[5];
 730           Double_t mcEt = TMath::Sqrt(mcEx*mcEx + mcEy*mcEy);
 731           if( TMath::Abs(mcEt)<1.e-4 ) break;
 732           Double_t mcSinPhi = mcEy / mcEt;
 733           Double_t mcDzDs   = mcEz / mcEt;
 734           Double_t mcQPt = mc.TPCPar()[6]/ mcEt;
 735           if( TMath::Abs(mcQPt)<1.e-4 ) break;
 736           Double_t mcPt = 1./TMath::Abs(mcQPt);
 737           if( mcPt<1. ) break;
 738           if( t.NHits() <  50 ) break;
 739           Double_t bz = fTracker->Slices()[0].Param().Bz();
 740           if( !p.TransportToXWithMaterial( mcX, bz ) ) break;
 741           if( p.GetCosPhi()*mcEx < 0 ){ // change direction
 742             mcSinPhi = -mcSinPhi;
 743             mcDzDs = -mcDzDs;
 744             mcQPt = -mcQPt;
 745           }
 746           const Double_t kCLight = 0.000299792458;
 747           Double_t k2QPt = 100;
 748           if( TMath::Abs(bz)>1.e-4 ) k2QPt= 1./(bz*kCLight);
 749           Double_t qPt = p.GetKappa()*k2QPt;
 750           Double_t pt = 100;
 751           if( TMath::Abs(qPt) >1.e-4 ) pt = 1./TMath::Abs(qPt);
 752
 753           fhResY->Fill( p.GetY() - mcY );
 754           fhResZ->Fill( p.GetZ() - mcZ );
 755           fhResSinPhi->Fill( p.GetSinPhi() - mcSinPhi );
 756           fhResDzDs->Fill( p.GetDzDs() - mcDzDs );
 757           fhResPt->Fill( ( pt - mcPt )/mcPt );
 758
 759           if( p.GetErr2Y()>0 ) fhPullY->Fill( (p.GetY() - mcY)/TMath::Sqrt(p.GetErr2Y()) );
 760           if( p.GetErr2Z()>0 ) fhPullZ->Fill( (p.GetZ() - mcZ)/TMath::Sqrt(p.GetErr2Z()) );
 761           if( p.GetErr2SinPhi()>0 ) fhPullSinPhi->Fill( (p.GetSinPhi() - mcSinPhi)/TMath::Sqrt(p.GetErr2SinPhi()) );
 762           if( p.GetErr2DzDs()>0 ) fhPullDzDs->Fill( (p.DzDs() - mcDzDs)/TMath::Sqrt(p.GetErr2DzDs()) );
 763           if( p.GetErr2Kappa()>0 ) fhPullQPt->Fill( (qPt - mcQPt)/TMath::Sqrt(p.GetErr2Kappa()*k2QPt*k2QPt) );
 764           break;
 765         }
 766       }
 767
 768       for (Int_t ipart=0; ipart<fNMCTracks; ipart++) {
 769         AliHLTTPCCAMCTrack &mc = fMCTracks[ipart];
 770         if( mc.Set()>0 ) fhGBEffVsP->Fill(mc.P(), ( mc.NReconstructed()>0 ?1 :0));
 771       }
 772
 773       if( traLabels ) delete[] traLabels;
 774       if( traPurity ) delete[] traPurity;
 775
 776       fStatGBNRecTot += nRecTot;
 777       fStatGBNRecOut += nRecOut;
 778       fStatGBNGhost  += nGhost;
 779       fStatGBNMCAll  += nMCAll;
 780       fStatGBNRecAll  += nRecAll;
 781       fStatGBNClonesAll  += nClonesAll;
 782       fStatGBNMCRef  += nMCRef;
 783       fStatGBNRecRef  += nRecRef;
 784       fStatGBNClonesRef  += nClonesRef;
 785   }
 786
 787
 788   // distribution of cluster errors
 789
 790   {
 791     Int_t nHits = fTracker->NHits();
 792     for( Int_t ih=0; ih<nHits; ih++ ){
 793       AliHLTTPCCAGBHit &hit = fTracker->Hits()[ih];
 794       fhHitErrY->Fill(hit.ErrY());
 795       fhHitErrZ->Fill(hit.ErrZ());
 796     }
 797   }
 798
 799   // cluster pulls
 800
 801   if( fDoClusterPulls && fNMCPoints>0 ) {
 802
 803     {
 804       for (Int_t ipart=0; ipart<fNMCTracks; ipart++) {
 805         AliHLTTPCCAMCTrack &mc = fMCTracks[ipart];
 806         mc.NMCPoints() = 0;
 807       }
 808       sort(fMCPoints, fMCPoints+fNMCPoints, AliHLTTPCCAMCPoint::Compare );
 809
 810       for( Int_t ip=0; ip<fNMCPoints; ip++ ){
 811         AliHLTTPCCAMCPoint &p = fMCPoints[ip];
 812         AliHLTTPCCAMCTrack &t = fMCTracks[p.TrackID()];
 813         if( t.NMCPoints()==0 ) t.FirstMCPointID() = ip;
 814         t.NMCPoints()++;
 815       }
 816     }
 817
 818     for( Int_t ih=0; ih<fNHits; ih++ ){
 819
 820       AliHLTTPCCAGBHit &hit = fTracker->Hits()[ih];
 821       AliHLTTPCCAHitLabel &l = fHitLabels[ih];
 822
 823       if( l.fLab[0]<0 || l.fLab[0]>=fNMCTracks
 824           || l.fLab[1]>=0 || l.fLab[2]>=0       ) continue;
 825
 826       Int_t lab = l.fLab[0];
 827
 828       AliHLTTPCCAMCTrack &track = fMCTracks[lab];
 829       //if( track.Pt()<1. ) continue;
 830       Int_t ip1=-1, ip2=-1;
 831       Double_t d1 = 1.e20, d2=1.e20;
 832       for( Int_t ip=0; ip<track.NMCPoints(); ip++ ){
 833         AliHLTTPCCAMCPoint &p = fMCPoints[track.FirstMCPointID() + ip];
 834         if( p.ISlice() != hit.ISlice() ) continue;
 835         Double_t dx = p.Sx()-hit.X();
 836         Double_t dy = p.Sy()-hit.Y();
 837         Double_t dz = p.Sz()-hit.Z();
 838         Double_t d = dx*dx + dy*dy + dz*dz;
 839         if( p.Sx()< hit.X() ){
 840           if( d<d1 ){
 841             d1 = d;
 842             ip1 = ip;
 843           }
 844         }else{
 845           if( d<d2 ){
 846             d2 = d;
 847             ip2 = ip;
 848           }
 849         }
 850       }
 851
 852       if( ip1<0 || ip2<0 ) continue;
 853
 854       AliHLTTPCCAMCPoint &p1 = fMCPoints[track.FirstMCPointID() + ip1];
 855       AliHLTTPCCAMCPoint &p2 = fMCPoints[track.FirstMCPointID() + ip2];
 856       Double_t dx = p2.Sx() - p1.Sx();
 857       Double_t dy = p2.Sy() - p1.Sy();
 858       Double_t dz = p2.Sz() - p1.Sz();
 859       if( TMath::Abs(dx)>1.e-8 && TMath::Abs(p1.Sx()-hit.X())<2. && TMath::Abs(p2.Sx()-hit.X())<2.  ){
 860         Double_t sx = hit.X();
 861         Double_t sy = p1.Sy() + dy/dx*(sx-p1.Sx());
 862         Double_t sz = p1.Sz() + dz/dx*(sx-p1.Sx());
 863
 864         Float_t errY, errZ;
 865         {
 866           AliHLTTPCCATrackParam t;
 867           t.Z() = sz;
 868           t.SinPhi() = dy/TMath::Sqrt(dx*dx+dy*dy);
 869           t.CosPhi() = dx/TMath::Sqrt(dx*dx+dy*dy);
 870           t.DzDs() = dz/TMath::Sqrt(dx*dx+dy*dy);
 871           fTracker->GetErrors2(hit,t,errY, errZ );
 872           errY = TMath::Sqrt(errY);
 873           errZ = TMath::Sqrt(errZ);
 874         }
 875
 876         fhHitResY->Fill((hit.Y()-sy));
 877         fhHitResZ->Fill((hit.Z()-sz));
 878         fhHitPullY->Fill((hit.Y()-sy)/errY);
 879         fhHitPullZ->Fill((hit.Z()-sz)/errZ);
 880         if( track.Pt()>=1. ){
 881           fhHitResY1->Fill((hit.Y()-sy));
 882           fhHitResZ1->Fill((hit.Z()-sz));
 883           fhHitPullY1->Fill((hit.Y()-sy)/errY);
 884           fhHitPullZ1->Fill((hit.Z()-sz)/errZ);
 885         }
 886       }
 887     }
 888   }
 889
 890   {
 891     cout<<"\nSlice tracker performance: \n"<<endl;
 892     cout<<" N tracks : "
 893         <<fStatNMCAll/fStatNEvents<<" mc all, "
 894         <<fStatNMCRef/fStatNEvents<<" mc ref, "
 895         <<fStatNRecTot/fStatNEvents<<" rec total, "
 896         <<fStatNRecAll/fStatNEvents<<" rec all, "
 897         <<fStatNClonesAll/fStatNEvents<<" clones all, "
 898         <<fStatNRecRef/fStatNEvents<<" rec ref, "
 899         <<fStatNClonesRef/fStatNEvents<<" clones ref, "
 900         <<fStatNRecOut/fStatNEvents<<" out, "
 901         <<fStatNGhost/fStatNEvents<<" ghost"<<endl;
 902
 903     Int_t nRecExtr = fStatNRecAll - fStatNRecRef;
 904     Int_t nMCExtr = fStatNMCAll - fStatNMCRef;
 905     Int_t nClonesExtr = fStatNClonesAll - fStatNClonesRef;
 906
 907     Double_t dRecTot = (fStatNRecTot>0 ) ? fStatNRecTot :1;
 908     Double_t dMCAll = (fStatNMCAll>0 ) ? fStatNMCAll :1;
 909     Double_t dMCRef = (fStatNMCRef>0 ) ? fStatNMCRef :1;
 910     Double_t dMCExtr = (nMCExtr>0 ) ? nMCExtr :1;
 911     Double_t dRecAll = (fStatNRecAll+fStatNClonesAll>0 ) ? fStatNRecAll+fStatNClonesAll :1;
 912     Double_t dRecRef = (fStatNRecRef+fStatNClonesRef>0 ) ? fStatNRecRef+fStatNClonesRef :1;
 913     Double_t dRecExtr = (nRecExtr+nClonesExtr>0 ) ? nRecExtr+nClonesExtr :1;
 914
 915     cout<<" EffRef = "<< fStatNRecRef/dMCRef
 916         <<", CloneRef = " << fStatNClonesRef/dRecRef <<endl;
 917     cout<<" EffExtra = "<<nRecExtr/dMCExtr
 918         <<", CloneExtra = " << nClonesExtr/dRecExtr<<endl;
 919     cout<<" EffAll = "<<fStatNRecAll/dMCAll
 920         <<", CloneAll = " << fStatNClonesAll/dRecAll<<endl;
 921     cout<<" Out = "<<fStatNRecOut/dRecTot
 922         <<", Ghost = "<<fStatNGhost/dRecTot<<endl;
 923     cout<<" Time = "<<fTracker->StatTime(0)/fTracker->StatNEvents()*1.e3<<" msec/event "<<endl;
 924     cout<<" Local timers = "
 925         <<fTracker->StatTime(1)/fTracker->StatNEvents()*1.e3<<" "
 926         <<fTracker->StatTime(2)/fTracker->StatNEvents()*1.e3<<" "
 927         <<fTracker->StatTime(3)/fTracker->StatNEvents()*1.e3<<" "
 928         <<fTracker->StatTime(4)/fTracker->StatNEvents()*1.e3<<" "
 929         <<fTracker->StatTime(5)/fTracker->StatNEvents()*1.e3<<"["
 930         <<fTracker->StatTime(6)/fTracker->StatNEvents()*1.e3<<"/"
 931         <<fTracker->StatTime(7)/fTracker->StatNEvents()*1.e3<<"] "
 932         <<fTracker->StatTime(8)/fTracker->StatNEvents()*1.e3<<" "
 933         <<" msec/event "<<endl;
 934   }
 935
 936   {
 937     cout<<"\nGlobal tracker performance: \n"<<endl;
 938     cout<<" N tracks : "
 939         <<fStatGBNMCAll/fStatNEvents<<" mc all, "
 940         <<fStatGBNMCRef/fStatNEvents<<" mc ref, "
 941         <<fStatGBNRecTot/fStatNEvents<<" rec total, "
 942         <<fStatGBNRecAll/fStatNEvents<<" rec all, "
 943         <<fStatGBNClonesAll/fStatNEvents<<" clones all, "
 944         <<fStatGBNRecRef/fStatNEvents<<" rec ref, "
 945         <<fStatGBNClonesRef/fStatNEvents<<" clones ref, "
 946         <<fStatGBNRecOut/fStatNEvents<<" out, "
 947         <<fStatGBNGhost/fStatNEvents<<" ghost"<<endl;
 948
 949     Int_t nRecExtr = fStatGBNRecAll - fStatGBNRecRef;
 950     Int_t nMCExtr = fStatGBNMCAll - fStatGBNMCRef;
 951     Int_t nClonesExtr = fStatGBNClonesAll - fStatGBNClonesRef;
 952
 953     Double_t dRecTot = (fStatGBNRecTot>0 ) ? fStatGBNRecTot :1;
 954     Double_t dMCAll = (fStatGBNMCAll>0 ) ? fStatGBNMCAll :1;
 955     Double_t dMCRef = (fStatGBNMCRef>0 ) ? fStatGBNMCRef :1;
 956     Double_t dMCExtr = (nMCExtr>0 ) ? nMCExtr :1;
 957     Double_t dRecAll = (fStatGBNRecAll+fStatGBNClonesAll>0 ) ? fStatGBNRecAll+fStatGBNClonesAll :1;
 958     Double_t dRecRef = (fStatGBNRecRef+fStatGBNClonesRef>0 ) ? fStatGBNRecRef+fStatGBNClonesRef :1;
 959     Double_t dRecExtr = (nRecExtr+nClonesExtr>0 ) ? nRecExtr+nClonesExtr :1;
 960
 961     cout<<" EffRef = "<< fStatGBNRecRef/dMCRef
 962         <<", CloneRef = " << fStatGBNClonesRef/dRecRef <<endl;
 963     cout<<" EffExtra = "<<nRecExtr/dMCExtr
 964         <<", CloneExtra = " << nClonesExtr/dRecExtr<<endl;
 965     cout<<" EffAll = "<<fStatGBNRecAll/dMCAll
 966         <<", CloneAll = " << fStatGBNClonesAll/dRecAll<<endl;
 967     cout<<" Out = "<<fStatGBNRecOut/dRecTot
 968         <<", Ghost = "<<fStatGBNGhost/dRecTot<<endl;
 969     cout<<" Time = "<<( fTracker->StatTime(0)+fTracker->StatTime(9) )/fTracker->StatNEvents()*1.e3<<" msec/event "<<endl;
 970     cout<<" Local timers: "<<endl;
 971     cout<<" slice tracker "<<fTracker->StatTime(0)/fTracker->StatNEvents()*1.e3<<": "
 972         <<fTracker->StatTime(1)/fTracker->StatNEvents()*1.e3<<" "
 973         <<fTracker->StatTime(2)/fTracker->StatNEvents()*1.e3<<" "
 974         <<fTracker->StatTime(3)/fTracker->StatNEvents()*1.e3<<" "
 975         <<fTracker->StatTime(4)/fTracker->StatNEvents()*1.e3<<" "
 976         <<fTracker->StatTime(5)/fTracker->StatNEvents()*1.e3<<"["
 977         <<fTracker->StatTime(6)/fTracker->StatNEvents()*1.e3<<"/"
 978         <<fTracker->StatTime(7)/fTracker->StatNEvents()*1.e3<<"] "
 979         <<fTracker->StatTime(8)/fTracker->StatNEvents()*1.e3
 980         <<" msec/event "<<endl;
 981     cout<<" GB merger "<<fTracker->StatTime(9)/fTracker->StatNEvents()*1.e3<<": "
 982         <<fTracker->StatTime(10)/fTracker->StatNEvents()*1.e3<<", "
 983         <<fTracker->StatTime(11)/fTracker->StatNEvents()*1.e3<<", "
 984         <<fTracker->StatTime(12)/fTracker->StatNEvents()*1.e3<<" "
 985         <<" msec/event "<<endl;
 986   }
 987
 988   WriteHistos();
 989 }
 990
 991 void AliHLTTPCCAPerformance::WriteMCEvent( ostream &out )
 992 {
 993   out<<fNMCTracks<<endl;
 994   for( Int_t it=0; it<fNMCTracks; it++ ){
 995     AliHLTTPCCAMCTrack &t = fMCTracks[it];
 996     out<<it<<" ";
 997     out<<t.PDG()<<endl;
 998     for( Int_t i=0; i<7; i++ ) out<<t.Par()[i]<<" ";
 999     out<<endl<<"    ";
1000     for( Int_t i=0; i<7; i++ ) out<<t.TPCPar()[i]<<" ";
1001     out<<endl<<"    ";
1002     out<< t.P()<<" ";
1003     out<< t.Pt()<<" ";
1004     out<< t.NMCPoints()<<" ";
1005     out<< t.FirstMCPointID()<<" ";
1006     out<< t.NHits()<<" ";
1007     out<< t.NReconstructed()<<" ";
1008     out<< t.Set()<<" ";
1009     out<< t.NTurns()<<endl;
1010   }
1011
1012   out<<fNHits<<endl;
1013   for( Int_t ih=0; ih<fNHits; ih++ ){
1014     AliHLTTPCCAHitLabel &l = fHitLabels[ih];
1015     out<<l.fLab[0]<<" "<<l.fLab[1]<<" "<<l.fLab[2]<<endl;
1016   }
1017 }
1018
1019 void AliHLTTPCCAPerformance::WriteMCPoints( ostream &out )
1020 {
1021   out<<fNMCPoints<<endl;
1022   for( Int_t ip=0; ip<fNMCPoints; ip++ ){
1023     AliHLTTPCCAMCPoint &p = fMCPoints[ip];
1024     out<< p.X()<<" ";
1025     out<< p.Y()<<" ";
1026     out<< p.Z()<<" ";
1027     out<< p.Sx()<<" ";
1028     out<< p.Sy()<<" ";
1029     out<< p.Sz()<<" ";
1030     out<< p.Time()<<" ";
1031     out<< p.ISlice()<<" ";
1032     out<< p.TrackID()<<endl;
1033   }
1034 }
1035
1036 void AliHLTTPCCAPerformance::ReadMCEvent( istream &in )
1037 {
1038   StartEvent();
1039   if( fMCTracks ) delete[] fMCTracks;
1040   fMCTracks = 0;
1041   fNMCTracks = 0;
1042   if( fHitLabels ) delete[] fHitLabels;
1043   fHitLabels = 0;
1044   fNHits = 0;
1045   if( fMCPoints ) delete[] fMCPoints;
1046   fMCPoints = 0;
1047   fNMCPoints = 0;
1048
1049   in>>fNMCTracks;
1050   fMCTracks = new AliHLTTPCCAMCTrack[fNMCTracks];
1051   for( Int_t it=0; it<fNMCTracks; it++ ){
1052     AliHLTTPCCAMCTrack &t = fMCTracks[it];
1053     Int_t j;
1054     in>>j;
1055     in>> t.PDG();
1056     for( Int_t i=0; i<7; i++ ) in>>t.Par()[i];
1057     for( Int_t i=0; i<7; i++ ) in>>t.TPCPar()[i];
1058     in>> t.P();
1059     in>> t.Pt();
1060     in>> t.NHits();
1061     in>> t.NMCPoints();
1062     in>> t.FirstMCPointID();
1063     in>> t.NReconstructed();
1064     in>> t.Set();
1065     in>> t.NTurns();
1066   }
1067
1068   in>>fNHits;
1069   fHitLabels = new AliHLTTPCCAHitLabel[fNHits];
1070   for( Int_t ih=0; ih<fNHits; ih++ ){
1071     AliHLTTPCCAHitLabel &l = fHitLabels[ih];
1072     in>>l.fLab[0]>>l.fLab[1]>>l.fLab[2];
1073   }
1074 }
1075
1076 void AliHLTTPCCAPerformance::ReadMCPoints( istream &in )
1077 {
1078   if( fMCPoints ) delete[] fMCPoints;
1079   fMCPoints = 0;
1080   fNMCPoints = 0;
1081
1082   in>>fNMCPoints;
1083   fMCPoints = new AliHLTTPCCAMCPoint[fNMCPoints];
1084   for( Int_t ip=0; ip<fNMCPoints; ip++ ){
1085     AliHLTTPCCAMCPoint &p = fMCPoints[ip];
1086     in>> p.X();
1087     in>> p.Y();
1088     in>> p.Z();
1089     in>> p.Sx();
1090     in>> p.Sy();
1091     in>> p.Sz();
1092     in>> p.Time();
1093     in>> p.ISlice();
1094     in>> p.TrackID();
1095   }
1096 }