]>
Commit | Line | Data |
---|---|---|
6eb67451 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
c68184b6 | 14 | **************************************************************************/ |
15 | ||
88cb7938 | 16 | /* $Id$ */ |
c68184b6 | 17 | |
6eb69b9a | 18 | /************************************************************************** |
6eb67451 | 19 | * * |
20 | * This class builds AliTPCtrack objects from generated tracks to feed * | |
21 | * ITS tracking (V2). The AliTPCtrack is built from its first hit in * | |
22 | * the TPC. The track is assigned a Kalman-like covariance matrix * | |
23 | * depending on its pT and pseudorapidity and track parameters are * | |
70521312 | 24 | * smeared according to this covariance matrix. * |
6eb69b9a | 25 | * Output file contains sorted tracks, ready for matching with ITS. * |
6eb67451 | 26 | * * |
e130146c | 27 | * For details: * |
b2bca9d4 | 28 | * Alice Internal Note (submitted - get it from andrea.dainese@pd.infn.it)* |
e130146c | 29 | * http://www.pd.infn.it/alipd/talks/soft/adIII02/TPCtrackingParam.htm * |
30 | * * | |
31 | * Test macro is: AliBarrelRec_TPCparam.C * | |
6eb67451 | 32 | * * |
b2bca9d4 | 33 | * 2002/10/01: Introduction of the database for pp collisions (B=0.4 T) * |
34 | * - Everything done separately for pions, kaons, protons, electrons and * | |
35 | * muons. * | |
36 | * - Now (only for pp) the tracks are built from the AliTrackReferences * | |
37 | * which contain the position and momentum of all tracks at R = 83 cm; * | |
38 | * This eliminates the loss of tracks due the dead zone of the TPC * | |
39 | * where the 1st hit is not created. * | |
40 | * - In AliBarrelRec_TPCparam.C there many possible ways of determining * | |
41 | * the z coordinate of the primary vertex in pp events (from pixels, * | |
42 | * from ITS tracks, smearing according to resolution given by tracks. * | |
43 | * * | |
6eb69b9a | 44 | * 2002/04/28: Major upgrade of the class * |
45 | * - Covariance matrices and pulls are now separeted for pions, kaons * | |
46 | * and electrons. * | |
47 | * - A parameterization of dE/dx in the TPC has been included and it is * | |
48 | * used to assign a mass to each track according to a rough dE/dx PID. * | |
49 | * - All the "numbers" have been moved to the file with the DataBase, * | |
50 | * they are read as objects of the class AliTPCkineGrid, and assigned * | |
51 | * to data memebers of the class AliTPCtrackerParam. * | |
52 | * - All the code necessary to create a BataBase has been included in * | |
53 | * class (see the macro AliTPCtrackingParamDB.C for the details). * | |
54 | * * | |
55 | * Origin: Andrea Dainese, Padova - e-mail: andrea.dainese@pd.infn.it * | |
6eb67451 | 56 | * * |
57 | **************************************************************************/ | |
b2bca9d4 | 58 | //------- Root headers -------- |
6eb69b9a | 59 | #include <TChain.h> |
60 | #include <TF1.h> | |
61 | #include <TGraph.h> | |
62 | #include <TGraphErrors.h> | |
63 | #include <TLegend.h> | |
64 | #include <TLine.h> | |
65 | #include <TMatrixD.h> | |
66 | #include <TNtuple.h> | |
67 | #include <TSystem.h> | |
b91c8ab5 | 68 | #include <TH2F.h> |
b2bca9d4 | 69 | //------ AliRoot headers ------ |
6eb67451 | 70 | #include "alles.h" |
6eb69b9a | 71 | #include "AliGausCorr.h" |
6eb67451 | 72 | #include "AliKalmanTrack.h" |
6eb69b9a | 73 | #include "AliMagF.h" |
70521312 | 74 | #include "AliMagFCM.h" |
88cb7938 | 75 | #include "AliRunLoader.h" |
76 | #include "AliTPCLoader.h" | |
6eb69b9a | 77 | #include "AliTPCkineGrid.h" |
78 | #include "AliTPCtrack.h" | |
b2bca9d4 | 79 | #include "AliTrackReference.h" |
6eb69b9a | 80 | #include "AliTPCtrackerParam.h" |
5d12ce38 | 81 | #include "AliMC.h" |
b2bca9d4 | 82 | //----------------------------- |
6eb69b9a | 83 | |
84 | Double_t RegFunc(Double_t *x,Double_t *par) { | |
85 | // This is the function used to regularize the covariance matrix | |
b2bca9d4 | 86 | Double_t value = par[0]+par[1]/TMath::Power(x[0],par[2]); |
6eb69b9a | 87 | return value; |
88 | } | |
70521312 | 89 | |
6eb69b9a | 90 | // structure for DB building |
91 | typedef struct { | |
92 | Int_t pdg; | |
93 | Int_t bin; | |
94 | Double_t r; | |
95 | Double_t p; | |
96 | Double_t pt; | |
97 | Double_t cosl; | |
98 | Double_t eta; | |
99 | Double_t dpt; | |
100 | Double_t dP0,dP1,dP2,dP3,dP4; | |
101 | Double_t c00,c10,c11,c20,c21,c22,c30,c31,c32,c33,c40,c41,c42,c43,c44; | |
102 | Double_t dEdx; | |
103 | } COMPTRACK; | |
104 | // cov matrix structure | |
105 | typedef struct { | |
106 | Double_t c00,c10,c11,c20,c21,c22,c30,c31,c32,c33,c40,c41,c42,c43,c44; | |
107 | } COVMATRIX; | |
70521312 | 108 | |
6eb67451 | 109 | ClassImp(AliTPCtrackerParam) |
110 | ||
6eb69b9a | 111 | //----------------------------------------------------------------------------- |
b2bca9d4 | 112 | AliTPCtrackerParam::AliTPCtrackerParam(const Int_t coll,const Double_t Bz, |
88cb7938 | 113 | const Int_t n, const char* evfoldname): |
114 | fEvFolderName(evfoldname) { | |
6eb69b9a | 115 | //----------------------------------------------------------------------------- |
e130146c | 116 | // This is the class conctructor |
6eb69b9a | 117 | //----------------------------------------------------------------------------- |
e130146c | 118 | |
b2bca9d4 | 119 | fNevents = n; // events to be processed |
6eb69b9a | 120 | fBz = Bz; // value of the z component of L3 field (Tesla) |
b2bca9d4 | 121 | fColl = coll; // collision code (0: PbPb6000; 1: pp) |
6eb69b9a | 122 | fSelAndSmear = kTRUE; // by default selection and smearing are done |
123 | ||
124 | if(fBz!=0.4) { | |
125 | cerr<<"AliTPCtrackerParam::AliTPCtrackerParam: Invalid field!\n"; | |
126 | cerr<<" Available: 0.4"<<endl; | |
127 | } | |
b2bca9d4 | 128 | if(fColl!=0 && fColl!=1) { |
6eb69b9a | 129 | cerr<<"AliTPCtrackerParam::AliTPCtrackerParam: Invalid collision!\n"; |
b2bca9d4 | 130 | cerr<<" Available: 0 -> PbPb6000"<<endl; |
131 | cerr<<" 1 -> pp"<<endl; | |
6eb69b9a | 132 | } |
133 | ||
134 | fDBfileName = gSystem->Getenv("ALICE_ROOT"); | |
135 | fDBfileName.Append("/TPC/CovMatrixDB_"); | |
b2bca9d4 | 136 | //fDBfileName = "CovMatrixDB_"; |
6eb69b9a | 137 | if(fColl==0) fDBfileName.Append("PbPb6000"); |
b2bca9d4 | 138 | if(fColl==1) fDBfileName.Append("pp"); |
6eb69b9a | 139 | if(fBz==0.4) fDBfileName.Append("_B0.4T.root"); |
6eb67451 | 140 | } |
6eb69b9a | 141 | //----------------------------------------------------------------------------- |
b2bca9d4 | 142 | AliTPCtrackerParam::~AliTPCtrackerParam() {} |
6eb69b9a | 143 | //----------------------------------------------------------------------------- |
b2bca9d4 | 144 | Int_t AliTPCtrackerParam::BuildTPCtracks(const TFile *inp, TFile *out) { |
6eb69b9a | 145 | //----------------------------------------------------------------------------- |
e130146c | 146 | // This function creates the TPC parameterized tracks |
6eb69b9a | 147 | //----------------------------------------------------------------------------- |
148 | ||
88cb7938 | 149 | Error("BuildTPCtracks","in and out parameters ignored. new io"); |
150 | ||
151 | /********************************************/ | |
152 | AliRunLoader* rl = AliRunLoader::GetRunLoader(fEvFolderName); | |
153 | if (rl == 0x0) | |
154 | { | |
155 | Error("BuildTPCtracks","Can not get Run Loader from event folder named %s.", | |
156 | fEvFolderName.Data()); | |
157 | return 2; | |
158 | } | |
159 | AliLoader* tpcloader = rl->GetLoader("TPCLoader"); | |
160 | if (tpcloader == 0x0) | |
161 | { | |
162 | Error("BuildTPCtracks","Can not get TPC Loader from Run Loader."); | |
163 | return 3; | |
164 | } | |
165 | ||
166 | /********************************************/ | |
167 | ||
6eb69b9a | 168 | TFile *fileDB=0; |
169 | TTree *covTreePi[50]; | |
170 | TTree *covTreeKa[50]; | |
b2bca9d4 | 171 | TTree *covTreePr[50]; |
6eb69b9a | 172 | TTree *covTreeEl[50]; |
b2bca9d4 | 173 | TTree *covTreeMu[50]; |
6eb69b9a | 174 | |
175 | if(fSelAndSmear) { | |
176 | cerr<<"+++\n+++ Reading DataBase from:\n+++ "<< | |
177 | fDBfileName.Data()<<"\n+++\n"; | |
178 | // Read paramters from DB file | |
179 | if(!ReadAllData(fDBfileName.Data())) { | |
c15d4c37 | 180 | cerr<<"AliTPCtrackerParam::BuildTPCtracks: \ |
6eb69b9a | 181 | Could not read data from DB\n\n"; return 1; |
182 | } | |
b2bca9d4 | 183 | |
6eb69b9a | 184 | // Read the trees with regularized cov. matrices from DB |
185 | TString str; | |
186 | fileDB = TFile::Open(fDBfileName.Data()); | |
187 | Int_t nBinsPi = fDBgridPi.GetTotBins(); | |
188 | for(Int_t l=0; l<nBinsPi; l++) { | |
189 | str = "/CovMatrices/Pions/CovTreePi_bin"; | |
190 | str+=l; | |
191 | covTreePi[l] = (TTree*)fileDB->Get(str.Data()); | |
192 | } | |
193 | Int_t nBinsKa = fDBgridKa.GetTotBins(); | |
194 | for(Int_t l=0; l<nBinsKa; l++) { | |
195 | str = "/CovMatrices/Kaons/CovTreeKa_bin"; | |
196 | str+=l; | |
197 | covTreeKa[l] = (TTree*)fileDB->Get(str.Data()); | |
198 | } | |
b2bca9d4 | 199 | Int_t nBinsPr = fDBgridPr.GetTotBins(); |
200 | for(Int_t l=0; l<nBinsPr; l++) { | |
201 | if(fColl==0) { | |
202 | str = "/CovMatrices/Pions/CovTreePi_bin"; | |
203 | } else { | |
204 | str = "/CovMatrices/Protons/CovTreePr_bin"; | |
205 | } | |
206 | str+=l; | |
207 | covTreePr[l] = (TTree*)fileDB->Get(str.Data()); | |
208 | } | |
6eb69b9a | 209 | Int_t nBinsEl = fDBgridEl.GetTotBins(); |
210 | for(Int_t l=0; l<nBinsEl; l++) { | |
211 | str = "/CovMatrices/Electrons/CovTreeEl_bin"; | |
212 | str+=l; | |
213 | covTreeEl[l] = (TTree*)fileDB->Get(str.Data()); | |
214 | } | |
b2bca9d4 | 215 | Int_t nBinsMu = fDBgridMu.GetTotBins(); |
216 | for(Int_t l=0; l<nBinsMu; l++) { | |
217 | if(fColl==0) { | |
218 | str = "/CovMatrices/Pions/CovTreePi_bin"; | |
219 | } else { | |
220 | str = "/CovMatrices/Muons/CovTreeMu_bin"; | |
221 | } | |
222 | str+=l; | |
223 | covTreeMu[l] = (TTree*)fileDB->Get(str.Data()); | |
224 | } | |
e130146c | 225 | |
6eb69b9a | 226 | } else cerr<<"\n ! Creating ALL TRUE tracks at TPC 1st hit !\n\n"; |
6eb67451 | 227 | |
228 | TFile *infile=(TFile*)inp; | |
6eb69b9a | 229 | infile->cd(); |
6eb67451 | 230 | |
6eb67451 | 231 | // Get gAlice object from file |
88cb7938 | 232 | |
233 | rl->LoadgAlice(); | |
234 | rl->LoadHeader(); | |
235 | tpcloader->LoadHits("read"); | |
236 | ||
237 | if(!(gAlice=rl->GetAliRun())) { | |
238 | cerr<<"Can not get gAlice from Run Loader !\n"; | |
6eb67451 | 239 | return 1; |
240 | } | |
241 | ||
6eb69b9a | 242 | // Check if value in the galice file is equal to selected one (fBz) |
b2bca9d4 | 243 | AliMagF *fiel = (AliMagF*)gAlice->Field(); |
70521312 | 244 | Double_t fieval=(Double_t)fiel->SolenoidField()/10.; |
245 | printf("Magnetic field is %6.2f Tesla\n",fieval); | |
246 | if(fBz!=fieval) { | |
247 | cerr<<"AliTPCtrackerParam::BuildTPCtracks: Invalid field!"<<endl; | |
248 | cerr<<"Field selected is: "<<fBz<<" T\n"; | |
249 | cerr<<"Field found on file is: "<<fieval<<" T\n"; | |
250 | return 0; | |
251 | } | |
252 | ||
6eb69b9a | 253 | // Get TPC detector |
254 | AliTPC *TPC=(AliTPC*)gAlice->GetDetector("TPC"); | |
255 | Int_t ver = TPC->IsVersion(); | |
256 | cerr<<"+++ TPC version "<<ver<<" has been found !\n"; | |
88cb7938 | 257 | |
258 | rl->CdGAFile(); | |
6eb69b9a | 259 | AliTPCParam *digp=(AliTPCParam*)infile->Get("75x40_100x60"); |
260 | if(digp){ | |
261 | delete digp; | |
262 | digp = new AliTPCParamSR(); | |
263 | } | |
264 | else digp=(AliTPCParam*)infile->Get("75x40_100x60_150x60"); | |
265 | ||
266 | if(!digp) { cerr<<"TPC parameters have not been found !\n"; return 1; } | |
267 | TPC->SetParam(digp); | |
268 | ||
269 | // Set the conversion constant between curvature and Pt | |
70521312 | 270 | AliKalmanTrack::SetConvConst(100/0.299792458/fBz); |
6eb67451 | 271 | |
b2bca9d4 | 272 | TParticle *Part=0; |
273 | AliTPCseedGeant *seed=0; | |
274 | AliTPCtrack *tpctrack=0; | |
275 | Double_t sPt,sEta; | |
276 | Int_t cl=0,bin,label,pdg,charge; | |
277 | Int_t tracks; | |
278 | Int_t nParticles,nSeeds,arrentr; | |
6eb69b9a | 279 | Char_t tname[100]; |
280 | //Int_t nSel=0,nAcc=0; | |
6eb67451 | 281 | |
b2bca9d4 | 282 | |
6eb67451 | 283 | // loop over first n events in file |
b2bca9d4 | 284 | for(Int_t evt=0; evt<fNevents; evt++){ |
6eb67451 | 285 | cerr<<"+++\n+++ Processing event "<<evt<<"\n+++\n"; |
b2bca9d4 | 286 | tracks=0; |
6eb67451 | 287 | |
6eb67451 | 288 | // tree for TPC tracks |
6eb67451 | 289 | sprintf(tname,"TreeT_TPC_%d",evt); |
290 | TTree *tracktree = new TTree(tname,"Tree with TPC tracks"); | |
291 | tracktree->Branch("tracks","AliTPCtrack",&tpctrack,20000,0); | |
292 | ||
293 | // array for TPC tracks | |
b2bca9d4 | 294 | TObjArray tArray(20000); |
295 | ||
296 | // array for TPC seeds with geant info | |
297 | TObjArray sArray(20000); | |
6eb67451 | 298 | |
6eb69b9a | 299 | // get the particles stack |
b2bca9d4 | 300 | nParticles = (Int_t)gAlice->GetEvent(evt); |
301 | ||
302 | Bool_t *done = new Bool_t[nParticles]; | |
303 | Int_t *pdgCodes = new Int_t[nParticles]; | |
304 | Double_t *ptkine = new Double_t[nParticles]; | |
305 | Double_t *pzkine = new Double_t[nParticles]; | |
6eb67451 | 306 | |
6eb69b9a | 307 | // loop on particles and store pdg codes |
308 | for(Int_t l=0; l<nParticles; l++) { | |
5d12ce38 | 309 | Part = (TParticle*)gAlice->GetMCApp()->Particle(l); |
6eb69b9a | 310 | pdgCodes[l] = Part->GetPdgCode(); |
b2bca9d4 | 311 | ptkine[l] = Part->Pt(); |
312 | pzkine[l] = Part->Pz(); | |
6eb69b9a | 313 | done[l] = kFALSE; |
7a09f434 | 314 | } |
6eb69b9a | 315 | cerr<<"+++\n+++ Number of particles in event "<<evt<<": "<<nParticles<< |
b2bca9d4 | 316 | "\n+++\n"; |
317 | ||
318 | cerr<<"\n ********** MAKING SEEDS *****************"<<endl<<endl; | |
319 | ||
320 | ||
321 | // Create the seeds for the TPC tracks at the inner radius of TPC | |
322 | if(fColl==0) { | |
323 | // Get TreeH with hits | |
88cb7938 | 324 | TTree *TH = tpcloader->TreeH(); |
b2bca9d4 | 325 | MakeSeedsFromHits(TPC,TH,sArray); |
326 | } else { | |
327 | // Get TreeTR with track references | |
88cb7938 | 328 | TTree *TTR = rl->TreeTR(); |
b2bca9d4 | 329 | MakeSeedsFromRefs(TTR,sArray); |
330 | } | |
331 | ||
332 | ||
333 | nSeeds = sArray.GetEntries(); | |
334 | cerr<<"\n\n+++\n+++ Number of seeds: "<<nSeeds<<"\n+++\n"; | |
335 | ||
336 | ||
337 | cerr<<"\n ********** BUILDING TRACKS **************"<<endl<<endl; | |
338 | ||
339 | // loop over entries in sArray | |
340 | for(Int_t l=0; l<nSeeds; l++) { | |
341 | //if(l%1000==0) cerr<<" --- Processing seed " | |
342 | // <<l<<" of "<<nSeeds<<" ---\r"; | |
343 | ||
344 | seed = (AliTPCseedGeant*)sArray.At(l); | |
345 | ||
346 | // this is TEMPORARY: only for reconstruction of pp production for charm | |
347 | if(fColl==1) cl = CheckLabel(seed,nParticles,ptkine,pzkine); | |
348 | if(cl) continue; | |
349 | ||
350 | // Get track label | |
351 | label = seed->GetLabel(); | |
352 | ||
353 | // check if this track has already been processed | |
354 | if(done[label]) continue; | |
355 | // PDG code & electric charge | |
356 | pdg = pdgCodes[label]; | |
357 | if(pdg>200 || pdg==-11 || pdg==-13) { charge=1; } | |
358 | else if(pdg<-200 || pdg==11 || pdg==13) { charge=-1; } | |
359 | else continue; | |
360 | pdg = TMath::Abs(pdg); | |
361 | if(pdg>3000) pdg=211; | |
362 | ||
363 | if(fSelAndSmear) SetParticle(pdg); | |
364 | ||
365 | ||
366 | sPt = seed->GetPt(); | |
367 | sEta = seed->GetEta(); | |
368 | ||
369 | // Apply selection according to TPC efficiency | |
370 | //if(TMath::Abs(pdg)==211) nAcc++; | |
371 | if(fSelAndSmear && !SelectedTrack(sPt,sEta)) continue; | |
372 | //if(TMath::Abs(pdg)==211) nSel++; | |
373 | ||
374 | // create AliTPCtrack object | |
375 | BuildTrack(seed,charge); | |
6eb69b9a | 376 | |
b2bca9d4 | 377 | if(fSelAndSmear) { |
378 | bin = fDBgrid->GetBin(sPt,sEta); | |
379 | switch (pdg) { | |
380 | case 211: | |
381 | fCovTree = covTreePi[bin]; | |
382 | break; | |
383 | case 321: | |
384 | fCovTree = covTreeKa[bin]; | |
385 | break; | |
386 | case 2212: | |
387 | fCovTree = covTreePr[bin]; | |
388 | break; | |
389 | case 11: | |
390 | fCovTree = covTreeEl[bin]; | |
391 | break; | |
392 | case 13: | |
393 | fCovTree = covTreeMu[bin]; | |
394 | break; | |
395 | } | |
396 | // deal with covariance matrix and smearing of parameters | |
397 | CookTrack(sPt,sEta); | |
398 | ||
399 | // assign the track a dE/dx and make a rough PID | |
400 | CookdEdx(sPt,sEta); | |
6eb67451 | 401 | } |
b2bca9d4 | 402 | |
403 | // put track in array | |
404 | AliTPCtrack *iotrack = new AliTPCtrack(fTrack); | |
405 | iotrack->SetLabel(label); | |
406 | tArray.AddLast(iotrack); | |
407 | // Mark track as "done" and register the pdg code | |
408 | done[label] = kTRUE; | |
409 | tracks++; | |
6eb67451 | 410 | |
b2bca9d4 | 411 | } // loop over entries in sArray |
412 | ||
6eb67451 | 413 | |
6eb67451 | 414 | // sort array with TPC tracks (decreasing pT) |
b2bca9d4 | 415 | tArray.Sort(); |
6eb67451 | 416 | |
b2bca9d4 | 417 | arrentr = tArray.GetEntriesFast(); |
6eb67451 | 418 | for(Int_t l=0; l<arrentr; l++) { |
b2bca9d4 | 419 | tpctrack=(AliTPCtrack*)tArray.UncheckedAt(l); |
6eb67451 | 420 | tracktree->Fill(); |
421 | } | |
422 | ||
b2bca9d4 | 423 | |
6eb67451 | 424 | // write the tree with tracks in the output file |
425 | out->cd(); | |
b2bca9d4 | 426 | tracktree->Write(); |
427 | ||
6eb67451 | 428 | delete tracktree; |
c15d4c37 | 429 | delete [] done; |
430 | delete [] pdgCodes; | |
b2bca9d4 | 431 | delete [] ptkine; |
432 | delete [] pzkine; | |
6eb67451 | 433 | |
434 | printf("\n\n+++\n+++ Number of TPC tracks: %d\n+++\n",tracks); | |
6eb69b9a | 435 | //cerr<<"Average Eff: "<<(Float_t)nSel/nAcc<<endl; |
6eb67451 | 436 | |
b2bca9d4 | 437 | sArray.Delete(); |
438 | tArray.Delete(); | |
439 | ||
440 | infile->cd(); | |
6eb67451 | 441 | } // loop on events |
442 | ||
6eb69b9a | 443 | if(fileDB) fileDB->Close(); |
444 | ||
6eb67451 | 445 | return 0; |
446 | } | |
6eb69b9a | 447 | //----------------------------------------------------------------------------- |
448 | void AliTPCtrackerParam::AnalyzedEdx(const Char_t *outName,Int_t pdg) { | |
449 | //----------------------------------------------------------------------------- | |
450 | // This function computes the dE/dx for pions, kaons, protons and electrons, | |
451 | // in the [pT,eta] bins. | |
452 | // Input file is CovMatrix_AllEvts.root. | |
453 | //----------------------------------------------------------------------------- | |
454 | ||
455 | gStyle->SetOptStat(0); | |
456 | gStyle->SetOptFit(10001); | |
457 | ||
458 | Char_t *part="PIONS"; | |
459 | Double_t ymax=500.; | |
460 | ||
b2bca9d4 | 461 | /* |
6eb69b9a | 462 | // create a chain with compared tracks |
b2bca9d4 | 463 | TChain *cmptrkchain = new ("cmptrktree"); |
464 | cmptrkchain.Add("CovMatrix_AllEvts.root"); | |
465 | //cmptrkchain.Add("CovMatrix_AllEvts_1.root"); | |
466 | //cmptrkchain.Add("CovMatrix_AllEvts_2.root"); | |
467 | //cmptrkchain.Add("CovMatrix_AllEvts_3.root"); | |
468 | */ | |
469 | ||
470 | ||
471 | TFile *infile = TFile::Open("CovMatrix_AllEvts.root"); | |
472 | TTree *cmptrktree = (TTree*)infile->Get("cmptrktree"); | |
6eb69b9a | 473 | |
474 | COMPTRACK cmptrk; | |
b2bca9d4 | 475 | cmptrktree->SetBranchAddress("comptracks",&cmptrk); |
476 | Int_t entries = (Int_t)cmptrktree->GetEntries(); | |
6eb69b9a | 477 | cerr<<" Number of entries: "<<entries<<endl; |
478 | ||
b2bca9d4 | 479 | InitializeKineGrid("DB"); |
480 | InitializeKineGrid("dEdx"); | |
6eb69b9a | 481 | |
482 | switch(pdg) { | |
483 | case 211: | |
484 | part = "PIONS"; | |
485 | ymax = 100.; | |
486 | break; | |
487 | case 321: | |
488 | part = "KAONS"; | |
489 | ymax = 300.; | |
490 | break; | |
491 | case 2212: | |
492 | part = "PROTONS"; | |
493 | ymax = 500.; | |
494 | break; | |
495 | case 11: | |
496 | part = "ELECTRONS"; | |
497 | ymax = 100.; | |
498 | break; | |
b2bca9d4 | 499 | case 13: |
500 | part = "MUONS"; | |
501 | ymax = 100.; | |
502 | break; | |
6eb69b9a | 503 | } |
504 | ||
505 | SetParticle(pdg); | |
506 | ||
507 | const Int_t nTotBins = fDBgrid->GetTotBins(); | |
508 | ||
509 | cerr<<" Fit bins: "<<nTotBins<<endl; | |
510 | ||
511 | Int_t bin=0; | |
b2bca9d4 | 512 | Int_t *n = new Int_t[nTotBins]; |
513 | Double_t *p = new Double_t[nTotBins]; | |
514 | Double_t *ep = new Double_t[nTotBins]; | |
515 | Double_t *mean = new Double_t[nTotBins]; | |
516 | Double_t *sigma = new Double_t[nTotBins]; | |
6eb69b9a | 517 | |
518 | for(Int_t l=0; l<nTotBins; l++) { | |
519 | n[l] = 1; // set to 1 to avoid divisions by 0 | |
520 | p[l] = mean[l] = sigma[l] = ep[l] = 0.; | |
521 | } | |
522 | ||
523 | // loop on chain entries for the mean | |
524 | for(Int_t l=0; l<entries; l++) { | |
b2bca9d4 | 525 | cmptrktree->GetEvent(l); |
6eb69b9a | 526 | if(TMath::Abs(cmptrk.pdg)!=pdg) continue; |
527 | bin = fDBgrid->GetBin(cmptrk.pt,cmptrk.eta); | |
528 | p[bin] += cmptrk.p; | |
529 | mean[bin] += cmptrk.dEdx; | |
530 | n[bin]++; | |
531 | } // loop on chain entries | |
532 | ||
533 | for(Int_t l=0; l<nTotBins; l++) { | |
534 | p[l] /= n[l]; | |
535 | mean[l] /= n[l]; | |
536 | n[l] = 1; // set to 1 to avoid divisions by 0 | |
537 | } | |
538 | ||
539 | // loop on chain entries for the sigma | |
540 | for(Int_t l=0; l<entries; l++) { | |
b2bca9d4 | 541 | cmptrktree->GetEvent(l); |
6eb69b9a | 542 | if(TMath::Abs(cmptrk.pdg)!=pdg) continue; |
543 | bin = fDBgrid->GetBin(cmptrk.pt,cmptrk.eta); | |
544 | if(cmptrk.p<1. && TMath::Abs(cmptrk.p-p[bin])>0.025) continue; | |
545 | n[bin]++; | |
546 | sigma[bin] += (cmptrk.dEdx-mean[bin])*(cmptrk.dEdx-mean[bin]); | |
547 | } // loop on chain entries | |
548 | ||
549 | for(Int_t l=0; l<nTotBins; l++) { | |
550 | sigma[l] = TMath::Sqrt(sigma[l]/n[l]); | |
551 | } | |
552 | ||
553 | // create the canvas | |
554 | TCanvas *canv = new TCanvas("canv","dEdx",0,0,900,700); | |
555 | ||
556 | // create the graph for dEdx vs p | |
557 | TGraphErrors *gr = new TGraphErrors(nTotBins,p,mean,ep,sigma); | |
558 | TString title(" : dE/dx vs momentum"); title.Prepend(part); | |
559 | TH2F *frame = new TH2F("frame1",title.Data(),2,0.1,50,2,0,ymax); | |
560 | frame->SetXTitle("p [GeV/c]"); | |
561 | frame->SetYTitle("dE/dx [a.u.]"); | |
562 | canv->SetLogx(); | |
563 | frame->Draw(); | |
564 | gr->Draw("P"); | |
565 | ||
566 | switch(pdg) { | |
567 | case 211: | |
568 | for(Int_t i=0; i<nTotBins; i++) { | |
569 | fdEdxMeanPi.SetParam(i,mean[i]); | |
570 | fdEdxRMSPi.SetParam(i,sigma[i]); | |
571 | } | |
572 | break; | |
573 | case 321: | |
574 | for(Int_t i=0; i<nTotBins; i++) { | |
575 | fdEdxMeanKa.SetParam(i,mean[i]); | |
576 | fdEdxRMSKa.SetParam(i,sigma[i]); | |
577 | } | |
578 | break; | |
579 | case 2212: | |
580 | for(Int_t i=0; i<nTotBins; i++) { | |
581 | fdEdxMeanPr.SetParam(i,mean[i]); | |
582 | fdEdxRMSPr.SetParam(i,sigma[i]); | |
583 | } | |
584 | break; | |
585 | case 11: | |
586 | for(Int_t i=0; i<nTotBins; i++) { | |
587 | fdEdxMeanEl.SetParam(i,mean[i]); | |
588 | fdEdxRMSEl.SetParam(i,sigma[i]); | |
589 | } | |
590 | break; | |
b2bca9d4 | 591 | case 13: |
592 | for(Int_t i=0; i<nTotBins; i++) { | |
593 | fdEdxMeanMu.SetParam(i,mean[i]); | |
594 | fdEdxRMSMu.SetParam(i,sigma[i]); | |
595 | } | |
596 | break; | |
6eb69b9a | 597 | } |
598 | ||
599 | // write results to file | |
600 | WritedEdx(outName,pdg); | |
601 | ||
c15d4c37 | 602 | delete [] n; |
603 | delete [] p; | |
604 | delete [] ep; | |
605 | delete [] mean; | |
606 | delete [] sigma; | |
6eb69b9a | 607 | |
608 | return; | |
609 | } | |
610 | //----------------------------------------------------------------------------- | |
611 | void AliTPCtrackerParam::AnalyzePulls(const Char_t *outName) { | |
612 | //----------------------------------------------------------------------------- | |
613 | // This function computes the pulls for pions, kaons and electrons, | |
614 | // in the [pT,eta] bins. | |
615 | // Input file is CovMatrix_AllEvts.root. | |
616 | // Output file is pulls.root. | |
617 | //----------------------------------------------------------------------------- | |
618 | ||
b2bca9d4 | 619 | /* |
6eb69b9a | 620 | // create a chain with compared tracks |
b2bca9d4 | 621 | TChain *cmptrkchain = new ("cmptrktree"); |
622 | cmptrkchain.Add("CovMatrix_AllEvts.root"); | |
623 | //cmptrkchain.Add("CovMatrix_AllEvts_1.root"); | |
624 | //cmptrkchain.Add("CovMatrix_AllEvts_2.root"); | |
625 | //cmptrkchain.Add("CovMatrix_AllEvts_3.root"); | |
626 | */ | |
627 | ||
628 | ||
629 | TFile *infile = TFile::Open("CovMatrix_AllEvts.root"); | |
630 | TTree *cmptrktree = (TTree*)infile->Get("cmptrktree"); | |
6eb69b9a | 631 | |
632 | COMPTRACK cmptrk; | |
b2bca9d4 | 633 | cmptrktree->SetBranchAddress("comptracks",&cmptrk); |
634 | Int_t entries = (Int_t)cmptrktree->GetEntries(); | |
6eb69b9a | 635 | cerr<<" Number of entries: "<<entries<<endl; |
636 | ||
637 | Int_t thisPdg=0; | |
638 | Char_t hname[100], htitle[100]; | |
639 | Int_t nTotBins,bin; | |
640 | ||
641 | AliTPCkineGrid pulls[5]; | |
642 | TH1F *hDum = new TH1F("name","title",100,-7.,7.); | |
643 | TF1 *g = new TF1("g","gaus"); | |
644 | ||
b2bca9d4 | 645 | InitializeKineGrid("pulls"); |
646 | InitializeKineGrid("DB"); | |
6eb69b9a | 647 | |
648 | ||
649 | ||
b2bca9d4 | 650 | // loop on the particles Pi,Ka,Pr,El,Mu |
651 | for(Int_t part=0; part<5; part++) { | |
6eb69b9a | 652 | |
653 | switch (part) { | |
654 | case 0: // pions | |
655 | thisPdg=211; | |
656 | cerr<<" Processing pions ...\n"; | |
657 | break; | |
658 | case 1: // kaons | |
659 | thisPdg=321; | |
660 | cerr<<" Processing kaons ...\n"; | |
661 | break; | |
b2bca9d4 | 662 | case 2: // protons |
663 | thisPdg=2212; | |
664 | cerr<<" Processing protons ...\n"; | |
665 | break; | |
666 | case 3: // electrons | |
6eb69b9a | 667 | thisPdg=11; |
668 | cerr<<" Processing electrons ...\n"; | |
669 | break; | |
b2bca9d4 | 670 | case 4: // muons |
671 | thisPdg=13; | |
672 | cerr<<" Processing muons ...\n"; | |
673 | break; | |
6eb69b9a | 674 | } |
675 | ||
676 | SetParticle(thisPdg); | |
6eb67451 | 677 | |
6eb69b9a | 678 | for(Int_t i=0;i<5;i++) { |
679 | pulls[i].~AliTPCkineGrid(); | |
680 | new(&pulls[i]) AliTPCkineGrid(*(fPulls+i)); | |
681 | } | |
682 | nTotBins = fDBgrid->GetTotBins(); | |
b2bca9d4 | 683 | cerr<<"nTotBins = "<<nTotBins<<endl; |
684 | ||
6eb69b9a | 685 | // create histograms for the all the bins |
686 | TH1F *hPulls0_=NULL; | |
687 | TH1F *hPulls1_=NULL; | |
688 | TH1F *hPulls2_=NULL; | |
689 | TH1F *hPulls3_=NULL; | |
690 | TH1F *hPulls4_=NULL; | |
691 | ||
692 | hPulls0_ = new TH1F[nTotBins]; | |
693 | hPulls1_ = new TH1F[nTotBins]; | |
694 | hPulls2_ = new TH1F[nTotBins]; | |
695 | hPulls3_ = new TH1F[nTotBins]; | |
696 | hPulls4_ = new TH1F[nTotBins]; | |
697 | ||
698 | ||
699 | for(Int_t i=0; i<nTotBins; i++) { | |
700 | sprintf(hname,"hPulls0_%d",i); | |
701 | sprintf(htitle,"P0 pulls for bin %d",i); | |
702 | hDum->SetName(hname); hDum->SetTitle(htitle); | |
703 | hPulls0_[i] = *hDum; | |
704 | sprintf(hname,"hPulls1_%d",i); | |
705 | sprintf(htitle,"P1 pulls for bin %d",i); | |
706 | hDum->SetName(hname); hDum->SetTitle(htitle); | |
707 | hPulls1_[i] = *hDum; | |
708 | sprintf(hname,"hPulls2_%d",i); | |
709 | sprintf(htitle,"P2 pulls for bin %d",i); | |
710 | hDum->SetName(hname); hDum->SetTitle(htitle); | |
711 | hPulls2_[i] = *hDum; | |
712 | sprintf(hname,"hPulls3_%d",i); | |
713 | sprintf(htitle,"P3 pulls for bin %d",i); | |
714 | hDum->SetName(hname); hDum->SetTitle(htitle); | |
715 | hPulls3_[i] = *hDum; | |
716 | sprintf(hname,"hPulls4_%d",i); | |
717 | sprintf(htitle,"P4 pulls for bin %d",i); | |
718 | hDum->SetName(hname); hDum->SetTitle(htitle); | |
719 | hPulls4_[i] = *hDum; | |
720 | } | |
721 | ||
722 | // loop on chain entries | |
723 | for(Int_t i=0; i<entries; i++) { | |
b2bca9d4 | 724 | cmptrktree->GetEvent(i); |
6eb69b9a | 725 | if(TMath::Abs(cmptrk.pdg)!=thisPdg) continue; |
726 | // fill histograms with the pulls | |
b2bca9d4 | 727 | bin = fDBgrid->GetBin(cmptrk.pt,cmptrk.eta); |
728 | //cerr<<" pt "<<cmptrk.pt<<" eta "<<cmptrk.eta<<" bin "<<bin<<endl; | |
6eb69b9a | 729 | hPulls0_[bin].Fill(cmptrk.dP0/TMath::Sqrt(cmptrk.c00)); |
730 | hPulls1_[bin].Fill(cmptrk.dP1/TMath::Sqrt(cmptrk.c11)); | |
731 | hPulls2_[bin].Fill(cmptrk.dP2/TMath::Sqrt(cmptrk.c22)); | |
732 | hPulls3_[bin].Fill(cmptrk.dP3/TMath::Sqrt(cmptrk.c33)); | |
733 | hPulls4_[bin].Fill(cmptrk.dP4/TMath::Sqrt(cmptrk.c44)); | |
734 | } // loop on chain entries | |
735 | ||
736 | // compute the sigma of the distributions | |
737 | for(Int_t i=0; i<nTotBins; i++) { | |
b2bca9d4 | 738 | if(hPulls0_[i].GetEntries()>10) { |
6eb69b9a | 739 | g->SetRange(-3.*hPulls0_[i].GetRMS(),3.*hPulls0_[i].GetRMS()); |
740 | hPulls0_[i].Fit("g","R,Q,N"); | |
741 | pulls[0].SetParam(i,g->GetParameter(2)); | |
742 | } else pulls[0].SetParam(i,-1.); | |
b2bca9d4 | 743 | if(hPulls1_[i].GetEntries()>10) { |
6eb69b9a | 744 | g->SetRange(-3.*hPulls1_[i].GetRMS(),3.*hPulls1_[i].GetRMS()); |
745 | hPulls1_[i].Fit("g","R,Q,N"); | |
746 | pulls[1].SetParam(i,g->GetParameter(2)); | |
747 | } else pulls[1].SetParam(i,-1.); | |
b2bca9d4 | 748 | if(hPulls2_[i].GetEntries()>10) { |
6eb69b9a | 749 | g->SetRange(-3.*hPulls2_[i].GetRMS(),3.*hPulls2_[i].GetRMS()); |
750 | hPulls2_[i].Fit("g","R,Q,N"); | |
751 | pulls[2].SetParam(i,g->GetParameter(2)); | |
752 | } else pulls[2].SetParam(i,-1.); | |
b2bca9d4 | 753 | if(hPulls3_[i].GetEntries()>10) { |
6eb69b9a | 754 | g->SetRange(-3.*hPulls3_[i].GetRMS(),3.*hPulls3_[i].GetRMS()); |
755 | hPulls3_[i].Fit("g","R,Q,N"); | |
756 | pulls[3].SetParam(i,g->GetParameter(2)); | |
757 | } else pulls[3].SetParam(i,-1.); | |
b2bca9d4 | 758 | if(hPulls4_[i].GetEntries()>10) { |
6eb69b9a | 759 | g->SetRange(-3.*hPulls4_[i].GetRMS(),3.*hPulls4_[i].GetRMS()); |
760 | hPulls4_[i].Fit("g","R,Q,N"); | |
761 | pulls[4].SetParam(i,g->GetParameter(2)); | |
762 | } else pulls[4].SetParam(i,-1.); | |
763 | } // loop on bins | |
764 | ||
765 | ||
766 | switch (part) { | |
767 | case 0: // pions | |
768 | for(Int_t i=0;i<5;i++) { | |
769 | fPullsPi[i].~AliTPCkineGrid(); | |
770 | new(&fPullsPi[i]) AliTPCkineGrid(pulls[i]); | |
771 | } | |
772 | break; | |
773 | case 1: // kaons | |
774 | for(Int_t i=0;i<5;i++) { | |
775 | fPullsKa[i].~AliTPCkineGrid(); | |
776 | new(&fPullsKa[i]) AliTPCkineGrid(pulls[i]); | |
777 | } | |
778 | break; | |
b2bca9d4 | 779 | case 2: // protons |
780 | for(Int_t i=0;i<5;i++) { | |
781 | fPullsPr[i].~AliTPCkineGrid(); | |
782 | new(&fPullsPr[i]) AliTPCkineGrid(pulls[i]); | |
783 | } | |
784 | break; | |
785 | case 3: // electrons | |
6eb69b9a | 786 | for(Int_t i=0;i<5;i++) { |
787 | fPullsEl[i].~AliTPCkineGrid(); | |
788 | new(&fPullsEl[i]) AliTPCkineGrid(pulls[i]); | |
789 | } | |
790 | break; | |
b2bca9d4 | 791 | case 4: // muons |
792 | for(Int_t i=0;i<5;i++) { | |
793 | fPullsMu[i].~AliTPCkineGrid(); | |
794 | new(&fPullsMu[i]) AliTPCkineGrid(pulls[i]); | |
795 | //cerr<<" mu pulls "<<i<<" "<<fPullsMu[i].GetParam(0)<<endl; | |
796 | } | |
797 | break; | |
6eb69b9a | 798 | } |
799 | ||
800 | delete [] hPulls0_; | |
801 | delete [] hPulls1_; | |
802 | delete [] hPulls2_; | |
803 | delete [] hPulls3_; | |
804 | delete [] hPulls4_; | |
805 | ||
806 | } // loop on particle species | |
807 | ||
808 | // write pulls to file | |
809 | WritePulls(outName); | |
810 | ||
811 | ||
812 | return; | |
813 | } | |
814 | //----------------------------------------------------------------------------- | |
b2bca9d4 | 815 | void AliTPCtrackerParam::AnalyzeResolutions(Int_t pdg) { |
816 | //----------------------------------------------------------------------------- | |
817 | // This function computes the resolutions: | |
818 | // - dy | |
819 | // - dC | |
820 | // - dPt/Pt | |
821 | // as a function of Pt | |
822 | // Input file is CovMatrix_AllEvts.root. | |
823 | //----------------------------------------------------------------------------- | |
824 | ||
825 | /* | |
826 | // create a chain with compared tracks | |
827 | TChain *cmptrkchain = new ("cmptrktree"); | |
828 | cmptrkchain.Add("CovMatrix_AllEvts.root"); | |
829 | //cmptrkchain.Add("CovMatrix_AllEvts_1.root"); | |
830 | //cmptrkchain.Add("CovMatrix_AllEvts_2.root"); | |
831 | //cmptrkchain.Add("CovMatrix_AllEvts_3.root"); | |
832 | */ | |
833 | ||
834 | ||
835 | TFile *infile = TFile::Open("CovMatrix_AllEvts.root"); | |
836 | TTree *cmptrktree = (TTree*)infile->Get("cmptrktree"); | |
837 | ||
838 | COMPTRACK cmptrk; | |
839 | cmptrktree->SetBranchAddress("comptracks",&cmptrk); | |
840 | Int_t entries = (Int_t)cmptrktree->GetEntries(); | |
841 | cerr<<" Number of entries: "<<entries<<endl; | |
842 | ||
843 | ||
844 | Int_t bin = 0; | |
845 | ||
846 | InitializeKineGrid("DB"); | |
847 | InitializeKineGrid("eff"); | |
848 | ||
849 | SetParticle(pdg); | |
850 | ||
851 | const Int_t nPtBins = fEff->GetPointsPt(); | |
852 | cerr<<"nPtBins = "<<nPtBins<<endl; | |
853 | Double_t *dP0 = new Double_t[nPtBins]; | |
854 | Double_t *dP4 = new Double_t[nPtBins]; | |
855 | Double_t *dPtToPt = new Double_t[nPtBins]; | |
856 | Double_t *pt = new Double_t[nPtBins]; | |
857 | fEff->GetArrayPt(pt); | |
858 | ||
859 | ||
860 | TH1F *hDumP0 = new TH1F("nameP0","dy",100,-0.3,0.3); | |
861 | TH1F *hDumP4 = new TH1F("nameP4","dC",100,-0.0005,0.0005); | |
862 | TH1F *hDumPt = new TH1F("namePt","dp_{T}/p_{T}",100,-0.5,0.5); | |
863 | ||
864 | TF1 *g = new TF1("g","gaus"); | |
865 | ||
866 | // create histograms for the all the bins | |
867 | TH1F *hP0_=NULL; | |
868 | TH1F *hP4_=NULL; | |
869 | TH1F *hPt_=NULL; | |
870 | ||
871 | hP0_ = new TH1F[nPtBins]; | |
872 | hP4_ = new TH1F[nPtBins]; | |
873 | hPt_ = new TH1F[nPtBins]; | |
874 | ||
875 | for(Int_t i=0; i<nPtBins; i++) { | |
876 | hP0_[i] = *hDumP0; | |
877 | hP4_[i] = *hDumP4; | |
878 | hPt_[i] = *hDumPt; | |
879 | } | |
880 | ||
881 | // loop on chain entries | |
882 | for(Int_t i=0; i<entries; i++) { | |
883 | cmptrktree->GetEvent(i); | |
884 | if(TMath::Abs(cmptrk.pdg)!=pdg) continue; | |
885 | // fill histograms with the residuals | |
886 | bin = (Int_t)fDBgrid->GetBin(cmptrk.pt,cmptrk.eta)/fDBgrid->GetBinsEta(); | |
887 | //cerr<<" pt "<<cmptrk.pt<<" eta "<<cmptrk.eta<<" bin "<<bin<<endl; | |
888 | hP0_[bin].Fill(cmptrk.dP0); | |
889 | hP4_[bin].Fill(cmptrk.dP4); | |
890 | hPt_[bin].Fill(cmptrk.dpt/cmptrk.pt); | |
891 | } // loop on chain entries | |
892 | ||
893 | ||
894 | TCanvas *cP0res = new TCanvas("cP0res","cP0res",0,0,1200,700); | |
895 | cP0res->Divide(5,2); | |
896 | TCanvas *cP4res = new TCanvas("cP4res","cP4res",0,0,1200,700); | |
897 | cP4res->Divide(5,2); | |
898 | TCanvas *cPtres = new TCanvas("cPtres","cPtres",0,0,1200,700); | |
899 | cPtres->Divide(5,2); | |
900 | ||
901 | // Draw histograms | |
902 | for(Int_t i=0; i<nPtBins; i++) { | |
903 | cP0res->cd(i+1); hP0_[i].Draw(); | |
904 | cP4res->cd(i+1); hP4_[i].Draw(); | |
905 | cPtres->cd(i+1); hPt_[i].Draw(); | |
906 | } | |
907 | ||
908 | ||
909 | // compute the sigma of the distributions | |
910 | for(Int_t i=0; i<nPtBins; i++) { | |
911 | if(hP0_[i].GetEntries()>10) { | |
912 | g->SetRange(-3.*hP0_[i].GetRMS(),3.*hP0_[i].GetRMS()); | |
913 | hP0_[i].Fit("g","R,Q,N"); | |
914 | dP0[i] = g->GetParameter(2); | |
915 | } else dP0[i] = 0.; | |
916 | if(hP4_[i].GetEntries()>10) { | |
917 | g->SetRange(-3.*hP4_[i].GetRMS(),3.*hP4_[i].GetRMS()); | |
918 | hP4_[i].Fit("g","R,Q,N"); | |
919 | dP4[i] = g->GetParameter(2); | |
920 | } else dP4[i] = 0.; | |
921 | if(hPt_[i].GetEntries()>10) { | |
922 | g->SetRange(-3.*hPt_[i].GetRMS(),3.*hPt_[i].GetRMS()); | |
923 | hPt_[i].Fit("g","R,Q,N"); | |
924 | dPtToPt[i] = 100.*g->GetParameter(2); | |
925 | } else dPtToPt[i] = 0.; | |
926 | } // loop on bins | |
927 | ||
928 | ||
929 | TGraph *grdP0 = new TGraph(nPtBins,pt,dP0); | |
930 | TGraph *grdP4 = new TGraph(nPtBins,pt,dP4); | |
931 | TGraph *grdPtToPt = new TGraph(nPtBins,pt,dPtToPt); | |
932 | ||
933 | grdP0->SetMarkerStyle(20); grdP0->SetMarkerColor(2); grdP0->SetMarkerSize(1.5); | |
934 | grdP4->SetMarkerStyle(21); grdP4->SetMarkerColor(3); grdP4->SetMarkerSize(1.5); | |
935 | grdPtToPt->SetMarkerStyle(22); grdPtToPt->SetMarkerColor(4); grdPtToPt->SetMarkerSize(1.5); | |
936 | ||
937 | // Plot Results | |
938 | gStyle->SetOptStat(0); | |
939 | TCanvas *c1 = new TCanvas("c1","dP0",0,0,900,700); | |
940 | c1->SetLogx(); | |
941 | c1->SetGridx(); | |
942 | c1->SetGridy(); | |
943 | ||
944 | TH2F *frame1 = new TH2F("frame1","y resolution VS p_{T} in TPC",2,0.1,30,2,0,0.1); | |
945 | frame1->SetXTitle("p_{T} [GeV/c]"); | |
946 | frame1->SetYTitle("#sigma(y) [cm]"); | |
947 | frame1->Draw(); | |
948 | grdP0->Draw("P"); | |
949 | ||
950 | ||
951 | TCanvas *c2 = new TCanvas("c2","dP4",0,0,900,700); | |
952 | c2->SetLogx(); | |
953 | c2->SetGridx(); | |
954 | c2->SetGridy(); | |
955 | ||
956 | TH2F *frame2 = new TH2F("frame2","C resolution VS p_{T} in TPC",2,0.1,30,2,0,0.0001); | |
957 | frame2->SetXTitle("p_{T} [GeV/c]"); | |
958 | frame2->SetYTitle("#sigma(C) [1/cm]"); | |
959 | frame2->Draw(); | |
960 | grdP4->Draw("P"); | |
961 | ||
962 | TCanvas *c3 = new TCanvas("c3","dPtToPt",0,0,900,700); | |
963 | c3->SetLogx(); | |
964 | c3->SetLogy(); | |
965 | c3->SetGridx(); | |
966 | c3->SetGridy(); | |
967 | ||
968 | TH2F *frame3 = new TH2F("frame3","Relative p_{T} resolution VS p_{T} in TPC",2,0.1,30,2,0.1,30.); | |
969 | frame3->SetXTitle("p_{T} [GeV/c]"); | |
970 | frame3->SetYTitle("dp_{T}/p_{T} (%)"); | |
971 | frame3->Draw(); | |
972 | grdPtToPt->Draw("P"); | |
973 | ||
974 | ||
975 | delete [] dP0; | |
976 | delete [] dP4; | |
977 | delete [] dPtToPt; | |
978 | delete [] pt; | |
979 | ||
980 | ||
981 | delete [] hP0_; | |
982 | delete [] hP4_; | |
983 | delete [] hPt_; | |
984 | ||
985 | return; | |
986 | } | |
987 | //----------------------------------------------------------------------------- | |
988 | void AliTPCtrackerParam::BuildTrack(AliTPCseedGeant *s,Int_t ch) { | |
6eb69b9a | 989 | //----------------------------------------------------------------------------- |
e130146c | 990 | // This function uses GEANT info to set true track parameters |
6eb69b9a | 991 | //----------------------------------------------------------------------------- |
b2bca9d4 | 992 | Double_t xref = s->GetXL(); |
6eb67451 | 993 | Double_t xx[5],cc[15]; |
994 | cc[0]=cc[2]=cc[5]=cc[9]=cc[14]=10.; | |
995 | cc[1]=cc[3]=cc[4]=cc[6]=cc[7]=cc[8]=cc[10]=cc[11]=cc[12]=cc[13]=0.; | |
996 | ||
997 | // Magnetic field | |
e130146c | 998 | TVector3 bfield(0.,0.,fBz); |
6eb67451 | 999 | |
1000 | ||
1001 | // radius [cm] of track projection in (x,y) | |
b2bca9d4 | 1002 | Double_t rho = s->GetPt()*100./0.299792458/bfield.Z(); |
6eb67451 | 1003 | // center of track projection in local reference frame |
b2bca9d4 | 1004 | TVector3 sMom,sPos; |
6eb67451 | 1005 | |
1006 | ||
b2bca9d4 | 1007 | // position (local) and momentum (local) at the seed |
6eb67451 | 1008 | // in the bending plane (z=0) |
b2bca9d4 | 1009 | sPos.SetXYZ(s->GetXL(),s->GetYL(),0.); |
1010 | sMom.SetXYZ(s->GetPx()*TMath::Cos(s->GetAlpha())+s->GetPy()*TMath::Sin(s->GetAlpha()),-s->GetPx()*TMath::Sin(s->GetAlpha())+s->GetPy()*TMath::Cos(s->GetAlpha()),0.); | |
1011 | TVector3 vrho = sMom.Cross(bfield); | |
e130146c | 1012 | vrho *= ch; |
1013 | vrho.SetMag(rho); | |
6eb67451 | 1014 | |
b2bca9d4 | 1015 | TVector3 vcenter = sPos+vrho; |
6eb67451 | 1016 | |
e130146c | 1017 | Double_t x0 = vcenter.X(); |
6eb67451 | 1018 | |
1019 | // fX = xref X-coordinate of this track (reference plane) | |
1020 | // fAlpha = Alpha Rotation angle the local (TPC sector) | |
1021 | // fP0 = YL Y-coordinate of a track | |
1022 | // fP1 = ZG Z-coordinate of a track | |
1023 | // fP2 = C*x0 x0 is center x in rotated frame | |
1024 | // fP3 = Tgl tangent of the track momentum dip angle | |
1025 | // fP4 = C track curvature | |
b2bca9d4 | 1026 | xx[0] = s->GetYL(); |
1027 | xx[1] = s->GetZL(); | |
1028 | xx[3] = s->GetPz()/s->GetPt(); | |
6eb67451 | 1029 | xx[4] = -ch/rho; |
1030 | xx[2] = xx[4]*x0; | |
1031 | ||
1032 | // create the object AliTPCtrack | |
b2bca9d4 | 1033 | AliTPCtrack track(0,xx,cc,xref,s->GetAlpha()); |
6eb69b9a | 1034 | new(&fTrack) AliTPCtrack(track); |
6eb67451 | 1035 | |
6eb69b9a | 1036 | return; |
6eb67451 | 1037 | } |
6eb69b9a | 1038 | //----------------------------------------------------------------------------- |
b2bca9d4 | 1039 | Int_t AliTPCtrackerParam::CheckLabel(AliTPCseedGeant *s,Int_t nPart, |
1040 | Double_t *ptkine,Double_t *pzkine) const { | |
1041 | //----------------------------------------------------------------------------- | |
1042 | // This function checks if the label of the seed has been correctly | |
1043 | // assigned (to do only for pp charm production with AliRoot v3-08-02) | |
1044 | //----------------------------------------------------------------------------- | |
1045 | ||
1046 | Int_t sLabel = s->GetLabel(); | |
1047 | Double_t sPt = s->GetPt(); | |
1048 | Double_t sPz = s->GetPz(); | |
1049 | ||
1050 | // check if the label is correct (comparing momentum) | |
1051 | if(sLabel<nPart && | |
1052 | TMath::Abs(sPt-ptkine[sLabel])* | |
1053 | TMath::Abs(sPz-pzkine[sLabel])<0.001) return 0; | |
1054 | ||
1055 | if((sLabel-30)>=nPart) return 1; | |
1056 | ||
1057 | Double_t diff=0,mindiff=1000.; | |
1058 | Int_t bestLabel=0; | |
1059 | ||
1060 | for(Int_t i=sLabel-30; i<sLabel; i++) { | |
1061 | if(i<0 || i>=nPart) continue; | |
1062 | diff = TMath::Abs(sPt-ptkine[i])*TMath::Abs(sPz-pzkine[i]); | |
1063 | if(diff<mindiff) { mindiff = diff; bestLabel = i; } | |
1064 | } | |
1065 | ||
1066 | if(mindiff>0.001) return 1; | |
1067 | s->SetLabel(bestLabel); | |
1068 | ||
1069 | return 0; | |
1070 | } | |
1071 | //----------------------------------------------------------------------------- | |
6eb69b9a | 1072 | void AliTPCtrackerParam::CompareTPCtracks( |
c15d4c37 | 1073 | const Char_t* galiceName, |
1074 | const Char_t* trkGeaName, | |
1075 | const Char_t* trkKalName, | |
1076 | const Char_t* covmatName, | |
b2bca9d4 | 1077 | const Char_t* tpceffasciiName, |
1078 | const Char_t* tpceffrootName) { | |
6eb69b9a | 1079 | //----------------------------------------------------------------------------- |
1080 | // This function compares tracks from TPC Kalman Filter V2 with | |
b2bca9d4 | 1081 | // geant tracks at TPC 1st hit. It gives: |
6eb69b9a | 1082 | // - a tree with Kalman cov. matrix elements, resolutions, dEdx |
1083 | // - the efficiencies as a function of particle type, pT, eta | |
1084 | //----------------------------------------------------------------------------- | |
1085 | ||
1086 | TFile *kalFile = TFile::Open(trkKalName); | |
1087 | TFile *geaFile = TFile::Open(trkGeaName); | |
1088 | TFile *galiceFile = TFile::Open(galiceName); | |
1089 | ||
b2bca9d4 | 1090 | // get the AliRun object |
6eb69b9a | 1091 | AliRun *gAlice = (AliRun*)galiceFile->Get("gAlice"); |
6eb69b9a | 1092 | |
6eb67451 | 1093 | |
b2bca9d4 | 1094 | // create the tree for comparison results |
1095 | COMPTRACK cmptrk; | |
1096 | TTree *cmptrktree = new TTree("cmptrktree","results of track comparison"); | |
1097 | cmptrktree->Branch("comptracks",&cmptrk,"pdg/I:bin:r/D:p:pt:cosl:eta:dpt:dP0:dP1:dP2:dP3:dP4:c00:c10:c11:c20:c21:c22:c30:c31:c32:c33:c40:c41:c42:c43:c44:dEdx"); | |
1098 | ||
1099 | InitializeKineGrid("eff"); | |
1100 | InitializeKineGrid("DB"); | |
1101 | Int_t effBins = fEffPi.GetTotPoints(); | |
1102 | Int_t effBinsPt = fEffPi.GetPointsPt(); | |
1103 | Double_t *pt = new Double_t[effBinsPt]; | |
1104 | fEffPi.GetArrayPt(pt); | |
6eb69b9a | 1105 | |
1106 | TParticle *Part; | |
b2bca9d4 | 1107 | Double_t ptgener; |
1108 | Bool_t usethis; | |
6eb69b9a | 1109 | Int_t label; |
1110 | Double_t cc[15],dAlpha; | |
6eb69b9a | 1111 | Int_t pi=0,ka=0,mu=0,el=0,pr=0; |
b2bca9d4 | 1112 | Int_t *geaPi = new Int_t[effBins]; |
1113 | Int_t *geaKa = new Int_t[effBins]; | |
1114 | Int_t *geaPr = new Int_t[effBins]; | |
1115 | Int_t *geaEl = new Int_t[effBins]; | |
1116 | Int_t *geaMu = new Int_t[effBins]; | |
1117 | Int_t *kalPi = new Int_t[effBins]; | |
1118 | Int_t *kalKa = new Int_t[effBins]; | |
1119 | Int_t *kalPr = new Int_t[effBins]; | |
1120 | Int_t *kalEl = new Int_t[effBins]; | |
1121 | Int_t *kalMu = new Int_t[effBins]; | |
1122 | Float_t *effPi = new Float_t[effBins]; | |
1123 | Float_t *effKa = new Float_t[effBins]; | |
1124 | Float_t *effPr = new Float_t[effBins]; | |
1125 | Float_t *effEl = new Float_t[effBins]; | |
1126 | Float_t *effMu = new Float_t[effBins]; | |
6eb69b9a | 1127 | Int_t bin; |
b2bca9d4 | 1128 | for(Int_t j=0; j<effBins; j++) { |
6eb69b9a | 1129 | geaPi[j]=geaKa[j]=geaPr[j]=geaEl[j]=geaMu[j]=0; |
1130 | kalPi[j]=kalKa[j]=kalPr[j]=kalEl[j]=kalMu[j]=0; | |
1131 | effPi[j]=effKa[j]=effPr[j]=effEl[j]=effMu[j]=-1.; | |
1132 | } | |
e130146c | 1133 | |
b2bca9d4 | 1134 | Char_t tname[100]; |
1135 | ||
1136 | // loop on events in file | |
1137 | for(Int_t evt=0; evt<fNevents; evt++) { | |
1138 | cerr<<"\n --- Reading tracks for event "<<evt<<" ---\n\n"; | |
1139 | sprintf(tname,"TreeT_TPC_%d",evt); | |
1140 | ||
1141 | // particles from TreeK | |
1142 | const Int_t nparticles = gAlice->GetEvent(evt); | |
1143 | ||
1144 | Int_t *kalLab = new Int_t[nparticles]; | |
1145 | for(Int_t i=0; i<nparticles; i++) kalLab[i] = -1; | |
1146 | ||
1147 | ||
1148 | // tracks from Kalman | |
1149 | TTree *kaltree=(TTree*)kalFile->Get(tname); | |
1150 | if(!kaltree) continue; | |
1151 | AliTPCtrack *kaltrack=new AliTPCtrack; | |
1152 | kaltree->SetBranchAddress("tracks",&kaltrack); | |
1153 | Int_t kalEntries = (Int_t)kaltree->GetEntries(); | |
1154 | ||
1155 | // tracks from 1st hit | |
1156 | TTree *geatree=(TTree*)geaFile->Get(tname); | |
1157 | if(!geatree) continue; | |
1158 | AliTPCtrack *geatrack=new AliTPCtrack; | |
1159 | geatree->SetBranchAddress("tracks",&geatrack); | |
1160 | Int_t geaEntries = (Int_t)geatree->GetEntries(); | |
1161 | ||
1162 | cerr<<"+++\n+++ Number of tracks: TPC Kalman = "<<kalEntries<<endl<<"+++ TPC 1st hit = "<<geaEntries<<endl<<"+++\n"; | |
1163 | ||
1164 | // set pointers for TPC tracks: | |
1165 | // the entry number of the track labelled l is stored in kalLab[l] | |
1166 | Int_t fake=0,mult=0; | |
1167 | for (Int_t j=0; j<kalEntries; j++) { | |
1168 | kaltree->GetEvent(j); | |
1169 | if(kaltrack->GetLabel()>=0) { | |
1170 | if(kalLab[kaltrack->GetLabel()]!=-1) mult++; | |
1171 | kalLab[kaltrack->GetLabel()] = j; | |
1172 | } | |
1173 | else fake++; | |
1174 | } | |
1175 | ||
1176 | cerr<<"+++ Number of fake tracks in TPC Kalman: "<<fake<<endl; | |
1177 | cerr<<"+++ Number of multiply found tracks in TPC Kalman: "<<mult<<endl; | |
6eb69b9a | 1178 | |
b2bca9d4 | 1179 | /* |
1180 | // Read the labels of the seeds | |
1181 | char sname[100]; | |
1182 | Int_t sLabel,ncol; | |
1183 | Bool_t *hasSeed = new Bool_t[nparticles]; | |
1184 | for(Int_t i=0; i<nparticles; i++) hasSeed[i] = kFALSE; | |
1185 | sprintf(sname,"seedLabels.%d.dat",evt); | |
1186 | FILE *seedFile = fopen(sname,"r"); | |
1187 | while(1) { | |
1188 | ncol = fscanf(seedFile,"%d",&sLabel); | |
1189 | if(ncol<1) break; | |
1190 | if(sLabel>0) hasSeed[sLabel]=kTRUE; | |
1191 | } | |
1192 | fclose(seedFile); | |
1193 | */ | |
1194 | ||
1195 | cerr<<"Doing track comparison...\n"; | |
1196 | // loop on tracks at 1st hit | |
1197 | for(Int_t j=0; j<geaEntries; j++) { | |
1198 | geatree->GetEvent(j); | |
1199 | ||
1200 | label = geatrack->GetLabel(); | |
5d12ce38 | 1201 | Part = (TParticle*)gAlice->GetMCApp()->Particle(label); |
b2bca9d4 | 1202 | |
1203 | // use only injected tracks with fixed values of pT | |
1204 | ptgener = Part->Pt(); | |
1205 | usethis = kFALSE; | |
1206 | for(Int_t l=0; l<fEffPi.GetPointsPt(); l++) { | |
1207 | if(TMath::Abs(ptgener-pt[l])<0.01) usethis = kTRUE; | |
1208 | } | |
1209 | if(!usethis) continue; | |
1210 | ||
1211 | // check if it has the seed | |
1212 | //if(!hasSeed[label]) continue; | |
1213 | ||
1214 | /* | |
1215 | // check if track is entirely contained in a TPC sector | |
1216 | Bool_t out = kFALSE; | |
1217 | for(Int_t l=0; l<10; l++) { | |
1218 | Double_t x = 85. + (250.-85.)*(Double_t)l/9.; | |
1219 | //cerr<<"x "<<x<<" X "<<geatrack->GetX()<<endl; | |
1220 | Double_t y = geatrack->GetY() + ( | |
1221 | TMath::Sqrt(1-(geatrack->GetC()*geatrack->GetX()-geatrack->GetEta())* | |
1222 | (geatrack->GetC()*geatrack->GetX()-geatrack->GetEta())) | |
1223 | -TMath::Sqrt(1-(geatrack->GetC()*x-geatrack->GetEta())* | |
1224 | (geatrack->GetC()*x-geatrack->GetEta())) | |
1225 | )/geatrack->GetC(); | |
1226 | y = TMath::Abs(y); | |
1227 | //cerr<<"y "<<y<<" Y "<<geatrack->GetY()<<endl; | |
1228 | if(y > 0.8*x*TMath::Tan(TMath::Pi()/18)) { out = kTRUE; break; } | |
1229 | } | |
1230 | if(out) continue; | |
1231 | */ | |
1232 | ||
1233 | cmptrk.pdg = Part->GetPdgCode(); | |
1234 | cmptrk.eta = Part->Eta(); | |
1235 | cmptrk.r = TMath::Sqrt(Part->Vx()*Part->Vx()+Part->Vy()*Part->Vy()); | |
1236 | ||
1237 | cmptrk.pt = 1/TMath::Abs(geatrack->Get1Pt()); | |
1238 | cmptrk.cosl = TMath::Cos(TMath::ATan(geatrack->GetTgl())); | |
1239 | cmptrk.p = cmptrk.pt/cmptrk.cosl; | |
1240 | ||
1241 | ||
1242 | bin = fDBgridPi.GetBin(cmptrk.pt,cmptrk.eta); | |
1243 | cmptrk.bin = bin; | |
1244 | if(abs(cmptrk.pdg)==211) geaPi[bin]++; | |
1245 | if(abs(cmptrk.pdg)==321) geaKa[bin]++; | |
1246 | if(abs(cmptrk.pdg)==2212) geaPr[bin]++; | |
1247 | if(abs(cmptrk.pdg)==11) geaEl[bin]++; | |
1248 | if(abs(cmptrk.pdg)==13) geaMu[bin]++; | |
1249 | ||
1250 | ||
1251 | // check if there is the corresponding track in the TPC kalman and get it | |
1252 | if(kalLab[label]==-1) continue; | |
1253 | kaltree->GetEvent(kalLab[label]); | |
1254 | ||
1255 | // and go on only if it has xref = 84.57 cm (inner pad row) | |
1256 | if(kaltrack->GetX()>90.) continue; | |
1257 | ||
1258 | if(abs(cmptrk.pdg)==211) { kalPi[bin]++; pi++; } | |
1259 | if(abs(cmptrk.pdg)==321) { kalKa[bin]++; ka++; } | |
1260 | if(abs(cmptrk.pdg)==2212) { kalPr[bin]++; pr++; } | |
1261 | if(abs(cmptrk.pdg)==11) { kalEl[bin]++; el++; } | |
1262 | if(abs(cmptrk.pdg)==13) { kalMu[bin]++; mu++; } | |
1263 | ||
1264 | kaltrack->PropagateTo(geatrack->GetX(),1.e-9,0.); | |
1265 | ||
1266 | cmptrk.dEdx = kaltrack->GetdEdx(); | |
1267 | ||
1268 | // compute errors on parameters | |
1269 | dAlpha = kaltrack->GetAlpha()-geatrack->GetAlpha(); | |
1270 | if(TMath::Abs(dAlpha)>0.1) { cerr<<" ! WRONG SECTOR !\n"; continue; } | |
1271 | ||
1272 | cmptrk.dP0 = kaltrack->GetY()-geatrack->GetY(); | |
1273 | cmptrk.dP1 = kaltrack->GetZ()-geatrack->GetZ(); | |
1274 | cmptrk.dP2 = kaltrack->GetEta()-geatrack->GetEta(); | |
1275 | cmptrk.dP3 = kaltrack->GetTgl()-geatrack->GetTgl(); | |
1276 | cmptrk.dP4 = kaltrack->GetC()-geatrack->GetC(); | |
1277 | cmptrk.dpt = 1/kaltrack->Get1Pt()-1/geatrack->Get1Pt(); | |
1278 | ||
1279 | // get covariance matrix | |
1280 | // beware: lines 3 and 4 in the matrix are inverted! | |
1281 | kaltrack->GetCovariance(cc); | |
1282 | ||
1283 | cmptrk.c00 = cc[0]; | |
1284 | cmptrk.c10 = cc[1]; | |
1285 | cmptrk.c11 = cc[2]; | |
1286 | cmptrk.c20 = cc[3]; | |
1287 | cmptrk.c21 = cc[4]; | |
1288 | cmptrk.c22 = cc[5]; | |
1289 | cmptrk.c30 = cc[10]; | |
1290 | cmptrk.c31 = cc[11]; | |
1291 | cmptrk.c32 = cc[12]; | |
1292 | cmptrk.c33 = cc[14]; | |
1293 | cmptrk.c40 = cc[6]; | |
1294 | cmptrk.c41 = cc[7]; | |
1295 | cmptrk.c42 = cc[8]; | |
1296 | cmptrk.c43 = cc[13]; | |
1297 | cmptrk.c44 = cc[9]; | |
1298 | ||
1299 | // fill tree | |
1300 | cmptrktree->Fill(); | |
1301 | ||
1302 | } // loop on tracks at TPC 1st hit | |
1303 | ||
1304 | delete [] kalLab; | |
1305 | //delete [] hasSeed; | |
1306 | ||
1307 | } // end loop on events in file | |
1308 | ||
6eb69b9a | 1309 | |
1310 | cerr<<"+++\n+++ Number of compared tracks: "<<pi+ka+el+mu+pr<<endl; | |
1311 | cerr<<"+++ Pions: "<<pi<<", Kaons: "<<ka<<", Protons : "<<pr<<", Electrons: "<<el<<", Muons: "<<mu<<endl; | |
1312 | cerr<<"+++"<<endl; | |
1313 | ||
1314 | // Write tree to file | |
1315 | TFile *outfile = new TFile(covmatName,"recreate"); | |
1316 | cmptrktree->Write(); | |
1317 | outfile->Close(); | |
1318 | ||
1319 | ||
b2bca9d4 | 1320 | // Write efficiencies to ascii file |
1321 | FILE *effFile = fopen(tpceffasciiName,"w"); | |
6eb69b9a | 1322 | //fprintf(effFile,"%d\n",kalEntries); |
b2bca9d4 | 1323 | for(Int_t j=0; j<effBins; j++) { |
1324 | if(geaPi[j]>=100) effPi[j]=(Float_t)kalPi[j]/geaPi[j]; | |
1325 | if(geaKa[j]>=100) effKa[j]=(Float_t)kalKa[j]/geaKa[j]; | |
1326 | if(geaPr[j]>=100) effPr[j]=(Float_t)kalPr[j]/geaPr[j]; | |
1327 | if(geaEl[j]>=500) effEl[j]=(Float_t)kalEl[j]/geaEl[j]; | |
1328 | if(geaMu[j]>=100) effMu[j]=(Float_t)kalMu[j]/geaMu[j]; | |
6eb69b9a | 1329 | fprintf(effFile,"%f %f %f %f %f\n",effPi[j],effKa[j],effPr[j],effEl[j],effMu[j]); |
1330 | } | |
6eb67451 | 1331 | |
b2bca9d4 | 1332 | for(Int_t j=0; j<effBins; j++) { |
6eb69b9a | 1333 | fprintf(effFile,"%d %d %d %d %d\n",geaPi[j],geaKa[j],geaPr[j],geaEl[j],geaMu[j]); |
1334 | } | |
b2bca9d4 | 1335 | for(Int_t j=0; j<effBins; j++) { |
6eb69b9a | 1336 | fprintf(effFile,"%d %d %d %d %d\n",kalPi[j],kalKa[j],kalPr[j],kalEl[j],kalMu[j]); |
1337 | } | |
1338 | fclose(effFile); | |
6eb67451 | 1339 | |
b2bca9d4 | 1340 | // Write efficiencies to root file |
1341 | for(Int_t j=0; j<effBins; j++) { | |
1342 | fEffPi.SetParam(j,(Double_t)effPi[j]); | |
1343 | fEffKa.SetParam(j,(Double_t)effKa[j]); | |
1344 | fEffPr.SetParam(j,(Double_t)effPr[j]); | |
1345 | fEffEl.SetParam(j,(Double_t)effEl[j]); | |
1346 | fEffMu.SetParam(j,(Double_t)effMu[j]); | |
1347 | } | |
1348 | WriteEffs(tpceffrootName); | |
1349 | ||
6eb69b9a | 1350 | // delete AliRun object |
1351 | delete gAlice; gAlice=0; | |
1352 | ||
1353 | // close all input files | |
1354 | kalFile->Close(); | |
1355 | geaFile->Close(); | |
1356 | galiceFile->Close(); | |
6eb67451 | 1357 | |
b2bca9d4 | 1358 | delete [] pt; |
1359 | delete [] effPi; | |
1360 | delete [] effKa; | |
1361 | delete [] effPr; | |
1362 | delete [] effEl; | |
1363 | delete [] effMu; | |
1364 | delete [] geaPi; | |
1365 | delete [] geaKa; | |
1366 | delete [] geaPr; | |
1367 | delete [] geaEl; | |
1368 | delete [] geaMu; | |
1369 | delete [] kalPi; | |
1370 | delete [] kalKa; | |
1371 | delete [] kalPr; | |
1372 | delete [] kalEl; | |
1373 | delete [] kalMu; | |
c15d4c37 | 1374 | |
6eb69b9a | 1375 | return; |
6eb67451 | 1376 | } |
6eb69b9a | 1377 | //----------------------------------------------------------------------------- |
1378 | void AliTPCtrackerParam::CookdEdx(Double_t pt,Double_t eta) { | |
1379 | //----------------------------------------------------------------------------- | |
1380 | // This function assigns the track a dE/dx and makes a rough PID | |
1381 | //----------------------------------------------------------------------------- | |
6eb67451 | 1382 | |
6eb69b9a | 1383 | Double_t mean = fdEdxMean->GetValueAt(pt,eta); |
1384 | Double_t rms = fdEdxRMS->GetValueAt(pt,eta); | |
1385 | ||
1386 | Double_t dEdx = gRandom->Gaus(mean,rms); | |
1387 | ||
1388 | fTrack.SetdEdx(dEdx); | |
1389 | ||
1390 | AliTPCtrackParam t(fTrack); | |
1391 | ||
1392 | //Very rough PID | |
1393 | Double_t p = TMath::Sqrt(1.+t.GetTgl()*t.GetTgl())*pt; | |
1394 | ||
1395 | if (p<0.6) { | |
1396 | if (dEdx < 39.+ 12./(p+0.25)/(p+0.25)) { | |
1397 | t.AssignMass(0.13957); new(&fTrack) AliTPCtrack(t); return; | |
6eb67451 | 1398 | } |
6eb69b9a | 1399 | if (dEdx < 39.+ 12./p/p) { |
1400 | t.AssignMass(0.49368); new(&fTrack) AliTPCtrack(t); return; | |
1401 | } | |
1402 | t.AssignMass(0.93827); new(&fTrack) AliTPCtrack(t); return; | |
1403 | } | |
6eb67451 | 1404 | |
6eb69b9a | 1405 | if (p<1.2) { |
1406 | if (dEdx < 39.+ 12./(p+0.25)/(p+0.25)) { | |
1407 | t.AssignMass(0.13957); new(&fTrack) AliTPCtrack(t); return; | |
6eb67451 | 1408 | } |
6eb69b9a | 1409 | t.AssignMass(0.93827); new(&fTrack) AliTPCtrack(t); return; |
6eb67451 | 1410 | } |
1411 | ||
6eb69b9a | 1412 | t.AssignMass(0.13957); new(&fTrack) AliTPCtrack(t); return; |
6eb67451 | 1413 | } |
6eb69b9a | 1414 | //----------------------------------------------------------------------------- |
1415 | void AliTPCtrackerParam::CookTrack(Double_t pt,Double_t eta) { | |
1416 | //----------------------------------------------------------------------------- | |
1417 | // This function deals with covariance matrix and smearing | |
1418 | //----------------------------------------------------------------------------- | |
1419 | ||
1420 | COVMATRIX covmat; | |
1421 | Double_t p,cosl; | |
b2bca9d4 | 1422 | Double_t trkKine[1],trkRegPar[3]; |
6eb69b9a | 1423 | Double_t xref,alpha,xx[5],xxsm[5],cc[15]; |
1424 | Int_t treeEntries; | |
1425 | ||
1426 | fCovTree->SetBranchAddress("matrix",&covmat); | |
1427 | ||
1428 | // get random entry from the tree | |
1429 | treeEntries = (Int_t)fCovTree->GetEntries(); | |
1430 | fCovTree->GetEvent(gRandom->Integer(treeEntries)); | |
1431 | ||
1432 | // get P and Cosl from track | |
1433 | cosl = TMath::Cos(TMath::ATan(fTrack.GetTgl())); | |
1434 | p = 1./TMath::Abs(fTrack.Get1Pt())/cosl; | |
1435 | ||
1436 | trkKine[0] = p; | |
6eb69b9a | 1437 | |
1438 | // get covariance matrix from regularized matrix | |
b2bca9d4 | 1439 | for(Int_t l=0;l<3;l++) trkRegPar[l] = (*fRegPar)(0,l); |
6eb69b9a | 1440 | cc[0] = covmat.c00*RegFunc(trkKine,trkRegPar); |
1441 | cc[1] = covmat.c10; | |
b2bca9d4 | 1442 | for(Int_t l=0;l<3;l++) trkRegPar[l] = (*fRegPar)(1,l); |
6eb69b9a | 1443 | cc[2] = covmat.c11*RegFunc(trkKine,trkRegPar); |
b2bca9d4 | 1444 | for(Int_t l=0;l<3;l++) trkRegPar[l] = (*fRegPar)(2,l); |
6eb69b9a | 1445 | cc[3] = covmat.c20*RegFunc(trkKine,trkRegPar); |
1446 | cc[4] = covmat.c21; | |
b2bca9d4 | 1447 | for(Int_t l=0;l<3;l++) trkRegPar[l] = (*fRegPar)(3,l); |
6eb69b9a | 1448 | cc[5] = covmat.c22*RegFunc(trkKine,trkRegPar); |
1449 | cc[6] = covmat.c30; | |
b2bca9d4 | 1450 | for(Int_t l=0;l<3;l++) trkRegPar[l] = (*fRegPar)(4,l); |
6eb69b9a | 1451 | cc[7] = covmat.c31*RegFunc(trkKine,trkRegPar); |
1452 | cc[8] = covmat.c32; | |
b2bca9d4 | 1453 | for(Int_t l=0;l<3;l++) trkRegPar[l] = (*fRegPar)(5,l); |
6eb69b9a | 1454 | cc[9] = covmat.c33*RegFunc(trkKine,trkRegPar); |
b2bca9d4 | 1455 | for(Int_t l=0;l<3;l++) trkRegPar[l] = (*fRegPar)(6,l); |
6eb69b9a | 1456 | cc[10]= covmat.c40*RegFunc(trkKine,trkRegPar); |
1457 | cc[11]= covmat.c41; | |
b2bca9d4 | 1458 | for(Int_t l=0;l<3;l++) trkRegPar[l] = (*fRegPar)(7,l); |
6eb69b9a | 1459 | cc[12]= covmat.c42*RegFunc(trkKine,trkRegPar); |
1460 | cc[13]= covmat.c43; | |
b2bca9d4 | 1461 | for(Int_t l=0;l<3;l++) trkRegPar[l] = (*fRegPar)(8,l); |
6eb69b9a | 1462 | cc[14]= covmat.c44*RegFunc(trkKine,trkRegPar); |
1463 | ||
1464 | TMatrixD covMatSmear(5,5); | |
1465 | ||
1466 | covMatSmear = GetSmearingMatrix(cc,pt,eta); | |
1467 | ||
1468 | // get track original parameters | |
1469 | xref =fTrack.GetX(); | |
1470 | alpha=fTrack.GetAlpha(); | |
1471 | xx[0]=fTrack.GetY(); | |
1472 | xx[1]=fTrack.GetZ(); | |
1473 | xx[2]=fTrack.GetX()*fTrack.GetC()-fTrack.GetSnp(); | |
1474 | xx[3]=fTrack.GetTgl(); | |
1475 | xx[4]=fTrack.GetC(); | |
1476 | ||
1477 | // use smearing matrix to smear the original parameters | |
1478 | SmearTrack(xx,xxsm,covMatSmear); | |
1479 | ||
1480 | AliTPCtrack track(0,xxsm,cc,xref,alpha); | |
1481 | new(&fTrack) AliTPCtrack(track); | |
1482 | ||
1483 | return; | |
1484 | } | |
1485 | //----------------------------------------------------------------------------- | |
c15d4c37 | 1486 | void AliTPCtrackerParam::DrawEffs(const Char_t* inName,Int_t pdg) { |
6eb69b9a | 1487 | //----------------------------------------------------------------------------- |
1488 | // This function draws the TPC efficiencies in the [pT,eta] bins | |
1489 | //----------------------------------------------------------------------------- | |
1490 | ||
1491 | ReadEffs(inName); | |
1492 | SetParticle(pdg); | |
1493 | ||
1494 | const Int_t n = fEff->GetPointsPt(); | |
b2bca9d4 | 1495 | Double_t *effsA = new Double_t[n]; |
1496 | Double_t *effsB = new Double_t[n]; | |
1497 | Double_t *effsC = new Double_t[n]; | |
1498 | Double_t *pt = new Double_t[n]; | |
c15d4c37 | 1499 | |
6eb69b9a | 1500 | fEff->GetArrayPt(pt); |
1501 | for(Int_t i=0;i<n;i++) { | |
1502 | effsA[i] = fEff->GetParam(i,0); | |
1503 | effsB[i] = fEff->GetParam(i,1); | |
1504 | effsC[i] = fEff->GetParam(i,2); | |
1505 | } | |
1506 | ||
1507 | TGraph *grA = new TGraph(n,pt,effsA); | |
1508 | TGraph *grB = new TGraph(n,pt,effsB); | |
1509 | TGraph *grC = new TGraph(n,pt,effsC); | |
1510 | ||
1511 | grA->SetMarkerStyle(20); grA->SetMarkerColor(2); grA->SetMarkerSize(1.5); | |
1512 | grB->SetMarkerStyle(21); grB->SetMarkerColor(3); grB->SetMarkerSize(1.5); | |
1513 | grC->SetMarkerStyle(22); grC->SetMarkerColor(4); grC->SetMarkerSize(1.5); | |
1514 | ||
1515 | TString title("Distribution of the TPC efficiencies"); | |
1516 | switch (pdg) { | |
1517 | case 211: | |
1518 | title.Prepend("PIONS - "); | |
1519 | break; | |
1520 | case 321: | |
1521 | title.Prepend("KAONS - "); | |
1522 | break; | |
1523 | case 2212: | |
1524 | title.Prepend("PROTONS - "); | |
1525 | break; | |
1526 | case 11: | |
1527 | title.Prepend("ELECTRONS - "); | |
1528 | break; | |
1529 | case 13: | |
1530 | title.Prepend("MUONS - "); | |
1531 | break; | |
1532 | } | |
1533 | ||
1534 | ||
1535 | gStyle->SetOptStat(0); | |
1536 | TCanvas *c = new TCanvas("c","effs",0,0,900,700); | |
1537 | c->SetLogx(); | |
1538 | c->SetGridx(); | |
1539 | c->SetGridy(); | |
1540 | ||
1541 | TH2F *frame = new TH2F("frame",title.Data(),2,0.1,30,2,0,1); | |
1542 | frame->SetXTitle("p_{T} [GeV/c]"); | |
1543 | frame->Draw(); | |
1544 | grA->Draw("P"); | |
1545 | grB->Draw("P"); | |
1546 | grC->Draw("P"); | |
1547 | ||
1548 | TLegend *leg = new TLegend(0.2,0.2,0.4,0.4); | |
1549 | leg->AddEntry(grA,"|#eta|<0.3","p"); | |
1550 | leg->AddEntry(grB,"0.3<|#eta|<0.6","p"); | |
1551 | leg->AddEntry(grC,"0.6<|#eta|<0.9","p"); | |
1552 | leg->SetFillColor(0); | |
1553 | leg->Draw(); | |
1554 | ||
c15d4c37 | 1555 | delete [] effsA; |
1556 | delete [] effsB; | |
1557 | delete [] effsC; | |
1558 | delete [] pt; | |
6eb69b9a | 1559 | |
1560 | return; | |
1561 | } | |
1562 | //----------------------------------------------------------------------------- | |
c15d4c37 | 1563 | void AliTPCtrackerParam::DrawPulls(const Char_t* inName,Int_t pdg, |
1564 | Int_t par) { | |
6eb69b9a | 1565 | //----------------------------------------------------------------------------- |
1566 | // This function draws the pulls in the [pT,eta] bins | |
1567 | //----------------------------------------------------------------------------- | |
1568 | ||
1569 | ReadPulls(inName); | |
1570 | SetParticle(pdg); | |
1571 | ||
1572 | const Int_t n = (fPulls+par)->GetPointsPt(); | |
b2bca9d4 | 1573 | Double_t *pullsA = new Double_t[n]; |
1574 | Double_t *pullsB = new Double_t[n]; | |
1575 | Double_t *pullsC = new Double_t[n]; | |
1576 | Double_t *pt = new Double_t[n]; | |
6eb69b9a | 1577 | (fPulls+par)->GetArrayPt(pt); |
1578 | for(Int_t i=0;i<n;i++) { | |
1579 | pullsA[i] = (fPulls+par)->GetParam(i,0); | |
1580 | pullsB[i] = (fPulls+par)->GetParam(i,1); | |
1581 | pullsC[i] = (fPulls+par)->GetParam(i,2); | |
1582 | } | |
1583 | ||
1584 | TGraph *grA = new TGraph(n,pt,pullsA); | |
1585 | TGraph *grB = new TGraph(n,pt,pullsB); | |
1586 | TGraph *grC = new TGraph(n,pt,pullsC); | |
1587 | ||
1588 | grA->SetMarkerStyle(20); grA->SetMarkerColor(2); grA->SetMarkerSize(1.5); | |
1589 | grB->SetMarkerStyle(21); grB->SetMarkerColor(3); grB->SetMarkerSize(1.5); | |
1590 | grC->SetMarkerStyle(22); grC->SetMarkerColor(4); grC->SetMarkerSize(1.5); | |
1591 | ||
1592 | TString title("Distribution of the pulls: "); | |
1593 | switch (pdg) { | |
1594 | case 211: | |
1595 | title.Prepend("PIONS - "); | |
1596 | break; | |
1597 | case 321: | |
1598 | title.Prepend("KAONS - "); | |
1599 | break; | |
b2bca9d4 | 1600 | case 2212: |
1601 | title.Prepend("PROTONS - "); | |
1602 | break; | |
6eb69b9a | 1603 | case 11: |
1604 | title.Prepend("ELECTRONS - "); | |
1605 | break; | |
b2bca9d4 | 1606 | case 13: |
1607 | title.Prepend("MUONS - "); | |
1608 | break; | |
6eb69b9a | 1609 | } |
1610 | switch (par) { | |
1611 | case 0: | |
1612 | title.Append("y"); | |
1613 | break; | |
1614 | case 1: | |
1615 | title.Append("z"); | |
1616 | break; | |
1617 | case 2: | |
1618 | title.Append(" #eta"); | |
1619 | break; | |
1620 | case 3: | |
1621 | title.Append("tg #lambda"); | |
1622 | break; | |
1623 | case 4: | |
1624 | title.Append("C"); | |
1625 | break; | |
1626 | } | |
1627 | ||
1628 | gStyle->SetOptStat(0); | |
1629 | TCanvas *c = new TCanvas("c","pulls",0,0,900,700); | |
1630 | c->SetLogx(); | |
1631 | c->SetGridx(); | |
1632 | c->SetGridy(); | |
1633 | ||
1634 | TH2F *frame = new TH2F("frame",title.Data(),2,0.1,30,2,0,2); | |
1635 | frame->SetXTitle("p_{T} [GeV/c]"); | |
1636 | frame->Draw(); | |
1637 | grA->Draw("P"); | |
1638 | grB->Draw("P"); | |
1639 | grC->Draw("P"); | |
1640 | ||
1641 | TLegend *leg = new TLegend(0.2,0.2,0.4,0.4); | |
1642 | leg->AddEntry(grA,"|#eta|<0.3","p"); | |
1643 | leg->AddEntry(grB,"0.3<|#eta|<0.6","p"); | |
1644 | leg->AddEntry(grC,"0.6<|#eta|<0.9","p"); | |
1645 | leg->SetFillColor(0); | |
1646 | leg->Draw(); | |
1647 | ||
c15d4c37 | 1648 | delete [] pullsA; |
1649 | delete [] pullsB; | |
1650 | delete [] pullsC; | |
1651 | delete [] pt; | |
6eb67451 | 1652 | |
6eb69b9a | 1653 | return; |
1654 | } | |
1655 | //----------------------------------------------------------------------------- | |
6eb69b9a | 1656 | TMatrixD AliTPCtrackerParam::GetSmearingMatrix(Double_t* cc,Double_t pt, |
1657 | Double_t eta) const { | |
1658 | //----------------------------------------------------------------------------- | |
e130146c | 1659 | // This function stretches the covariance matrix according to the pulls |
6eb69b9a | 1660 | //----------------------------------------------------------------------------- |
1661 | ||
e130146c | 1662 | TMatrixD covMat(5,5); |
1663 | ||
1664 | covMat(0,0)=cc[0]; | |
1665 | covMat(1,0)=cc[1]; covMat(0,1)=covMat(1,0); | |
1666 | covMat(1,1)=cc[2]; | |
1667 | covMat(2,0)=cc[3]; covMat(0,2)=covMat(2,0); | |
1668 | covMat(2,1)=cc[4]; covMat(1,2)=covMat(2,1); | |
1669 | covMat(2,2)=cc[5]; | |
1670 | covMat(3,0)=cc[6]; covMat(0,3)=covMat(3,0); | |
1671 | covMat(3,1)=cc[7]; covMat(1,3)=covMat(3,1); | |
1672 | covMat(3,2)=cc[8]; covMat(2,3)=covMat(3,2); | |
1673 | covMat(3,3)=cc[9]; | |
1674 | covMat(4,0)=cc[10]; covMat(0,4)=covMat(4,0); | |
1675 | covMat(4,1)=cc[11]; covMat(1,4)=covMat(4,1); | |
1676 | covMat(4,2)=cc[12]; covMat(2,4)=covMat(4,2); | |
1677 | covMat(4,3)=cc[13]; covMat(3,4)=covMat(4,3); | |
1678 | covMat(4,4)=cc[14]; | |
1679 | ||
b2bca9d4 | 1680 | |
e130146c | 1681 | TMatrixD stretchMat(5,5); |
6eb67451 | 1682 | for(Int_t k=0;k<5;k++) { |
1683 | for(Int_t l=0;l<5;l++) { | |
e130146c | 1684 | stretchMat(k,l)=0.; |
6eb67451 | 1685 | } |
1686 | } | |
1687 | ||
b2bca9d4 | 1688 | for(Int_t i=0;i<5;i++) { |
1689 | stretchMat(i,i) = (fPulls+i)->GetValueAt(pt,eta); | |
1690 | if(stretchMat(i,i)==0.) stretchMat(i,i) = 1.; | |
1691 | } | |
6eb69b9a | 1692 | |
1693 | TMatrixD mat(stretchMat,TMatrixD::kMult,covMat); | |
1694 | TMatrixD covMatSmear(mat,TMatrixD::kMult,stretchMat); | |
1695 | ||
1696 | return covMatSmear; | |
1697 | } | |
1698 | //----------------------------------------------------------------------------- | |
b2bca9d4 | 1699 | void AliTPCtrackerParam::InitializeKineGrid(Option_t* which) { |
6eb69b9a | 1700 | //----------------------------------------------------------------------------- |
1701 | // This function initializes ([pt,eta] points) the data members AliTPCkineGrid | |
1702 | // which = "DB" -> initialize fDBgrid... members | |
1703 | // "eff" -> initialize fEff... members | |
1704 | // "pulls" -> initialize fPulls... members | |
1705 | // "dEdx" -> initialize fdEdx... members | |
6eb69b9a | 1706 | //----------------------------------------------------------------------------- |
1707 | ||
6eb69b9a | 1708 | const char *DB = strstr(which,"DB"); |
1709 | const char *eff = strstr(which,"eff"); | |
1710 | const char *pulls = strstr(which,"pulls"); | |
1711 | const char *dEdx = strstr(which,"dEdx"); | |
1712 | ||
1713 | ||
b2bca9d4 | 1714 | Int_t nEta=0, nPt=0; |
6eb69b9a | 1715 | |
1716 | Double_t etaPoints[2] = {0.3,0.6}; | |
1717 | Double_t etaBins[3] = {0.15,0.45,0.75}; | |
b2bca9d4 | 1718 | |
1719 | Double_t ptPoints[9] = {0.4,0.6,0.8,1.2,1.7,3.,5.,8.,15.}; | |
1720 | Double_t ptBins[10] = {0.3,0.5,0.7,1.,1.5,2.,4.,6.,10.,20.}; | |
1721 | ||
1722 | ||
1723 | Double_t *eta=0,*pt=0; | |
1724 | ||
1725 | if(DB) { | |
1726 | nEta = 2; | |
1727 | nPt = 9; | |
6eb69b9a | 1728 | eta = etaPoints; |
b2bca9d4 | 1729 | pt = ptPoints; |
1730 | } else { | |
1731 | nEta = 3; | |
1732 | nPt = 10; | |
6eb69b9a | 1733 | eta = etaBins; |
b2bca9d4 | 1734 | pt = ptBins; |
6eb69b9a | 1735 | } |
1736 | ||
1737 | AliTPCkineGrid *dummy=0; | |
1738 | ||
1739 | if(DB) { | |
b2bca9d4 | 1740 | dummy = new AliTPCkineGrid(nPt,nEta,pt,eta); |
6eb69b9a | 1741 | new(&fDBgridPi) AliTPCkineGrid(*dummy); |
6eb69b9a | 1742 | new(&fDBgridKa) AliTPCkineGrid(*dummy); |
b2bca9d4 | 1743 | new(&fDBgridPr) AliTPCkineGrid(*dummy); |
6eb69b9a | 1744 | new(&fDBgridEl) AliTPCkineGrid(*dummy); |
b2bca9d4 | 1745 | new(&fDBgridMu) AliTPCkineGrid(*dummy); |
6eb69b9a | 1746 | delete dummy; |
1747 | } | |
1748 | if(eff) { | |
b2bca9d4 | 1749 | dummy = new AliTPCkineGrid(nPt,nEta,pt,eta); |
6eb69b9a | 1750 | new(&fEffPi) AliTPCkineGrid(*dummy); |
6eb69b9a | 1751 | new(&fEffKa) AliTPCkineGrid(*dummy); |
6eb69b9a | 1752 | new(&fEffPr) AliTPCkineGrid(*dummy); |
6eb69b9a | 1753 | new(&fEffEl) AliTPCkineGrid(*dummy); |
6eb69b9a | 1754 | new(&fEffMu) AliTPCkineGrid(*dummy); |
1755 | delete dummy; | |
1756 | } | |
1757 | if(pulls) { | |
b2bca9d4 | 1758 | dummy = new AliTPCkineGrid(nPt,nEta,pt,eta); |
6eb69b9a | 1759 | for(Int_t i=0;i<5;i++) new(&fPullsPi[i]) AliTPCkineGrid(*dummy); |
6eb69b9a | 1760 | for(Int_t i=0;i<5;i++) new(&fPullsKa[i]) AliTPCkineGrid(*dummy); |
b2bca9d4 | 1761 | for(Int_t i=0;i<5;i++) new(&fPullsPr[i]) AliTPCkineGrid(*dummy); |
6eb69b9a | 1762 | for(Int_t i=0;i<5;i++) new(&fPullsEl[i]) AliTPCkineGrid(*dummy); |
b2bca9d4 | 1763 | for(Int_t i=0;i<5;i++) new(&fPullsMu[i]) AliTPCkineGrid(*dummy); |
6eb69b9a | 1764 | delete dummy; |
1765 | } | |
1766 | if(dEdx) { | |
b2bca9d4 | 1767 | dummy = new AliTPCkineGrid(nPt,nEta,pt,eta); |
6eb69b9a | 1768 | new(&fdEdxMeanPi) AliTPCkineGrid(*dummy); |
1769 | new(&fdEdxRMSPi) AliTPCkineGrid(*dummy); | |
6eb69b9a | 1770 | new(&fdEdxMeanKa) AliTPCkineGrid(*dummy); |
1771 | new(&fdEdxRMSKa) AliTPCkineGrid(*dummy); | |
6eb69b9a | 1772 | new(&fdEdxMeanPr) AliTPCkineGrid(*dummy); |
1773 | new(&fdEdxRMSPr) AliTPCkineGrid(*dummy); | |
6eb69b9a | 1774 | new(&fdEdxMeanEl) AliTPCkineGrid(*dummy); |
1775 | new(&fdEdxRMSEl) AliTPCkineGrid(*dummy); | |
b2bca9d4 | 1776 | new(&fdEdxMeanMu) AliTPCkineGrid(*dummy); |
1777 | new(&fdEdxRMSMu) AliTPCkineGrid(*dummy); | |
1778 | delete dummy; | |
6eb69b9a | 1779 | } |
1780 | ||
1781 | return; | |
1782 | } | |
1783 | //----------------------------------------------------------------------------- | |
1784 | void AliTPCtrackerParam::MakeDataBase() { | |
1785 | //----------------------------------------------------------------------------- | |
1786 | // This function creates the DB file and store in it: | |
1787 | // - TPC Efficiencies for pi,ka,pr,el,mu (AliTPCkineGrid class) | |
1788 | // - Pulls for pi,ka,el (AliTPCkineGrid class) | |
1789 | // - Regularization parameters for pi,ka,el (TMatrixD class) | |
1790 | // - dE/dx parameterization for pi,ka,pr,el (AliTPCkineGrid class) | |
1791 | // - Regularized cov. matrices for pi,ka,el (COVMATRIX structure) | |
1792 | //----------------------------------------------------------------------------- | |
1793 | ||
1794 | // define some file names | |
b2bca9d4 | 1795 | Char_t *effFile ="TPCeff.root"; |
1796 | Char_t *pullsFile ="pulls.root"; | |
1797 | Char_t *regPiFile ="regPi.root"; | |
1798 | Char_t *regKaFile ="regKa.root"; | |
1799 | Char_t *regPrFile ="regPr.root"; | |
1800 | Char_t *regElFile ="regEl.root"; | |
1801 | Char_t *regMuFile ="regMu.root"; | |
6eb69b9a | 1802 | Char_t *dEdxPiFile="dEdxPi.root"; |
1803 | Char_t *dEdxKaFile="dEdxKa.root"; | |
1804 | Char_t *dEdxPrFile="dEdxPr.root"; | |
1805 | Char_t *dEdxElFile="dEdxEl.root"; | |
b2bca9d4 | 1806 | Char_t *dEdxMuFile="dEdxMu.root"; |
1807 | Char_t *cmFile ="CovMatrix_AllEvts.root"; | |
1808 | /* | |
6eb69b9a | 1809 | Char_t *cmFile1 ="CovMatrix_AllEvts_1.root"; |
1810 | Char_t *cmFile2 ="CovMatrix_AllEvts_2.root"; | |
1811 | Char_t *cmFile3 ="CovMatrix_AllEvts_3.root"; | |
b2bca9d4 | 1812 | */ |
6eb69b9a | 1813 | |
1814 | // store the effieciencies | |
1815 | ReadEffs(effFile); | |
1816 | WriteEffs(fDBfileName.Data()); | |
1817 | ||
1818 | // store the pulls | |
1819 | ReadPulls(pullsFile); | |
1820 | WritePulls(fDBfileName.Data()); | |
1821 | ||
1822 | //store the regularization parameters | |
1823 | ReadRegParams(regPiFile,211); | |
1824 | WriteRegParams(fDBfileName.Data(),211); | |
1825 | ReadRegParams(regKaFile,321); | |
1826 | WriteRegParams(fDBfileName.Data(),321); | |
b2bca9d4 | 1827 | ReadRegParams(regPrFile,2212); |
1828 | WriteRegParams(fDBfileName.Data(),2212); | |
6eb69b9a | 1829 | ReadRegParams(regElFile,11); |
1830 | WriteRegParams(fDBfileName.Data(),11); | |
b2bca9d4 | 1831 | ReadRegParams(regMuFile,13); |
1832 | WriteRegParams(fDBfileName.Data(),13); | |
6eb69b9a | 1833 | |
1834 | // store the dEdx parameters | |
1835 | ReaddEdx(dEdxPiFile,211); | |
1836 | WritedEdx(fDBfileName.Data(),211); | |
1837 | ReaddEdx(dEdxKaFile,321); | |
1838 | WritedEdx(fDBfileName.Data(),321); | |
1839 | ReaddEdx(dEdxPrFile,2212); | |
1840 | WritedEdx(fDBfileName.Data(),2212); | |
1841 | ReaddEdx(dEdxElFile,11); | |
1842 | WritedEdx(fDBfileName.Data(),11); | |
b2bca9d4 | 1843 | ReaddEdx(dEdxMuFile,13); |
1844 | WritedEdx(fDBfileName.Data(),13); | |
6eb69b9a | 1845 | |
1846 | ||
1847 | // | |
1848 | // store the regularized covariance matrices | |
1849 | // | |
b2bca9d4 | 1850 | InitializeKineGrid("DB"); |
6eb69b9a | 1851 | |
1852 | const Int_t nBinsPi = fDBgridPi.GetTotBins(); | |
1853 | const Int_t nBinsKa = fDBgridKa.GetTotBins(); | |
b2bca9d4 | 1854 | const Int_t nBinsPr = fDBgridPr.GetTotBins(); |
6eb69b9a | 1855 | const Int_t nBinsEl = fDBgridEl.GetTotBins(); |
b2bca9d4 | 1856 | const Int_t nBinsMu = fDBgridMu.GetTotBins(); |
6eb69b9a | 1857 | |
1858 | ||
1859 | // create the trees for cov. matrices | |
1860 | // trees for pions | |
1861 | TTree *CovTreePi_ = NULL; | |
1862 | CovTreePi_ = new TTree[nBinsPi]; | |
1863 | // trees for kaons | |
1864 | TTree *CovTreeKa_ = NULL; | |
1865 | CovTreeKa_ = new TTree[nBinsKa]; | |
b2bca9d4 | 1866 | // trees for protons |
1867 | TTree *CovTreePr_ = NULL; | |
1868 | CovTreePr_ = new TTree[nBinsPr]; | |
6eb69b9a | 1869 | // trees for electrons |
1870 | TTree *CovTreeEl_ = NULL; | |
1871 | CovTreeEl_ = new TTree[nBinsEl]; | |
b2bca9d4 | 1872 | // trees for muons |
1873 | TTree *CovTreeMu_ = NULL; | |
1874 | CovTreeMu_ = new TTree[nBinsMu]; | |
6eb69b9a | 1875 | |
1876 | Char_t hname[100], htitle[100]; | |
1877 | COVMATRIX covmat; | |
1878 | ||
1879 | ||
1880 | for(Int_t i=0; i<nBinsPi; i++) { | |
1881 | sprintf(hname,"CovTreePi_bin%d",i); | |
1882 | sprintf(htitle,"Tree with cov matrix elements for bin %d",i); | |
1883 | CovTreePi_[i].SetName(hname); CovTreePi_[i].SetTitle(htitle); | |
1884 | CovTreePi_[i].Branch("matrix",&covmat,"c00/D:c10:c11:c20:c21:c22:c30:c31:c32:c33:c40:c41:c42:c43:c44",5000000); | |
1885 | } | |
1886 | for(Int_t i=0; i<nBinsKa; i++) { | |
1887 | sprintf(hname,"CovTreeKa_bin%d",i); | |
1888 | sprintf(htitle,"Tree with cov matrix elements for bin %d",i); | |
1889 | CovTreeKa_[i].SetName(hname); CovTreeKa_[i].SetTitle(htitle); | |
1890 | CovTreeKa_[i].Branch("matrix",&covmat,"c00/D:c10:c11:c20:c21:c22:c30:c31:c32:c33:c40:c41:c42:c43:c44",1000000); | |
1891 | } | |
b2bca9d4 | 1892 | for(Int_t i=0; i<nBinsPr; i++) { |
1893 | sprintf(hname,"CovTreePr_bin%d",i); | |
1894 | sprintf(htitle,"Tree with cov matrix elements for bin %d",i); | |
1895 | CovTreePr_[i].SetName(hname); CovTreePr_[i].SetTitle(htitle); | |
1896 | CovTreePr_[i].Branch("matrix",&covmat,"c00/D:c10:c11:c20:c21:c22:c30:c31:c32:c33:c40:c41:c42:c43:c44",1000000); | |
1897 | } | |
6eb69b9a | 1898 | for(Int_t i=0; i<nBinsEl; i++) { |
1899 | sprintf(hname,"CovTreeEl_bin%d",i); | |
1900 | sprintf(htitle,"Tree with cov matrix elements for bin %d",i); | |
1901 | CovTreeEl_[i].SetName(hname); CovTreeEl_[i].SetTitle(htitle); | |
1902 | CovTreeEl_[i].Branch("matrix",&covmat,"c00/D:c10:c11:c20:c21:c22:c30:c31:c32:c33:c40:c41:c42:c43:c44",1000000); | |
1903 | } | |
b2bca9d4 | 1904 | for(Int_t i=0; i<nBinsMu; i++) { |
1905 | sprintf(hname,"CovTreeMu_bin%d",i); | |
1906 | sprintf(htitle,"Tree with cov matrix elements for bin %d",i); | |
1907 | CovTreeMu_[i].SetName(hname); CovTreeMu_[i].SetTitle(htitle); | |
1908 | CovTreeMu_[i].Branch("matrix",&covmat,"c00/D:c10:c11:c20:c21:c22:c30:c31:c32:c33:c40:c41:c42:c43:c44",1000000); | |
1909 | } | |
1910 | ||
1911 | /* | |
6eb69b9a | 1912 | // create the chain with the compared tracks |
b2bca9d4 | 1913 | TChain *cmptrktree = new TChain("cmptrktree"); |
6eb69b9a | 1914 | cmptrkchain.Add(cmFile1); |
1915 | cmptrkchain.Add(cmFile2); | |
1916 | cmptrkchain.Add(cmFile3); | |
b2bca9d4 | 1917 | */ |
6eb69b9a | 1918 | |
b2bca9d4 | 1919 | TFile *infile = TFile::Open(cmFile); |
1920 | TTree *cmptrktree = (TTree*)infile->Get("cmptrktree"); | |
1921 | ||
6eb69b9a | 1922 | COMPTRACK cmptrk; |
b2bca9d4 | 1923 | cmptrktree->SetBranchAddress("comptracks",&cmptrk); |
1924 | Int_t entries = (Int_t)cmptrktree->GetEntries(); | |
6eb69b9a | 1925 | cerr<<" Number of entries: "<<entries<<endl; |
1926 | ||
1927 | Int_t trkPdg,trkBin; | |
b2bca9d4 | 1928 | Double_t trkKine[1],trkRegPar[3]; |
1929 | Int_t *nPerBinPi = new Int_t[nBinsPi]; | |
c15d4c37 | 1930 | for(Int_t k=0;k<nBinsPi;k++) nPerBinPi[k]=0; |
b2bca9d4 | 1931 | Int_t *nPerBinKa = new Int_t[nBinsKa]; |
c15d4c37 | 1932 | for(Int_t k=0;k<nBinsKa;k++) nPerBinKa[k]=0; |
b2bca9d4 | 1933 | Int_t *nPerBinMu = new Int_t[nBinsMu]; |
1934 | for(Int_t k=0;k<nBinsMu;k++) nPerBinMu[k]=0; | |
1935 | Int_t *nPerBinEl = new Int_t[nBinsEl]; | |
c15d4c37 | 1936 | for(Int_t k=0;k<nBinsEl;k++) nPerBinEl[k]=0; |
b2bca9d4 | 1937 | Int_t *nPerBinPr = new Int_t[nBinsPr]; |
1938 | for(Int_t k=0;k<nBinsPr;k++) nPerBinPr[k]=0; | |
6eb69b9a | 1939 | |
1940 | // loop on chain entries | |
1941 | for(Int_t l=0; l<entries; l++) { | |
1942 | if(l % 10000 == 0) cerr<<"--- Processing track "<<l<<" of "<<entries<<" ---"<<endl; | |
1943 | // get the event | |
b2bca9d4 | 1944 | cmptrktree->GetEvent(l); |
6eb69b9a | 1945 | // get the pdg code |
1946 | trkPdg = TMath::Abs(cmptrk.pdg); | |
b2bca9d4 | 1947 | // use only pions, kaons, protons, electrons, muons |
1948 | if(trkPdg!=211 && trkPdg!=321 && trkPdg!=2212 && trkPdg!=11 && trkPdg!=13) continue; | |
6eb69b9a | 1949 | SetParticle(trkPdg); |
1950 | trkBin = fDBgrid->GetBin(cmptrk.pt,cmptrk.eta); | |
1951 | //cerr<<cmptrk.pt<<" "<<cmptrk.eta<<" "<<trkBin<<endl; | |
1952 | ||
1953 | if(trkPdg==211 && nPerBinPi[trkBin]>=5000) continue; | |
1954 | if(trkPdg==321 && nPerBinKa[trkBin]>=5000) continue; | |
b2bca9d4 | 1955 | if(trkPdg==2212 && nPerBinPr[trkBin]>=5000) continue; |
6eb69b9a | 1956 | if(trkPdg==11 && nPerBinEl[trkBin]>=5000) continue; |
b2bca9d4 | 1957 | if(trkPdg==13 && nPerBinMu[trkBin]>=5000) continue; |
6eb69b9a | 1958 | |
1959 | trkKine[0] = cmptrk.p; | |
6eb69b9a | 1960 | |
1961 | // get regularized covariance matrix | |
b2bca9d4 | 1962 | for(Int_t k=0;k<3;k++) trkRegPar[k] = (*fRegPar)(0,k); |
6eb69b9a | 1963 | covmat.c00 = cmptrk.c00/RegFunc(trkKine,trkRegPar); |
1964 | covmat.c10 = cmptrk.c10; | |
b2bca9d4 | 1965 | for(Int_t k=0;k<3;k++) trkRegPar[k] = (*fRegPar)(1,k); |
6eb69b9a | 1966 | covmat.c11 = cmptrk.c11/RegFunc(trkKine,trkRegPar); |
b2bca9d4 | 1967 | for(Int_t k=0;k<3;k++) trkRegPar[k] = (*fRegPar)(2,k); |
6eb69b9a | 1968 | covmat.c20 = cmptrk.c20/RegFunc(trkKine,trkRegPar); |
1969 | covmat.c21 = cmptrk.c21; | |
b2bca9d4 | 1970 | for(Int_t k=0;k<3;k++) trkRegPar[k] = (*fRegPar)(3,k); |
6eb69b9a | 1971 | covmat.c22 = cmptrk.c22/RegFunc(trkKine,trkRegPar); |
1972 | covmat.c30 = cmptrk.c30; | |
b2bca9d4 | 1973 | for(Int_t k=0;k<3;k++) trkRegPar[k] = (*fRegPar)(4,k); |
6eb69b9a | 1974 | covmat.c31 = cmptrk.c31/RegFunc(trkKine,trkRegPar); |
1975 | covmat.c32 = cmptrk.c32; | |
b2bca9d4 | 1976 | for(Int_t k=0;k<3;k++) trkRegPar[k] = (*fRegPar)(5,k); |
6eb69b9a | 1977 | covmat.c33 = cmptrk.c33/RegFunc(trkKine,trkRegPar); |
b2bca9d4 | 1978 | for(Int_t k=0;k<3;k++) trkRegPar[k] = (*fRegPar)(6,k); |
6eb69b9a | 1979 | covmat.c40 = cmptrk.c40/RegFunc(trkKine,trkRegPar); |
1980 | covmat.c41 = cmptrk.c41; | |
b2bca9d4 | 1981 | for(Int_t k=0;k<3;k++) trkRegPar[k] = (*fRegPar)(7,k); |
6eb69b9a | 1982 | covmat.c42 = cmptrk.c42/RegFunc(trkKine,trkRegPar); |
1983 | covmat.c43 = cmptrk.c43; | |
b2bca9d4 | 1984 | for(Int_t k=0;k<3;k++) trkRegPar[k] = (*fRegPar)(8,k); |
6eb69b9a | 1985 | covmat.c44 = cmptrk.c44/RegFunc(trkKine,trkRegPar); |
1986 | ||
1987 | // fill the tree | |
1988 | switch (trkPdg) { | |
1989 | case 211: // pions | |
1990 | CovTreePi_[trkBin].Fill(); | |
1991 | nPerBinPi[trkBin]++; | |
1992 | break; | |
1993 | case 321: // kaons | |
1994 | CovTreeKa_[trkBin].Fill(); | |
1995 | nPerBinKa[trkBin]++; | |
1996 | break; | |
b2bca9d4 | 1997 | case 2212: // protons |
1998 | CovTreePr_[trkBin].Fill(); | |
1999 | nPerBinPr[trkBin]++; | |
2000 | break; | |
6eb69b9a | 2001 | case 11: // electrons |
2002 | CovTreeEl_[trkBin].Fill(); | |
2003 | nPerBinEl[trkBin]++; | |
2004 | break; | |
b2bca9d4 | 2005 | case 13: // muons |
2006 | CovTreeMu_[trkBin].Fill(); | |
2007 | nPerBinMu[trkBin]++; | |
2008 | break; | |
6eb69b9a | 2009 | } |
2010 | } // loop on chain entries | |
2011 | ||
2012 | // store all trees the DB file | |
2013 | TFile *DBfile = new TFile(fDBfileName.Data(),"update"); | |
2014 | DBfile->mkdir("CovMatrices"); | |
2015 | gDirectory->cd("/CovMatrices"); | |
2016 | gDirectory->mkdir("Pions"); | |
2017 | gDirectory->mkdir("Kaons"); | |
b2bca9d4 | 2018 | gDirectory->mkdir("Protons"); |
6eb69b9a | 2019 | gDirectory->mkdir("Electrons"); |
b2bca9d4 | 2020 | gDirectory->mkdir("Muons"); |
6eb69b9a | 2021 | // store pions |
2022 | gDirectory->cd("/CovMatrices/Pions"); | |
2023 | fDBgridPi.SetName("DBgridPi"); fDBgridPi.Write(); | |
2024 | for(Int_t i=0;i<nBinsPi;i++) CovTreePi_[i].Write(); | |
2025 | // store kaons | |
2026 | gDirectory->cd("/CovMatrices/Kaons"); | |
2027 | fDBgridKa.SetName("DBgridKa"); fDBgridKa.Write(); | |
2028 | for(Int_t i=0;i<nBinsKa;i++) CovTreeKa_[i].Write(); | |
b2bca9d4 | 2029 | // store kaons |
2030 | gDirectory->cd("/CovMatrices/Protons"); | |
2031 | fDBgridPr.SetName("DBgridPr"); fDBgridPr.Write(); | |
2032 | for(Int_t i=0;i<nBinsPr;i++) CovTreePr_[i].Write(); | |
6eb69b9a | 2033 | // store electrons |
2034 | gDirectory->cd("/CovMatrices/Electrons"); | |
2035 | fDBgridEl.SetName("DBgridEl"); fDBgridEl.Write(); | |
2036 | for(Int_t i=0;i<nBinsEl;i++) CovTreeEl_[i].Write(); | |
b2bca9d4 | 2037 | // store kaons |
2038 | gDirectory->cd("/CovMatrices/Muons"); | |
2039 | fDBgridMu.SetName("DBgridMu"); fDBgridMu.Write(); | |
2040 | for(Int_t i=0;i<nBinsMu;i++) CovTreeMu_[i].Write(); | |
6eb69b9a | 2041 | |
2042 | DBfile->Close(); | |
c15d4c37 | 2043 | delete [] nPerBinPi; |
2044 | delete [] nPerBinKa; | |
b2bca9d4 | 2045 | delete [] nPerBinPr; |
c15d4c37 | 2046 | delete [] nPerBinEl; |
b2bca9d4 | 2047 | delete [] nPerBinMu; |
2048 | ||
2049 | return; | |
2050 | } | |
2051 | //----------------------------------------------------------------------------- | |
2052 | void AliTPCtrackerParam::MakeSeedsFromHits(AliTPC *TPC,TTree *TH, | |
2053 | TObjArray &seedArray) const { | |
2054 | //----------------------------------------------------------------------------- | |
2055 | // This function makes the seeds for tracks from the 1st hits in the TPC | |
2056 | //----------------------------------------------------------------------------- | |
2057 | ||
2058 | Double_t xg,yg,zg,px,py,pz,pt; | |
2059 | Int_t label; | |
2060 | Int_t nTracks=(Int_t)TH->GetEntries(); | |
2061 | ||
2062 | cerr<<"+++\n+++ Number of \"primary tracks\"(entries in TreeH): "<<nTracks<< | |
2063 | "\n+++\n\n"; | |
2064 | ||
2065 | AliTPChit *tpcHit=0; | |
2066 | ||
2067 | // loop over entries in TreeH | |
2068 | for(Int_t l=0; l<nTracks; l++) { | |
2069 | if(l%1000==0) cerr<<" --- Processing primary track " | |
2070 | <<l<<" of "<<nTracks<<" ---\r"; | |
2071 | TPC->ResetHits(); | |
2072 | TH->GetEvent(l); | |
2073 | // Get FirstHit | |
2074 | tpcHit=(AliTPChit*)TPC->FirstHit(-1); | |
2075 | for( ; tpcHit; tpcHit=(AliTPChit*)TPC->NextHit() ) { | |
2076 | if(tpcHit->fQ !=0.) continue; | |
2077 | // Get particle momentum at hit | |
2078 | px=tpcHit->X(); py=tpcHit->Y(); pz=tpcHit->Z(); | |
2079 | ||
2080 | pt=TMath::Sqrt(px*px+py*py); | |
2081 | // reject hits with Pt<mag*0.45 GeV/c | |
2082 | if(pt<(fBz*0.45)) continue; | |
2083 | ||
2084 | // Get track label | |
2085 | label=tpcHit->Track(); | |
2086 | ||
2087 | if((tpcHit=(AliTPChit*)TPC->NextHit())==0) break; | |
2088 | if(tpcHit->fQ != 0.) continue; | |
2089 | // Get global coordinates of hit | |
2090 | xg=tpcHit->X(); yg=tpcHit->Y(); zg=tpcHit->Z(); | |
2091 | if(TMath::Sqrt(xg*xg+yg*yg)>90.) continue; | |
2092 | ||
2093 | // create the seed | |
2094 | AliTPCseedGeant *ioseed = new AliTPCseedGeant(xg,yg,zg,px,py,pz,label); | |
2095 | ||
2096 | // reject tracks which are not in the TPC acceptance | |
2097 | if(!ioseed->InTPCAcceptance()) { delete ioseed; continue; } | |
2098 | ||
2099 | // put seed in array | |
2100 | seedArray.AddLast(ioseed); | |
2101 | } | |
2102 | ||
2103 | } // loop over entries in TreeH | |
2104 | ||
2105 | return; | |
2106 | } | |
2107 | //----------------------------------------------------------------------------- | |
2108 | void AliTPCtrackerParam::MakeSeedsFromRefs(TTree *TTR, | |
2109 | TObjArray &seedArray) const { | |
2110 | //----------------------------------------------------------------------------- | |
2111 | // This function makes the seeds for tracks from the track references | |
2112 | //----------------------------------------------------------------------------- | |
2113 | ||
2114 | Double_t xg,yg,zg,px,py,pz,pt; | |
2115 | Int_t label; | |
2116 | Int_t nnn,k; | |
2117 | ||
2118 | TClonesArray *tkRefArray = new TClonesArray("AliTrackReference"); | |
2119 | ||
2120 | TBranch *b =(TBranch*)TTR->GetBranch("TPC"); | |
2121 | if(!b) {cerr<<"TPC branch of TreeTR not found"<<endl; return; } | |
2122 | b->SetAddress(&tkRefArray); | |
2123 | Int_t nTkRef = (Int_t)b->GetEntries(); | |
2124 | cerr<<"+++\n+++ Number of entries in TreeTR(TPC): "<<nTkRef<< | |
2125 | "\n+++\n\n"; | |
2126 | ||
2127 | // loop on track references | |
2128 | for(Int_t l=0; l<nTkRef; l++){ | |
2129 | if(l%1000==0) cerr<<" --- Processing primary track " | |
2130 | <<l<<" of "<<nTkRef<<" ---\r"; | |
2131 | if(!b->GetEvent(l)) continue; | |
2132 | nnn = tkRefArray->GetEntriesFast(); | |
2133 | ||
2134 | if(nnn <= 0) continue; | |
2135 | for(k=0; k<nnn; k++) { | |
2136 | AliTrackReference *tkref = (AliTrackReference*)tkRefArray->UncheckedAt(k); | |
2137 | ||
2138 | xg = tkref->X(); | |
2139 | yg = tkref->Y(); | |
2140 | zg = tkref->Z(); | |
2141 | px = tkref->Px(); | |
2142 | py = tkref->Py(); | |
2143 | pz = tkref->Pz(); | |
2144 | label = tkref->GetTrack(); | |
2145 | ||
2146 | pt=TMath::Sqrt(px*px+py*py); | |
2147 | // reject hits with Pt<mag*0.45 GeV/c | |
2148 | if(pt<(fBz*0.45)) continue; | |
2149 | ||
2150 | // create the seed | |
2151 | AliTPCseedGeant *ioseed = new AliTPCseedGeant(xg,yg,zg,px,py,pz,label); | |
2152 | ||
2153 | // reject if not at the inner part of TPC | |
2154 | if(TMath::Abs(ioseed->GetXL()-82.9701) > 0.1) { | |
2155 | delete ioseed; continue; | |
2156 | } | |
2157 | ||
2158 | // reject tracks which are not in the TPC acceptance | |
2159 | if(!ioseed->InTPCAcceptance()) { | |
2160 | delete ioseed; continue; | |
2161 | } | |
2162 | ||
2163 | // put seed in array | |
2164 | seedArray.AddLast(ioseed); | |
2165 | ||
2166 | } | |
2167 | ||
2168 | ||
2169 | } // loop on track references | |
2170 | ||
2171 | delete tkRefArray; | |
2172 | ||
2173 | ||
6eb69b9a | 2174 | return; |
2175 | } | |
2176 | //----------------------------------------------------------------------------- | |
c15d4c37 | 2177 | void AliTPCtrackerParam::MergeEvents(Int_t evFirst,Int_t evLast) { |
6eb69b9a | 2178 | //----------------------------------------------------------------------------- |
2179 | // This function: 1) merges the files from track comparison | |
2180 | // (beware: better no more than 100 events per file) | |
2181 | // 2) computes the average TPC efficiencies | |
2182 | //----------------------------------------------------------------------------- | |
2183 | ||
2184 | Char_t *outName="TPCeff.root"; | |
2185 | ||
b2bca9d4 | 2186 | // merge files with tracks |
2187 | cerr<<" ******** MERGING FILES **********\n\n"; | |
6eb67451 | 2188 | |
6eb69b9a | 2189 | // create the chain for the tree of compared tracks |
2190 | TChain ch1("cmptrktree"); | |
2191 | TChain ch2("cmptrktree"); | |
2192 | TChain ch3("cmptrktree"); | |
6eb67451 | 2193 | |
6eb69b9a | 2194 | for(Int_t j=evFirst; j<=evLast; j++) { |
2195 | cerr<<"Processing event "<<j<<endl; | |
6eb67451 | 2196 | |
6eb69b9a | 2197 | TString covName("CovMatrix."); |
6eb69b9a | 2198 | covName+=j; |
6eb69b9a | 2199 | covName.Append(".root"); |
6eb67451 | 2200 | |
b2bca9d4 | 2201 | if(gSystem->AccessPathName(covName.Data(),kFileExists)) continue; |
2202 | ||
6eb69b9a | 2203 | |
2204 | if(j<=100) ch1.Add(covName.Data()); | |
2205 | if(j>100 && j<=200) ch2.Add(covName.Data()); | |
2206 | if(j>200) ch3.Add(covName.Data()); | |
2207 | ||
6eb69b9a | 2208 | } // loop on events |
2209 | ||
2210 | // merge chain in one file | |
2211 | TFile *covOut=0; | |
2212 | covOut = new TFile("CovMatrix_AllEvts_1.root","recreate"); | |
2213 | ch1.Merge(covOut,1000000000); | |
2214 | covOut->Close(); | |
2215 | delete covOut; | |
2216 | covOut = new TFile("CovMatrix_AllEvts_2.root","recreate"); | |
2217 | ch2.Merge(covOut,1000000000); | |
2218 | covOut->Close(); | |
2219 | delete covOut; | |
2220 | covOut = new TFile("CovMatrix_AllEvts_3.root","recreate"); | |
2221 | ch3.Merge(covOut,1000000000); | |
2222 | covOut->Close(); | |
2223 | delete covOut; | |
2224 | ||
2225 | ||
2226 | // efficiencies | |
b2bca9d4 | 2227 | cerr<<" ***** EFFICIENCIES ******\n\n"; |
2228 | ||
2229 | ReadEffs("TPCeff.1.root"); | |
2230 | ||
2231 | Int_t n = fEffPi.GetTotPoints(); | |
2232 | Double_t *avEffPi = new Double_t[n]; | |
2233 | Double_t *avEffKa = new Double_t[n]; | |
2234 | Double_t *avEffPr = new Double_t[n]; | |
2235 | Double_t *avEffEl = new Double_t[n]; | |
2236 | Double_t *avEffMu = new Double_t[n]; | |
2237 | Int_t *evtsPi = new Int_t[n]; | |
2238 | Int_t *evtsKa = new Int_t[n]; | |
2239 | Int_t *evtsPr = new Int_t[n]; | |
2240 | Int_t *evtsEl = new Int_t[n]; | |
2241 | Int_t *evtsMu = new Int_t[n]; | |
2242 | ||
2243 | for(Int_t j=0; j<n; j++) { | |
2244 | avEffPi[j]=avEffKa[j]=avEffPr[j]=avEffEl[j]=avEffMu[j]=0.; | |
2245 | evtsPi[j]=evtsKa[j]=evtsPr[j]=evtsEl[j]=evtsMu[j]=0; | |
2246 | } | |
2247 | ||
2248 | for(Int_t j=evFirst; j<=evLast; j++) { | |
2249 | cerr<<"Processing event "<<j<<endl; | |
6eb69b9a | 2250 | |
b2bca9d4 | 2251 | TString effName("TPCeff."); |
2252 | effName+=j; | |
2253 | effName.Append(".root"); | |
2254 | ||
2255 | if(gSystem->AccessPathName(effName.Data(),kFileExists)) continue; | |
2256 | ||
2257 | ReadEffs(effName.Data()); | |
2258 | ||
2259 | for(Int_t k=0; k<n; k++) { | |
2260 | if(fEffPi.GetParam(k)>=0.) {avEffPi[k]+=fEffPi.GetParam(k); evtsPi[k]++;} | |
2261 | if(fEffKa.GetParam(k)>=0.) {avEffKa[k]+=fEffKa.GetParam(k); evtsKa[k]++;} | |
2262 | if(fEffPr.GetParam(k)>=0.) {avEffPr[k]+=fEffPr.GetParam(k); evtsPr[k]++;} | |
2263 | if(fEffEl.GetParam(k)>=0.) {avEffEl[k]+=fEffEl.GetParam(k); evtsEl[k]++;} | |
2264 | if(fEffMu.GetParam(k)>=0.) {avEffMu[k]+=fEffMu.GetParam(k); evtsMu[k]++;} | |
2265 | } | |
2266 | ||
2267 | } // loop on events | |
6eb69b9a | 2268 | |
b2bca9d4 | 2269 | // compute average efficiencies |
2270 | for(Int_t j=0; j<n; j++) { | |
6eb69b9a | 2271 | if(evtsPi[j]==0) evtsPi[j]++; |
2272 | fEffPi.SetParam(j,(Double_t)avEffPi[j]/evtsPi[j]); | |
b2bca9d4 | 2273 | if(evtsKa[j]==0) evtsKa[j]++; |
2274 | fEffKa.SetParam(j,(Double_t)avEffKa[j]/evtsKa[j]); | |
6eb69b9a | 2275 | if(evtsPr[j]==0) evtsPr[j]++; |
2276 | fEffPr.SetParam(j,(Double_t)avEffPr[j]/evtsPr[j]); | |
6eb69b9a | 2277 | if(evtsEl[j]==0) evtsEl[j]++; |
2278 | fEffEl.SetParam(j,(Double_t)avEffEl[j]/evtsEl[j]); | |
6eb69b9a | 2279 | if(evtsMu[j]==0) evtsMu[j]++; |
2280 | fEffMu.SetParam(j,(Double_t)avEffMu[j]/evtsMu[j]); | |
2281 | } | |
b2bca9d4 | 2282 | |
6eb69b9a | 2283 | // write efficiencies to a file |
2284 | WriteEffs(outName); | |
2285 | ||
b2bca9d4 | 2286 | delete [] avEffPi; |
2287 | delete [] avEffKa; | |
2288 | delete [] avEffPr; | |
2289 | delete [] avEffEl; | |
2290 | delete [] avEffMu; | |
2291 | delete [] evtsPi; | |
2292 | delete [] evtsKa; | |
2293 | delete [] evtsPr; | |
2294 | delete [] evtsEl; | |
2295 | delete [] evtsMu; | |
2296 | ||
6eb69b9a | 2297 | return; |
2298 | } | |
2299 | //----------------------------------------------------------------------------- | |
2300 | Int_t AliTPCtrackerParam::ReadAllData(const Char_t* inName) { | |
2301 | //----------------------------------------------------------------------------- | |
2302 | // This function reads all parameters from the DB | |
2303 | //----------------------------------------------------------------------------- | |
2304 | ||
2305 | if(!ReadEffs(inName)) return 0; | |
2306 | if(!ReadPulls(inName)) return 0; | |
2307 | if(!ReadRegParams(inName,211)) return 0; | |
2308 | if(!ReadRegParams(inName,321)) return 0; | |
b2bca9d4 | 2309 | if(!ReadRegParams(inName,2212)) return 0; |
6eb69b9a | 2310 | if(!ReadRegParams(inName,11)) return 0; |
b2bca9d4 | 2311 | if(!ReadRegParams(inName,13)) return 0; |
6eb69b9a | 2312 | if(!ReaddEdx(inName,211)) return 0; |
2313 | if(!ReaddEdx(inName,321)) return 0; | |
2314 | if(!ReaddEdx(inName,2212)) return 0; | |
2315 | if(!ReaddEdx(inName,11)) return 0; | |
b2bca9d4 | 2316 | if(!ReaddEdx(inName,13)) return 0; |
6eb69b9a | 2317 | if(!ReadDBgrid(inName)) return 0; |
2318 | ||
2319 | return 1; | |
6eb67451 | 2320 | } |
6eb69b9a | 2321 | //----------------------------------------------------------------------------- |
2322 | Int_t AliTPCtrackerParam::ReaddEdx(const Char_t* inName,Int_t pdg) { | |
2323 | //----------------------------------------------------------------------------- | |
2324 | // This function reads the dEdx parameters from the DB | |
2325 | //----------------------------------------------------------------------------- | |
2326 | ||
2327 | if(gSystem->AccessPathName(inName,kFileExists)) { | |
2328 | cerr<<"AliTPCtrackerParam::ReaddEdx: "<<inName<<" not found\n"; | |
2329 | return 0; } | |
2330 | ||
2331 | TFile *inFile = TFile::Open(inName); | |
2332 | switch (pdg) { | |
2333 | case 211: | |
2334 | fdEdxMean = (AliTPCkineGrid*)inFile->Get("/dEdx/Pions/dEdxMeanPi"); | |
2335 | fdEdxMeanPi.~AliTPCkineGrid(); | |
2336 | new(&fdEdxMeanPi) AliTPCkineGrid(*fdEdxMean); | |
2337 | fdEdxRMS = (AliTPCkineGrid*)inFile->Get("/dEdx/Pions/dEdxRMSPi"); | |
2338 | fdEdxRMSPi.~AliTPCkineGrid(); | |
2339 | new(&fdEdxRMSPi) AliTPCkineGrid(*fdEdxRMS); | |
2340 | break; | |
2341 | case 321: | |
2342 | fdEdxMean = (AliTPCkineGrid*)inFile->Get("/dEdx/Kaons/dEdxMeanKa"); | |
2343 | fdEdxMeanKa.~AliTPCkineGrid(); | |
2344 | new(&fdEdxMeanKa) AliTPCkineGrid(*fdEdxMean); | |
2345 | fdEdxRMS = (AliTPCkineGrid*)inFile->Get("/dEdx/Kaons/dEdxRMSKa"); | |
2346 | fdEdxRMSKa.~AliTPCkineGrid(); | |
2347 | new(&fdEdxRMSKa) AliTPCkineGrid(*fdEdxRMS); | |
2348 | break; | |
2349 | case 2212: | |
2350 | fdEdxMean = (AliTPCkineGrid*)inFile->Get("/dEdx/Protons/dEdxMeanPr"); | |
2351 | fdEdxMeanPr.~AliTPCkineGrid(); | |
2352 | new(&fdEdxMeanPr) AliTPCkineGrid(*fdEdxMean); | |
2353 | fdEdxRMS = (AliTPCkineGrid*)inFile->Get("/dEdx/Protons/dEdxRMSPr"); | |
2354 | fdEdxRMSPr.~AliTPCkineGrid(); | |
2355 | new(&fdEdxRMSPr) AliTPCkineGrid(*fdEdxRMS); | |
2356 | break; | |
2357 | case 11: | |
2358 | fdEdxMean = (AliTPCkineGrid*)inFile->Get("/dEdx/Electrons/dEdxMeanEl"); | |
2359 | fdEdxMeanEl.~AliTPCkineGrid(); | |
2360 | new(&fdEdxMeanEl) AliTPCkineGrid(*fdEdxMean); | |
2361 | fdEdxRMS = (AliTPCkineGrid*)inFile->Get("/dEdx/Electrons/dEdxRMSEl"); | |
2362 | fdEdxRMSEl.~AliTPCkineGrid(); | |
2363 | new(&fdEdxRMSEl) AliTPCkineGrid(*fdEdxRMS); | |
2364 | break; | |
b2bca9d4 | 2365 | case 13: |
2366 | if(fColl==0) { | |
2367 | fdEdxMean = (AliTPCkineGrid*)inFile->Get("/dEdx/Pions/dEdxMeanPi"); | |
2368 | fdEdxRMS = (AliTPCkineGrid*)inFile->Get("/dEdx/Pions/dEdxRMSPi"); | |
2369 | } else { | |
2370 | fdEdxMean = (AliTPCkineGrid*)inFile->Get("/dEdx/Muons/dEdxMeanMu"); | |
2371 | fdEdxRMS = (AliTPCkineGrid*)inFile->Get("/dEdx/Muons/dEdxRMSMu"); | |
2372 | } | |
2373 | fdEdxMeanMu.~AliTPCkineGrid(); | |
2374 | new(&fdEdxMeanMu) AliTPCkineGrid(*fdEdxMean); | |
2375 | fdEdxRMSMu.~AliTPCkineGrid(); | |
2376 | new(&fdEdxRMSMu) AliTPCkineGrid(*fdEdxRMS); | |
2377 | break; | |
6eb69b9a | 2378 | } |
2379 | inFile->Close(); | |
6eb67451 | 2380 | |
6eb69b9a | 2381 | return 1; |
2382 | } | |
2383 | //----------------------------------------------------------------------------- | |
2384 | Int_t AliTPCtrackerParam::ReadDBgrid(const Char_t* inName) { | |
2385 | //----------------------------------------------------------------------------- | |
2386 | // This function reads the kine grid from the DB | |
2387 | //----------------------------------------------------------------------------- | |
2388 | ||
2389 | if(gSystem->AccessPathName(inName,kFileExists)) { | |
2390 | cerr<<"AliTPCtrackerParam::ReadCovMatrices: "<<inName<<" not found\n"; | |
2391 | return 0; } | |
2392 | ||
2393 | TFile *inFile = TFile::Open(inName); | |
2394 | ||
2395 | // first read the DB grid for the different particles | |
2396 | fDBgrid = (AliTPCkineGrid*)inFile->Get("/CovMatrices/Pions/DBgridPi"); | |
2397 | fDBgridPi.~AliTPCkineGrid(); | |
2398 | new(&fDBgridPi) AliTPCkineGrid(*fDBgrid); | |
2399 | fDBgrid = (AliTPCkineGrid*)inFile->Get("/CovMatrices/Kaons/DBgridKa"); | |
2400 | fDBgridKa.~AliTPCkineGrid(); | |
2401 | new(&fDBgridKa) AliTPCkineGrid(*fDBgrid); | |
b2bca9d4 | 2402 | if(fColl==0) { |
2403 | fDBgrid = (AliTPCkineGrid*)inFile->Get("/CovMatrices/Pions/DBgridPi"); | |
2404 | } else { | |
2405 | fDBgrid = (AliTPCkineGrid*)inFile->Get("/CovMatrices/Protons/DBgridPr"); | |
2406 | } | |
2407 | fDBgridPr.~AliTPCkineGrid(); | |
2408 | new(&fDBgridPr) AliTPCkineGrid(*fDBgrid); | |
6eb69b9a | 2409 | fDBgrid = (AliTPCkineGrid*)inFile->Get("/CovMatrices/Electrons/DBgridEl"); |
b2bca9d4 | 2410 | fDBgridEl.~AliTPCkineGrid(); |
6eb69b9a | 2411 | new(&fDBgridEl) AliTPCkineGrid(*fDBgrid); |
b2bca9d4 | 2412 | if(fColl==0) { |
2413 | fDBgrid = (AliTPCkineGrid*)inFile->Get("/CovMatrices/Pions/DBgridPi"); | |
2414 | } else { | |
2415 | fDBgrid = (AliTPCkineGrid*)inFile->Get("/CovMatrices/Muons/DBgridMu"); | |
2416 | } | |
2417 | fDBgridMu.~AliTPCkineGrid(); | |
2418 | new(&fDBgridMu) AliTPCkineGrid(*fDBgrid); | |
6eb69b9a | 2419 | |
2420 | inFile->Close(); | |
2421 | ||
2422 | return 1; | |
2423 | } | |
2424 | //----------------------------------------------------------------------------- | |
2425 | Int_t AliTPCtrackerParam::ReadEffs(const Char_t* inName) { | |
2426 | //----------------------------------------------------------------------------- | |
2427 | // This function reads the TPC efficiencies from the DB | |
2428 | //----------------------------------------------------------------------------- | |
2429 | ||
2430 | if(gSystem->AccessPathName(inName,kFileExists)) { | |
2431 | cerr<<"AliTPCtrackerParam::ReadEffs: "<<inName<<" not found\n"; | |
2432 | return 0; } | |
2433 | ||
2434 | TFile *inFile = TFile::Open(inName); | |
2435 | ||
2436 | fEff = (AliTPCkineGrid*)inFile->Get("/Efficiencies/Pions/EffPi"); | |
2437 | fEffPi.~AliTPCkineGrid(); | |
2438 | new(&fEffPi) AliTPCkineGrid(*fEff); | |
2439 | fEff = (AliTPCkineGrid*)inFile->Get("/Efficiencies/Kaons/EffKa"); | |
2440 | fEffKa.~AliTPCkineGrid(); | |
2441 | new(&fEffKa) AliTPCkineGrid(*fEff); | |
2442 | fEff = (AliTPCkineGrid*)inFile->Get("/Efficiencies/Protons/EffPr"); | |
2443 | fEffPr.~AliTPCkineGrid(); | |
2444 | new(&fEffPr) AliTPCkineGrid(*fEff); | |
2445 | fEff = (AliTPCkineGrid*)inFile->Get("/Efficiencies/Electrons/EffEl"); | |
2446 | fEffEl.~AliTPCkineGrid(); | |
2447 | new(&fEffEl) AliTPCkineGrid(*fEff); | |
2448 | fEff = (AliTPCkineGrid*)inFile->Get("/Efficiencies/Muons/EffMu"); | |
2449 | fEffMu.~AliTPCkineGrid(); | |
2450 | new(&fEffMu) AliTPCkineGrid(*fEff); | |
2451 | ||
2452 | inFile->Close(); | |
2453 | ||
2454 | return 1; | |
2455 | } | |
2456 | //----------------------------------------------------------------------------- | |
2457 | Int_t AliTPCtrackerParam::ReadPulls(const Char_t* inName) { | |
2458 | //----------------------------------------------------------------------------- | |
2459 | // This function reads the pulls from the DB | |
2460 | //----------------------------------------------------------------------------- | |
2461 | ||
2462 | if(gSystem->AccessPathName(inName,kFileExists)) { | |
2463 | cerr<<"AliTPCtrackerParam::ReadPulls: "<<inName<<" not found\n"; | |
2464 | return 0; } | |
2465 | ||
2466 | TFile *inFile = TFile::Open(inName); | |
2467 | ||
2468 | for(Int_t i=0; i<5; i++) { | |
2469 | TString pi("/Pulls/Pions/PullsPi_"); pi+=i; | |
2470 | TString ka("/Pulls/Kaons/PullsKa_"); ka+=i; | |
b2bca9d4 | 2471 | TString pr("/Pulls/Protons/PullsPr_"); pr+=i; |
6eb69b9a | 2472 | TString el("/Pulls/Electrons/PullsEl_"); el+=i; |
b2bca9d4 | 2473 | TString mu("/Pulls/Muons/PullsMu_"); mu+=i; |
2474 | ||
6eb69b9a | 2475 | fPulls = (AliTPCkineGrid*)inFile->Get(pi.Data()); |
2476 | fPullsPi[i].~AliTPCkineGrid(); | |
2477 | new(&fPullsPi[i]) AliTPCkineGrid(*fPulls); | |
b2bca9d4 | 2478 | |
6eb69b9a | 2479 | fPulls = (AliTPCkineGrid*)inFile->Get(ka.Data()); |
2480 | fPullsKa[i].~AliTPCkineGrid(); | |
2481 | new(&fPullsKa[i]) AliTPCkineGrid(*fPulls); | |
b2bca9d4 | 2482 | |
2483 | if(fColl==0) { | |
2484 | fPulls = (AliTPCkineGrid*)inFile->Get(pi.Data()); | |
2485 | } else { | |
2486 | fPulls = (AliTPCkineGrid*)inFile->Get(pr.Data()); | |
2487 | } | |
2488 | fPullsPr[i].~AliTPCkineGrid(); | |
2489 | new(&fPullsPr[i]) AliTPCkineGrid(*fPulls); | |
2490 | ||
6eb69b9a | 2491 | fPulls = (AliTPCkineGrid*)inFile->Get(el.Data()); |
2492 | fPullsEl[i].~AliTPCkineGrid(); | |
2493 | new(&fPullsEl[i]) AliTPCkineGrid(*fPulls); | |
b2bca9d4 | 2494 | |
2495 | if(fColl==0) { | |
2496 | fPulls = (AliTPCkineGrid*)inFile->Get(pi.Data()); | |
2497 | } else { | |
2498 | fPulls = (AliTPCkineGrid*)inFile->Get(mu.Data()); | |
2499 | } | |
2500 | fPullsMu[i].~AliTPCkineGrid(); | |
2501 | new(&fPullsMu[i]) AliTPCkineGrid(*fPulls); | |
6eb69b9a | 2502 | } |
2503 | ||
2504 | inFile->Close(); | |
2505 | ||
2506 | return 1; | |
2507 | } | |
2508 | //----------------------------------------------------------------------------- | |
2509 | Int_t AliTPCtrackerParam::ReadRegParams(const Char_t* inName,Int_t pdg) { | |
2510 | //----------------------------------------------------------------------------- | |
2511 | // This function reads the regularization parameters from the DB | |
2512 | //----------------------------------------------------------------------------- | |
2513 | ||
2514 | if(gSystem->AccessPathName(inName,kFileExists)) { | |
2515 | cerr<<"AliTPCtrackerParam::ReadRegParams: "<<inName<<" not found\n"; | |
2516 | return 0; } | |
2517 | ||
2518 | TFile *inFile = TFile::Open(inName); | |
2519 | switch (pdg) { | |
2520 | case 211: | |
2521 | fRegPar = (TMatrixD*)inFile->Get("/RegParams/Pions/RegPions"); | |
2522 | new(&fRegParPi) TMatrixD(*fRegPar); | |
2523 | break; | |
2524 | case 321: | |
2525 | fRegPar = (TMatrixD*)inFile->Get("/RegParams/Kaons/RegKaons"); | |
2526 | new(&fRegParKa) TMatrixD(*fRegPar); | |
2527 | break; | |
b2bca9d4 | 2528 | case 2212: |
2529 | if(fColl==0) { | |
2530 | fRegPar = (TMatrixD*)inFile->Get("/RegParams/Pions/RegPions"); | |
2531 | } else { | |
2532 | fRegPar = (TMatrixD*)inFile->Get("/RegParams/Protons/RegProtons"); | |
2533 | } | |
2534 | new(&fRegParPr) TMatrixD(*fRegPar); | |
2535 | break; | |
6eb69b9a | 2536 | case 11: |
2537 | fRegPar = (TMatrixD*)inFile->Get("/RegParams/Electrons/RegElectrons"); | |
2538 | new(&fRegParEl) TMatrixD(*fRegPar); | |
2539 | break; | |
b2bca9d4 | 2540 | case 13: |
2541 | if(fColl==0) { | |
2542 | fRegPar = (TMatrixD*)inFile->Get("/RegParams/Pions/RegPions"); | |
2543 | } else { | |
2544 | fRegPar = (TMatrixD*)inFile->Get("/RegParams/Muons/RegMuons"); | |
2545 | } | |
2546 | new(&fRegParMu) TMatrixD(*fRegPar); | |
2547 | break; | |
6eb69b9a | 2548 | } |
2549 | inFile->Close(); | |
2550 | ||
2551 | return 1; | |
2552 | } | |
2553 | //----------------------------------------------------------------------------- | |
2554 | void AliTPCtrackerParam::RegularizeCovMatrix(const Char_t *outName,Int_t pdg) { | |
2555 | //----------------------------------------------------------------------------- | |
2556 | // This function regularizes the elements of the covariance matrix | |
2557 | // that show a momentum depence: | |
2558 | // c00, c11, c22, c33, c44, c20, c24, c40, c31 | |
2559 | // The regularization is done separately for pions, kaons, electrons: | |
2560 | // give "Pion","Kaon" or "Electron" as first argument. | |
2561 | //----------------------------------------------------------------------------- | |
2562 | ||
2563 | gStyle->SetOptStat(0); | |
2564 | gStyle->SetOptFit(10001); | |
2565 | ||
2566 | Int_t thisPdg=211; | |
2567 | Char_t *part="Pions - "; | |
2568 | ||
b2bca9d4 | 2569 | InitializeKineGrid("DB"); |
6eb69b9a | 2570 | SetParticle(pdg); |
2571 | const Int_t fitbins = fDBgrid->GetBinsPt(); | |
2572 | cerr<<" Fit bins: "<<fitbins<<endl; | |
2573 | ||
2574 | switch (pdg) { | |
2575 | case 211: // pions | |
2576 | thisPdg=211; | |
2577 | part="Pions - "; | |
2578 | cerr<<" Processing pions ...\n"; | |
2579 | break; | |
2580 | case 321: //kaons | |
2581 | thisPdg=321; | |
2582 | part="Kaons - "; | |
2583 | cerr<<" Processing kaons ...\n"; | |
2584 | break; | |
b2bca9d4 | 2585 | case 2212: //protons |
2586 | thisPdg=2212; | |
2587 | part="Protons - "; | |
2588 | cerr<<" Processing protons ...\n"; | |
2589 | break; | |
6eb69b9a | 2590 | case 11: // electrons |
2591 | thisPdg= 11; | |
2592 | part="Electrons - "; | |
2593 | cerr<<" Processing electrons ...\n"; | |
2594 | break; | |
b2bca9d4 | 2595 | case 13: // muons |
2596 | thisPdg= 13; | |
2597 | part="Muons - "; | |
2598 | cerr<<" Processing muons ...\n"; | |
2599 | break; | |
6eb69b9a | 2600 | } |
2601 | ||
b2bca9d4 | 2602 | |
2603 | /* | |
2604 | // create a chain with compared tracks | |
2605 | TChain *cmptrkchain = new ("cmptrktree"); | |
2606 | cmptrkchain.Add("CovMatrix_AllEvts.root"); | |
2607 | //cmptrkchain.Add("CovMatrix_AllEvts_1.root"); | |
2608 | //cmptrkchain.Add("CovMatrix_AllEvts_2.root"); | |
2609 | //cmptrkchain.Add("CovMatrix_AllEvts_3.root"); | |
2610 | */ | |
2611 | ||
2612 | ||
2613 | TFile *infile = TFile::Open("CovMatrix_AllEvts.root"); | |
2614 | TTree *cmptrktree = (TTree*)infile->Get("cmptrktree"); | |
6eb69b9a | 2615 | |
2616 | COMPTRACK cmptrk; | |
b2bca9d4 | 2617 | cmptrktree->SetBranchAddress("comptracks",&cmptrk); |
2618 | Int_t entries = (Int_t)cmptrktree->GetEntries(); | |
6eb69b9a | 2619 | |
2620 | ||
2621 | Int_t pbin; | |
b2bca9d4 | 2622 | Int_t *n = new Int_t[fitbins]; |
2623 | Int_t *n00 = new Int_t[fitbins]; | |
2624 | Int_t *n11 = new Int_t[fitbins]; | |
2625 | Int_t *n20 = new Int_t[fitbins]; | |
2626 | Int_t *n22 = new Int_t[fitbins]; | |
2627 | Int_t *n31 = new Int_t[fitbins]; | |
2628 | Int_t *n33 = new Int_t[fitbins]; | |
2629 | Int_t *n40 = new Int_t[fitbins]; | |
2630 | Int_t *n42 = new Int_t[fitbins]; | |
2631 | Int_t *n44 = new Int_t[fitbins]; | |
2632 | Double_t *p = new Double_t[fitbins]; | |
2633 | Double_t *ep = new Double_t[fitbins]; | |
2634 | Double_t *mean00 = new Double_t[fitbins]; | |
2635 | Double_t *mean11 = new Double_t[fitbins]; | |
2636 | Double_t *mean20 = new Double_t[fitbins]; | |
2637 | Double_t *mean22 = new Double_t[fitbins]; | |
2638 | Double_t *mean31 = new Double_t[fitbins]; | |
2639 | Double_t *mean33 = new Double_t[fitbins]; | |
2640 | Double_t *mean40 = new Double_t[fitbins]; | |
2641 | Double_t *mean42 = new Double_t[fitbins]; | |
2642 | Double_t *mean44 = new Double_t[fitbins]; | |
2643 | Double_t *sigma00 = new Double_t[fitbins]; | |
2644 | Double_t *sigma11 = new Double_t[fitbins]; | |
2645 | Double_t *sigma20 = new Double_t[fitbins]; | |
2646 | Double_t *sigma22 = new Double_t[fitbins]; | |
2647 | Double_t *sigma31 = new Double_t[fitbins]; | |
2648 | Double_t *sigma33 = new Double_t[fitbins]; | |
2649 | Double_t *sigma40 = new Double_t[fitbins]; | |
2650 | Double_t *sigma42 = new Double_t[fitbins]; | |
2651 | Double_t *sigma44 = new Double_t[fitbins]; | |
2652 | Double_t *rmean = new Double_t[fitbins]; | |
2653 | Double_t *rsigma = new Double_t[fitbins]; | |
2654 | Double_t fitpar[3]; | |
6eb69b9a | 2655 | |
2656 | for(Int_t l=0; l<fitbins; l++) { | |
2657 | n[l]=1; | |
2658 | n00[l]=n11[l]=n20[l]=n22[l]=n31[l]=n33[l]=n40[l]=n42[l]=n44[l]=1; | |
2659 | p[l ]=ep[l]=0.; | |
2660 | mean00[l]=mean11[l]=mean20[l]=mean22[l]=mean31[l]=mean33[l]=mean40[l]=mean42[l]=mean44[l]=0.; | |
2661 | sigma00[l]=sigma11[l]=sigma20[l]=sigma22[l]=sigma31[l]=sigma33[l]=sigma40[l]=sigma42[l]=sigma44[l]=0.; | |
2662 | } | |
2663 | ||
2664 | // loop on chain entries for mean | |
2665 | for(Int_t l=0; l<entries; l++) { | |
b2bca9d4 | 2666 | cmptrktree->GetEvent(l); |
6eb69b9a | 2667 | if(TMath::Abs(cmptrk.pdg)!=thisPdg) continue; |
2668 | pbin = (Int_t)fDBgrid->GetBin(cmptrk.pt,cmptrk.eta)/fDBgrid->GetBinsEta(); | |
2669 | n[pbin]++; | |
2670 | p[pbin]+=cmptrk.p; | |
2671 | mean00[pbin]+=cmptrk.c00; | |
2672 | mean11[pbin]+=cmptrk.c11; | |
2673 | mean20[pbin]+=cmptrk.c20; | |
2674 | mean22[pbin]+=cmptrk.c22; | |
2675 | mean31[pbin]+=cmptrk.c31; | |
2676 | mean33[pbin]+=cmptrk.c33; | |
2677 | mean40[pbin]+=cmptrk.c40; | |
2678 | mean42[pbin]+=cmptrk.c42; | |
2679 | mean44[pbin]+=cmptrk.c44; | |
2680 | } // loop on chain entries | |
2681 | ||
2682 | for(Int_t l=0; l<fitbins; l++) { | |
2683 | p[l]/=n[l]; | |
2684 | mean00[l]/=n[l]; | |
2685 | mean11[l]/=n[l]; | |
2686 | mean20[l]/=n[l]; | |
2687 | mean22[l]/=n[l]; | |
2688 | mean31[l]/=n[l]; | |
2689 | mean33[l]/=n[l]; | |
2690 | mean40[l]/=n[l]; | |
2691 | mean42[l]/=n[l]; | |
2692 | mean44[l]/=n[l]; | |
2693 | } | |
2694 | ||
2695 | // loop on chain entries for sigma | |
2696 | for(Int_t l=0; l<entries; l++) { | |
b2bca9d4 | 2697 | cmptrktree->GetEvent(l); |
6eb69b9a | 2698 | if(TMath::Abs(cmptrk.pdg)!=thisPdg) continue; |
2699 | pbin = (Int_t)fDBgrid->GetBin(cmptrk.pt,cmptrk.eta)/fDBgrid->GetBinsEta(); | |
2700 | if(TMath::Abs(cmptrk.c00-mean00[pbin])<0.4*mean00[pbin]) { n00[pbin]++; | |
2701 | sigma00[pbin]+=(cmptrk.c00-mean00[pbin])*(cmptrk.c00-mean00[pbin]); } | |
2702 | if(TMath::Abs(cmptrk.c11-mean11[pbin])<0.4*mean11[pbin]) { n11[pbin]++; | |
2703 | sigma11[pbin]+=(cmptrk.c11-mean11[pbin])*(cmptrk.c11-mean11[pbin]); } | |
2704 | if(TMath::Abs(cmptrk.c20-mean20[pbin])<0.4*mean20[pbin]) { n20[pbin]++; | |
2705 | sigma20[pbin]+=(cmptrk.c20-mean20[pbin])*(cmptrk.c20-mean20[pbin]); } | |
2706 | if(TMath::Abs(cmptrk.c22-mean22[pbin])<0.4*mean22[pbin]) { n22[pbin]++; | |
2707 | sigma22[pbin]+=(cmptrk.c22-mean22[pbin])*(cmptrk.c22-mean22[pbin]); } | |
2708 | if(TMath::Abs(cmptrk.c31-mean31[pbin])<-0.4*mean31[pbin]) { n31[pbin]++; | |
2709 | sigma31[pbin]+=(cmptrk.c31-mean31[pbin])*(cmptrk.c31-mean31[pbin]); } | |
2710 | if(TMath::Abs(cmptrk.c33-mean33[pbin])<0.4*mean33[pbin]) { n33[pbin]++; | |
2711 | sigma33[pbin]+=(cmptrk.c33-mean33[pbin])*(cmptrk.c33-mean33[pbin]); } | |
2712 | if(TMath::Abs(cmptrk.c40-mean40[pbin])<0.4*mean40[pbin]) { n40[pbin]++; | |
2713 | sigma40[pbin]+=(cmptrk.c40-mean40[pbin])*(cmptrk.c40-mean40[pbin]); } | |
2714 | if(TMath::Abs(cmptrk.c42-mean42[pbin])<0.4*mean42[pbin]) { n42[pbin]++; | |
2715 | sigma42[pbin]+=(cmptrk.c42-mean42[pbin])*(cmptrk.c42-mean42[pbin]); } | |
2716 | if(TMath::Abs(cmptrk.c44-mean44[pbin])<0.4*mean44[pbin]) { n44[pbin]++; | |
2717 | sigma44[pbin]+=(cmptrk.c44-mean44[pbin])*(cmptrk.c44-mean44[pbin]); } | |
2718 | } // loop on chain entries | |
2719 | ||
2720 | for(Int_t l=0; l<fitbins; l++) { | |
2721 | sigma00[l] = TMath::Sqrt(sigma00[l]/n00[l]); | |
2722 | sigma11[l] = TMath::Sqrt(sigma11[l]/n11[l]); | |
2723 | sigma20[l] = TMath::Sqrt(sigma20[l]/n20[l]); | |
2724 | sigma22[l] = TMath::Sqrt(sigma22[l]/n22[l]); | |
2725 | sigma31[l] = TMath::Sqrt(sigma31[l]/n31[l]); | |
2726 | sigma33[l] = TMath::Sqrt(sigma33[l]/n33[l]); | |
2727 | sigma40[l] = TMath::Sqrt(sigma40[l]/n40[l]); | |
2728 | sigma42[l] = TMath::Sqrt(sigma42[l]/n42[l]); | |
2729 | sigma44[l] = TMath::Sqrt(sigma44[l]/n44[l]); | |
2730 | } | |
2731 | ||
2732 | ||
2733 | // Fit function | |
b2bca9d4 | 2734 | TF1 *func = new TF1("RegFunc",RegFunc,0.23,50.,3); |
2735 | func->SetParNames("A_meas","A_scatt","B"); | |
6eb69b9a | 2736 | |
2737 | // line to draw on the plots | |
2738 | TLine *lin = new TLine(-1,1,1.69,1); | |
2739 | lin->SetLineStyle(2); | |
2740 | lin->SetLineWidth(2); | |
2741 | ||
2742 | // matrix used to store fit results | |
b2bca9d4 | 2743 | TMatrixD fitRes(9,3); |
6eb69b9a | 2744 | |
2745 | // --- c00 --- | |
2746 | ||
2747 | // create the canvas | |
2748 | TCanvas *canv00 = new TCanvas("canv00","c00",0,0,700,900); | |
2749 | canv00->Divide(1,2); | |
2750 | // create the graph for cov matrix | |
2751 | TGraphErrors *gr00 = new TGraphErrors(fitbins,p,mean00,ep,sigma00); | |
2752 | TString title00("C(y,y)"); title00.Prepend(part); | |
2753 | TH2F *frame00 = new TH2F("frame00",title00.Data(),2,0.1,50,2,0,5e-3); | |
2754 | frame00->SetXTitle("p [GeV/c]"); | |
2755 | canv00->cd(1); gPad->SetLogx(); | |
2756 | frame00->Draw(); | |
2757 | gr00->Draw("P"); | |
2758 | // Sets initial values for parameters | |
b2bca9d4 | 2759 | func->SetParameters(1.6e-3,1.9e-4,1.5); |
6eb69b9a | 2760 | // Fit points in range defined by function |
b2bca9d4 | 2761 | gr00->Fit("RegFunc","R,Q"); |
6eb69b9a | 2762 | func->GetParameters(fitpar); |
b2bca9d4 | 2763 | for(Int_t i=0; i<3; i++) fitRes(0,i)=fitpar[i]; |
6eb69b9a | 2764 | for(Int_t l=0; l<fitbins; l++) { |
b2bca9d4 | 2765 | rmean[l] = mean00[l]/RegFunc(&p[l],fitpar); |
2766 | rsigma[l] = sigma00[l]/RegFunc(&p[l],fitpar); | |
6eb69b9a | 2767 | } |
2768 | // create the graph the regularized cov. matrix | |
2769 | TGraphErrors *gr00reg = new TGraphErrors(fitbins,p,rmean,ep,rsigma); | |
b2bca9d4 | 2770 | TString regtitle00("C(y,y)/(A_meas+A_scatt/p^{B})"); |
6eb69b9a | 2771 | regtitle00.Prepend(part); |
2772 | TH2F *frame00reg = new TH2F("frame00reg",regtitle00.Data(),2,0.1,50,2,0,2); | |
2773 | frame00reg->SetXTitle("p [GeV/c]"); | |
2774 | canv00->cd(2); gPad->SetLogx(); | |
2775 | frame00reg->Draw(); | |
2776 | gr00reg->Draw("P"); | |
2777 | lin->Draw("same"); | |
2778 | ||
2779 | ||
2780 | // --- c11 --- | |
b2bca9d4 | 2781 | |
6eb69b9a | 2782 | // create the canvas |
2783 | TCanvas *canv11 = new TCanvas("canv11","c11",0,0,700,900); | |
2784 | canv11->Divide(1,2); | |
2785 | // create the graph for cov matrix | |
2786 | TGraphErrors *gr11 = new TGraphErrors(fitbins,p,mean11,ep,sigma11); | |
2787 | TString title11("C(z,z)"); title11.Prepend(part); | |
2788 | TH2F *frame11 = new TH2F("frame11",title11.Data(),2,0.1,50,2,0,6e-3); | |
2789 | frame11->SetXTitle("p [GeV/c]"); | |
2790 | canv11->cd(1); gPad->SetLogx(); | |
2791 | frame11->Draw(); | |
2792 | gr11->Draw("P"); | |
2793 | // Sets initial values for parameters | |
b2bca9d4 | 2794 | func->SetParameters(1.2e-3,8.1e-4,1.); |
6eb69b9a | 2795 | // Fit points in range defined by function |
b2bca9d4 | 2796 | gr11->Fit("RegFunc","R,Q"); |
6eb69b9a | 2797 | func->GetParameters(fitpar); |
b2bca9d4 | 2798 | for(Int_t i=0; i<3; i++) fitRes(1,i)=fitpar[i]; |
6eb69b9a | 2799 | for(Int_t l=0; l<fitbins; l++) { |
b2bca9d4 | 2800 | rmean[l] = mean11[l]/RegFunc(&p[l],fitpar); |
2801 | rsigma[l] = sigma11[l]/RegFunc(&p[l],fitpar); | |
6eb69b9a | 2802 | } |
2803 | // create the graph the regularized cov. matrix | |
2804 | TGraphErrors *gr11reg = new TGraphErrors(fitbins,p,rmean,ep,rsigma); | |
b2bca9d4 | 2805 | TString regtitle11("C(z,z)/(A_meas+A_scatt/p^{B})"); |
6eb69b9a | 2806 | regtitle11.Prepend(part); |
2807 | TH2F *frame11reg = new TH2F("frame11reg",regtitle11.Data(),2,0.1,50,2,0,2); | |
2808 | frame11reg->SetXTitle("p [GeV/c]"); | |
2809 | canv11->cd(2); gPad->SetLogx(); | |
2810 | frame11reg->Draw(); | |
2811 | gr11reg->Draw("P"); | |
2812 | lin->Draw("same"); | |
2813 | ||
2814 | ||
2815 | // --- c20 --- | |
2816 | ||
2817 | // create the canvas | |
2818 | TCanvas *canv20 = new TCanvas("canv20","c20",0,0,700,900); | |
2819 | canv20->Divide(1,2); | |
2820 | // create the graph for cov matrix | |
2821 | TGraphErrors *gr20 = new TGraphErrors(fitbins,p,mean20,ep,sigma20); | |
2822 | TString title20("C(#eta, y)"); title20.Prepend(part); | |
2823 | TH2F *frame20 = new TH2F("frame20",title20.Data(),2,0.1,50,2,0,2.5e-4); | |
2824 | frame20->SetXTitle("p [GeV/c]"); | |
2825 | canv20->cd(1); gPad->SetLogx(); | |
2826 | frame20->Draw(); | |
2827 | gr20->Draw("P"); | |
2828 | // Sets initial values for parameters | |
b2bca9d4 | 2829 | func->SetParameters(7.3e-5,1.2e-5,1.5); |
6eb69b9a | 2830 | // Fit points in range defined by function |
b2bca9d4 | 2831 | gr20->Fit("RegFunc","R,Q"); |
6eb69b9a | 2832 | func->GetParameters(fitpar); |
b2bca9d4 | 2833 | for(Int_t i=0; i<3; i++) fitRes(2,i)=fitpar[i]; |
6eb69b9a | 2834 | for(Int_t l=0; l<fitbins; l++) { |
b2bca9d4 | 2835 | rmean[l] = mean20[l]/RegFunc(&p[l],fitpar); |
2836 | rsigma[l] = sigma20[l]/RegFunc(&p[l],fitpar); | |
6eb69b9a | 2837 | } |
2838 | // create the graph the regularized cov. matrix | |
2839 | TGraphErrors *gr20reg = new TGraphErrors(fitbins,p,rmean,ep,rsigma); | |
b2bca9d4 | 2840 | TString regtitle20("C(#eta, y)/(A_meas+A_scatt/p^{B})"); |
6eb69b9a | 2841 | regtitle20.Prepend(part); |
2842 | TH2F *frame20reg = new TH2F("frame20reg",regtitle20.Data(),2,0.1,50,2,0,2); | |
2843 | frame20reg->SetXTitle("p [GeV/c]"); | |
2844 | canv20->cd(2); gPad->SetLogx(); | |
2845 | frame20reg->Draw(); | |
2846 | gr20reg->Draw("P"); | |
2847 | lin->Draw("same"); | |
2848 | ||
2849 | ||
2850 | // --- c22 --- | |
2851 | ||
2852 | // create the canvas | |
2853 | TCanvas *canv22 = new TCanvas("canv22","c22",0,0,700,900); | |
2854 | canv22->Divide(1,2); | |
2855 | // create the graph for cov matrix | |
2856 | TGraphErrors *gr22 = new TGraphErrors(fitbins,p,mean22,ep,sigma22); | |
2857 | TString title22("C(#eta, #eta)"); title22.Prepend(part); | |
2858 | TH2F *frame22 = new TH2F("frame22",title22.Data(),2,0.1,50,2,0,3e-5); | |
2859 | frame22->SetXTitle("p [GeV/c]"); | |
2860 | canv22->cd(1); gPad->SetLogx(); | |
2861 | frame22->Draw(); | |
2862 | gr22->Draw("P"); | |
2863 | // Sets initial values for parameters | |
b2bca9d4 | 2864 | func->SetParameters(5.2e-6,1.1e-6,2.); |
6eb69b9a | 2865 | // Fit points in range defined by function |
b2bca9d4 | 2866 | gr22->Fit("RegFunc","R,Q"); |
6eb69b9a | 2867 | func->GetParameters(fitpar); |
b2bca9d4 | 2868 | for(Int_t i=0; i<3; i++) fitRes(3,i)=fitpar[i]; |
6eb69b9a | 2869 | for(Int_t l=0; l<fitbins; l++) { |
b2bca9d4 | 2870 | rmean[l] = mean22[l]/RegFunc(&p[l],fitpar); |
2871 | rsigma[l] = sigma22[l]/RegFunc(&p[l],fitpar); | |
6eb69b9a | 2872 | } |
2873 | // create the graph the regularized cov. matrix | |
2874 | TGraphErrors *gr22reg = new TGraphErrors(fitbins,p,rmean,ep,rsigma); | |
b2bca9d4 | 2875 | TString regtitle22("C(#eta, #eta)/(A_meas+A_scatt/p^{B})"); |
6eb69b9a | 2876 | regtitle22.Prepend(part); |
2877 | TH2F *frame22reg = new TH2F("frame22reg",regtitle22.Data(),2,0.1,50,2,0,2); | |
2878 | frame22reg->SetXTitle("p [GeV/c]"); | |
2879 | canv22->cd(2); gPad->SetLogx(); | |
2880 | frame22reg->Draw(); | |
2881 | gr22reg->Draw("P"); | |
2882 | lin->Draw("same"); | |
2883 | ||
2884 | ||
2885 | // --- c31 --- | |
2886 | ||
2887 | // create the canvas | |
2888 | TCanvas *canv31 = new TCanvas("canv31","c31",0,0,700,900); | |
2889 | canv31->Divide(1,2); | |
2890 | // create the graph for cov matrix | |
2891 | TGraphErrors *gr31 = new TGraphErrors(fitbins,p,mean31,ep,sigma31); | |
2892 | TString title31("C(tg #lambda,z)"); title31.Prepend(part); | |
2893 | TH2F *frame31 = new TH2F("frame31",title31.Data(),2,0.1,50,2,-2e-4,0); | |
2894 | frame31->SetXTitle("p [GeV/c]"); | |
2895 | canv31->cd(1); gPad->SetLogx(); | |
2896 | frame31->Draw(); | |
2897 | gr31->Draw("P"); | |
2898 | // Sets initial values for parameters | |
b2bca9d4 | 2899 | func->SetParameters(-1.2e-5,-1.2e-5,1.5); |
6eb69b9a | 2900 | // Fit points in range defined by function |
b2bca9d4 | 2901 | gr31->Fit("RegFunc","R,Q"); |
6eb69b9a | 2902 | func->GetParameters(fitpar); |
b2bca9d4 | 2903 | for(Int_t i=0; i<3; i++) fitRes(4,i)=fitpar[i]; |
6eb69b9a | 2904 | for(Int_t l=0; l<fitbins; l++) { |
b2bca9d4 | 2905 | rmean[l] = mean31[l]/RegFunc(&p[l],fitpar); |
2906 | rsigma[l] = -sigma31[l]/RegFunc(&p[l],fitpar); | |
6eb69b9a | 2907 | } |
2908 | // create the graph the regularized cov. matrix | |
2909 | TGraphErrors *gr31reg = new TGraphErrors(fitbins,p,rmean,ep,rsigma); | |
b2bca9d4 | 2910 | TString regtitle31("C(tg #lambda,z)/(A_meas+A_scatt/p^{B})"); |
6eb69b9a | 2911 | regtitle31.Prepend(part); |
2912 | TH2F *frame31reg = new TH2F("frame31reg",regtitle31.Data(),2,0.1,50,2,0,2); | |
2913 | frame31reg->SetXTitle("p [GeV/c]"); | |
2914 | canv31->cd(2); gPad->SetLogx(); | |
2915 | frame31reg->Draw(); | |
2916 | gr31reg->Draw("P"); | |
2917 | lin->Draw("same"); | |
2918 | ||
2919 | ||
2920 | // --- c33 --- | |
2921 | ||
2922 | // create the canvas | |
2923 | TCanvas *canv33 = new TCanvas("canv33","c33",0,0,700,900); | |
2924 | canv33->Divide(1,2); | |
2925 | // create the graph for cov matrix | |
2926 | TGraphErrors *gr33 = new TGraphErrors(fitbins,p,mean33,ep,sigma33); | |
2927 | TString title33("C(tg #lambda,tg #lambda)"); title33.Prepend(part); | |
2928 | TH2F *frame33 = new TH2F("frame33",title33.Data(),2,0.1,50,2,0,1e-5); | |
2929 | frame33->SetXTitle("p [GeV/c]"); | |
2930 | canv33->cd(1); gPad->SetLogx(); | |
2931 | frame33->Draw(); | |
2932 | gr33->Draw("P"); | |
2933 | // Sets initial values for parameters | |
b2bca9d4 | 2934 | func->SetParameters(1.3e-7,4.6e-7,1.7); |
6eb69b9a | 2935 | // Fit points in range defined by function |
b2bca9d4 | 2936 | gr33->Fit("RegFunc","R,Q"); |
6eb69b9a | 2937 | func->GetParameters(fitpar); |
b2bca9d4 | 2938 | for(Int_t i=0; i<3; i++) fitRes(5,i)=fitpar[i]; |
6eb69b9a | 2939 | for(Int_t l=0; l<fitbins; l++) { |
b2bca9d4 | 2940 | rmean[l] = mean33[l]/RegFunc(&p[l],fitpar); |
2941 | rsigma[l] = sigma33[l]/RegFunc(&p[l],fitpar); | |
6eb69b9a | 2942 | } |
2943 | // create the graph the regularized cov. matrix | |
2944 | TGraphErrors *gr33reg = new TGraphErrors(fitbins,p,rmean,ep,rsigma); | |
b2bca9d4 | 2945 | TString regtitle33("C(tg #lambda,tg #lambda)/(A_meas+A_scatt/p^{B})"); |
6eb69b9a | 2946 | regtitle33.Prepend(part); |
2947 | TH2F *frame33reg = new TH2F("frame33reg",regtitle33.Data(),2,0.1,50,2,0,2); | |
2948 | frame33reg->SetXTitle("p [GeV/c]"); | |
2949 | canv33->cd(2); gPad->SetLogx(); | |
2950 | frame33reg->Draw(); | |
2951 | gr33reg->Draw("P"); | |
2952 | lin->Draw("same"); | |
2953 | ||
2954 | ||
2955 | // --- c40 --- | |
2956 | ||
2957 | // create the canvas | |
2958 | TCanvas *canv40 = new TCanvas("canv40","c40",0,0,700,900); | |
2959 | canv40->Divide(1,2); | |
2960 | // create the graph for cov matrix | |
2961 | TGraphErrors *gr40 = new TGraphErrors(fitbins,p,mean40,ep,sigma40); | |
2962 | TString title40("C(C,y)"); title40.Prepend(part); | |
2963 | TH2F *frame40 = new TH2F("frame40",title40.Data(),2,0.1,50,2,0,1e-6); | |
2964 | frame40->SetXTitle("p [GeV/c]"); | |
2965 | canv40->cd(1); gPad->SetLogx(); | |
2966 | frame40->Draw(); | |
2967 | gr40->Draw("P"); | |
2968 | // Sets initial values for parameters | |
b2bca9d4 | 2969 | func->SetParameters(4.e-7,4.4e-8,1.5); |
6eb69b9a | 2970 | // Fit points in range defined by function |
b2bca9d4 | 2971 | gr40->Fit("RegFunc","R,Q"); |
6eb69b9a | 2972 | func->GetParameters(fitpar); |
b2bca9d4 | 2973 | for(Int_t i=0; i<3; i++) fitRes(6,i)=fitpar[i]; |
6eb69b9a | 2974 | for(Int_t l=0; l<fitbins; l++) { |
b2bca9d4 | 2975 | rmean[l] = mean40[l]/RegFunc(&p[l],fitpar); |
2976 | rsigma[l] = sigma40[l]/RegFunc(&p[l],fitpar); | |
6eb69b9a | 2977 | } |
2978 | // create the graph the regularized cov. matrix | |
2979 | TGraphErrors *gr40reg = new TGraphErrors(fitbins,p,rmean,ep,rsigma); | |
b2bca9d4 | 2980 | TString regtitle40("C(C,y)/(A_meas+A_scatt/p^{B})"); |
6eb69b9a | 2981 | regtitle40.Prepend(part); |
2982 | TH2F *frame40reg = new TH2F("frame40reg",regtitle40.Data(),2,0.1,50,2,0,2); | |
2983 | frame40reg->SetXTitle("p [GeV/c]"); | |
2984 | canv40->cd(2); gPad->SetLogx(); | |
2985 | frame40reg->Draw(); | |
2986 | gr40reg->Draw("P"); | |
2987 | lin->Draw("same"); | |
2988 | ||
2989 | ||
2990 | // --- c42 --- | |
2991 | ||
2992 | // create the canvas | |
2993 | TCanvas *canv42 = new TCanvas("canv42","c42",0,0,700,900); | |
2994 | canv42->Divide(1,2); | |
2995 | // create the graph for cov matrix | |
2996 | TGraphErrors *gr42 = new TGraphErrors(fitbins,p,mean42,ep,sigma42); | |
2997 | TString title42("C(C, #eta)"); title42.Prepend(part); | |
2998 | TH2F *frame42 = new TH2F("frame42",title42.Data(),2,0.1,50,2,0,2.2e-7); | |
2999 | frame42->SetXTitle("p [GeV/c]"); | |
3000 | canv42->cd(1); gPad->SetLogx(); | |
3001 | frame42->Draw(); | |
3002 | gr42->Draw("P"); | |
3003 | // Sets initial values for parameters | |
b2bca9d4 | 3004 | func->SetParameters(3.e-8,8.2e-9,2.); |
6eb69b9a | 3005 | // Fit points in range defined by function |
b2bca9d4 | 3006 | gr42->Fit("RegFunc","R,Q"); |
6eb69b9a | 3007 | func->GetParameters(fitpar); |
b2bca9d4 | 3008 | for(Int_t i=0; i<3; i++) fitRes(7,i)=fitpar[i]; |
6eb69b9a | 3009 | for(Int_t l=0; l<fitbins; l++) { |
b2bca9d4 | 3010 | rmean[l] = mean42[l]/RegFunc(&p[l],fitpar); |
3011 | rsigma[l] = sigma42[l]/RegFunc(&p[l],fitpar); | |
6eb69b9a | 3012 | } |
3013 | // create the graph the regularized cov. matrix | |
3014 | TGraphErrors *gr42reg = new TGraphErrors(fitbins,p,rmean,ep,rsigma); | |
b2bca9d4 | 3015 | TString regtitle42("C(C, #eta)/(A_meas+A_scatt/p^{B})"); |
6eb69b9a | 3016 | regtitle42.Prepend(part); |
3017 | TH2F *frame42reg = new TH2F("frame42reg",regtitle42.Data(),2,0.1,50,2,0,2); | |
3018 | frame42reg->SetXTitle("p [GeV/c]"); | |
3019 | canv42->cd(2); gPad->SetLogx(); | |
3020 | frame42reg->Draw(); | |
3021 | gr42reg->Draw("P"); | |
3022 | lin->Draw("same"); | |
3023 | ||
3024 | ||
3025 | // --- c44 --- | |
3026 | ||
3027 | // create the canvas | |
3028 | TCanvas *canv44 = new TCanvas("canv44","c44",0,0,700,900); | |
3029 | canv44->Divide(1,2); | |
3030 | // create the graph for cov matrix | |
3031 | TGraphErrors *gr44 = new TGraphErrors(fitbins,p,mean44,ep,sigma44); | |
3032 | TString title44("C(C,C)"); title44.Prepend(part); | |
3033 | TH2F *frame44 = new TH2F("frame44",title44.Data(),2,0.1,50,2,0,2e-9); | |
3034 | frame44->SetXTitle("p [GeV/c]"); | |
3035 | canv44->cd(1); gPad->SetLogx(); | |
3036 | frame44->Draw(); | |
3037 | gr44->Draw("P"); | |
3038 | // Sets initial values for parameters | |
b2bca9d4 | 3039 | func->SetParameters(1.8e-10,5.8e-11,2.); |
6eb69b9a | 3040 | // Fit points in range defined by function |
b2bca9d4 | 3041 | gr44->Fit("RegFunc","R,Q"); |
6eb69b9a | 3042 | func->GetParameters(fitpar); |
b2bca9d4 | 3043 | for(Int_t i=0; i<3; i++) fitRes(8,i)=fitpar[i]; |
6eb69b9a | 3044 | for(Int_t l=0; l<fitbins; l++) { |
b2bca9d4 | 3045 | rmean[l] = mean44[l]/RegFunc(&p[l],fitpar); |
3046 | rsigma[l] = sigma44[l]/RegFunc(&p[l],fitpar); | |
6eb69b9a | 3047 | } |
3048 | // create the graph the regularized cov. matrix | |
3049 | TGraphErrors *gr44reg = new TGraphErrors(fitbins,p,rmean,ep,rsigma); | |
b2bca9d4 | 3050 | TString regtitle44("C(C,C)/(A_meas+A_scatt/p^{B})"); |
6eb69b9a | 3051 | regtitle44.Prepend(part); |
3052 | TH2F *frame44reg = new TH2F("frame44reg",regtitle44.Data(),2,0.1,50,2,0,2); | |
3053 | frame44reg->SetXTitle("p [GeV/c]"); | |
3054 | canv44->cd(2); gPad->SetLogx(); | |
3055 | frame44reg->Draw(); | |
3056 | gr44reg->Draw("P"); | |
3057 | lin->Draw("same"); | |
3058 | ||
3059 | ||
3060 | ||
3061 | ||
3062 | switch (pdg) { | |
3063 | case 211: | |
3064 | new(&fRegParPi) TMatrixD(fitRes); | |
3065 | break; | |
3066 | case 321: | |
3067 | new(&fRegParKa) TMatrixD(fitRes); | |
3068 | break; | |
b2bca9d4 | 3069 | case 2212: |
3070 | new(&fRegParPr) TMatrixD(fitRes); | |
3071 | break; | |
6eb69b9a | 3072 | case 11: |
3073 | new(&fRegParEl) TMatrixD(fitRes); | |
3074 | break; | |
b2bca9d4 | 3075 | case 13: |
3076 | new(&fRegParMu) TMatrixD(fitRes); | |
3077 | break; | |
6eb69b9a | 3078 | } |
3079 | ||
3080 | // write fit parameters to file | |
3081 | WriteRegParams(outName,pdg); | |
3082 | ||
c15d4c37 | 3083 | delete [] n; |
3084 | delete [] n00; | |
3085 | delete [] n11; | |
3086 | delete [] n20; | |
3087 | delete [] n22; | |
3088 | delete [] n31; | |
3089 | delete [] n33; | |
3090 | delete [] n40; | |
3091 | delete [] n42; | |
3092 | delete [] n44; | |
3093 | delete [] p; | |
3094 | delete [] ep; | |
3095 | delete [] mean00; | |
3096 | delete [] mean11; | |
3097 | delete [] mean20; | |
3098 | delete [] mean22; | |
3099 | delete [] mean31; | |
3100 | delete [] mean33; | |
3101 | delete [] mean40; | |
3102 | delete [] mean42; | |
3103 | delete [] mean44; | |
3104 | delete [] sigma00; | |
3105 | delete [] sigma11; | |
3106 | delete [] sigma20; | |
3107 | delete [] sigma22; | |
3108 | delete [] sigma31; | |
3109 | delete [] sigma33; | |
3110 | delete [] sigma40; | |
3111 | delete [] sigma42; | |
3112 | delete [] sigma44; | |
3113 | delete [] rmean; | |
3114 | delete [] rsigma; | |
3115 | ||
6eb69b9a | 3116 | return; |
3117 | } | |
3118 | //----------------------------------------------------------------------------- | |
3119 | Bool_t AliTPCtrackerParam::SelectedTrack(Double_t pt,Double_t eta) const { | |
3120 | //----------------------------------------------------------------------------- | |
3121 | // This function makes a selection according to TPC tracking efficiency | |
3122 | //----------------------------------------------------------------------------- | |
3123 | ||
3124 | Double_t eff=0.; | |
3125 | ||
3126 | eff = fEff->GetValueAt(pt,eta); | |
3127 | ||
3128 | if(gRandom->Rndm() < eff) return kTRUE; | |
3129 | ||
3130 | return kFALSE; | |
3131 | } | |
3132 | //----------------------------------------------------------------------------- | |
3133 | void AliTPCtrackerParam::SetParticle(Int_t pdg) { | |
3134 | //----------------------------------------------------------------------------- | |
3135 | // This function set efficiencies, pulls, etc... for the particle | |
3136 | // specie of the current particle | |
3137 | //----------------------------------------------------------------------------- | |
3138 | ||
3139 | switch (pdg) { | |
3140 | case 211: | |
3141 | fDBgrid = &fDBgridPi; | |
3142 | fEff = &fEffPi; | |
3143 | fPulls = fPullsPi; | |
3144 | fRegPar = &fRegParPi; | |
3145 | fdEdxMean = &fdEdxMeanPi; | |
3146 | fdEdxRMS = &fdEdxRMSPi; | |
3147 | break; | |
3148 | case 321: | |
3149 | fDBgrid = &fDBgridKa; | |
3150 | fEff = &fEffKa; | |
3151 | fPulls = fPullsKa; | |
3152 | fRegPar = &fRegParKa; | |
3153 | fdEdxMean = &fdEdxMeanKa; | |
3154 | fdEdxRMS = &fdEdxRMSKa; | |
3155 | break; | |
3156 | case 2212: | |
b2bca9d4 | 3157 | fDBgrid = &fDBgridPr; |
6eb69b9a | 3158 | fEff = &fEffPr; |
b2bca9d4 | 3159 | fPulls = fPullsPr; |
3160 | fRegPar = &fRegParPr; | |
6eb69b9a | 3161 | fdEdxMean = &fdEdxMeanPr; |
3162 | fdEdxRMS = &fdEdxRMSPr; | |
3163 | break; | |
3164 | case 11: | |
3165 | fDBgrid = &fDBgridEl; | |
3166 | fEff = &fEffEl; | |
3167 | fPulls = fPullsEl; | |
3168 | fRegPar = &fRegParEl; | |
3169 | fdEdxMean = &fdEdxMeanEl; | |
3170 | fdEdxRMS = &fdEdxRMSEl; | |
3171 | break; | |
3172 | case 13: | |
b2bca9d4 | 3173 | fDBgrid = &fDBgridMu; |
6eb69b9a | 3174 | fEff = &fEffMu; |
b2bca9d4 | 3175 | fPulls = fPullsMu; |
3176 | fRegPar = &fRegParMu; | |
3177 | fdEdxMean = &fdEdxMeanMu; | |
3178 | fdEdxRMS = &fdEdxRMSMu; | |
6eb69b9a | 3179 | break; |
3180 | } | |
3181 | ||
3182 | return; | |
3183 | } | |
3184 | //----------------------------------------------------------------------------- | |
3185 | void AliTPCtrackerParam::SmearTrack(Double_t *xx,Double_t *xxsm,TMatrixD cov) | |
3186 | const { | |
3187 | //----------------------------------------------------------------------------- | |
e130146c | 3188 | // This function smears track parameters according to streched cov. matrix |
6eb69b9a | 3189 | //----------------------------------------------------------------------------- |
70521312 | 3190 | AliGausCorr *corgen = new AliGausCorr(cov,5); |
3191 | TArrayD corr(5); | |
3192 | corgen->GetGaussN(corr); | |
3193 | delete corgen; | |
3194 | corgen = 0; | |
6eb67451 | 3195 | |
3196 | for(Int_t l=0;l<5;l++) { | |
3197 | xxsm[l] = xx[l]+corr[l]; | |
3198 | } | |
3199 | ||
3200 | return; | |
3201 | } | |
6eb69b9a | 3202 | //----------------------------------------------------------------------------- |
3203 | Int_t AliTPCtrackerParam::WritedEdx(const Char_t *outName,Int_t pdg) { | |
3204 | //----------------------------------------------------------------------------- | |
3205 | // This function writes the dEdx parameters to the DB | |
3206 | //----------------------------------------------------------------------------- | |
3207 | ||
3208 | Option_t *opt; | |
3209 | Char_t *dirName="Pions"; | |
3210 | Char_t *meanName="dEdxMeanPi"; | |
3211 | Char_t *RMSName="dEdxRMSPi"; | |
3212 | ||
3213 | SetParticle(pdg); | |
3214 | ||
3215 | if(gSystem->AccessPathName(outName,kFileExists)) { opt="recreate"; | |
3216 | } else { opt="update"; } | |
3217 | ||
3218 | switch (pdg) { | |
3219 | case 211: | |
3220 | dirName="Pions"; | |
3221 | meanName="dEdxMeanPi"; | |
3222 | RMSName="dEdxRMSPi"; | |
3223 | break; | |
3224 | case 321: | |
3225 | dirName="Kaons"; | |
3226 | meanName="dEdxMeanKa"; | |
3227 | RMSName="dEdxRMSKa"; | |
3228 | break; | |
3229 | case 2212: | |
3230 | dirName="Protons"; | |
3231 | meanName="dEdxMeanPr"; | |
3232 | RMSName="dEdxRMSPr"; | |
3233 | break; | |
3234 | case 11: | |
3235 | dirName="Electrons"; | |
3236 | meanName="dEdxMeanEl"; | |
3237 | RMSName="dEdxRMSEl"; | |
3238 | break; | |
b2bca9d4 | 3239 | case 13: |
3240 | dirName="Muons"; | |
3241 | meanName="dEdxMeanMu"; | |
3242 | RMSName="dEdxRMSMu"; | |
3243 | break; | |
6eb69b9a | 3244 | } |
6eb67451 | 3245 | |
6eb69b9a | 3246 | TFile *outFile = new TFile(outName,opt); |
3247 | if(!gDirectory->cd("/dEdx")) { | |
3248 | outFile->mkdir("dEdx"); | |
3249 | gDirectory->cd("/dEdx"); | |
3250 | } | |
3251 | TDirectory *dir2 = gDirectory->mkdir(dirName); | |
3252 | dir2->cd(); | |
3253 | fdEdxMean->SetName(meanName); fdEdxMean->Write(); | |
3254 | fdEdxRMS->SetName(RMSName); fdEdxRMS->Write(); | |
6eb67451 | 3255 | |
6eb69b9a | 3256 | outFile->Close(); |
3257 | delete outFile; | |
6eb67451 | 3258 | |
6eb67451 | 3259 | |
6eb69b9a | 3260 | return 1; |
3261 | } | |
3262 | //----------------------------------------------------------------------------- | |
3263 | Int_t AliTPCtrackerParam::WriteEffs(const Char_t *outName) { | |
3264 | //----------------------------------------------------------------------------- | |
3265 | // This function writes the TPC efficiencies to the DB | |
3266 | //----------------------------------------------------------------------------- | |
3267 | ||
3268 | ||
3269 | Option_t *opt; | |
3270 | ||
3271 | if(gSystem->AccessPathName(outName,kFileExists)) { opt="recreate"; | |
3272 | } else { opt="update"; } | |
3273 | ||
3274 | TFile *outFile = new TFile(outName,opt); | |
3275 | ||
3276 | outFile->mkdir("Efficiencies"); | |
3277 | gDirectory->cd("/Efficiencies"); | |
3278 | gDirectory->mkdir("Pions"); | |
3279 | gDirectory->mkdir("Kaons"); | |
3280 | gDirectory->mkdir("Protons"); | |
3281 | gDirectory->mkdir("Electrons"); | |
3282 | gDirectory->mkdir("Muons"); | |
3283 | ||
3284 | gDirectory->cd("/Efficiencies/Pions"); | |
3285 | fEffPi.SetName("EffPi"); | |
3286 | fEffPi.Write(); | |
3287 | gDirectory->cd("/Efficiencies/Kaons"); | |
3288 | fEffKa.SetName("EffKa"); | |
3289 | fEffKa.Write(); | |
3290 | gDirectory->cd("/Efficiencies/Protons"); | |
3291 | fEffPr.SetName("EffPr"); | |
3292 | fEffPr.Write(); | |
3293 | gDirectory->cd("/Efficiencies/Electrons"); | |
3294 | fEffEl.SetName("EffEl"); | |
3295 | fEffEl.Write(); | |
3296 | gDirectory->cd("/Efficiencies/Muons"); | |
3297 | fEffMu.SetName("EffMu"); | |
3298 | fEffMu.Write(); | |
3299 | ||
3300 | outFile->Close(); | |
3301 | ||
3302 | delete outFile; | |
3303 | ||
3304 | return 1; | |
3305 | } | |
3306 | //----------------------------------------------------------------------------- | |
3307 | Int_t AliTPCtrackerParam::WritePulls(const Char_t *outName) { | |
3308 | //----------------------------------------------------------------------------- | |
3309 | // This function writes the pulls to the DB | |
3310 | //----------------------------------------------------------------------------- | |
3311 | ||
3312 | Option_t *opt; | |
3313 | ||
3314 | if(gSystem->AccessPathName(outName,kFileExists)) { opt="recreate"; | |
3315 | } else { opt="update"; } | |
3316 | ||
3317 | TFile *outFile = new TFile(outName,opt); | |
3318 | ||
3319 | outFile->mkdir("Pulls"); | |
3320 | gDirectory->cd("/Pulls"); | |
3321 | gDirectory->mkdir("Pions"); | |
3322 | gDirectory->mkdir("Kaons"); | |
b2bca9d4 | 3323 | gDirectory->mkdir("Protons"); |
6eb69b9a | 3324 | gDirectory->mkdir("Electrons"); |
b2bca9d4 | 3325 | gDirectory->mkdir("Muons"); |
6eb69b9a | 3326 | |
3327 | for(Int_t i=0;i<5;i++) { | |
3328 | TString pi("PullsPi_"); pi+=i; | |
3329 | TString ka("PullsKa_"); ka+=i; | |
b2bca9d4 | 3330 | TString pr("PullsPr_"); pr+=i; |
6eb69b9a | 3331 | TString el("PullsEl_"); el+=i; |
b2bca9d4 | 3332 | TString mu("PullsMu_"); mu+=i; |
6eb69b9a | 3333 | fPullsPi[i].SetName(pi.Data()); |
3334 | fPullsKa[i].SetName(ka.Data()); | |
b2bca9d4 | 3335 | fPullsPr[i].SetName(pr.Data()); |
6eb69b9a | 3336 | fPullsEl[i].SetName(el.Data()); |
b2bca9d4 | 3337 | fPullsMu[i].SetName(mu.Data()); |
6eb69b9a | 3338 | gDirectory->cd("/Pulls/Pions"); |
3339 | fPullsPi[i].Write(); | |
3340 | gDirectory->cd("/Pulls/Kaons"); | |
3341 | fPullsKa[i].Write(); | |
b2bca9d4 | 3342 | gDirectory->cd("/Pulls/Protons"); |
3343 | fPullsPr[i].Write(); | |
6eb69b9a | 3344 | gDirectory->cd("/Pulls/Electrons"); |
3345 | fPullsEl[i].Write(); | |
b2bca9d4 | 3346 | gDirectory->cd("/Pulls/Muons"); |
3347 | fPullsMu[i].Write(); | |
6eb69b9a | 3348 | } |
3349 | outFile->Close(); | |
3350 | delete outFile; | |
6eb67451 | 3351 | |
6eb69b9a | 3352 | return 1; |
3353 | } | |
3354 | //----------------------------------------------------------------------------- | |
3355 | Int_t AliTPCtrackerParam::WriteRegParams(const Char_t *outName,Int_t pdg) { | |
3356 | //----------------------------------------------------------------------------- | |
3357 | // This function writes the regularization parameters to the DB | |
3358 | //----------------------------------------------------------------------------- | |
3359 | ||
3360 | Option_t *opt; | |
3361 | Char_t *dirName="Pions"; | |
3362 | Char_t *keyName="RegPions"; | |
3363 | ||
3364 | SetParticle(pdg); | |
3365 | ||
3366 | if(gSystem->AccessPathName(outName,kFileExists)) { opt="recreate"; | |
3367 | } else { opt="update"; } | |
3368 | ||
3369 | switch (pdg) { | |
3370 | case 211: | |
3371 | dirName="Pions"; | |
3372 | keyName="RegPions"; | |
3373 | break; | |
3374 | case 321: | |
3375 | dirName="Kaons"; | |
3376 | keyName="RegKaons"; | |
3377 | break; | |
b2bca9d4 | 3378 | case 2212: |
3379 | dirName="Protons"; | |
3380 | keyName="RegProtons"; | |
3381 | break; | |
6eb69b9a | 3382 | case 11: |
3383 | dirName="Electrons"; | |
3384 | keyName="RegElectrons"; | |
3385 | break; | |
b2bca9d4 | 3386 | case 13: |
3387 | dirName="Muons"; | |
3388 | keyName="RegMuons"; | |
3389 | break; | |
6eb69b9a | 3390 | } |
70521312 | 3391 | |
6eb69b9a | 3392 | TFile *outFile = new TFile(outName,opt); |
3393 | if(!gDirectory->cd("/RegParams")) { | |
3394 | outFile->mkdir("RegParams"); | |
3395 | gDirectory->cd("/RegParams"); | |
6eb67451 | 3396 | } |
6eb69b9a | 3397 | TDirectory *dir2 = gDirectory->mkdir(dirName); |
3398 | dir2->cd(); | |
3399 | fRegPar->Write(keyName); | |
6eb67451 | 3400 | |
6eb69b9a | 3401 | outFile->Close(); |
3402 | delete outFile; | |
70521312 | 3403 | |
3404 | ||
6eb69b9a | 3405 | return 1; |
6eb67451 | 3406 | } |
b2bca9d4 | 3407 | //----------------------------------------------------------------------------- |
3408 | //***************************************************************************** | |
3409 | //----------------------------------------------------------------------------- | |
3410 | ||
3411 | ||
3412 | ||
3413 | ||
3414 | ||
3415 | ||
3416 | ||
3417 | ||
3418 | ||
3419 | ||
3420 | ||
3421 | ||
3422 | ||
3423 | ||
6eb67451 | 3424 | |
3425 | ||
3426 | ||
3427 | ||
3428 | ||
3429 | ||
3430 | ||
3431 | ||
3432 |