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[15000];
  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==15000) {
  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      //Int_t labTPC=labITS;  //provvisoria
 127           Double_t phistate=iotrack->GetStatePhi();
 128           Double_t tgl=iotrack->GetStateTgl();
 129           Double_t Zstate=iotrack->GetStateZ();
 130           Double_t Dr=iotrack->GetStateD();
 131           Double_t C=iotrack->GetStateC();
 132           Double_t px=iotrack->GetPx();
 133           Double_t py=iotrack->GetPy();
 134           Double_t pz=iotrack->GetPz();
 135           Double_t pt=TMath::Sqrt(px*px+py*py);
 136           Int_t c = iotrack->GetCharge();
 137           Double_t x=iotrack->GetX();
 138           Double_t y=iotrack->GetY();
 139           Double_t z= iotrack->GetZ();
 140         //  Double_t Dz=z;   //non e' vero bisogna levare vertice
 141          // Double_t Dtot= TMath::Sqrt(Dr*Dr+Dz*Dz);
 142
 143      // cout<<" track label = "<<label<<"\n";
 144      // cout<<" phi z D tanl C = "<<phistate<<" "<<Zstate<<" "<<Dr<<" "<<tgl<<" "<<C<<"\n";
 145
 146   //  Int_t ilab=TMath::Abs(iotrack->GetLabel());
 147     Int_t flaglab=0;
 148     Int_t iii=0;
 149    Double_t ptg=0.,pxg=0.,pyg=0.,pzg=0.;
 150         Double_t xo=0., yo=0., zo=0.;
 151
 152     Int_t mingood=0, maxgood=0, jj=0;
 153     if(nev==evNumber1) mingood=0;
 154     else
 155     {for(jj=evNumber1; jj<nev; jj++) mingood += ntrackevtpc(jj);}
 156      for(jj=evNumber1; jj<=nev; jj++) maxgood+= ntrackevtpc(jj);
 157
 158    Int_t ilab=TMath::Abs(labTPC);
 159    // for (iii=0;iii<ngood;iii++) {
 160    for(iii=mingood; iii<maxgood; iii++){
 161          //cout<<" ilab, gt[iii].lab = "<<ilab<<" "<<gt[iii].lab<<"\n"; getchar();
 162       if (ilab==gt[iii].lab) {
 163         flaglab=1;
 164         ptg=gt[iii].ptg;
 165         pxg=gt[iii].pxg;
 166         pyg=gt[iii].pyg;
 167         pzg=gt[iii].pzg;
 168         xo=gt[iii].x;
 169         yo=gt[iii].y;
 170         zo=gt[iii].z;
 171         break;
 172       }
 173     }
 174
 175     if (!flaglab) continue;
 176     //cout<<" j =  " <<j<<"\n";  getchar();
 177    TVector dataOut(10);
 178     Int_t kkk=0;
 179
 180     dataOut(kkk) = ptg; kkk++; dataOut(kkk)=labITS; kkk++; dataOut(kkk)=labTPC; kkk++;
 181
 182       ///////////////////////////////
 183       Double_t difpt= (pt-ptg)/ptg*100.;
 184       //cout<<" difpt = "<<difpt<<"\n"; getchar();
 185       dataOut(kkk)=difpt; kkk++;
 186       Double_t lambdag=TMath::ATan(pzg/ptg);
 187       Double_t   lam=TMath::ATan(tgl);
 188       Double_t diflam = (lam - lambdag)*1000.;
 189       dataOut(kkk) = diflam; kkk++;
 190       Double_t phig=TMath::ATan2(pyg,pxg);  if(phig<0) phig=2.*TMath::Pi()+phig;
 191      // Double_t phi=phivertex;
 192        Double_t phi=TMath::ACos(px/pt);
 193      Double_t duepi=2.*TMath::Pi();
 194      if(phi>duepi) phi-=duepi;
 195      if(phi<0.) phi+=duepi;
 196       Double_t signC=0.;
 197       if(c>0) signC=1.; else signC=-1.;
 198           Double_t Dz=z-zo;   // vertex subtraction
 199           Double_t Dtot= TMath::Sqrt(Dr*Dr+Dz*Dz);
 200       Double_t difphi = (phi - phig)*1000.;
 201       dataOut(kkk)=difphi; kkk++;
 202       dataOut(kkk)=Dtot*1.e4; kkk++;
 203       dataOut(kkk)=Dr*1.e4; kkk++;
 204       dataOut(kkk)=Dz*1.e4; kkk++;
 205       dataOut(kkk)=signC; kkk++;
 206       Int_t r;
 207       for (r=0; r<10; r++) { out<<dataOut(r)<<" ";}
 208       out<<"\n";
 209       kkk=0;
 210
 211     delete iotrack;
 212    }
 213   //out.close();
 214    }   // event loop
 215   out.close();
 216   in.close();
 217   file->Close();
 218 }
 219