--- /dev/null
+#include "iostream.h"
+
+void AliITSComparisonV1
+(const char *foundfile = "itstracks.root", const char *truefile = "galice_tracks_0.root", Int_t evnum = 0) {
+
+ // Make sure that ALICE objects are loaded
+ if (gClassTable->GetID("AliRun") < 0) {
+ gROOT->LoadMacro("loadlibs.C");
+ loadlibs();
+ }
+ else {
+ delete gAlice;
+ gAlice=0;
+ }
+
+ // Load the tree with track data
+ TFile *fileTrueTracks = new TFile(truefile);
+ TTree *treeTrueTracks = (TTree*)fileTrueTracks->Get("Tracks");
+ TFile *fileFoundTracks = new TFile(foundfile);
+ TTree *treeFoundTracks = (TTree*)fileFoundTracks->Get(Form("TreeT%d", evnum));
+
+ Int_t nTrue = (Int_t)treeTrueTracks->GetEntries();
+ Int_t *cnt = new Int_t[nTrue];
+ for (Int_t i = 0; i < nTrue; i++) cnt[i] = 0;
+
+ Int_t label, tpc_ok, entry, pdg_code;
+ Double_t px, py, pz, pt;
+ treeTrueTracks->SetBranchAddress("label", &label);
+ treeTrueTracks->SetBranchAddress("entry", &entry);
+ treeTrueTracks->SetBranchAddress("tpc_ok", &tpc_ok);
+ treeTrueTracks->SetBranchAddress("pdg_code", &pdg_code);
+ treeTrueTracks->SetBranchAddress("px", &px);
+ treeTrueTracks->SetBranchAddress("py", &py);
+ treeTrueTracks->SetBranchAddress("pz", &pz);
+ treeTrueTracks->SetBranchAddress("pt", &pt);
+
+ AliITSIOTrack *iotrack = 0;
+ treeFoundTracks->SetBranchAddress("ITStracks", &iotrack);
+
+ const Int_t npos = 36, nneg = 31;
+ Int_t pos[npos] = {2212, -11, -13, 211, 321, 3222, 213, 323, 10323, 3224,
+ 2224, 2214, -1114, -3112, -3312, 3224, -3114, -3314, 411,
+ 431, 413, 433, -15, 4232, 4222, 4322, 4422, 4412, 4432,
+ 4224, 4214, 4324, 4424, 4414, 4434, 4444};
+ Int_t neg[nneg] = {2212, 11, 13, -211, -321, 3112, -213, -323, -10323, 3114,
+ 1114, -2224, -2214, 33112, -3222, 3114, 3314, 3334, -3224,
+ -411, -431, -413, -433, 15, -4422, -4432, -4214, -4324,
+ -4424, -4434, -444};
+
+ // Evaluation tree definition
+ Int_t labITS, labTPC, signC;
+ Double_t difpt, diflambda, difphi, Dz, Dr, Dtot, ptg;
+ TTree *treeEvaluation = new TTree("Eval", "Parameters for tracking evaluation");
+ treeEvaluation->Branch("labITS" , &labITS , "labITS/I" );
+ treeEvaluation->Branch("labTPC" , &labTPC , "labTPC/I" );
+ treeEvaluation->Branch("difpt" , &difpt , "difpt/D" );
+ treeEvaluation->Branch("diflambda", &diflambda, "diflambda/D");
+ treeEvaluation->Branch("difphi" , &difphi , "difphi/D" );
+ treeEvaluation->Branch("Dz" , &Dz , "Dz/D" );
+ treeEvaluation->Branch("Dr" , &Dr , "Dr/D" );
+ treeEvaluation->Branch("Dtot" , &Dtot , "Dtot/D" );
+ treeEvaluation->Branch("ptg" , &ptg , "ptg/D" );
+ treeEvaluation->Branch("signC" , &signC , "signC/I" );
+
+ // Make comparison
+ Double_t *var = 0;
+ Bool_t isGood;
+ Int_t nOK, trueEntry;
+ Double_t found_px, found_py, found_pz, found_pt;
+ Double_t found_tgl, found_lambda, found_phi, found_curv;
+ Double_t true_lambda, true_phi, true_px, true_py, true_pz, true_pt;
+ Double_t duepi = 2.*TMath::Pi();
+ Bool_t ispos, isneg;
+ Int_t countpos = 0, countneg = 0, found_charge;
+ for (Int_t i = 0; i < treeFoundTracks->GetEntries(); i++) {
+ treeFoundTracks->GetEntry(i);
+ labITS = iotrack->GetLabel();
+ labTPC = iotrack->GetTPCLabel();
+ isGood = (labITS >= 0);
+ nOK = treeTrueTracks->Draw("px:py:pz", Form("label==%d && tpc_ok==1", abs(labITS)), "goff");
+ if (!nOK) {
+ cerr << "ITS label not found among findable tracks:";
+ cerr << " labITS = " << labITS;
+ cerr << " labTPC = " << labTPC;
+ cerr << endl;
+ isGood = kFALSE;
+ }
+ if (nOK > 1) {
+ cerr << "More than 1 tracks with label " << labITS << ": taking the first" << endl;
+ }
+ true_px = *treeTrueTracks->GetV1();
+ true_py = *treeTrueTracks->GetV2();
+ true_pz = *treeTrueTracks->GetV3();
+ true_pt = TMath::Sqrt(true_px * true_px + true_py * true_py);
+ true_phi = TMath::ATan2(true_py, true_px);
+ if(true_phi < 0) true_phi += duepi;
+ true_phi *= 1000.0;
+ true_lambda = TMath::ATan(true_pz / true_pt) * 1000.0;
+ ptg = true_pt;
+
+ // checks if two found good tracks have the same label in ITS
+ treeTrueTracks->Draw("entry", Form("label==%d && tpc_ok==1", abs(labITS)), "goff");
+ trueEntry = (Int_t)*treeTrueTracks->GetV1();
+ if (isGood && cnt[trueEntry] == 0)
+ cnt[trueEntry] = 1;
+ else if (isGood) {
+ cout << "Entry " << trueEntry << " already analyzed!" << endl;
+ continue;
+ }
+
+ // charge matching
+ found_charge = (Int_t)iotrack->GetCharge();
+ ispos = isneg = kFALSE;
+ for (Int_t j = 0; j < npos; j++) ispos = (pdg_code == pos[j]);
+ for (Int_t j = 0; j < nneg; j++) isneg = (pdg_code == neg[j]);
+ if (ispos) countpos++;
+ if (isneg) countneg++;
+
+ // pt resolution (%)
+ found_px = iotrack->GetPx();
+ found_py = iotrack->GetPy();
+ found_pz = iotrack->GetPz();
+ found_pt = TMath::Sqrt(found_px*found_px + found_py*found_py);
+ difpt = ((2. * found_pt - true_pt) / true_pt) * 100.;
+
+ // lambda (mrad)
+ found_tgl = iotrack->GetStateTgl();
+ found_lambda = TMath::ATan(found_tgl) * 1000.0;
+ diflambda = found_lambda - true_lambda;
+// cout << "lambda " << found_lambda << " " << true_lambda << " " << diflambda << endl;
+
+ // phi (mrad)
+ found_phi = TMath::ACos(found_px / found_pt);
+ if(found_phi > duepi) found_phi -= duepi;
+ if(found_phi < 0.) found_phi += duepi;
+ found_phi *= 1000.0;
+ difphi = found_phi - true_phi;
+// cout << "phi " << found_phi << " " << true_phi << " " << difphi << endl;
+
+ // impact parameters (microns)
+ Dr = iotrack->GetStateD() * 1.e4;
+ Dz = iotrack->GetDz() * 1.e4;
+ Dtot = sqrt(Dr*Dr + Dz*Dz);
+
+ // curvature sign
+ found_curv = iotrack->GetStateC();
+ signC = (found_curv > 0.) ? 1 : -1;
+
+ // fill tree
+ treeEvaluation->Fill();
+ }
+
+
+ TFile *outFile = new TFile("AliITSComparisonV1.root", "recreate");
+ treeEvaluation->Write("Eval", TObject::kOverwrite);
+ outFile->Close();
+}
+
#include <iostream.h>
#include <fstream.h>
-void AliITSPlotTracksV1(){
+void AliITSPlotTracksV1()
+{
+ TFile *file = new TFile("AliITSComparisonV1.root");
+ TTree *tree = (TTree*)file->Get("Eval");
+
+ if (!tree) {
+ cerr << "Evaluation tree not found!" << endl;
+ return;
+ }
+
+ // histogram definition - I (efficiency)
+ TH1D *hFindables = new TH1D("hFindables", "Findable tracks", 10, 0, 2);
+ TH1D *hGood = new TH1D("hGood", "Good found tracks", 10, 0, 2);
+ TH1D *hFake = new TH1D("hFake", "Fake found tracks", 10, 0, 2);
+ TH1D *hRatioG = new TH1D("hRatioG", "", 10, 0, 2);
+ TH1D *hRatioF = new TH1D("hRatioF", "", 10, 0, 2);
+ hGood->Sumw2();
+ hFake->Sumw2();
+ hFindables->Sumw2();
+ hRatioG->SetLineColor(kBlue);
+ hRatioG->SetLineWidth(2);
+ hRatioF->SetLineColor(kRed);
+ hRatioF->SetLineWidth(2);
+
+ // histograms definition - II (resolution)
+ TH1D *hPhi = new TH1D("hPhi", "#Phi resolution", 50, -15., 15.);
+ TH1D *hLambda = new TH1D("hLambda", "#lambda resolution", 50, -15., 15.);
+ TH1D *hPt = new TH1D("hPt", "Relative Pt resolution", 40, -10., 10.);
+ TH1D *hDtot = new TH1D("hDtot", "Total impact parameter distribution", 100, 0., 2000.);
+ TH1D *hDr = new TH1D("hDr", "Transv. impact parameter distribution", 50, -1000., 1000.);
+ TH1D *hDz = new TH1D("hDz", "Long. impact parameter distribution", 50, -1000., 1000.);
+ hPhi->SetFillColor(4);
+ hLambda->SetFillColor(4);
+ hPt->SetFillColor(2);
+ hDtot->SetFillColor(6);
+ hDr->SetFillColor(kGreen);
+ hDz->SetFillColor(kBlue);
+
+ // Evaluation tree settings
+ Int_t labITS, labTPC, signC;
+ Int_t i, j, tot = (Int_t)tree->GetEntries();
+ Double_t difpt, diflambda, difphi, Dz, Dr, Dtot, ptg;
+ tree->SetBranchAddress("labITS" , &labITS );
+ tree->SetBranchAddress("labTPC" , &labTPC );
+ tree->SetBranchAddress("difpt" , &difpt );
+ tree->SetBranchAddress("diflambda", &diflambda);
+ tree->SetBranchAddress("difphi" , &difphi );
+ tree->SetBranchAddress("Dz" , &Dz );
+ tree->SetBranchAddress("Dr" , &Dr );
+ tree->SetBranchAddress("Dtot" , &Dtot );
+ tree->SetBranchAddress("ptg" , &ptg );
+ tree->SetBranchAddress("signC" , &signC );
+
+ // Filling the histogram of findable tracks (w.r. to momentum)
+ for(i = 0; i < tot; i++) {
+ tree->GetEntry(i);
+ hFindables->Fill(ptg);
+ }
+
+ // Filling the evaluation histograms
+ Int_t neglabs = 0;
+ for(i = 0; i < tot; i++) {
+ tree->GetEntry(i);
+// if(signC<0) continue;
+ if (labITS < 0) neglabs++;
+ if (labITS >= 0) {
+ hGood->Fill(ptg);
+ hPt->Fill(difpt);
+ hLambda->Fill(diflambda);
+ hPhi->Fill(difphi);
+ hDtot->Fill(Dtot);
+ hDr->Fill(Dr);
+ hDz->Fill(Dz);
+ }
+ else {
+ hFake->Fill(ptg);
+ neglabs++;
+ }
+ }
- ifstream in ("AliITSTra.out");
-
- TVector DataOut(10);
-
- ///////////////////////////////// Histograms definition ///////////////////////////////////////////
-
-
- TH1F *hp=new TH1F("hp","PHI resolution",50,-15.,15.); hp->SetFillColor(4);
- TH1F *hl=new TH1F("hl","LAMBDA resolution",50,-15.,15.); hl->SetFillColor(4);
-
- // TH1F *hp=new TH1F("hp","PHI resolution",50,-5.,5.); hp->SetFillColor(4);
- //TH1F *hl=new TH1F("hl","LAMBDA resolution",50,-5.,5.); hl->SetFillColor(4);
-
- // TH1F *hpt=new TH1F("hpt","Relative Pt resolution",50,-5.,5.);
- TH1F *hpt=new TH1F("hpt","Relative Pt resolution",40,-10.,10.);
- hpt->SetFillColor(2);
- TH1F *hd=new TH1F("hd","Impact parameter distribution ",100,0,2000);
- hd->SetFillColor(6);
-
- TH1F *hdr=new TH1F("hdr","Dr ",50,-1000,1000);
- hdr->SetFillColor(kGreen);
- TH1F *hdz=new TH1F("hdz","Dz ",50,-1000,1000);
- hdz->SetFillColor(kBlue);
-
- /*
- TH1F *hdr=new TH1F("hdr","Dr ",50,-100,100);
- hdr->SetFillColor(kGreen);
- TH1F *hdz=new TH1F("hdz","Dz ",50,-1000,1000);
- hdz->SetFillColor(kBlue);
- */
-
- TH1F *hgood=new TH1F("hgood","Good tracks",10,0,2);
- hgood->Sumw2(); // aggiunto il 15-01-01
- TH1F *hfound=new TH1F("hfound","Found tracks",10,0,2);
- hfound->Sumw2(); // aggiunto il 15-01-01
- TH1F *hfake=new TH1F("hfake","Fake tracks",10,0,2);
- hfake->Sumw2();
- TH1F *hg=new TH1F("hg","",10,0,2); //efficiency for good tracks
- hg->SetLineColor(4); hg->SetLineWidth(2);
- TH1F *hf=new TH1F("hf","Efficiency for fake tracks",10,0,2);
- /*hf->SetFillColor(1); hf->SetFillStyle(3013);*/ hf->SetLineColor(4); hf->SetLineWidth(2);
-
- ///////////////////////////////////////////////////////////////////////////////////////////////////
-
- ifstream in1 ("AliITSTrag.out");
-
- Double_t ptg;
- for(;;) {
- in1 >> ptg;
- if( in1.eof() ) break;
- hgood->Fill(ptg);
- }
- in1.close();
-
- Int_t neglabs=0;
- for (;;){
- for (int r=0; r<10; r++) in>>DataOut(r);
- if( in.eof() ) break;
-
- Double_t ptg=DataOut(0); Double_t labITS=DataOut(1); Double_t labTPC=DataOut(2); Double_t ptperc=DataOut(3);
- Double_t deltalam=DataOut(4); Double_t deltaphi=DataOut(5);
- Double_t Dtot=DataOut(6); Double_t Dr=DataOut(7); Double_t Dz=DataOut(8);
- Double_t signC=DataOut(9);
-
- // if(signC<0) continue;
-
- if(labITS<0) neglabs++;
- if(labITS>=0) hfound->Fill(ptg); else { hfake->Fill(ptg);}
-
- if(labITS>=0 ) {
- hpt->Fill(ptperc);
- hl->Fill(deltalam);
- hp->Fill(deltaphi);
- hd->Fill(Dtot);
- hdr->Fill(Dr);
- hdz->Fill(Dz);
- }
- }
-
-
- in.close();
- Stat_t ngood=hgood->GetEntries(); cerr<<"Good tracks "<<ngood<<endl;
- Stat_t nfound=hfound->GetEntries(); cerr<<"Found tracks "<<nfound<<endl;
- Stat_t nfake=hfake->GetEntries(); cerr<<"Fake tracks "<<nfake<<endl;
- gStyle->SetOptStat(111110);
- gStyle->SetOptFit(1);
- TCanvas *c1=new TCanvas("c1","",0,0,700,700);
- TPad *p1=new TPad("p1","",0,0.5,0.5,1); p1->Draw(); hp->SetXTitle("(mrad)");
- p1->cd(); hp->Draw(); hp->Fit("gaus"); c1->cd();
- TPad *p2=new TPad("p2","",0.5,0.5,1,1); p2->Draw(); hl->SetXTitle("(mrad)");
- p2->cd(); hl->Draw(); hl->Fit("gaus"); c1->cd();
- TPad *p3=new TPad("p3","",0,0,0.5,0.5); p3->Draw(); hpt->SetXTitle("(%)");
- p3->cd(); hpt->Draw(); hpt->Fit("gaus"); c1->cd();
- TPad *p4=new TPad("p4","",0.5,0,1,0.5); p4->Draw(); hd->SetXTitle("(micron)");
- p4->cd(); hd->Draw(); c1->cd();
-
- TCanvas *c3=new TCanvas("c3","",200,200,800,500);
- hfound->Print("all"); // aggiunto il 16-01-01
- hgood->Print("all"); // aggiunto il 16-01-01
- TPad *p7=new TPad("p7","",0,0,0.333,1); p7->Draw(); p7->cd(); hfound->Draw(); c3->cd();
- TPad *p8=new TPad("p8","",0.333,0,0.666,1); p8->Draw(); p8->cd(); hfake->Draw(); c3->cd();
- TPad *p9=new TPad("p9","",0.666,0,1,1); p9->Draw(); p9->cd(); hgood->Draw(); c3->cd();
-
- TCanvas *c4=new TCanvas("c4","",300,300,800,500);
- hg->Divide(hfound,hgood,1.,1.); //,"b");
- hf->Divide(hfake,hgood,1.,1.); //,"b");
- hg->SetMaximum(1.4);
- hg->SetYTitle("Tracking efficiency");
- hg->SetXTitle("Pt (GeV/c)");
- hg->Print("all"); // aggiunto il 16-01-01
- hg->Draw(); // to not plot the erro bars hg->Draw("histo");
- TLine *line1 = new TLine(0,1.0,2,1.0); line1->SetLineStyle(4);
- line1->Draw("same");
- TLine *line2 = new TLine(0,0.9,2,0.9); line2->SetLineStyle(4);
- line2->Draw("histosame");
-
-
- hf->SetFillColor(1);
- hf->SetFillStyle(3013);
- hf->SetLineColor(2);
- hf->SetLineWidth(2);
- hf->Draw("same"); // to not plot the error bars hf->Draw("histosame");
-
- TText *text = new TText(0.461176,0.248448,"Fake tracks");
- text->SetTextSize(0.05);
- text->Draw();
- text = new TText(0.453919,1.11408,"Good tracks");
- text->SetTextSize(0.05);
- text->Draw();
-
- TCanvas *c2=new TCanvas("c2","",100,100,700,400);
- TPad *p5=new TPad("p5","",0,0,0.5,1); p5->Draw(); hdr->SetXTitle("(micron)");
- p5->cd(); hdr->Draw(); hdr->Fit("gaus"); c2->cd();
- TPad *p6=new TPad("p6","",0.5,0,1,1); p6->Draw(); hdz->SetXTitle("(micron)");
- p6->cd(); hdz->Draw(); hdz->Fit("gaus"); c2->cd();
-
- cout<<"neglabs = "<<neglabs<<"\n"; // provvisoria
+ // Drawing
+ cerr << "Findable tracks : " << hFindables->GetEntries() << endl;
+ cerr << "Good found tracks : " << hGood->GetEntries() << endl;
+ cerr << "Fake found tracks : " << hFake->GetEntries() << endl;
+
+ gStyle->SetOptStat(111110);
+ gStyle->SetOptFit(1);
+
+ TCanvas *c1 = new TCanvas("c1","Parameter resolutions",0,0,700,700);
+ c1->Divide(2, 2, 0.001, 0.001);
+ c1->cd(1); hPhi->SetXTitle("(mrad)"); hPhi->Draw(); hPhi->Fit("gaus", "Q");
+ c1->cd(2); hLambda->SetXTitle("(mrad)"); hLambda->Draw(); hLambda->Fit("gaus", "Q");
+ c1->cd(3); hPt->SetXTitle("(%)"); hPt->Draw(); hPt->Fit("gaus", "Q");
+ c1->cd(4); hDtot->SetXTitle("(micron)"); hDtot->Draw();
+ c1->Update();
+
+ TCanvas *c2 = new TCanvas("c2","Impact parameters resolution",100,100,700,400);
+ c2->Divide(2, 1, 0.001, 0.001);
+ c2->cd(1); hDr->SetXTitle("(micron)"); hDr->Draw(); hDr->Fit("gaus", "Q");
+ c2->cd(2); hDz->SetXTitle("(micron)"); hDz->Draw(); hDz->Fit("gaus", "Q");
+ c2->Update();
+
+ TCanvas *c3 = new TCanvas("c3","Momentum distributions",200,200,800,500);
+ c3->Divide(3, 1, 0.001, 0.001);
+ c3->cd(1); hGood->Draw();
+ c3->cd(2); hFake->Draw();
+ c3->cd(3); hFindables->Draw();
+ c3->Update();
+
+ TCanvas *c4 = new TCanvas("c4","Tracking efficiency",300,300,800,500);
+ hRatioG->Divide(hGood, hFindables, 1., 1.);
+ hRatioF->Divide(hFake, hFindables, 1., 1.);
+ hRatioG->SetMaximum(1.4);
+ hRatioG->SetYTitle("Tracking efficiency");
+ hRatioG->SetXTitle("Pt (GeV/c)");
+ hRatioG->Draw(); // to not plot the erro bars hg->Draw("histo");
+ hRatioF->SetFillColor(1);
+ hRatioF->SetFillStyle(3013);
+ hRatioF->SetLineColor(2);
+ hRatioF->SetLineWidth(2);
+ hRatioF->Draw("same"); // to not plot the error bars hRatioF->Draw("histosame");
+ // a line to mark the best efficiency
+ TLine *line1 = new TLine(0,1.0,2,1.0); line1->SetLineStyle(4);
+ line1->Draw("same");
+ TLine *line2 = new TLine(0,0.9,2,0.9); line2->SetLineStyle(4);
+ line2->Draw("histosame");
+ // a text explanation
+ TText *text = new TText(0.461176,0.248448,"Fake tracks");
+ text->SetTextSize(0.05);
+ text->Draw();
+ text = new TText(0.453919,1.11408,"Good tracks");
+ text->SetTextSize(0.05);
+ text->Draw();
+ c4->Update();
+
+ cout << "neglabs = " << neglabs << endl; // provvisoria
}
--- /dev/null
+Int_t AliITSStoreFindableTracks
+(Int_t nMinClusters = 5, const Text_t *evname = "galice", Int_t evnum = 0)
+{
+ gSystem->SetIncludePath("-I- -I$ALICE_ROOT/ITS -I$ALICE_ROOT/STEER");
+
+ gROOT->LoadMacro("$ALICE_ROOT/ITS/AliITSStoreFindableTracksCompiled.C+");
+ AliITSStoreFindableTracksCompiled(nMinClusters, evname, evnum);
+}
--- /dev/null
+#include <fstream.h>
+
+#include <TClassTable.h>
+#include <TFile.h>
+#include <TTree.h>
+#include <TString.h>
+#include <TClonesArray.h>
+#include <TParticle.h>
+
+#include "AliRun.h"
+#include "AliITS.h"
+#include "AliITSgeom.h"
+#include "AliITSRecPoint.h"
+
+Int_t AliITSStoreFindableTracksCompiled
+(Int_t nMinClusters = 5, const Text_t *evname = "galice", Int_t evnum = 0)
+{
+ // Make sure that ALICE objects are loaded
+ if (gAlice) {
+ delete gAlice;
+ gAlice = 0;
+ }
+
+ // Define the names of all involved files
+ TString strEventFile(evname);
+ TString strOutputFile(evname);
+ strEventFile.Append(".root");
+ strOutputFile.Append("_tracks_");
+ strOutputFile += evnum;
+ strOutputFile.Append(".root");
+
+ // Connect the Root Galice file containing Geometry, Kine and Hits
+ TFile *fileEvent = (TFile*)gROOT->GetListOfFiles()->FindObject(strEventFile);
+ if (!fileEvent) fileEvent = new TFile(strEventFile,"UPDATE");
+
+ // Get AliRun object from file
+ gAlice = (AliRun*)fileEvent->Get("gAlice");
+ if (gAlice) cout << "OK, found an AliRun object in file" << endl;
+
+ // Get ITS related objects and data
+ AliITS* ITS =(AliITS *)gAlice->GetDetector("ITS");
+ if (!ITS) {
+ cerr << "ITS object not found!" << endl;
+ return 1;
+ }
+ AliITSgeom *geometry = ITS->GetITSgeom();
+ if (!geometry) {
+ cerr << "ITS geometry object not found!" << endl;
+ return 2;
+ }
+
+ // Count the number of modules per layer
+ Int_t nLadders, nDetectors, mod_min[6], mod_max[6];
+ for(Int_t i = 0; i < 6; i++) {
+ nLadders = geometry->GetNladders(i + 1);
+ nDetectors = geometry->GetNdetectors(i + 1);
+ mod_min[i] = geometry->GetModuleIndex(i + 1, 1, 1);
+ mod_max[i] = geometry->GetModuleIndex(i + 1, nLadders, nDetectors);
+ }
+
+ // Load event and ITS recpoints
+ Int_t nParticles = gAlice->GetEvent(evnum);
+ cout << "Event number: " << evnum << endl;
+ cout << "# particles : " << nParticles <<endl;
+ if (nParticles <= 0) {
+ cerr << "Can't have <= 0 particles!" << endl;
+ return 3;
+ }
+ AliITSRecPoint *recp = 0;
+ TClonesArray *recPoints = ITS->RecPoints();
+ TObjArray *particles = gAlice->Particles();
+ Int_t nTracks = gAlice->GetNtrack(); //FCA correction
+ Bool_t *hitITSLayer[6];
+ for (Int_t i = 0; i < 6; i++) {
+ hitITSLayer[i] = new Bool_t[nTracks];
+ for (Int_t j = 0; j < nTracks; j++) hitITSLayer[i][j] = kFALSE;
+ }
+
+ // Load recpoints in event
+ TTree *TR = gAlice->TreeR();
+ if (!TR) {
+ cerr << "TreeR object not found!" << endl;
+ return 4;
+ }
+
+ // Scan recpoints and define findable tracks
+ Int_t nModules = (Int_t)TR->GetEntries(), nPoints = 0;
+ cout << "Found " << nModules;
+ cout << " entries in the TreeR (must be one per module!)" << endl;
+ for (Int_t layer = 1; layer <= 6; layer++) {
+ for (Int_t mod = mod_min[layer - 1]; mod <= mod_max[layer - 1]; mod++) {
+ ITS->ResetRecPoints();
+ TR->GetEntry(mod);
+ nPoints = recPoints->GetEntries();
+ if(!nPoints) {
+ cout << "Module " << mod << " is empty..." << endl;
+ continue;
+ }
+ for (Int_t point = 0; point < nPoints; point++) {
+ recp = (AliITSRecPoint*)recPoints->UncheckedAt(point);
+ for (Int_t it = 0; it < 3; it++) {
+ Int_t track = recp->GetLabel(it);
+ if(track < 0) continue;
+ if(track > nTracks) {
+ cout << "Found track index " << track;
+ cout << " whilw gAlice->GetNtrack() = " << nTracks << endl;
+ continue;
+ }
+ hitITSLayer[layer - 1][track] = kTRUE;
+ } // loop over recpoint labels
+ } //loop over points
+ } //loop over modules
+ } //loop over layers
+
+ // Scan the file of tracks in TPC to retrieve the findable TPC tracks
+ TString strLabelsTPC;
+ Int_t label, pdg_code, nFindablesTPC = 0;
+ Double_t dummy;
+ ifstream tpc("good_tracks_tpc");
+ while (tpc >> label >> pdg_code) {
+ for (Int_t i = 0; i < 6; i++) tpc >> dummy;
+ nFindablesTPC++;
+ strLabelsTPC.Append(Form("[%d]", label));
+ }
+
+ // Define the TTree with tracks data by means of a set of variables
+ Int_t nFindablesITS = 0, nFindablesITSTPC = 0;
+ Int_t nhits, tpc_ok, entry = 0;
+ Double_t vx, vy, vz;
+ Double_t px, py, pz, pt;
+
+ TTree *tree = new TTree("Tracks", "Findable tracks in ITS");
+
+ tree->Branch("vx", &vx, "vx/D");
+ tree->Branch("vy", &vy, "vy/D");
+ tree->Branch("vz", &vz, "vz/D");
+ tree->Branch("px", &px, "px/D");
+ tree->Branch("py", &py, "py/D");
+ tree->Branch("pz", &pz, "pz/D");
+ tree->Branch("pt", &pt, "pt/D");
+ tree->Branch("label", &label, "label/I");
+ tree->Branch("entry", &entry, "entry/I");
+ tree->Branch("pdg_code", &pdg_code, "pdg_code/I");
+ tree->Branch("nhits", &nhits, "nhits/I");
+ tree->Branch("tpc_ok", &tpc_ok, "tpc_ok/I");
+
+ // Fill the tree
+ cout << endl;
+ TParticle *p = 0;
+ for (Int_t i = 0; i < nTracks; i++) {
+ p = gAlice->Particle(i);
+ px = p->Px();
+ py = p->Py();
+ pz = p->Pz();
+ pt = p->Pt();
+ vx = p->Vx();
+ vy = p->Vy();
+ vz = p->Vz();
+ nhits = 0;
+ for (Int_t j = 0; j < 6; j++) if (hitITSLayer[j][i]) nhits++;
+ if (nhits < nMinClusters) continue;
+ cout << "Track " << i << " stored\r" << flush;
+ tpc_ok = (strLabelsTPC.Contains(Form("[%d]", i)));
+ pdg_code = p->GetPdgCode();
+ label = i;
+ nFindablesITS++;
+ if (tpc_ok) nFindablesITSTPC++;
+ tree->Fill();
+ entry++;
+ }
+
+ // Save into a file
+ TFile *fileOutput = new TFile(strOutputFile, "recreate");
+ tree->Write();
+ fileOutput->Close();
+
+ cout << "# findable tracks in TPC : " << nFindablesTPC << endl;
+ cout << "# findable tracks in ITS : " << nFindablesITS << endl;
+ cout << "# findable tracks in ITS+TPC : " << nFindablesITSTPC << endl;
+
+ return 0;
+}
/*
$Log$
+Revision 1.22 2002/10/22 14:45:36 alibrary
+Introducing Riostream.h
+
Revision 1.21 2002/02/05 09:12:26 hristov
Small mods for gcc 3.02
kkprov->SetConvConst(100/0.299792458/0.2/fFieldFactor);
TFile *cf=TFile::Open("AliTPCclusters.root");
- AliTPCParam *digp= (AliTPCParam*)cf->Get("75x40_100x60");
+ AliTPCParam *digp= (AliTPCParam*)cf->Get("75x40_100x60_150x60");
if (!digp) { cerr<<"TPC parameters have not been found !\n"; getchar();}
cf->cd();
#include "iostream.h"
#endif
+Bool_t TPCSortTracks(Int_t event = 0)
+{
+ TFile *fileTracks = TFile::Open("AliTPCtracks.root");
+ TFile *fileClusters = TFile::Open("AliTPCclusters.root");
+ TFile *fileEvent = TFile::Open("galice.root");
+
+ // get TPC parameterization
+ AliTPCParam *param=(AliTPCParam *)fileEvent->Get("75x40_100x60_150x60");
+ if (!param) {
+ cerr << "(TPCSortTracks) ERROR: TPC parameters have not been found!" << endl;
+ return kFALSE;
+ }
+
+ // read and sort tracks
+ Int_t i;
+ TSortedList tracks_list;
+ AliTPCtrack *iotrack = 0;
+ TTree *tracktree = (TTree*)fileTracks->Get(Form("TreeT_TPC_%d", event));
+ Int_t nentr = (Int_t)tracktree->GetEntries();
+ for (i = 0; i < nentr; i++) {
+ iotrack = new AliTPCtrack;
+ tracktree->SetBranchAddress("tracks", &iotrack);
+ tracktree->GetEvent(i);
+ tracks_list.Add(iotrack);
+ }
+ delete tracktree;
+
+ // assign to each track its GEANT label
+ fileClusters->cd();
+ AliTPCtracker *tracker = new AliTPCtracker(param, event);
+ tracker->LoadInnerSectors();
+ tracker->LoadOuterSectors();
+ TListIter iter(&tracks_list);
+ for (i = 0; i < nentr; i++) {
+ iotrack = (AliTPCtrack*)iter.Next();
+ if (!iotrack) {
+ cerr << "(TPCSortTracks) WARNING: Track no. " << i << " is NULL!!!" << endl;
+ continue;
+ }
+ tracker->CookLabel(iotrack, 0.1);
+ }
+ delete tracker;
+
+ // create the new TTree of TPC tracks sorted w.r. to Pt
+ tracktree = new TTree(Form("TreeT_TPC_%d", event),"Tree with TPC tracks sorted w.r to pt");
+ tracktree->Branch("tracks", "AliTPCtrack", &iotrack, 32000, 0);
+ iter.Reset();
+ for (i = 0; i < nentr; i++) {
+ iotrack = (AliTPCtrack*)iter.Next();
+ if (!iotrack) {
+ cerr << "(TPCSortTracks) WARNING: Track no. " << i << " is NULL!!!" << endl;
+ continue;
+ }
+ tracktree->Fill();
+ }
+
+ // save the new tree into new file
+ TFile *fileOutput = TFile::Open("AliTPCtracksSorted.root","recreate");
+ tracktree->Write();
+ fileOutput->Close();
+ fileEvent->Close();
+ fileClusters->Close();
+ fileTracks->Close();
+
+ return kTRUE;
+}
+
+
void AliITSTrackingV1(Int_t evNumber1=0,Int_t evNumber2=0, Int_t min_t=-1, Int_t max_t=0,Bool_t flagvert=1, Bool_t realmass=0, const char *filename="galice.root") {
-
///////////////// Dynamically link some shared libs ////////////////////////////////
if (gClassTable->GetID("AliRun") < 0) {
delete gAlice;
gAlice=0;
}
+
+ cout << "Sorting TPC tracks w.r. to transverse momentum...";
+ Bool_t success_sorting = TPCSortTracks();
+ if (success_sorting) {
+ cout << "DONE!" << endl;
+ }
+ else {
+ cout << "Some error occurred..." << endl;
+ return 1;
+ }
// Connect the Root Galice file containing Geometry, Kine and Hits
TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
./AliITSDeleteOldFiles.sh
# run the hit generation
-aliroot -q -b "$ALICE_ROOT/macros/grun.C($1)"
+aliroot -q -b "$ALICE_ROOT/macros/grun.C"
# digitize TPC
-aliroot -q -b "$ALICE_ROOT/TPC/AliTPCHits2Digits.C($1)"
+aliroot -q -b "$ALICE_ROOT/TPC/AliTPCHits2Digits.C"
-# prepare TPC tracks for matching with the ITS
-aliroot -q -b "$ALICE_ROOT/ITS/AliTPCTracking4ITS.C($1)"
+# find clusters in TPC
+aliroot -q -b "$ALICE_ROOT/TPC/AliTPCFindClusters.C"
+
+# find tracks in TPC
+aliroot -q -b "$ALICE_ROOT/TPC/AliTPCFindTracks.C"
+
+# do TPC tracking comparison
+aliroot -q -b "$ALICE_ROOT/TPC/AliTPCComparison.C"
# create summable digits for the ITS
aliroot -q -b "$ALICE_ROOT/ITS/AliITSHits2SDigits.C"
aliroot -q -b "$ALICE_ROOT/ITS/AliITSSDigits2Digits.C"
# create reconstructed points for the ITS
-aliroot -q -b "$ALICE_ROOT/ITS/AliITSDigits2RecPoints.C(0,$[$1-1])"
-
-# prepare for tracking
-aliroot -q -b "$ALICE_ROOT/ITS/AliITSTracksV1.C(0,$[$1-1])"
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSDigits2RecPoints.C"
# do the tracking
-aliroot -q -b "$ALICE_ROOT/ITS/AliITSTrackingV1.C(0,$[$1-1])"
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSTrackingV1.C"
+
+# prepare for comparison
+# aliroot -q -b "$ALICE_ROOT/ITS/AliITSTracksV1.C"
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSStoreFindableTracks.C"
-# do the comparison
-aliroot -q -b "$ALICE_ROOT/ITS/AliITSComparisonV1.C(0,$[$1-1])"
+# do ITS tracking comparison
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSComparisonV1.C"
#
# after all of the above you can run AliITSPlotTracksV1.C macro under
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff