]>
Commit | Line | Data |
---|---|---|
75ec0f41 | 1 | // |
2 | // Script to test the dN/dEta analysis using the dNdEtaAnalysis and | |
3 | // dNdEtaCorrection classes. Note that there is a cut on the events, | |
4 | // so the measurement will be biassed. | |
5 | // | |
6 | ||
7 | testAnalysis(Char_t* dataDir, Int_t nRuns=20) { | |
8 | ||
9 | Char_t str[256]; | |
10 | ||
11 | gSystem->Load("../esdTrackCuts/libESDtrackQuality.so"); | |
12 | gSystem->Load("libdNdEta.so"); | |
13 | ||
14 | // ######################################################## | |
15 | // selection of esd tracks | |
16 | ESDtrackQualityCuts* esdTrackCuts = new ESDtrackQualityCuts(); | |
17 | esdTrackCuts->DefineHistograms(1); | |
18 | ||
19 | esdTrackCuts->SetMinNClustersTPC(50); | |
20 | esdTrackCuts->SetMaxChi2PerClusterTPC(3.5); | |
21 | esdTrackCuts->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); | |
22 | esdTrackCuts->SetRequireTPCRefit(kTRUE); | |
23 | ||
24 | esdTrackCuts->SetMinNsigmaToVertex(3); | |
25 | esdTrackCuts->SetAcceptKingDaughters(kFALSE); | |
26 | ||
27 | AliLog::SetClassDebugLevel("ESDtrackQualityCuts",1); | |
28 | ||
29 | // ######################################################## | |
30 | // definition of dNdEta objects | |
31 | ||
32 | dNdEtaCorrection* dNdEtaMap = new dNdEtaCorrection(); | |
33 | ||
34 | dNdEtaMap->LoadHistograms("correction_map.root","dndeta_correction"); | |
35 | dNdEtaMap->RemoveEdges(2,0,2); | |
36 | ||
37 | dNdEtaAnalysis* analyse = new dNdEtaAnalysis("dndeta"); | |
38 | ||
39 | // ######################################################## | |
40 | // get the data dir | |
41 | Char_t execDir[256]; | |
42 | sprintf(execDir,gSystem->pwd()); | |
43 | TSystemDirectory* baseDir = new TSystemDirectory(".",dataDir); | |
44 | TList* dirList = baseDir->GetListOfFiles(); | |
45 | Int_t nDirs = dirList->GetEntries(); | |
46 | // go back to the dir where this script is executed | |
47 | gSystem->cd(execDir); | |
48 | ||
49 | ||
50 | // ######################################################## | |
51 | // definition of used pointers | |
52 | TFile* esdFile; | |
53 | TTree* esdTree; | |
54 | TBranch* esdBranch; | |
55 | ||
56 | AliESD* esd =0; | |
57 | ||
58 | // ######################################################## | |
59 | // loop over runs | |
60 | Int_t nRunCounter = 0; | |
61 | for (Int_t r=1; r<=nDirs; r++) { | |
62 | ||
63 | TSystemFile* presentDir = (TSystemFile*)dirList->At(r); | |
64 | if (!presentDir->IsDirectory()) | |
65 | continue; | |
66 | // check that the files are there | |
ceb5d1b5 | 67 | TString currentDataDir; |
68 | currentDataDir.Form("%s/%s",dataDir, presentDir->GetName()); | |
69 | if ((!gSystem->Which(currentDataDir.Data(),"galice.root")) || | |
70 | (!gSystem->Which(currentDataDir.Data(),"AliESDs.root"))) | |
75ec0f41 | 71 | continue; |
72 | ||
73 | if (nRunCounter++ >= nRuns) | |
74 | break; | |
75 | ||
76 | cout << "run #" << nRunCounter << endl; | |
77 | ||
78 | // ######################################################### | |
79 | // setup galice and runloader | |
80 | if (gAlice) { | |
81 | delete gAlice->GetRunLoader(); | |
82 | delete gAlice; | |
83 | gAlice=0; | |
84 | } | |
85 | ||
ceb5d1b5 | 86 | sprintf(str,"%s/galice.root",currentDataDir.Data()); |
75ec0f41 | 87 | AliRunLoader* runLoader = AliRunLoader::Open(str); |
88 | ||
89 | runLoader->LoadgAlice(); | |
90 | gAlice = runLoader->GetAliRun(); | |
91 | runLoader->LoadKinematics(); | |
92 | runLoader->LoadHeader(); | |
93 | ||
94 | // ######################################################### | |
95 | // open esd file and get the tree | |
96 | ||
ceb5d1b5 | 97 | sprintf(str,"%s/AliESDs.root",currentDataDir.Data()); |
75ec0f41 | 98 | // close it first to avoid memory leak |
99 | if (esdFile) | |
100 | if (esdFile->IsOpen()) | |
101 | esdFile->Close(); | |
102 | ||
103 | esdFile = TFile::Open(str); | |
104 | esdTree = (TTree*)esdFile->Get("esdTree"); | |
105 | if (!esdTree) | |
106 | continue; | |
107 | esdBranch = esdTree->GetBranch("ESD"); | |
108 | esdBranch->SetAddress(&esd); | |
109 | if (!esdBranch) | |
110 | continue; | |
111 | ||
112 | // ######################################################## | |
113 | // Magnetic field | |
114 | AliTracker::SetFieldMap(gAlice->Field(),kTRUE); // kTRUE means uniform magnetic field | |
115 | //AliKalmanTrack::SetConvConst(1000/0.299792458/5.); | |
116 | ||
117 | // ######################################################## | |
118 | // getting number of events | |
119 | Int_t nEvents = (Int_t)runLoader->GetNumberOfEvents(); | |
120 | Int_t nESDEvents = esdBranch->GetEntries(); | |
121 | ||
122 | if (nEvents!=nESDEvents) { | |
123 | cout << " Different number of events from runloader and esdtree!!!" << nEvents << " / " << nESDEvents << endl; | |
124 | return; | |
125 | } | |
126 | // ######################################################## | |
127 | // loop over number of events | |
128 | cout << " looping over events..." << endl; | |
129 | for(Int_t i=1; i<nEvents; i++) { | |
130 | ||
131 | esdBranch->GetEntry(i); | |
132 | runLoader->GetEvent(i); | |
133 | ||
134 | ||
135 | // ######################################################## | |
136 | // get the EDS vertex | |
137 | AliESDVertex* vtxESD = esd->GetVertex(); | |
138 | ||
139 | Double_t vtx[3]; | |
140 | Double_t vtx_res[3]; | |
141 | vtxESD->GetXYZ(vtx); | |
142 | ||
143 | vtx_res[0] = vtxESD->GetXRes(); | |
144 | vtx_res[1] = vtxESD->GetYRes(); | |
145 | vtx_res[2] = vtxESD->GetZRes(); | |
146 | ||
147 | // we need a good vertex | |
148 | // => there will be a bias on the measurement, since this cuts away some events | |
149 | if (strcmp(vtxESD->GetName(),"default")==0) | |
150 | continue; | |
151 | if (vtx_res[2]==0 || vtx_res[2]>0.1) | |
152 | continue; | |
153 | ||
154 | // ######################################################## | |
155 | // loop over esd tracks | |
156 | Int_t nTracks = esd->GetNumberOfTracks(); | |
157 | for (Int_t t=0; t<nTracks; t++) { | |
158 | AliESDtrack* esdTrack = esd->GetTrack(t); | |
159 | ||
160 | if (!esdTrackCuts->AcceptTrack(esdTrack)) | |
161 | continue; | |
162 | ||
163 | AliTPCtrack* tpcTrack = new AliTPCtrack(*esdTrack); | |
164 | ||
165 | if (tpcTrack->GetAlpha()==0) { | |
166 | cout << " Warning esd track alpha = 0" << endl; | |
167 | continue; | |
168 | } | |
169 | ||
170 | Float_t theta = tpcTrack->Theta(); | |
171 | Float_t eta = -TMath::Log(TMath::Tan(theta/2.)); | |
172 | Float_t correction = dNdEtaMap->GetCorrection(vtx[2], eta); | |
173 | ||
174 | dNdEtaMap->FillMeas(vtx[2], eta); | |
175 | ||
176 | analyse ->FillTrack(vtx[2], eta, correction); | |
177 | ||
178 | } // end of track loop | |
179 | analyse->FillEvent(vtx[2]); | |
180 | ||
181 | } // end of event loop | |
182 | } // end of run loop | |
183 | ||
184 | analyse->Finish(); | |
185 | ||
186 | TFile* fout = new TFile("out.root","RECREATE"); | |
187 | ||
188 | esdTrackCuts->SaveHistograms("esd_tracks_cuts"); | |
189 | dNdEtaMap->SaveHistograms(); | |
190 | analyse ->SaveHistograms(); | |
191 | ||
192 | fout->Write(); | |
193 | fout->Close(); | |
194 | ||
195 | } |