New AliITSclustererV2 class: update of V2 test macros
authormasera <masera@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 26 Feb 2003 13:43:41 +0000 (13:43 +0000)
committermasera <masera@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 26 Feb 2003 13:43:41 +0000 (13:43 +0000)
ITS/AliITSComparisonV2.C
ITS/AliITSFindClustersV2.C
ITS/AliITStestV2.C

index 52a2a70..9d47628 100644 (file)
@@ -12,6 +12,7 @@
   #include "TFile.h"
   #include "TTree.h"
   #include "TH1.h"
+  #include "TH2.h"
   #include "TObjArray.h"
   #include "TStyle.h"
   #include "TCanvas.h"
@@ -42,6 +43,7 @@ Int_t AliITSComparisonV2(Int_t event=0) {
      if (!tracktree) {cerr<<"Can't get a tree with ITS tracks !\n"; return 4;}
      TBranch *tbranch=tracktree->GetBranch("tracks");
      nentr=(Int_t)tracktree->GetEntries();
+
      for (Int_t i=0; i<nentr; i++) {
         AliITStrackV2 *iotrack=new AliITStrackV2;
         tbranch->SetAddress(&iotrack);
@@ -88,7 +90,6 @@ Int_t AliITSComparisonV2(Int_t event=0) {
       } else cerr<<"Can not open file (good_tracks_its) !\n";
       out.close();
    }
-   cerr<<"Number of good tracks : "<<ngood<<endl;
 
    TH1F *hp=new TH1F("hp","PHI resolution",50,-20.,20.); hp->SetFillColor(4);
    TH1F *hl=new TH1F("hl","LAMBDA resolution",50,-20,20);hl->SetFillColor(4);
@@ -97,7 +98,6 @@ Int_t AliITSComparisonV2(Int_t event=0) {
    TH1F *hmpt=new TH1F("hmpt","Transverse impact parameter",30,-300,300); 
    hmpt->SetFillColor(6);
    TH1F *hz=new TH1F("hz","Longitudinal impact parameter",30,-300,300); 
-   //hmpt->SetFillColor(6);
 
    AliITStrackV2 *trk=(AliITStrackV2*)tarray.UncheckedAt(0);
    Double_t pmin=0.1*(100/0.299792458/0.2/trk->GetConvConst());
@@ -111,34 +111,57 @@ Int_t AliITSComparisonV2(Int_t event=0) {
    TH1F *hf=new TH1F("hf","Efficiency for fake tracks",30,pmin,pmax);
    hf->SetFillColor(1); hf->SetFillStyle(3013); hf->SetLineWidth(2);
 
-   TH1F *hptw=new TH1F("hptw","Weghted pt",30,pmax,pmin);
+   //TH1F *hptw=new TH1F("hptw","Weghted pt",30,pmax,pmin);
+
+   TH1F *he =new TH1F("he","dE/dX for pions with 0.4<p<0.5 GeV/c",50,0.,100.);
+   TH2F *hep=new TH2F("hep","dE/dX vs momentum",50,0.,2.,50,0.,200.);
+   hep->SetMarkerStyle(8);
+   hep->SetMarkerSize(0.4);
 
-   while (ngood--) {
-      Int_t lab=gt[ngood].lab, tlab=-1;
-      Double_t pxg=gt[ngood].px, pyg=gt[ngood].py, pzg=gt[ngood].pz;
+
+   Int_t notf[MAX], nnotf=0;
+   Int_t fake[MAX], nfake=0;
+   Int_t mult[MAX], numb[MAX], nmult=0;
+
+   Int_t ng;
+   for (ng=0; ng<ngood; ng++) {
+      Int_t lab=gt[ng].lab, tlab=-1;
+      Double_t pxg=gt[ng].px, pyg=gt[ng].py, pzg=gt[ng].pz;
       Double_t ptg=TMath::Sqrt(pxg*pxg+pyg*pyg);
 
-      if (ptg<pmin) continue;
 
-      hgood->Fill(ptg);
+      if (ptg>pmin) hgood->Fill(ptg);
 
       AliITStrackV2 *track=0;
+      Int_t cnt=0;
       Int_t j;
       for (j=0; j<nentr; j++) {
-          track=(AliITStrackV2*)tarray.UncheckedAt(j);
-          tlab=track->GetLabel();
-          if (lab==TMath::Abs(tlab)) break;
+          AliITStrackV2 *trk=(AliITStrackV2*)tarray.UncheckedAt(j);
+          Int_t lbl=trk->GetLabel();
+          if (lab==TMath::Abs(lbl)) {
+           if (cnt==0) {track=trk; tlab=lbl;}
+            cnt++;
+          }
       }
-      if (j==nentr) {
-       cerr<<"Track "<<lab<<" was not found !\n";
+      if (cnt==0) {
+        notf[nnotf++]=lab;
         continue;
+      } else if (cnt>1){
+        mult[nmult]=lab;
+        numb[nmult]=cnt; nmult++;        
       }
+
+      if (ptg>pmin) {
+        if (lab==tlab) hfound->Fill(ptg);
+        else {
+          fake[nfake++]=lab;
+          hfake->Fill(ptg); 
+        }
+      }
+
       track->PropagateTo(3.,0.0028,65.19);
       track->PropagateToVertex();
 
-      if (lab==tlab) hfound->Fill(ptg);
-      else { hfake->Fill(ptg); cerr<<lab<<" fake\n";}
-
       Double_t xv,par[5]; track->GetExternalParameters(xv,par);
       Float_t phi=TMath::ASin(par[2]) + track->GetAlpha();
       if (phi<-TMath::Pi()) phi+=2*TMath::Pi();
@@ -155,20 +178,45 @@ Int_t AliITSComparisonV2(Int_t event=0) {
       Double_t d=10000*track->GetD();
       hmpt->Fill(d);
 
-      hptw->Fill(ptg,TMath::Abs(d));
+      //hptw->Fill(ptg,TMath::Abs(d));
 
       Double_t z=10000*track->GetZ();
       hz->Fill(z);
 
-//if (TMath::Abs(gt[ngood].code)==11 && ptg>4.)
       hpt->Fill((pt_1 - 1/ptg)/(1/ptg)*100.);
 
+      Float_t mom=1./(pt_1*TMath::Cos(lam));
+      Float_t dedx=track->GetdEdx();
+      hep->Fill(mom,dedx,1.);
+      if (TMath::Abs(gt[ng].code)==211)
+         if (mom>0.4 && mom<0.5) he->Fill(dedx,1.);
+   }
+
+
+   cout<<"\nList of Not found tracks :\n";
+   for (ng=0; ng<nnotf; ng++){
+     cout<<notf[ng]<<"\t";
+     if ((ng%9)==8) cout<<"\n";
+   }
+   cout<<"\n\nList of fake  tracks :\n";
+   for (ng=0; ng<nfake; ng++){
+     cout<<fake[ng]<<"\t";
+     if ((ng%9)==8) cout<<"\n";
+   }
+   cout<<"\n\nList of multiple found tracks :\n";
+   for (ng=0; ng<nmult; ng++) {
+       cout<<"id.   "<<mult[ng]
+           <<"     found - "<<numb[ng]<<"times\n";
    }
+   cout<<endl;
 
-   Stat_t ng=hgood->GetEntries(); cerr<<"Good tracks "<<ng<<endl;
-   Stat_t nf=hfound->GetEntries();
-   if (ng!=0) 
-      cerr<<"Integral efficiency is about "<<nf/ng*100.<<" %\n";
+   cerr<<"Number of good tracks : "<<ngood<<endl;
+   cerr<<"Number of found tracks : "<<nentr<<endl;
+
+   ng=(Int_t)hgood->GetEntries(); //cerr<<"Good tracks "<<ng<<endl;
+   if (ng!=0)  
+   cerr<<"Integral efficiency is about "<<hfound->GetEntries()/ng*100.<<" %\n";
+   cerr<<endl;
 
    gStyle->SetOptStat(111110);
    gStyle->SetOptFit(1);
@@ -227,6 +275,22 @@ Int_t AliITSComparisonV2(Int_t event=0) {
    text->SetTextSize(0.05);
    text->Draw();
 
+
+
+   TCanvas *c2=new TCanvas("c2","",320,32,530,590);
+
+   TPad *p6=new TPad("p6","",0.,0.,1.,.5); p6->Draw();
+   p6->cd(); p6->SetFillColor(42); p6->SetFrameFillColor(10); 
+   he->SetFillColor(2); he->SetFillStyle(3005);  
+   he->SetXTitle("Arbitrary Units"); 
+   he->Fit("gaus"); c2->cd();
+
+   TPad *p7=new TPad("p7","",0.,0.5,1.,1.); p7->Draw(); 
+   p7->cd(); p7->SetFillColor(42); p7->SetFrameFillColor(10);
+   hep->SetXTitle("p (Gev/c)"); hep->SetYTitle("dE/dX (Arb. Units)"); 
+   hep->Draw(); c1->cd();
+
+
    return 0;
 }
 
index 0703baa..9521a59 100644 (file)
+/****************************************************************************
+ *           Origin: I.Belikov, CERN, Jouri.Belikov@cern.ch                 *
+ ****************************************************************************/
+
 #ifndef __CINT__
-  #include <Riostream.h>
+  #include <iostream.h>
 
   #include "AliRun.h"
   #include "AliITS.h"
   #include "AliITSgeom.h"
-  #include "AliITSRecPoint.h"
-  #include "AliITSclusterV2.h"
+  #include "AliITSclustererV2.h"
 
   #include "TFile.h"
-  #include "TTree.h"
-  #include "TParticle.h"
+  #include "TStopwatch.h"
 #endif
 
-Int_t AliITSFindClustersV2(Char_t SlowOrFast='f') {
-/****************************************************************
- *  This macro converts AliITSRecPoint(s) to AliITSclusterV2(s) *
- ****************************************************************/
-   cerr<<"AliITSRecPoint(s) -> AliITSclusterV2(s)...\n";
+Int_t AliITSFindClustersV2(Char_t SlowOrFast='s',Int_t eventn=1) {
 
-   if (gAlice) {delete gAlice; gAlice=0;}
+   cerr<<"Looking for clusters...\n";
 
    TFile *in=TFile::Open("galice.root");
-   if (!in->IsOpen()) { cerr<<"Can't open galice.root !\n"; return 1; }
-
+   if (!in->IsOpen()) {cerr<<"Can't open galice.root !\n"; return 2;}
    if (!(gAlice=(AliRun*)in->Get("gAlice"))) {
       cerr<<"Can't find gAlice !\n";
       return 2;
    }
-   gAlice->GetEvent(0);
-
    AliITS *ITS  = (AliITS*)gAlice->GetModule("ITS");
    if (!ITS) { cerr<<"Can't find the ITS !\n"; return 3; }
+
    AliITSgeom *geom=ITS->GetITSgeom();
+
    TFile *out=TFile::Open("AliITSclustersV2.root","new");
    if (!out->IsOpen()) {
-      cerr<<"Delete old AliITSclustersV2.root !\n"; 
-      return 4;
-   }
+      cerr<<"Delete old AliITSclustersV2.root !\n"; return 1;}
    geom->Write();
 
-   TClonesArray *clusters=new TClonesArray("AliITSclusterV2",10000);
-   TTree *cTree=new TTree("TreeC_ITS_0","ITS clusters");
-   cTree->Branch("Clusters",&clusters);
-
-   TTree *pTree=gAlice->TreeR();
-   if (!pTree) { cerr<<"Can't get TreeR !\n"; return 5; }
-   TBranch *branch = 0;
-   if (SlowOrFast=='f') {
-     branch = pTree->GetBranch("ITSRecPointsF");
-   }
-   else {
-     branch = pTree->GetBranch("ITSRecPoints");
-   }
-   if (!branch) { cerr<<"Can't get ITSRecPoints branch !\n"; return 6; }
-   TClonesArray *points=new TClonesArray("AliITSRecPoint",10000);
-   branch->SetAddress(&points);
-
-   TClonesArray &cl=*clusters;
-   Int_t nclusters=0;
-   Int_t nentr=(Int_t)branch->GetEntries();
-
-   cerr<<"Number of entries: "<<nentr<<endl;
-
-   //   Float_t lp[5]; Int_t lab[6]; //Why can't it be inside a loop ?
-   Float_t * lp = new Float_t[5]; 
-   Int_t * lab = new Int_t[6];
-
-   for (Int_t i=0; i<nentr; i++) {
-       points->Clear();
-       branch->GetEvent(i);
-       Int_t ncl=points->GetEntriesFast(); if (ncl==0){cTree->Fill();continue;}
-       Int_t lay,lad,det; geom->GetModuleId(i,lay,lad,det);
-       Float_t x,y,zshift; geom->GetTrans(lay,lad,det,x,y,zshift); 
-       Double_t rot[9];    geom->GetRotMatrix(lay,lad,det,rot);
-       Double_t yshift = x*rot[0] + y*rot[1];
-       Int_t ndet=(lad-1)*geom->GetNdetectors(lay) + (det-1);
-       nclusters+=ncl;
-
-       Float_t kmip=1; // ADC->mip normalization factor for the SDD and SSD 
-       if(lay==4 || lay==3){kmip=280.;};
-       if(lay==6 || lay==5){kmip=38.;};
-
-       for (Int_t j=0; j<ncl; j++) {
-          AliITSRecPoint *p=(AliITSRecPoint*)points->UncheckedAt(j);
-          //Float_t lp[5];
-          lp[0]=-p->GetX()-yshift; if (lay==1) lp[0]=-lp[0];
-          lp[1]=p->GetZ()+zshift;
-          lp[2]=p->GetSigmaX2();
-          lp[3]=p->GetSigmaZ2();
-          lp[4]=p->GetQ(); lp[4]/=kmip;
-          //Int_t lab[6]; 
-          lab[0]=p->GetLabel(0);lab[1]=p->GetLabel(1);lab[2]=p->GetLabel(2);
-          lab[3]=ndet;
-
-          Int_t label=lab[0];
-          if (label>=0) {
-             TParticle *part=(TParticle*)gAlice->Particle(label);
-             label=-3;
-             while (part->P() < 0.005) {
-                Int_t m=part->GetFirstMother();
-                if (m<0) {cerr<<"Primary momentum: "<<part->P()<<endl; break;}
-                label=m;
-                part=(TParticle*)gAlice->Particle(label);
-             }
-             if      (lab[1]<0) lab[1]=label;
-             else if (lab[2]<0) lab[2]=label;
-             else cerr<<"No empty labels !\n";
-         }
-
-          new(cl[j]) AliITSclusterV2(lab,lp);
-       }
-       cTree->Fill(); clusters->Delete();
-       points->Delete();
+   TStopwatch timer;
+   AliITSclustererV2 clusterer(geom);
+   for (Int_t i=0; i<eventn; i++) {
+       cerr<<"Processing event number: "<<i<<endl;
+       gAlice->GetEvent(i);
+       //ITS->MakeTreeC(); //To make the V1 cluster finders happy
+       clusterer.SetEvent(i);
+       if (SlowOrFast=='s') clusterer.Digits2Clusters(in,out);
+       else                 clusterer.Hits2Clusters(in,out);
    }
-   cTree->Write();
-
-   cerr<<"Number of clusters: "<<nclusters<<endl;
-
-   delete [] lp;
-   delete [] lab;
-
-   delete cTree; delete clusters; delete points;
+   timer.Stop(); timer.Print();
 
    delete gAlice; gAlice=0;
-
-   in->Close();
    out->Close();
+   in->Close();
 
    return 0;
-
 }
 
 
 
-
-
-
-
-
-
-
-
-
-
-
index 4da897e..09a3c58 100644 (file)
@@ -16,17 +16,17 @@ Int_t AliITStestV2(Char_t SlowOrFast='s') {
    in->Close();
 
    if (SlowOrFast=='f') {
-      cerr<<"Fast AliITSRecPoint(s) !\n";
-      gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSHits2FastRecPoints.C");
-      AliITSHits2FastRecPoints();
+      //cerr<<"Fast AliITSRecPoint(s) !\n";
+      //gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSHits2FastRecPoints.C");
+      //AliITSHits2FastRecPoints();
    } else {
       cerr<<"Slow AliITSRecPoint(s) !\n";
       gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSHits2SDigits.C");
       AliITSHits2SDigits();
       gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSSDigits2Digits.C");
       AliITSSDigits2Digits();
-      gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSDigits2RecPoints.C");
-      AliITSDigits2RecPoints();
+      //gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSDigits2RecPoints.C");
+      //AliITSDigits2RecPoints();
    }
    gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSFindClustersV2.C");
    if (rc=AliITSFindClustersV2(SlowOrFast)) return rc;
@@ -36,7 +36,7 @@ Int_t AliITStestV2(Char_t SlowOrFast='s') {
 
    gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSComparisonV2.C");
    if (rc=AliITSComparisonV2()) return rc;
-
+/*
    gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliV0FindVertices.C");
    if (rc=AliV0FindVertices()) return rc;
 
@@ -49,6 +49,6 @@ Int_t AliITStestV2(Char_t SlowOrFast='s') {
    gSystem->SetIncludePath("-I$ROOTSYS/include -I$ALICE_ROOT/include -I$ALICE_ROOT/ITS -I$ALICE_ROOT/TPC -I$ALICE_ROOT/CONTAINERS");
    gROOT->ProcessLine(".L $(ALICE_ROOT)/ITS/AliCascadeComparison.C+");
    if (rc=AliCascadeComparison()) return rc;
-
+*/
    return rc;
 }