ITS/AliITSComparisonV1.C

   1 #include "iostream.h"
   2
   3
   4   void AliITSComparisonV1(Int_t evNumber1=0,Int_t evNumber2=0) {
   5
   6   const char *filename="itstracks.root";
   7
   8   ///////////////// Dynamically link some shared libs ////////////////////////////////
   9
  10   if (gClassTable->GetID("AliRun") < 0) {
  11     gROOT->LoadMacro("loadlibs.C");
  12     loadlibs();
  13   } else {
  14     delete gAlice;
  15     gAlice=0;
  16   }
  17
  18 // Connect the Root Galice file containing Geometry, Kine and Hits
  19    TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
  20    if (!file) file = new TFile(filename);
  21
  22    ifstream in("itsgood_tracks");
  23 //
  24 //   Loop over events
  25 //
  26
  27    char tname[30];
  28
  29
  30   struct GoodTrack {
  31    Int_t fEventN;
  32     Int_t lab,code;
  33     Float_t px,py,pz,x,y,z,pxg,pyg,pzg,ptg;
  34     Bool_t flag;
  35   };
  36     Int_t i;
  37   GoodTrack gt[50000];
  38   Int_t ngood=0;
  39   cerr<<"Reading itsgood tracks...\n";
  40   while (in>>gt[ngood].fEventN>>gt[ngood].lab>>gt[ngood].code
  41           >>gt[ngood].px >>gt[ngood].py>>gt[ngood].pz
  42           >>gt[ngood].x  >>gt[ngood].y >>gt[ngood].z
  43           >>gt[ngood].pxg  >>gt[ngood].pyg >>gt[ngood].pzg
  44           >>gt[ngood].ptg >>gt[ngood].flag) {
  45     ngood++;
  46     cerr<<ngood<<'\r';
  47     if (ngood==50000) {
  48       cerr<<"Too many good tracks !\n";
  49       break;
  50     }
  51   }
  52   if (!in.eof()) cerr<<"Read error (itsgood_tracks) !\n";
  53    ofstream out1 ("AliITSTrag.out");
  54   Int_t countpos=0,countneg=0;
  55
  56
  57   //for (i=0; i<ngood; i++) out1 << gt[i].ptg << "\n";
  58   for (i=0; i<ngood; i++) {
  59     out1 << gt[i].ptg << "\n";
  60     Int_t codpar=gt[i].code;
  61     if(codpar==2212 || codpar==-11 || codpar==-13 || codpar==211 || codpar==321 || codpar==3222
  62          || codpar==213 || codpar==323 || codpar==10323 || codpar==3224 || codpar==2224 || codpar==2214
  63          || codpar==-1114 || codpar==-3112 || codpar==-3312 || codpar==3224 || codpar==-3114 || codpar==-3314
  64          || codpar==411 || codpar==431 || codpar==413 || codpar==433 || codpar==-15 || codpar==4232
  65          || codpar==4222 || codpar==4322 || codpar==4422 || codpar==4412 || codpar==4432 || codpar==4224
  66          ||codpar==4214 || codpar==4324 || codpar==4424 || codpar==4414 || codpar==4434 || codpar==4444)
  67     countpos++;
  68          if(codpar==-2212 || codpar==11 || codpar==13 || codpar==-211 || codpar==-321 || codpar==3112
  69          || codpar==-213 || codpar==-323 || codpar==-10323 || codpar==3114 || codpar==1114 || codpar==-2224
  70          || codpar==-2214 || codpar==33112 || codpar==-3222 || codpar==3114 || codpar==3314 || codpar==3334
  71          || codpar==-3224 || codpar==-411 || codpar==-431 || codpar==-413 || codpar==-433 || codpar==15
  72          || codpar==-4422 || codpar==-4432 || codpar==-4214 || codpar==-4324 || codpar==-4424 || codpar==-4434
  73          || codpar==-444)
  74          countneg++;
  75   }
  76   out1.close();
  77    ofstream out ("AliITSTra.out");
  78 ///   definition of nm of good track for each event
  79    TVector ntrackevtpc(evNumber2-evNumber1 +1);
  80    Int_t nchange=-1;
  81    Int_t nmev=-100;
  82
  83
  84     Int_t i;
  85     for( i =0; i<ngood ; i++){
  86     if(gt[i].fEventN != nmev ){
  87     nmev=gt[i].fEventN;
  88     nchange++;
  89     }
  90     if(nmev == gt[i].fEventN) ntrackevtpc(nchange)++;
  91     }
  92 //////////////////////////////////////////////////////////////
  93
  94    for (Int_t nev=evNumber1; nev<= evNumber2; nev++) {
  95
  96
  97    sprintf(tname,"TreeT%d",nev);
  98    TTree *tracktree=(TTree*)file->Get(tname);
  99    TBranch *tbranch=tracktree->GetBranch("ITStracks");
 100    cout<<" nev = "<<nev<<"\n";
 101         //cout<<" open the file \n";
 102
 103    Int_t nentr=tracktree->GetEntries();
 104
 105    TObjArray tarray(nentr);
 106   // AliITSIOTrack *iotrack=0;
 107    printf("nentr %d\n",nentr);
 108
 109    for (Int_t i=0; i<nentr; i++) {
 110       AliITSIOTrack *iotrack=new AliITSIOTrack;
 111       // tarray.AddAt(new AliITSIOTrack,i);
 112       // iotrack=(AliITSiotrack*)tarray.UncheckedAt(i);
 113        tbranch->SetAddress(&iotrack);
 114        tracktree->GetEvent(i);
 115         tarray.AddLast(iotrack);
 116    }
 117         //file->Close();
 118
 119           AliITSIOTrack *iotrack;
 120    for (Int_t i=0; i<nentr; i++) {
 121      AliITSIOTrack *iotrack=new AliITSIOTrack;
 122          iotrack=(AliITSIOTrack*)tarray.UncheckedAt(i);
 123          if(!iotrack) continue;
 124      Int_t labITS=iotrack->GetLabel();
 125      Int_t labTPC=iotrack->GetTPCLabel();
 126           Double_t phistate=iotrack->GetStatePhi();
 127           Double_t tgl=iotrack->GetStateTgl();
 128           Double_t Zstate=iotrack->GetStateZ();
 129           Double_t Dr=iotrack->GetStateD();
 130           Double_t C=iotrack->GetStateC();
 131           Double_t px=iotrack->GetPx();
 132           Double_t py=iotrack->GetPy();
 133           Double_t pz=iotrack->GetPz();
 134           Double_t pt=TMath::Sqrt(px*px+py*py);
 135           Int_t c = iotrack->GetCharge();
 136           Double_t x=iotrack->GetX();
 137           Double_t y=iotrack->GetY();
 138           Double_t z= iotrack->GetZ();
 139           Double_t Dz=iotrack->GetDz();
 140
 141
 142      // cout<<" track label = "<<label<<"\n";
 143      // cout<<" phi z D tanl C = "<<phistate<<" "<<Zstate<<" "<<Dr<<" "<<tgl<<" "<<C<<"\n";
 144
 145   //  Int_t ilab=TMath::Abs(iotrack->GetLabel());
 146     Int_t flaglab=0;
 147     Int_t iii=0;
 148    Double_t ptg=0.,pxg=0.,pyg=0.,pzg=0.;
 149         Double_t xo=0., yo=0., zo=0.;
 150
 151     Int_t mingood=0, maxgood=0, jj=0;
 152     if(nev==evNumber1) mingood=0;
 153     else
 154     {for(jj=evNumber1; jj<nev; jj++) mingood += ntrackevtpc(jj);}
 155      for(jj=evNumber1; jj<=nev; jj++) maxgood+= ntrackevtpc(jj);
 156
 157    Int_t ilab=TMath::Abs(labTPC);
 158    // for (iii=0;iii<ngood;iii++) {
 159    for(iii=mingood; iii<maxgood; iii++){
 160          //cout<<" ilab, gt[iii].lab = "<<ilab<<" "<<gt[iii].lab<<"\n"; getchar();
 161       if (ilab==gt[iii].lab) {
 162         flaglab=1;
 163         ptg=gt[iii].ptg;
 164         pxg=gt[iii].pxg;
 165         pyg=gt[iii].pyg;
 166         pzg=gt[iii].pzg;
 167         xo=gt[iii].x;
 168         yo=gt[iii].y;
 169         zo=gt[iii].z;
 170         break;
 171       }
 172     }
 173
 174     if (!flaglab) continue;
 175     //cout<<" j =  " <<j<<"\n";  getchar();
 176    TVector dataOut(10);
 177     Int_t kkk=0;
 178
 179     dataOut(kkk) = ptg; kkk++; dataOut(kkk)=labITS; kkk++; dataOut(kkk)=labTPC; kkk++;
 180
 181       ///////////////////////////////
 182       Double_t difpt= (pt-ptg)/ptg*100.;
 183       //cout<<" difpt = "<<difpt<<"\n"; getchar();
 184       dataOut(kkk)=difpt; kkk++;
 185       Double_t lambdag=TMath::ATan(pzg/ptg);
 186       Double_t   lam=TMath::ATan(tgl);
 187       Double_t diflam = (lam - lambdag)*1000.;
 188       dataOut(kkk) = diflam; kkk++;
 189       Double_t phig=TMath::ATan2(pyg,pxg);  if(phig<0) phig=2.*TMath::Pi()+phig;
 190      // Double_t phi=phivertex;
 191        Double_t phi=TMath::ACos(px/pt);
 192      Double_t duepi=2.*TMath::Pi();
 193      if(phi>duepi) phi-=duepi;
 194      if(phi<0.) phi+=duepi;
 195       Double_t signC=0.;
 196       if(c>0) signC=1.; else signC=-1.;
 197           Double_t Dtot= TMath::Sqrt(Dr*Dr+Dz*Dz);
 198       Double_t difphi = (phi - phig)*1000.;
 199       dataOut(kkk)=difphi; kkk++;
 200       dataOut(kkk)=Dtot*1.e4; kkk++;
 201       dataOut(kkk)=Dr*1.e4; kkk++;
 202       dataOut(kkk)=Dz*1.e4; kkk++;
 203       dataOut(kkk)=signC; kkk++;
 204       Int_t r;
 205       for (r=0; r<10; r++) { out<<dataOut(r)<<" ";}
 206       out<<"\n";
 207       kkk=0;
 208
 209     delete iotrack;
 210    }
 211   //out.close();
 212    }   // event loop
 213   out.close();
 214   in.close();
 215   file->Close();
 216 }
 217