]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TPC/AliTPCcalibGainMult.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Adding additional variable
[u/mrichter/AliRoot.git] / TPC / AliTPCcalibGainMult.cxx
CommitLineData
f72219cb 1
2/**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16
17/*
18
19
20Basic calibration and QA class for the TPC gain calibration based on tracks from BEAM events.
21
22
23Send comments etc. to: A.Kalweit@gsi.de, marian.ivanov@cern.ch
24*/
25
26
27#include "Riostream.h"
28#include "TChain.h"
29#include "TTree.h"
30#include "TH1F.h"
31#include "TH2F.h"
32#include "TList.h"
33#include "TMath.h"
34#include "TCanvas.h"
35#include "TFile.h"
36#include "TF1.h"
6feb400f 37#include "TVectorF.h"
f72219cb 38#include "TProfile.h"
39
40#include "AliTPCcalibDB.h"
41#include "AliTPCclusterMI.h"
42#include "AliTPCClusterParam.h"
43#include "AliTPCseed.h"
44#include "AliESDVertex.h"
45#include "AliESDEvent.h"
46#include "AliESDfriend.h"
47#include "AliESDInputHandler.h"
48#include "AliAnalysisManager.h"
49#include "AliTPCParam.h"
50
51#include "AliComplexCluster.h"
52#include "AliTPCclusterMI.h"
53
54#include "AliLog.h"
55
56#include "AliTPCcalibGainMult.h"
57
58#include "TTreeStream.h"
6feb400f 59#include "TDatabasePDG.h"
60#include "AliKFParticle.h"
61#include "AliKFVertex.h"
62#include "AliESDv0.h"
63#include "AliESDkink.h"
64#include "AliRecoParam.h"
65#include "AliTracker.h"
66#include "AliTPCTransform.h"
67#include "AliTPCROC.h"
68#include "TROOT.h"
f72219cb 69
70ClassImp(AliTPCcalibGainMult)
71
72
73AliTPCcalibGainMult::AliTPCcalibGainMult()
74 :AliTPCcalibBase(),
75 fMIP(0),
76 fLowerTrunc(0),
77 fUpperTrunc(0),
78 fUseMax(kFALSE),
79 fHistNTracks(0),
80 fHistClusterShape(0),
81 fHistQA(0),
82 fHistGainSector(0),
83 fHistPadEqual(0),
6feb400f 84 fHistGainMult(0),
85 fPIDMatrix(0),
86 fHistdEdxMap(0),
87 fHistdEdxMax(0),
88 fHistdEdxTot(0),
89 fdEdxTree(0),
90 fBBParam(0)
f72219cb 91{
92 //
93 // Empty default cosntructor
94 //
95 AliInfo("Default Constructor");
96}
97
98
99AliTPCcalibGainMult::AliTPCcalibGainMult(const Text_t *name, const Text_t *title)
100 :AliTPCcalibBase(),
101 fMIP(0),
102 fLowerTrunc(0),
103 fUpperTrunc(0),
104 fUseMax(kFALSE),
105 fHistNTracks(0),
106 fHistClusterShape(0),
107 fHistQA(0),
108 fHistGainSector(0),
109 fHistPadEqual(0),
6feb400f 110 fHistGainMult(0),
111 fPIDMatrix(0),
112 fHistdEdxMap(0),
113 fHistdEdxMax(0),
114 fHistdEdxTot(0),
115 fdEdxTree(0),
116 fBBParam(0)
f72219cb 117{
118 //
119 //
120 //
121 SetName(name);
122 SetTitle(title);
123 //
124 fMIP = 50.;
125 fLowerTrunc = 0.02; // IMPORTANT CHANGE --> REMOVE HARDWIRED TRUNCATION FROM TRACKER
126 fUpperTrunc = 0.6;
127 fUseMax = kTRUE; // IMPORTANT CHANGE FOR PbPb; standard: kFALSE;
128 //
129 fHistNTracks = new TH1F("ntracks","Number of Tracks per Event; number of tracks per event; number of tracks",1001,-0.5,1000.5);
130 fHistClusterShape = new TH1F("fHistClusterShape","cluster shape; rms meas. / rms exp.;",300,0,3);
131 fHistQA = new TH3F("fHistQA","dEdx; momentum p (GeV); TPC signal (a.u.); pad",500,0.1,20.,500,0.,500,6,-0.5,5.5);
6feb400f 132 AliTPCcalibBase::BinLogX(fHistQA);
f72219cb 133 //
134 //
135 // MIP, sect, pad (short,med,long,full,oroc), run, ncl
136 Int_t binsGainSec[5] = { 145, 72, 4, 10000000, 65};
137 Double_t xminGainSec[5] = { 10., -0.5, -0.5, -0.5, -0.5};
138 Double_t xmaxGainSec[5] = {300., 71.5, 3.5, 9999999.5, 64.5};
139 TString axisNameSec[5]={"Q","sector","pad type","run", "ncl"};
140 TString axisTitleSec[5]={"Q (a.u)","sector","pad type","run","ncl"};
141 //
142 fHistGainSector = new THnSparseF("fHistGainSector","0:MIP, 1:sect, 2:pad, 3:run, 4:ncl", 5, binsGainSec, xminGainSec, xmaxGainSec);
143 for (Int_t iaxis=0; iaxis<5;iaxis++){
144 fHistGainSector->GetAxis(iaxis)->SetName(axisNameSec[iaxis]);
145 fHistGainSector->GetAxis(iaxis)->SetTitle(axisTitleSec[iaxis]);
146 }
147 //
148 //
149 //
bf07be7e 150 Int_t binsPadEqual[6] = { 400, 400, 4, 20, 50, 100};
151 Double_t xminPadEqual[6] = { 0.0, 0.0, -0.5, 0, -250, 0};
152 Double_t xmaxPadEqual[6] = { 2.0, 2.0, 3.5, 13000, 250, 3};
f72219cb 153 TString axisNamePadEqual[6] = {"dEdxRatioMax","dEdxRatioTot","padType","mult","driftlength", "1_pt"};
154 TString axisTitlePadEqual[6] = {"dEdx_padRegion/mean_dEdx Qmax", "dEdx_padRegion/mean_dEdx Qtot","padType","mult","driftlength", "1/pt"};
155 //
156 fHistPadEqual = new THnSparseF("fHistPadEqual","0:dEdx_pad/dEdx_mean, 1:pad, 2:mult, 3:drift, 4:1/pt", 6, binsPadEqual, xminPadEqual, xmaxPadEqual);
157 for (Int_t iaxis=0; iaxis<6;iaxis++){
158 fHistPadEqual->GetAxis(iaxis)->SetName(axisNamePadEqual[iaxis]);
159 fHistPadEqual->GetAxis(iaxis)->SetTitle(axisTitlePadEqual[iaxis]);
160 }
161 //
162 //
163 // MIP Qmax, MIP Qtot, z, pad, vtx. contribut., ncl
164 Int_t binsGainMult[6] = { 145, 145, 25, 4, 100, 80};
165 Double_t xminGainMult[6] = { 10., 10., 0, -0.5, 0, -0.5};
012c4694 166 Double_t xmaxGainMult[6] = {300., 300., 250, 3.5, 13000, 159.5};
f72219cb 167 TString axisNameMult[6]={"Qmax","Qtot","drift","padtype""multiplicity","ncl"};
168 TString axisTitleMult[6]={"Qmax (a.u)","Qtot (a.u.)","driftlenght l (cm)","Pad Type","multiplicity","ncl"};
169 //
170 fHistGainMult = new THnSparseF("fHistGainMult","MIP Qmax, MIP Qtot, z, type, vtx. contribut.", 6, binsGainMult, xminGainMult, xmaxGainMult);
171 for (Int_t iaxis=0; iaxis<6;iaxis++){
172 fHistGainMult->GetAxis(iaxis)->SetName(axisNameMult[iaxis]);
173 fHistGainMult->GetAxis(iaxis)->SetTitle(axisTitleMult[iaxis]);
174 }
175 //
f72219cb 176 //
6feb400f 177 // dedx maps - bigger granulatity in phi -
178 // to construct the dedx sector/phi map
179 Int_t binsGainMap[4] = { 100, 90, 10, 6};
180 Double_t xminGainMap[4] = { 0.3, -TMath::Pi(), 0, 0};
181 Double_t xmaxGainMap[4] = { 2, TMath::Pi(), 1, 6};
182 TString axisNameMap[4] = {"Q_Qexp","phi", "1/Qexp","Pad Type"};
183 TString axisTitleMap[4] = {"Q/Q_{exp}","#phi (a.u.)","1/Q_{exp}","Pad Type"};
184 //
185 fHistdEdxMap = new THnSparseF("fHistdEdxMap","fHistdEdxMap", 4, binsGainMap, xminGainMap, xmaxGainMap);
186 for (Int_t iaxis=0; iaxis<4;iaxis++){
187 fHistdEdxMap->GetAxis(iaxis)->SetName(axisNameMap[iaxis]);
188 fHistdEdxMap->GetAxis(iaxis)->SetTitle(axisTitleMap[iaxis]);
189 }
190 //
191 //
192 //
193 // dedx maps
194 Int_t binsGainMax[6] = { 100, 10, 10, 10, 5, 3};
195 Double_t xminGainMax[6] = { 0.5, 0, 0, 0, 0, 0};
196 Double_t xmaxGainMax[6] = { 1.5, 1, 1.0, 1.0, 3000, 3};
197 TString axisNameMax[6] = {"Q_Qexp","1/Qexp", "phi","theta","mult", "Pad Type"};
198 TString axisTitleMax[6] = {"Q/Q_{exp}","1/Qexp", "#phi","#theta","mult","Pad Type"};
199 //
200 fHistdEdxMax = new THnSparseF("fHistdEdxMax","fHistdEdxMax", 6, binsGainMax, xminGainMax, xmaxGainMax);
201 fHistdEdxTot = new THnSparseF("fHistdEdxTot","fHistdEdxTot", 6, binsGainMax, xminGainMax, xmaxGainMax);
202 for (Int_t iaxis=0; iaxis<6;iaxis++){
203 fHistdEdxMax->GetAxis(iaxis)->SetName(axisNameMax[iaxis]);
204 fHistdEdxMax->GetAxis(iaxis)->SetTitle(axisTitleMax[iaxis]);
205 fHistdEdxTot->GetAxis(iaxis)->SetName(axisNameMax[iaxis]);
206 fHistdEdxTot->GetAxis(iaxis)->SetTitle(axisTitleMax[iaxis]);
207 }
208 //
209 AliInfo("Non Default Constructor");
f72219cb 210}
211
212
213AliTPCcalibGainMult::~AliTPCcalibGainMult(){
214 //
6feb400f 215 // Destructor
f72219cb 216 //
b9ab8e40 217 delete fHistNTracks; // histogram showing number of ESD tracks per event
218 delete fHistClusterShape; // histogram to check the cluster shape
219 delete fHistQA; // dE/dx histogram showing the final spectrum
220 //
221 delete fHistGainSector; // histogram which shows MIP peak for each of the 3x36 sectors (pad region)
222 delete fHistPadEqual; // histogram for the equalization of the gain in the different pad regions -> pass0
223 delete fHistGainMult; // histogram which shows decrease of MIP signal as a function
6feb400f 224 //
225 delete fHistdEdxMap;
226 delete fHistdEdxMax;
227 delete fHistdEdxTot;
228 delete fdEdxTree;
f72219cb 229}
230
231
232
233void AliTPCcalibGainMult::Process(AliESDEvent *event) {
234 //
6feb400f 235 // Main function of the class
236 // 1. Select Identified particles - for identified particles the flag in the PID matrix is stored
237 // 1.a) ProcessV0s - select Electron (gamma coversion) and pion canditates (from K0s)
238 // 1.b) ProcessTOF - select - Proton, kaon and pions candidates
239 // AS THE TOF not calibrated yet in Pass0 - we are calibrating the TOF T0 in this function
240 // 1.c) ProcessCosmic - select cosmic mumn candidates - too few entries - not significant for the calibration
241 // 1.d) ProcessKinks - select Kaon and pion candidates. From our experience (see Kink debug streamer), the angular cut for kink daughter is not sufficient - big contamination - delta rays, hadronic interaction (proton)
242 // - NOT USED for the
243 //
244 // 2. Loop over tracks
245 // 2.a DumpTrack() - for identified particles dump the track and dEdx information into the tree (for later fitting)
246 // 3. Actual fitting for the moment macro
247
f72219cb 248 //
6feb400f 249 // Criteria for the track selection
f72219cb 250 //
6feb400f 251 const Int_t kMinNCL=80; // minimal number of cluster - tracks accepted for the dedx calibration
252 const Double_t kMaxEta=0.8; // maximal eta fo the track to be accepted
253 const Double_t kMaxDCAR=10; // maximal DCA R of the track
254 const Double_t kMaxDCAZ=5; // maximal DCA Z of the track
255 const Double_t kMIPPt=0.45; // MIP pt
256
f72219cb 257 if (!event) {
258 Printf("ERROR: ESD not available");
259 return;
260 }
6feb400f 261 fCurrentEvent=event ;
262 fMagF = event->GetMagneticField();
263 Int_t ntracks=event->GetNumberOfTracks();
f72219cb 264 AliESDfriend *esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
265 if (!esdFriend) {
6feb400f 266 //Printf("ERROR: esdFriend not available");
267 delete fPIDMatrix;
268 return;
269 }
270 if (!(esdFriend->TestSkipBit())) fPIDMatrix= new TMatrixD(ntracks,5);
271 fHistNTracks->Fill(ntracks);
272 ProcessCosmic(event); // usually not enogh statistic
273
274 if (esdFriend->TestSkipBit()) {
275 return;
f72219cb 276 }
6feb400f 277 //
278 ProcessV0s(event); //
279 ProcessTOF(event); //
280 ProcessKinks(event); // not relyable
281 DumpHPT(event); //
f72219cb 282 UInt_t runNumber = event->GetRunNumber();
012c4694 283 Int_t nContributors = event->GetNumberOfTracks();
f72219cb 284 //
285 // track loop
286 //
287 for (Int_t i=0;i<ntracks;++i) {
288 //
289 AliESDtrack *track = event->GetTrack(i);
290 if (!track) continue;
291 //
292 AliExternalTrackParam * trackIn = (AliExternalTrackParam *)track->GetInnerParam();
293 if (!trackIn) continue;
294
295 // calculate necessary track parameters
296 Double_t meanP = trackIn->GetP();
297 Int_t ncls = track->GetTPCNcls();
298
6feb400f 299 if (ncls < kMinNCL) continue;
f72219cb 300 // exclude tracks which do not look like primaries or are simply too short or on wrong sectors
6feb400f 301 if (TMath::Abs(trackIn->Eta()) > kMaxEta) continue;
f72219cb 302
f72219cb 303 UInt_t status = track->GetStatus();
304 if ((status&AliESDtrack::kTPCrefit)==0) continue;
305 //if (track->GetNcls(0) < 3) continue; // ITS clusters
306 Float_t dca[2], cov[3];
307 track->GetImpactParameters(dca,cov);
308 Float_t primVtxDCA = TMath::Sqrt(dca[0]*dca[0]);
6feb400f 309 if (primVtxDCA > kMaxDCAR || primVtxDCA < 0.00001) continue;
310 if (TMath::Abs(dca[1]) > kMaxDCAZ) continue;
f72219cb 311 //
f72219cb 312 //
6feb400f 313 // fill Alexander QA histogram
f72219cb 314 //
f72219cb 315 if (primVtxDCA < 3 && track->GetNcls(0) > 3 && track->GetKinkIndex(0) == 0 && ncls > 100) fHistQA->Fill(meanP, track->GetTPCsignal(), 5);
316
317 // Get seeds
318 AliESDfriendTrack *friendTrack = esdFriend->GetTrack(i);
319 if (!friendTrack) continue;
320 TObject *calibObject;
321 AliTPCseed *seed = 0;
322 for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) {
323 if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
324 }
6feb400f 325 if (seed) DumpTrack(track, friendTrack,seed,i); // MI implementation for the identified particles
326 //
327 if (seed) { // seed the container with track parameters and the clusters
328 //
329 const AliExternalTrackParam * trackOut = friendTrack->GetTPCOut(); // tack at the outer radius of TPC
f72219cb 330 if (!trackIn) continue;
331 if (!trackOut) continue;
332 Double_t meanDrift = 250 - 0.5*TMath::Abs(trackIn->GetZ() + trackOut->GetZ());
333 //
334 for (Int_t irow =0; irow<160;irow++) {
335 AliTPCTrackerPoint * point = seed->GetTrackPoint(irow);
336 if (point==0) continue;
337 AliTPCclusterMI * cl = seed->GetClusterPointer(irow);
338 if (cl==0) continue;
339 //
340 Float_t rsigmay = TMath::Sqrt(point->GetSigmaY());
341 fHistClusterShape->Fill(rsigmay);
342 }
343 //
344 Int_t row0 = 0;
345 Int_t row1 = 160;
346 //
347 Double_t signalShortMax = seed->CookdEdxAnalytical(fLowerTrunc,fUpperTrunc,1,0,62);
348 Double_t signalMedMax = seed->CookdEdxAnalytical(fLowerTrunc,fUpperTrunc,1,63,126);
349 Double_t signalLongMax = seed->CookdEdxAnalytical(fLowerTrunc,fUpperTrunc,1,127,159);
350 Double_t signalMax = seed->CookdEdxAnalytical(fLowerTrunc,fUpperTrunc,1,row0,row1);
351 Double_t signalArrayMax[4] = {signalShortMax, signalMedMax, signalLongMax, signalMax};
352 //
353 Double_t signalShortTot = seed->CookdEdxAnalytical(fLowerTrunc,fUpperTrunc,0,0,62);
354 Double_t signalMedTot = seed->CookdEdxAnalytical(fLowerTrunc,fUpperTrunc,0,63,126);
355 Double_t signalLongTot = seed->CookdEdxAnalytical(fLowerTrunc,fUpperTrunc,0,127,159);
356 Double_t signalTot = seed->CookdEdxAnalytical(fLowerTrunc,fUpperTrunc,0,row0,row1);
357 Double_t signalArrayTot[4] = {signalShortTot, signalMedTot, signalLongTot, signalTot};
358 //
359 Double_t mipSignalShort = fUseMax ? signalShortMax : signalShortTot;
360 Double_t mipSignalMed = fUseMax ? signalMedMax : signalMedTot;
361 Double_t mipSignalLong = fUseMax ? signalLongMax : signalLongTot;
362 Double_t mipSignalOroc = seed->CookdEdxAnalytical(fLowerTrunc,fUpperTrunc,fUseMax,63,159);
363 Double_t signal = fUseMax ? signalMax : signalTot;
364 //
365 fHistQA->Fill(meanP, mipSignalShort, 0);
366 fHistQA->Fill(meanP, mipSignalMed, 1);
367 fHistQA->Fill(meanP, mipSignalLong, 2);
368 fHistQA->Fill(meanP, signal, 3);
369 fHistQA->Fill(meanP, mipSignalOroc, 4);
370 //
371 // "dEdxRatioMax","dEdxRatioTot","padType","mult","driftlength", "1_pt"
372 Float_t meanMax = (1/3.)*(signalArrayMax[0] + signalArrayMax[1] + signalArrayMax[2]);
b4ffdb90 373 Float_t meanTot = (1/3.)*(signalArrayTot[0] + signalArrayTot[1] + signalArrayTot[2]);
4a4c6be2 374 if (meanMax < 1e-5 || meanTot < 1e-5) continue;
f72219cb 375 for(Int_t ipad = 0; ipad < 4; ipad ++) {
376 Double_t vecPadEqual[6] = {signalArrayMax[ipad]/meanMax, signalArrayTot[ipad]/meanTot, ipad, nContributors, meanDrift, track->OneOverPt()};
377 fHistPadEqual->Fill(vecPadEqual);
378 }
379 //
6feb400f 380 // if (meanP > 0.4 && meanP < 0.55) {
381 if ( TMath::Abs(meanP-kMIPPt)<0.05 ) {
f72219cb 382 Double_t vecMult[6] = {seed->CookdEdxAnalytical(fLowerTrunc,fUpperTrunc,1,row0,row1),
383 seed->CookdEdxAnalytical(fLowerTrunc,fUpperTrunc,0,row0,row1),
384 meanDrift,
385 3,
386 nContributors,
387 ncls};
388 //
389 fHistGainMult->Fill(vecMult);
390 vecMult[0]=mipSignalShort; vecMult[1]=mipSignalShort; vecMult[3]=0;
391 fHistGainMult->Fill(vecMult);
392 vecMult[0]=mipSignalMed; vecMult[1]=mipSignalMed; vecMult[3]=1;
393 fHistGainMult->Fill(vecMult);
394 vecMult[0]=mipSignalLong; vecMult[1]=mipSignalLong; vecMult[3]=2;
395 fHistGainMult->Fill(vecMult);
396 //
397 }
398 //
399 //
6feb400f 400 if ( TMath::Abs(meanP-kMIPPt)>0.05 ) continue; // only MIP pions
401 //if (meanP > 0.5 || meanP < 0.4) continue; // only MIP pions
f72219cb 402 //
403 // for each track, we look at the three different pad regions, split it into tracklets, check if the sector does not change and fill the histogram
404 //
405 Bool_t isNotSplit[3] = {kTRUE, kTRUE, kTRUE}; // short, medium, long (true if the track is not split between two chambers)
406 //
407 Double_t sector[4] = {-1, -1, -1, -1}; // sector number short, medium, long, all
408 Int_t ncl[3] = {0,0,0};
409 //
6feb400f 410
f72219cb 411 for (Int_t irow=0; irow < 159; irow++){
412 Int_t padRegion = 0;
413 if (irow > 62) padRegion = 1;
414 if (irow > 126) padRegion = 2;
415 //
416 AliTPCclusterMI* cluster = seed->GetClusterPointer(irow);
417 if (!cluster) continue;
418 if (sector[padRegion] == -1) {
419 sector[padRegion] = cluster->GetDetector();
420 continue;
421 }
422 if (sector[padRegion] != -1 && sector[padRegion] != cluster->GetDetector()) isNotSplit[padRegion] = kFALSE;
423 ncl[padRegion]++;
424 }
425 //
426 // MIP, sect, pad, run
427 //
428 Double_t vecMip[5] = {mipSignalShort, mipSignalMed, mipSignalLong, signal, mipSignalOroc};
429 //
430 for(Int_t ipad = 0; ipad < 3; ipad++) {
6feb400f 431 // AK. - run Number To be removed - not needed
f72219cb 432 Double_t vecGainSec[5] = {vecMip[ipad], sector[ipad], ipad, runNumber, ncl[ipad]};
433 if (isNotSplit[ipad]) fHistGainSector->Fill(vecGainSec);
434 }
435 }
436
437 }
6feb400f 438
439 delete fPIDMatrix;
f72219cb 440}
441
442
a8ef8a9c 443void AliTPCcalibGainMult::MakeLookup(THnSparse * /*hist*/, Char_t * /*outputFile*/) {
444 //
445 // Not yet implemented
446 //
4a4c6be2 447}
f72219cb 448
449
450void AliTPCcalibGainMult::Analyze() {
6feb400f 451 //
452 // Not implemented
453 //
f72219cb 454
455 return;
456
457}
458
459
460Long64_t AliTPCcalibGainMult::Merge(TCollection *li) {
6feb400f 461 //
462 // merging of the component
463 //
f72219cb 464
465 TIterator* iter = li->MakeIterator();
466 AliTPCcalibGainMult* cal = 0;
467
468 while ((cal = (AliTPCcalibGainMult*)iter->Next())) {
469 if (!cal->InheritsFrom(AliTPCcalibGainMult::Class())) {
470 Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
471 return -1;
472 }
473
474 if (cal->GetHistNTracks()) fHistNTracks->Add(cal->GetHistNTracks());
475 if (cal->GetHistClusterShape()) fHistClusterShape->Add(cal->GetHistClusterShape());
476 if (cal->GetHistQA()) fHistQA->Add(cal->GetHistQA());
477 if (cal->GetHistGainSector()) fHistGainSector->Add(cal->GetHistGainSector());
478 if (cal->GetHistPadEqual()) fHistPadEqual->Add(cal->GetHistPadEqual());
479 if (cal->GetHistGainMult()) fHistGainMult->Add(cal->GetHistGainMult());
6feb400f 480 if (cal->fHistdEdxMap){
481 if (fHistdEdxMap) fHistdEdxMap->Add(cal->fHistdEdxMap);
482 }
483 if (cal->fHistdEdxMax){
484 if (fHistdEdxMax) fHistdEdxMax->Add(cal->fHistdEdxMax);
485 }
486 if (cal->fHistdEdxTot){
487 if (fHistdEdxTot) fHistdEdxTot->Add(cal->fHistdEdxTot);
488 }
489 //
490 // Originally we tireied to write the tree to the same file as other calibration components
491 // We failed in merging => therefore this optio was disabled
492 //
493 // if (cal->fdEdxTree && cal->fdEdxTree->GetEntries()>0) {
494 // if (fdEdxTree) {
495 // const Int_t kMax=100000;
496 // Int_t entriesSum = (Int_t)fdEdxTree->GetEntries();
497 // Int_t entriesCurrent = (Int_t)cal->fdEdxTree->GetEntries();
498 // Int_t entriesCp=TMath::Min((Int_t) entriesCurrent*(kMax*entriesSum),entriesCurrent);
499 // // cal->fdEdxTree->SetBranchStatus("track*",0);
500 // // cal->fdEdxTree->SetBranchStatus("vertex*",0);
501 // // cal->fdEdxTree->SetBranchStatus("tpcOut*",0);
502 // // cal->fdEdxTree->SetBranchStatus("vec*",0);
503 // // fdEdxTree->SetBranchStatus("track*",0);
504 // // fdEdxTree->SetBranchStatus("vertex*",0);
505 // // fdEdxTree->SetBranchStatus("tpcOut*",0);
506 // // fdEdxTree->SetBranchStatus("vec*",0);
507 // fdEdxTree->Print();
508 // fdEdxTree->Dump();
509 // fdEdxTree->GetEntry(0);
510 // for (Int_t i=0; i<entriesCurrent; i++){
511 // cal->fdEdxTree->CopyAddresses(fdEdxTree);
512 // cal->fdEdxTree->GetEntry(i);
513 // fdEdxTree->Fill();
514 // }
515 // TObjArray *brarray = cal->fdEdxTree->GetListOfBranches();
516 // for (Int_t i=0; i<brarray->GetEntries(); i++) {TBranch * br = (TBranch *)brarray->At(i); br->SetAddress(0); }
517 // }
518 // else{
519 // fdEdxTree = (TTree*)(cal->fdEdxTree->Clone());
520 // TObjArray *brarray = fdEdxTree->GetListOfBranches();
521 // for (Int_t i=0; i<brarray->GetEntries(); i++) {TBranch * br = (TBranch *)brarray->At(i); br->SetAddress(0);}
522 // }
523 //}
524
f72219cb 525 }
526
527 return 0;
528
529}
530
531
532
f72219cb 533
534
535void AliTPCcalibGainMult::UpdateGainMap() {
536 //
537 // read in the old gain map and scale it appropriately...
538 //
539 /*
540 gSystem->Load("libANALYSIS");
541 gSystem->Load("libTPCcalib");
542 //
543 TFile jj("Run0_999999999_v1_s0.root");
544 AliTPCCalPad * pad = AliCDBEntry->GetObject()->Clone();
545 TFile hh("output.root");
546 AliTPCcalibGainMult * gain = calibTracksGain;
547 TH2D * histGainSec = gain->GetHistGainSector()->Projection(0,1);
548 //
549 TObjArray arr;
550 histGainSec->FitSlicesY(0, 0, -1, 0, "QNR", &arr);
551 TH1D * meanGainSec = arr->At(1);
552 Double_t gainsIROC[36];
553 Double_t gainsOROC[36];
554 Double_t gains[72];
555 //
556 for(Int_t isec = 1; isec < meanGainSec->GetNbinsX() + 1; isec++) {
557 cout << isec << " " << meanGainSec->GetXaxis()->GetBinCenter(isec) << " " <<meanGainSec->GetBinContent(isec) << endl;
558 gains[isec-1] = meanGainSec->GetBinContent(isec);
559 if (isec < 37) {
560 gainsIROC[isec-1] = meanGainSec->GetBinContent(isec);
561 } else {
562 gainsOROC[isec - 36 -1] = meanGainSec->GetBinContent(isec);
563 }
564 }
565 Double_t meanIroc = TMath::Mean(36, gainsIROC);
566 Double_t meanOroc = TMath::Mean(36, gainsIROC);
567 for(Int_t i = 0; i < 36; i++) gains[i] /= meanIroc;
568 for(Int_t i = 36; i < 72; i++) gains[i] /= meanOroc;
569 //
570 for(Int_t i = 0; i < 72; i++) {
571 AliTPCCalROC * chamber = pad->GetCalROC(i);
572 chamber->Multiply(gains[i]);
573 cout << i << " "<< chamber->GetMean() << endl;
574 }
575 //
576 // update the OCDB
577 //
578 AliCDBMetaData *metaData= new AliCDBMetaData();
579 metaData->SetObjectClassName("AliTPCCalPad");
580 metaData->SetResponsible("Alexander Kalweit");
581 metaData->SetBeamPeriod(1);
582 metaData->SetAliRootVersion("04-19-05"); //root version
583 metaData->SetComment("New gain map for 1600V OROC gain increase and equalization. Valid for runs starting after technical stop beginning of September.");
584 AliCDBId id1("TPC/Calib/GainFactorDedx", 131541, AliCDBRunRange::Infinity()); // important: new gain runs here..
585 AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage("local:///d/alice05/akalweit/projects/OCDBupdate/HighGain_2010-09-03/OCDB/");
586 gStorage->Put(pad, id1, metaData);
587 */
588
589}
590
591void AliTPCcalibGainMult::UpdateClusterParam() {
592 //
593 //
594 //
595 /*
596 gSystem->Load("libANALYSIS");
597 gSystem->Load("libTPCcalib");
598 //
599 TFile ff("OldClsParam.root");
600 AliTPCClusterParam * param = AliCDBEntry->GetObject()->Clone();
601
602 TFile hh("output.root");
603 AliTPCcalibGainMult * gain = calibGainMult;
604 TH2D * histGainSec = gain->GetHistGainSector()->Projection(0,2);
605 TObjArray arr;
606 histGainSec->FitSlicesY(0, 0, -1, 0, "QNR", &arr);
607 histGainSec->Draw("colz");
608 TH1D * fitVal = arr.At(1);
609 fitVal->Draw("same");
610 param->GetQnormCorrMatrix()->Print();
611 param->GetQnormCorrMatrix()(0,5) *= fitVal->GetBinContent(1)/fitVal->GetBinContent(1); // short pads Qtot
612 param->GetQnormCorrMatrix()(1,5) *= fitVal->GetBinContent(2)/fitVal->GetBinContent(1); // med pads Qtot
613 param->GetQnormCorrMatrix()(2,5) *= fitVal->GetBinContent(3)/fitVal->GetBinContent(1); // long pads Qtot
614 //
615 param->GetQnormCorrMatrix()(3,5) *= fitVal->GetBinContent(1)/fitVal->GetBinContent(1); // short pads Qmax -> scaling assumed
616 param->GetQnormCorrMatrix()(4,5) *= fitVal->GetBinContent(2)/fitVal->GetBinContent(1); // med pads Qmax -> scaling assumed
617 param->GetQnormCorrMatrix()(5,5) *= fitVal->GetBinContent(3)/fitVal->GetBinContent(1); // long pads Qmax -> scaling assumed
618 //
619 TFile jj("newClusterParam.root","RECREATE");
620 param->Write();
621 param->GetQnormCorrMatrix()->Print();
622 //
623 // update the OCDB
624 //
625 AliCDBMetaData *metaData= new AliCDBMetaData();
626 metaData->SetObjectClassName("AliTPCClusterParam");
627 metaData->SetResponsible("Alexander Kalweit");
628 metaData->SetBeamPeriod(1);
629 metaData->SetAliRootVersion("04-19-04"); //root version
630 metaData->SetComment("1600V OROC / hard thres. / new algorithm");
631 AliCDBId id1("TPC/Calib/ClusterParam", 0, AliCDBRunRange::Infinity()); // important: new gain runs here..
632 AliCDBStorage * gStorage = AliCDBManager::Instance()->GetStorage("local:///lustre/alice/akalweit/baseline/CalibrationEntries/OldThres_NewAlgo_PP");
633 gStorage->Put(param, id1, metaData);
634 */
635
636
637}
638
6feb400f 639
640void AliTPCcalibGainMult::DumpTrack(AliESDtrack * track, AliESDfriendTrack *ftrack, AliTPCseed * seed, Int_t index){
641 //
642 // dump interesting tracks
643 // 1. track at MIP region
644 // 2. highly ionizing protons
645 // pidType: 0 - unselected
646 // 1 - TOF
647 // 2 - V0
648 // 4 - Cosmic
649 // or of value
650 //
651 const Int_t kMax=10000;
652 const Int_t kMinRows=80;
653 const Double_t kDCAcut=30;
654 //
655 // Bethe Bloch paramerization
656 //
657 Double_t kp1= 0.0851148;
658 Double_t kp2= 9.25771;
659 Double_t kp3= 2.6558e-05;
660 Double_t kp4= 2.32742;
661 Double_t kp5= 1.83039;
662 if (fBBParam){
663 kp1=(*fBBParam)[0];
664 kp2=(*fBBParam)[1];
665 kp3=(*fBBParam)[2];
666 kp4=(*fBBParam)[3];
667 kp5=(*fBBParam)[4];
668 }
669 //
f4162f93 670 static const AliTPCROC *roc = AliTPCROC::Instance();
cb4201d4 671 static const TDatabasePDG *pdg = TDatabasePDG::Instance();
6feb400f 672
673 Int_t nclITS = track->GetNcls(0);
674 Int_t ncl = track->GetTPCncls();
675 Double_t ncl21 = track->GetTPCClusterInfo(3,1);
676 Double_t ncl20 = track->GetTPCClusterInfo(3,0);
677 //
678 if (!seed) return;
679 if (ncl21<kMinRows) return;
680 static Int_t counter=0;
681 static Int_t counterHPT=0;
682 //
683 static TH1F *hisBB=new TH1F("hisBB","hisBB",20,0.1,1.00); // bethe bloch histogram =
684 // used to cover more homegenously differnt dEdx regions
685 static Double_t massPi = pdg->GetParticle("pi-")->Mass(); //
686 static Double_t massK = pdg->GetParticle("K-")->Mass();
687 static Double_t massP = pdg->GetParticle("proton")->Mass();
688 static Double_t massE = pdg->GetParticle("e-")->Mass();
689 static Double_t massMuon = pdg->GetParticle("mu-")->Mass();
690 static Double_t radius0= roc->GetPadRowRadiiLow(roc->GetNRows(0)/2);
691 static Double_t radius1= roc->GetPadRowRadiiUp(30);
692 static Double_t radius2= roc->GetPadRowRadiiUp(roc->GetNRows(36)-15);
693
694 AliESDVertex *vertex= (AliESDVertex *)fCurrentEvent->GetPrimaryVertex();
695 //
696 // Estimate current MIP position -
697 //
698 static Double_t mipArray[kMax]; // mip array
699 static Int_t counterMIP0=0;
700 static Double_t medianMIP0=100000; // current MIP median position - estimated after some mimnimum number of MIP tracks
701
702 if (TMath::Abs(track->GetP()-0.5)<0.1&&track->GetTPCsignal()/medianMIP0<1.5){
703 mipArray[counterMIP0%kMax]= track->GetTPCsignal();
704 counterMIP0++;
705 if (counterMIP0>10) medianMIP0=TMath::Median(counterMIP0%kMax, mipArray);
706 }
707 // the PID as defiend from the external sources
708 //
709 Int_t isElectron = TMath::Nint((*fPIDMatrix)(index,0));
710 Int_t isMuon = TMath::Nint((*fPIDMatrix)(index,1));
711 Int_t isPion = TMath::Nint((*fPIDMatrix)(index,2));
712 Int_t isKaon = TMath::Nint((*fPIDMatrix)(index,3));
713 Int_t isProton = TMath::Nint((*fPIDMatrix)(index,4));
714 Float_t dca[2];
715 track->GetImpactParameters(dca[0],dca[1]);
716 //
717 if ( (isMuon==0 && isElectron==0) && (TMath::Sqrt(dca[0]*dca[0]+dca[1]*dca[1])>kDCAcut) ) return;
718 Double_t normdEdx= track->GetTPCsignal()/(medianMIP0); // TPC signal normalized to the MIP
719 //
720 AliExternalTrackParam * trackIn = (AliExternalTrackParam *)track->GetInnerParam();
721 AliExternalTrackParam * trackOut = (AliExternalTrackParam *)track->GetOuterParam();
722 AliExternalTrackParam * tpcOut = (AliExternalTrackParam *)ftrack->GetTPCOut();
723 if (!trackIn) return;
724 if (!trackOut) return;
725 if (!tpcOut) return;
726 if (trackIn->GetZ()*trackOut->GetZ()<0) return; // remove crossing tracks
727 //
728 // calculate local and global angle
729 Int_t side = (trackIn->GetZ()>0)? 1:-1;
730 Double_t tgl=trackIn->GetTgl();
731 Double_t gangle[3]={0,0,0};
732 Double_t langle[3]={0,0,0};
733 Double_t length[3]={0,0,0};
734 Double_t pxpypz[3]={0,0,0};
735 Double_t bz=fMagF;
736 trackIn->GetXYZAt(radius0,bz,pxpypz); // get the global position at the middle of the IROC
737 gangle[0]=TMath::ATan2(pxpypz[1],pxpypz[0]); // global angle IROC
738 trackIn->GetXYZAt(radius1,bz,pxpypz); // get the global position at the middle of the OROC - medium pads
739 gangle[1]=TMath::ATan2(pxpypz[1],pxpypz[0]); // global angle OROC
740 trackOut->GetXYZAt(radius2,bz,pxpypz); // get the global position at the middle of OROC - long pads
741 gangle[2]=TMath::ATan2(pxpypz[1],pxpypz[0]);
742 //
743 trackIn->GetPxPyPzAt(radius0,bz,pxpypz); //get momentum vector
744 langle[0]=TMath::ATan2(pxpypz[1],pxpypz[0])-gangle[0]; //local angle between padrow and track IROC
745 trackIn->GetPxPyPzAt(radius1,bz,pxpypz);
746 langle[1]=TMath::ATan2(pxpypz[1],pxpypz[0])-gangle[1];
747 trackOut->GetPxPyPzAt(radius2,bz,pxpypz); // OROC medium
748 langle[2]=TMath::ATan2(pxpypz[1],pxpypz[0])-gangle[2];
749 for (Int_t i=0; i<3; i++){
750 if (langle[i]>TMath::Pi()) langle[i]-=TMath::TwoPi();
751 if (langle[i]<-TMath::Pi()) langle[i]+=TMath::TwoPi();
752 length[i]=TMath::Sqrt(1+langle[i]*langle[i]+tgl*tgl); // the tracklet length
753 }
754 //
755 // Select the kaons and Protons which are "isolated" in TPC dedx curve
756 //
757 //
758 Double_t dedxP = AliExternalTrackParam::BetheBlochAleph(track->GetInnerParam()->GetP()/massP,kp1,kp2,kp3,kp4,kp5);
759 Double_t dedxK = AliExternalTrackParam::BetheBlochAleph(track->GetInnerParam()->GetP()/massK,kp1,kp2,kp3,kp4,kp5);
760 if (dedxP>2 || dedxK>2){
761 if (track->GetP()<1.2 && normdEdx>1.8&&counterMIP0>10){ // not enough from TOF and V0s triggered by high dedx
762 // signing the Proton and kaon - using the "bitmask" bit 1 and 2 is dedicated for V0s and TOF selected
763 if ( TMath::Abs(normdEdx/dedxP-1)<0.3) isProton+=4;
764 if ( TMath::Abs(normdEdx/dedxK-1)<0.3) isKaon+=4;
765 if (normdEdx/dedxK>1.3) isProton+=8;
766 if (normdEdx/dedxP<0.7) isKaon+=8;
767 }
768 }
769 //
770 //
771 //
772 Double_t mass = 0;
773 Bool_t isHighPt = ((TMath::Power(1/track->Pt(),4)*gRandom->Rndm())<0.005); // rnadomly selected HPT tracks
774 // there are selected for the QA of the obtained calibration
775 Bool_t isMIP = TMath::Abs(track->GetInnerParam()->P()-0.4)<0.005&&(counter<kMax); //
776 // REMINDER - it is not exactly MIP - we select the regtion where the Kaon and Electrons are well separated
777
778 if (isElectron>0) mass = massE;
779 if (isProton>0) mass = massP;
780 if (isKaon>0) mass = massK;
781 if (isMuon>0) mass = massMuon;
782 if (isPion>0) mass = massPi;
783 if (isHighPt) mass = massPi; //assign mass of pions
784 if (isMIP&&track->GetTPCsignal()/medianMIP0<1.5) mass = massPi; //assign mass of pions
785 if (mass==0) return;
786 //
787 // calculate expected dEdx
788 Double_t dedxDef= 0;
789 Double_t dedxDefPion= 0,dedxDefProton=0, dedxDefKaon=0;
790 Double_t pin=trackIn->GetP();
791 Double_t pout=trackOut->GetP();
792 Double_t p=(pin+pout)*0.5; // momenta as the mean between tpc momenta at inner and outer wall of the TPC
793 if (mass>0) dedxDef = AliExternalTrackParam::BetheBlochAleph(p/mass,kp1,kp2,kp3,kp4,kp5);
794 dedxDefPion = AliExternalTrackParam::BetheBlochAleph(p/massPi,kp1,kp2,kp3,kp4,kp5);
795 dedxDefProton = AliExternalTrackParam::BetheBlochAleph(p/massP,kp1,kp2,kp3,kp4,kp5);
796 dedxDefKaon = AliExternalTrackParam::BetheBlochAleph(p/massK,kp1,kp2,kp3,kp4,kp5);
797 //
798 // dEdx Truncated mean vectros with differnt tuncation
799 // 11 truncations array - 0-10 - 0~50% 11=100%
800 // 3 Regions - IROC,OROC0, OROC1
801 // 2 Q - total charge and max charge
802 // Log - Logarithmic mean used
803 // Up/Dwon - Upper half or lower half of truncation used
804 // RMS - rms of the distribction (otherwise truncated mean)
805 // M2 suffix - second moment ( truncated)
806 TVectorF truncUp(11);
807 TVectorF truncDown(11);
808 TVectorF vecAllMax(11);
809 TVectorF vecIROCMax(11);
810 TVectorF vecOROCMax(11);
811 TVectorF vecOROC0Max(11);
812 TVectorF vecOROC1Max(11);
813 //
814 TVectorF vecAllTot(11);
815 TVectorF vecIROCTot(11);
816 TVectorF vecOROCTot(11);
817 TVectorF vecOROC0Tot(11);
818 TVectorF vecOROC1Tot(11);
819 //
820 TVectorF vecAllTotLog(11);
821 TVectorF vecIROCTotLog(11);
822 TVectorF vecOROCTotLog(11);
823 TVectorF vecOROC0TotLog(11);
824 TVectorF vecOROC1TotLog(11);
825 //
826 TVectorF vecAllTotUp(11);
827 TVectorF vecIROCTotUp(11);
828 TVectorF vecOROCTotUp(11);
829 TVectorF vecOROC0TotUp(11);
830 TVectorF vecOROC1TotUp(11);
831 //
832 TVectorF vecAllTotDown(11);
833 TVectorF vecIROCTotDown(11);
834 TVectorF vecOROCTotDown(11);
835 TVectorF vecOROC0TotDown(11);
836 TVectorF vecOROC1TotDown(11);
837
838 TVectorF vecAllTotRMS(11);
839 TVectorF vecIROCTotRMS(11);
840 TVectorF vecOROCTotRMS(11);
841 TVectorF vecOROC0TotRMS(11);
842 TVectorF vecOROC1TotRMS(11);
843 //
844 TVectorF vecAllTotM2(11);
845 TVectorF vecIROCTotM2(11);
846 TVectorF vecOROCTotM2(11);
847 TVectorF vecOROC0TotM2(11);
848 TVectorF vecOROC1TotM2(11);
849 //
850 TVectorF vecAllTotMS(11);
851 TVectorF vecIROCTotMS(11);
852 TVectorF vecOROCTotMS(11);
853 TVectorF vecOROC0TotMS(11);
854 TVectorF vecOROC1TotMS(11);
855 //
856 // Differnt number of clusters definitions - in separate regions of the TPC
857 // 20 - ratio - found/findabel
858 // 21 - number of clusters used for given dEdx calculation
859 //
860 // suffix - 3 or 4 - number of padrows before and after given row to define findable row
861 //
862 Double_t ncl20All = seed->CookdEdxAnalytical(0.0,1, 1 ,0,159,3);
863 Double_t ncl20IROC = seed->CookdEdxAnalytical(0.,1, 1 ,0,63,3);
864 Double_t ncl20OROC = seed->CookdEdxAnalytical(0.,1, 1 ,64,159,3);
865 Double_t ncl20OROC0= seed->CookdEdxAnalytical(0.,1, 1 ,64,128,3);
866 Double_t ncl20OROC1= seed->CookdEdxAnalytical(0.,1, 1 ,129,159,3);
867 //
868 Double_t ncl20All4 = seed->CookdEdxAnalytical(0.0,1, 1 ,0,159,3,4);
869 Double_t ncl20IROC4 = seed->CookdEdxAnalytical(0.,1, 1 ,0,63,3,4);
870 Double_t ncl20OROC4 = seed->CookdEdxAnalytical(0.,1, 1 ,64,159,3,4);
871 Double_t ncl20OROC04= seed->CookdEdxAnalytical(0.,1, 1 ,64,128,3,4);
872 Double_t ncl20OROC14= seed->CookdEdxAnalytical(0.,1, 1 ,129,159,3,4);
873 //
874 Double_t ncl20All3 = seed->CookdEdxAnalytical(0.0,1, 1 ,0,159,3,3);
875 Double_t ncl20IROC3 = seed->CookdEdxAnalytical(0.,1, 1 ,0,63,3,3);
876 Double_t ncl20OROC3 = seed->CookdEdxAnalytical(0.,1, 1 ,64,159,3,3);
877 Double_t ncl20OROC03= seed->CookdEdxAnalytical(0.,1, 1 ,64,128,3,3);
878 Double_t ncl20OROC13= seed->CookdEdxAnalytical(0.,1, 1 ,129,159,3,3);
879 //
880 Double_t ncl21All = seed->CookdEdxAnalytical(0.0,1, 1 ,0,159,2);
881 Double_t ncl21IROC = seed->CookdEdxAnalytical(0.,1, 1 ,0,63,2);
882 Double_t ncl21OROC = seed->CookdEdxAnalytical(0.,1, 1 ,64,159,2);
883 Double_t ncl21OROC0= seed->CookdEdxAnalytical(0.,1, 1 ,64,128,2);
884 Double_t ncl21OROC1= seed->CookdEdxAnalytical(0.,1, 1 ,129,159,2);
885 // calculate truncated dEdx - mean rms M2 ...
886 Int_t ifrac=0;
887 for (Int_t ifracDown=0; ifracDown<1; ifracDown++){
888 for (Int_t ifracUp=0; ifracUp<11; ifracUp++){
889 Double_t fracDown = 0.0+Double_t(ifracDown)*0.05;
890 Double_t fracUp = 0.5+Double_t(ifracUp)*0.05;
891 vecAllMax[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 1 ,0,159,0);
892 vecIROCMax[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 1 ,0,63,0);
893 vecOROCMax[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 1 ,64,159,0);
894 vecOROC0Max[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 1 ,64,128,0);
895 vecOROC1Max[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 1 ,129,159,0);
896 //
897 vecAllTot[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,159,0);
898 vecIROCTot[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,63,0);
899 vecOROCTot[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,159,0);
900 vecOROC0Tot[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,128,0);
901 vecOROC1Tot[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,129,159,0);
902 //
903 vecAllTotLog[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,159,0,2,1);
904 vecIROCTotLog[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,63,0,2,1);
905 vecOROCTotLog[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,159,0,2,1);
906 vecOROC0TotLog[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,128,0,2,1);
907 vecOROC1TotLog[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,129,159,0,2,1);
908 //
909 vecAllTotUp[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,159,4,2,1);
910 vecIROCTotUp[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,63,4,2,1);
911 vecOROCTotUp[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,159,4,2,1);
912 vecOROC0TotUp[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,128,4,2,1);
913 vecOROC1TotUp[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,129,159,4,2,1);
914 //
915 vecAllTotDown[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,159,5,2,1);
916 vecIROCTotDown[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,63,5,2,1);
917 vecOROCTotDown[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,159,5,2,1);
918 vecOROC0TotDown[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,128,5,2,1);
919 vecOROC1TotDown[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,129,159,5,2,1);
920 //
921 vecAllTotRMS[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,159,1,2,0);
922 vecIROCTotRMS[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,63,1,2,0);
923 vecOROCTotRMS[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,159,1,2,0);
924 vecOROC0TotRMS[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,128,1,2,0);
925 vecOROC1TotRMS[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,129,159,1,2,0);
926 //
927 vecAllTotM2[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,159,6,2,1);
928 vecIROCTotM2[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,63,6,2,1);
929 vecOROCTotM2[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,159,6,2,1);
930 vecOROC0TotM2[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,128,6,2,1);
931 vecOROC1TotM2[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,129,159,6,2,1);
932 //
933 vecAllTotMS[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,159,8,2,1);
934 vecIROCTotMS[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,0,63,8,2,1);
935 vecOROCTotMS[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,159,8,2,1);
936 vecOROC0TotMS[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,64,128,8,2,1);
937 vecOROC1TotMS[ifrac]= seed->CookdEdxAnalytical(fracDown,fracUp, 0 ,129,159,8,2,1);
938 truncUp[ifrac]=fracUp;
939 truncDown[ifrac]=fracDown;
940 ifrac++;
941 }
942 }
943 //
944 // Fill histograms
945 //
946 if ((isKaon||isProton||isPion||isElectron||isMIP||isMuon&&(!isHighPt)) && dedxDef>0) {
947 //
948 Int_t ncont = vertex->GetNContributors();
949 for (Int_t ipad=0; ipad<3; ipad++){
950 // histogram with enahanced phi granularity - to make gain phi maps
951 Double_t xxx[4]={0,gangle[ipad],1./dedxDef,ipad*2+((side>0)?0:1)};
952 Double_t nclR=0;
953 if (ipad==0) {xxx[0]=vecIROCTot[4]/medianMIP0; nclR=ncl21IROC/63.;}
954 if (ipad==1) {xxx[0]=vecOROC0Tot[4]/medianMIP0;nclR=ncl21OROC0/63.;}
955 if (ipad==2) {xxx[0]=vecOROC1Tot[4]/medianMIP0;nclR=ncl21OROC1/32.;}
956 xxx[0]/=dedxDef;
957 if (xxx[0]>0) xxx[0]=1./xxx[0];
958 if (TMath::Abs(langle[ipad])<0.25&&nclR>0.4) fHistdEdxMap->Fill(xxx);
959 }
960 for (Int_t ipad=0; ipad<3; ipad++){
961 //
962 // this are histogram to define overall main gain correction
963 // Maybe dead end - we can not put all info which can be used into the THnSparse
964 // It is keeped there for educational point of view
965 //
966 Double_t xxx[6]={0,1./dedxDef, TMath::Abs(langle[ipad]), TMath::Abs(tgl), ncont, ipad };
967 if (ipad==0) {xxx[0]=vecIROCTot[4]/medianMIP0;}
968 if (ipad==1) {xxx[0]=vecOROC0Tot[4]/medianMIP0;}
969 if (ipad==2) {xxx[0]=vecOROC1Tot[4]/medianMIP0;}
970 xxx[0]/=dedxDef;
971 if (xxx[0]>0) xxx[0]=1./xxx[0];
972 if (xxx[0]>0) fHistdEdxTot->Fill(xxx);
973 if (ipad==0) {xxx[0]=vecIROCMax[4]/medianMIP0;}
974 if (ipad==1) {xxx[0]=vecOROC0Max[4]/medianMIP0;}
975 if (ipad==2) {xxx[0]=vecOROC1Max[4]/medianMIP0;}
976 xxx[0]=dedxDef;
977 if (xxx[0]>0) xxx[0]=1./xxx[0];
978 fHistdEdxMax->Fill(xxx);
979 }
980 }
981 //
982 // Downscale selected tracks before filling the tree
983 //
984 Bool_t isSelected = kFALSE;
985 //
986 if (isKaon||isProton||isPion||isElectron||isMIP||isMuon) isSelected=kTRUE;
987 isHighPt = kFALSE;
988 if (!isSelected) isHighPt = ((TMath::Power(1/track->Pt(),4)*gRandom->Rndm())<0.005);
989 if (counter>kMax && ((1/track->Pt()*gRandom->Rndm())>kMax/counter)) return;
990 isSelected|=isHighPt;
991 //
992 //
993 //
994 // Equalize statistic in BB bins - special care about pions
995 Int_t entriesBB = (Int_t)hisBB->GetEntries();
996 if ((isElectron==0 &&isMuon==0 && p<2.) && entriesBB>20 &&dedxDef>0.01){
997 Int_t bin = hisBB->GetXaxis()->FindBin(1./dedxDef);
998 Double_t cont = hisBB->GetBinContent(bin);
999 Double_t mean =(entriesBB)/20.;
1000 if ((isPion>0) && gRandom->Rndm()*cont > 0.1*mean) return;
1001 if ((isPion==0) && gRandom->Rndm()*cont > 0.25*mean) return;
1002 }
1003 if (!isSelected) return;
1004 if (dedxDef>0.01) hisBB->Fill(1./dedxDef);
1005 //
1006 if (isHighPt) counterHPT++;
1007 counter++;
1008 //
1009 TTreeSRedirector * pcstream = GetDebugStreamer();
1010 Double_t ptrel0 = AliTPCcalibDB::GetPTRelative(fTime,fRun,0);
1011 Double_t ptrel1 = AliTPCcalibDB::GetPTRelative(fTime,fRun,1);
1012 Int_t sectorIn = Int_t(18+9*(trackIn->GetAlpha()/TMath::Pi()))%18;
1013 Int_t sectorOut = Int_t(18+9*(trackOut->GetAlpha()/TMath::Pi()))%18;
1014 //
1015 if (pcstream){
1016 (*pcstream)<<"dump"<<
1017 "vertex.="<<vertex<<
1018 "bz="<<fMagF<<
1019 "ptrel0="<<ptrel0<<
1020 "ptrel1="<<ptrel1<<
1021 "sectorIn="<<sectorIn<<
1022 "sectorOut="<<sectorOut<<
1023 "side="<<side<<
1024 // pid type
1025 "isMuon="<<isMuon<<
1026 "isProton="<<isProton<<
1027 "isKaon="<<isKaon<<
1028 "isPion="<<isPion<<
1029 "isElectron="<<isElectron<<
1030 "isMIP="<<isMIP<<
1031 "isHighPt="<<isHighPt<<
1032 "mass="<<mass<<
1033 "dedxDef="<<dedxDef<<
1034 "dedxDefPion="<<dedxDefPion<<
1035 "dedxDefKaon="<<dedxDefKaon<<
1036 "dedxDefProton="<<dedxDefProton<<
1037 //
1038 "nclITS="<<nclITS<<
1039 "ncl="<<ncl<<
1040 "ncl21="<<ncl21<<
1041 "ncl20="<<ncl20<<
1042 //
1043 "ncl20All="<<ncl20All<<
1044 "ncl20OROC="<<ncl20OROC<<
1045 "ncl20IROC="<<ncl20IROC<<
1046 "ncl20OROC0="<<ncl20OROC0<<
1047 "ncl20OROC1="<<ncl20OROC1<<
1048 //
1049 "ncl20All4="<<ncl20All4<<
1050 "ncl20OROC4="<<ncl20OROC4<<
1051 "ncl20IROC4="<<ncl20IROC4<<
1052 "ncl20OROC04="<<ncl20OROC04<<
1053 "ncl20OROC14="<<ncl20OROC14<<
1054 //
1055 "ncl20All3="<<ncl20All3<<
1056 "ncl20OROC3="<<ncl20OROC3<<
1057 "ncl20IROC3="<<ncl20IROC3<<
1058 "ncl20OROC03="<<ncl20OROC03<<
1059 "ncl20OROC13="<<ncl20OROC13<<
1060 //
1061 "ncl21All="<<ncl21All<<
1062 "ncl21OROC="<<ncl21OROC<<
1063 "ncl21IROC="<<ncl21IROC<<
1064 "ncl21OROC0="<<ncl21OROC0<<
1065 "ncl21OROC1="<<ncl21OROC1<<
1066 //track properties
1067 "langle0="<<langle[0]<<
1068 "langle1="<<langle[1]<<
1069 "langle2="<<langle[2]<<
1070 "gangle0="<<gangle[0]<< //global angle phi IROC
1071 "gangle1="<<gangle[1]<< // OROC medium
1072 "gangle2="<<gangle[2]<< // OROC long
1073 "L0="<<length[0]<<
1074 "L1="<<length[1]<<
1075 "L2="<<length[2]<<
1076 "p="<<p<<
1077 "pin="<<pin<<
1078 "pout="<<pout<<
1079 "tgl="<<tgl<<
1080 "track.="<<track<<
1081 "trackIn.="<<trackIn<<
1082 "trackOut.="<<trackOut<<
1083 "tpcOut.="<<tpcOut<<
4ebbe4be 1084 "seed.="<<seed<<
6feb400f 1085 "medianMIP0="<<medianMIP0<< // median MIP position as estimated from the array of (not only) "MIPS"
1086 //dedx
1087 "truncUp.="<<&truncUp<<
1088 "truncDown.="<<&truncDown<<
1089 "vecAllMax.="<<&vecAllMax<<
1090 "vecIROCMax.="<<&vecIROCMax<<
1091 "vecOROCMax.="<<&vecOROCMax<<
1092 "vecOROC0Max.="<<&vecOROC0Max<<
1093 "vecOROC1Max.="<<&vecOROC1Max<<
1094 //
1095 "vecAllTot.="<<&vecAllTot<<
1096 "vecIROCTot.="<<&vecIROCTot<<
1097 "vecOROCTot.="<<&vecOROCTot<<
1098 "vecOROC0Tot.="<<&vecOROC0Tot<<
1099 "vecOROC1Tot.="<<&vecOROC1Tot<<
1100 //
1101 "vecAllTotLog.="<<&vecAllTotLog<<
1102 "vecIROCTotLog.="<<&vecIROCTotLog<<
1103 "vecOROCTotLog.="<<&vecOROCTotLog<<
1104 "vecOROC0TotLog.="<<&vecOROC0TotLog<<
1105 "vecOROC1TotLog.="<<&vecOROC1TotLog<<
1106 //
1107 "vecAllTotUp.="<<&vecAllTotUp<<
1108 "vecIROCTotUp.="<<&vecIROCTotUp<<
1109 "vecOROCTotUp.="<<&vecOROCTotUp<<
1110 "vecOROC0TotUp.="<<&vecOROC0TotUp<<
1111 "vecOROC1TotUp.="<<&vecOROC1TotUp<<
1112 //
1113 "vecAllTotDown.="<<&vecAllTotDown<<
1114 "vecIROCTotDown.="<<&vecIROCTotDown<<
1115 "vecOROCTotDown.="<<&vecOROCTotDown<<
1116 "vecOROC0TotDown.="<<&vecOROC0TotDown<<
1117 "vecOROC1TotDown.="<<&vecOROC1TotDown<<
1118 //
1119 "vecAllTotRMS.="<<&vecAllTotRMS<<
1120 "vecIROCTotRMS.="<<&vecIROCTotRMS<<
1121 "vecOROCTotRMS.="<<&vecOROCTotRMS<<
1122 "vecOROC0TotRMS.="<<&vecOROC0TotRMS<<
1123 "vecOROC1TotRMS.="<<&vecOROC1TotRMS<<
1124 //
1125 "vecAllTotM2.="<<&vecAllTotM2<<
1126 "vecIROCTotM2.="<<&vecIROCTotM2<<
1127 "vecOROCTotM2.="<<&vecOROCTotM2<<
1128 "vecOROC0TotM2.="<<&vecOROC0TotM2<<
1129 "vecOROC1TotM2.="<<&vecOROC1TotM2<<
1130 //
1131 "vecAllTotMS.="<<&vecAllTotMS<<
1132 "vecIROCTotMS.="<<&vecIROCTotMS<<
1133 "vecOROCTotMS.="<<&vecOROCTotMS<<
1134 "vecOROC0TotMS.="<<&vecOROC0TotMS<<
1135 "vecOROC1TotMS.="<<&vecOROC1TotMS<<
1136 "\n";
1137 }
1138}
1139
1140
1141
1142
1143void AliTPCcalibGainMult::ProcessV0s(AliESDEvent * event){
1144 //
1145 // Select the K0s and gamma - and sign daughter products
1146 //
1147 TTreeSRedirector * pcstream = GetDebugStreamer();
1148 AliKFParticle::SetField(event->GetMagneticField());
1149 AliESDfriend *esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
1150 if (!esdFriend) {
1151 //Printf("ERROR: esdFriend not available");
1152 return;
1153 }
1154 if (esdFriend->TestSkipBit()) return;
1155 //
1156 //
cb4201d4 1157 static const TDatabasePDG *pdg = TDatabasePDG::Instance();
6feb400f 1158 const Double_t kChi2Cut=5;
1159 const Double_t kMinR=2;
1160 const Int_t kMinNcl=110;
1161 const Double_t kMinREl=5;
1162 const Double_t kMaxREl=70;
1163 //
1164 Int_t nv0 = event->GetNumberOfV0s();
1165 AliESDVertex *vertex= (AliESDVertex *)event->GetPrimaryVertex();
1166 AliKFVertex kfvertex=*vertex;
1167 //
1168 for (Int_t iv0=0;iv0<nv0;iv0++){
1169 AliESDv0 *v0 = event->GetV0(iv0);
1170 if (!v0) continue;
1171 if (v0->GetOnFlyStatus()<0.5) continue;
1172 if (v0->GetPindex()<0) continue;
1173 if (v0->GetNindex()<0) continue;
1174 if (TMath::Max(v0->GetPindex(), v0->GetNindex())>event->GetNumberOfTracks()) continue;
1175 //
1176 //
1177 AliExternalTrackParam pp=(v0->GetParamP()->GetSign()>0) ? (*(v0->GetParamP())):(*(v0->GetParamN()));
1178 AliExternalTrackParam pn=(v0->GetParamP()->GetSign()>0) ? (*(v0->GetParamN())):(*(v0->GetParamP()));
1179 AliKFParticle kfp1( pp, 211 );
1180 AliKFParticle kfp2( pn, -211 );
1181 //
1182 AliKFParticle *v0KFK0 = new AliKFParticle(kfp1,kfp2);
1183 AliKFParticle *v0KFK0CV = new AliKFParticle(*v0KFK0);
1184 v0KFK0CV->SetProductionVertex(kfvertex);
1185 v0KFK0CV->TransportToProductionVertex();
1186 Double_t chi2K0 = v0KFK0CV->GetChi2();
1187 if (chi2K0>kChi2Cut) continue;
1188 if (v0->GetRr()<kMinR) continue;
1189 Bool_t isOKC=TMath::Max(v0->GetCausalityP()[0],v0->GetCausalityP()[1])<0.7&&TMath::Min(v0->GetCausalityP()[2],v0->GetCausalityP()[3])>0.2;
1190 //
1191 Double_t effMass22=v0->GetEffMass(2,2);
1192 Double_t effMass42=v0->GetEffMass(4,2);
1193 Double_t effMass24=v0->GetEffMass(2,4);
1194 Double_t effMass00=v0->GetEffMass(0,0);
1195 AliKFParticle *v0KFK0CVM = new AliKFParticle(*v0KFK0CV);
1196 v0KFK0CVM->SetMassConstraint(pdg->GetParticle("K_S0")->Mass());
1197 Bool_t isV0= kFALSE;
1198 //
1199 Double_t d22 = TMath::Abs(effMass22-pdg->GetParticle("K_S0")->Mass());
1200 Double_t d42 = TMath::Abs(effMass42-pdg->GetParticle("Lambda0")->Mass());
1201 Double_t d24 = TMath::Abs(effMass24-pdg->GetParticle("Lambda0")->Mass());
1202 Double_t d00 = TMath::Abs(effMass00);
1203 //
1204 Bool_t isKaon = d22<0.01 && d22< 0.3 * TMath::Min(TMath::Min(d42,d24),d00);
1205 Bool_t isLambda = d42<0.01 && d42< 0.3 * TMath::Min(TMath::Min(d22,d24),d00);
1206 Bool_t isAntiLambda= d24<0.01 && d24< 0.3 * TMath::Min(TMath::Min(d22,d42),d00);
1207 Bool_t isGamma = d00<0.02 && d00< 0.3 * TMath::Min(TMath::Min(d42,d24),d22);
1208 //
1209 if (isGamma && (isKaon||isLambda||isAntiLambda)) continue;
1210 if (isLambda && (isKaon||isGamma||isAntiLambda)) continue;
1211 if (isKaon && (isLambda||isGamma||isAntiLambda)) continue;
1212 Double_t sign= v0->GetParamP()->GetSign()* v0->GetParamN()->GetSign();
1213 if (sign>0) continue;
1214 isV0=isKaon||isLambda||isAntiLambda||isGamma;
1215 if (!(isV0)) continue;
1216 if (isGamma&&v0->GetRr()<kMinREl) continue;
1217 if (isGamma&&v0->GetRr()>kMaxREl) continue;
1218 if (!isOKC) continue;
1219 //
1220 Int_t pindex = (v0->GetParamP()->GetSign()>0) ? v0->GetPindex() : v0->GetNindex();
1221 Int_t nindex = (v0->GetParamP()->GetSign()>0) ? v0->GetNindex() : v0->GetPindex();
1222 AliESDtrack * trackP = event->GetTrack(pindex);
1223 AliESDtrack * trackN = event->GetTrack(nindex);
1224 if (!trackN) continue;
1225 if (!trackP) continue;
1226 Int_t nclP= (Int_t)trackP->GetTPCClusterInfo(2,1);
1227 Int_t nclN= (Int_t)trackN->GetTPCClusterInfo(2,1);
1228 if (TMath::Min(nclP,nclN)<kMinNcl) continue;
1229 Double_t eta = TMath::Max(TMath::Abs(trackP->Eta()), TMath::Abs(trackN->Eta()));
1230 if (TMath::Abs(eta)>1) continue;
1231 //
1232 //
1233 AliESDfriendTrack *friendTrackP = esdFriend->GetTrack(pindex);
1234 AliESDfriendTrack *friendTrackN = esdFriend->GetTrack(nindex);
1235 if (!friendTrackP) continue;
1236 if (!friendTrackN) continue;
1237 TObject *calibObject;
1238 AliTPCseed *seedP = 0;
1239 AliTPCseed *seedN = 0;
1240 for (Int_t l=0;(calibObject=friendTrackP->GetCalibObject(l));++l) {
1241 if ((seedP=dynamic_cast<AliTPCseed*>(calibObject))) break;
1242 }
1243 for (Int_t l=0;(calibObject=friendTrackN->GetCalibObject(l));++l) {
1244 if ((seedN=dynamic_cast<AliTPCseed*>(calibObject))) break;
1245 }
1246 if (isGamma){
1247 if ( TMath::Abs((trackP->GetTPCsignal()/(trackN->GetTPCsignal()+0.0001)-1)>0.3)) continue;
1248 }
1249 if (isGamma) (*fPIDMatrix)(pindex, 0)+=2;
1250 if (isGamma) (*fPIDMatrix)(nindex, 0)+=2;
1251 //
1252 if (isKaon) (*fPIDMatrix)(pindex, 2)+=2;
1253 if (isKaon) (*fPIDMatrix)(nindex, 2)+=2;
1254 //
1255 //
1256 if (pcstream){
1257 (*pcstream)<<"v0s"<<
1258 "isGamma="<<isGamma<<
1259 "isKaon="<<isKaon<<
1260 "isLambda="<<isLambda<<
1261 "isAntiLambda="<<isAntiLambda<<
1262 "chi2="<<chi2K0<<
1263 "trackP.="<<trackP<<
1264 "trackN.="<<trackN<<
1265 "v0.="<<v0<<
1266 "\n";
1267 }
1268 }
1269}
1270
1271
1272
1273
1274void AliTPCcalibGainMult::ProcessCosmic(const AliESDEvent * event) {
1275 //
1276 // Find cosmic pairs trigger by random trigger
1277 //
1278 //
1279 AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;
1280 AliTPCParam *param = AliTPCcalibDB::Instance()->GetParameters();
1281
1282 AliESDVertex *vertexSPD = (AliESDVertex *)event->GetPrimaryVertexSPD();
1283 AliESDVertex *vertexTPC = (AliESDVertex *)event->GetPrimaryVertexTPC();
1284 AliESDfriend *esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
1285 const Double_t kMinPt=4;
1286 const Double_t kMinPtMax=0.8;
1287 const Double_t kMinNcl=159*0.5;
1288 const Double_t kMaxDelta[5]={2,600,0.02,0.02,0.1};
1289 Int_t ntracks=event->GetNumberOfTracks();
1290 const Double_t kMaxImpact=80;
1291 // Float_t dcaTPC[2]={0,0};
1292 // Float_t covTPC[3]={0,0,0};
1293
1294 UInt_t specie = event->GetEventSpecie(); // skip laser events
1295 if (specie==AliRecoParam::kCalib) return;
1296
1297
1298 for (Int_t itrack0=0;itrack0<ntracks;itrack0++) {
1299 AliESDtrack *track0 = event->GetTrack(itrack0);
1300 if (!track0) continue;
1301 if (!track0->IsOn(AliESDtrack::kTPCrefit)) continue;
1302
1303 if (TMath::Abs(AliTracker::GetBz())>1&&track0->Pt()<kMinPt) continue;
1304 if (track0->GetTPCncls()<kMinNcl) continue;
1305 if (TMath::Abs(track0->GetY())<2*kMaxDelta[0]) continue;
1306 if (TMath::Abs(track0->GetY())>kMaxImpact) continue;
1307 if (track0->GetKinkIndex(0)>0) continue;
1308 const Double_t * par0=track0->GetParameter(); //track param at rhe DCA
1309 //rm primaries
1310 //
1311 for (Int_t itrack1=itrack0+1;itrack1<ntracks;itrack1++) {
1312 AliESDtrack *track1 = event->GetTrack(itrack1);
1313 if (!track1) continue;
1314 if (!track1->IsOn(AliESDtrack::kTPCrefit)) continue;
1315 if (track1->GetKinkIndex(0)>0) continue;
1316 if (TMath::Abs(AliTracker::GetBz())>1&&track1->Pt()<kMinPt) continue;
1317 if (track1->GetTPCncls()<kMinNcl) continue;
1318 if (TMath::Abs(AliTracker::GetBz())>1&&TMath::Max(track1->Pt(), track0->Pt())<kMinPtMax) continue;
1319 if (TMath::Abs(track1->GetY())<2*kMaxDelta[0]) continue;
1320 if (TMath::Abs(track1->GetY())>kMaxImpact) continue;
1321 //
1322 const Double_t* par1=track1->GetParameter(); //track param at rhe DCA
1323 //
1324 Bool_t isPair=kTRUE;
1325 for (Int_t ipar=0; ipar<5; ipar++){
1326 if (ipar==4&&TMath::Abs(AliTracker::GetBz())<1) continue; // 1/pt not defined for B field off
1327 if (TMath::Abs(TMath::Abs(par0[ipar])-TMath::Abs(par1[ipar]))>kMaxDelta[ipar]) isPair=kFALSE;
1328 }
1329 if (!isPair) continue;
1330 if (TMath::Abs(TMath::Abs(track0->GetAlpha()-track1->GetAlpha())-TMath::Pi())>kMaxDelta[2]) isPair=kFALSE;
1331 //delta with correct sign
1332 if (TMath::Abs(par0[0]+par1[0])>kMaxDelta[0]) isPair=kFALSE; //delta y opposite sign
1333 if (TMath::Abs(par0[3]+par1[3])>kMaxDelta[3]) isPair=kFALSE; //delta tgl opposite sign
1334 if (TMath::Abs(AliTracker::GetBz())>1 && TMath::Abs(par0[4]+par1[4])>kMaxDelta[4]) isPair=kFALSE; //delta 1/pt opposite sign
1335 if (!isPair) continue;
1336 TString filename(AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile()->GetName());
1337 Int_t eventNumber = event->GetEventNumberInFile();
1338 Bool_t hasFriend=(esdFriend) ? (esdFriend->GetTrack(itrack0)!=0):0;
1339 Bool_t hasITS=(track0->GetNcls(0)+track1->GetNcls(0)>4);
1340 printf("DUMPHPTCosmic:%s|%f|%d|%d|%d\n",filename.Data(),(TMath::Min(track0->Pt(),track1->Pt())), eventNumber,hasFriend,hasITS);
1341 //
1342 //
1343 TTreeSRedirector * pcstream = GetDebugStreamer();
1344 Int_t ntracksSPD = vertexSPD->GetNContributors();
1345 Int_t ntracksTPC = vertexTPC->GetNContributors();
1346 //
1347 if (pcstream){
1348 (*pcstream)<<"cosmicPairsAll"<<
1349 "run="<<fRun<< // run number
1350 "event="<<fEvent<< // event number
1351 "time="<<fTime<< // time stamp of event
1352 "trigger="<<fTrigger<< // trigger
1353 "triggerClass="<<&fTriggerClass<< // trigger
1354 "bz="<<fMagF<< // magnetic field
1355 //
1356 "nSPD="<<ntracksSPD<<
1357 "nTPC="<<ntracksTPC<<
1358 "vSPD.="<<vertexSPD<< //primary vertex -SPD
1359 "vTPC.="<<vertexTPC<< //primary vertex -TPC
1360 "t0.="<<track0<< //track0
1361 "t1.="<<track1<< //track1
1362 "\n";
1363 }
1364 //
1365 AliESDfriendTrack *friendTrack0 = esdFriend->GetTrack(itrack0);
1366 if (!friendTrack0) continue;
1367 AliESDfriendTrack *friendTrack1 = esdFriend->GetTrack(itrack1);
1368 if (!friendTrack1) continue;
1369 TObject *calibObject;
1370 AliTPCseed *seed0 = 0;
1371 AliTPCseed *seed1 = 0;
1372 //
1373 for (Int_t l=0;(calibObject=friendTrack0->GetCalibObject(l));++l) {
1374 if ((seed0=dynamic_cast<AliTPCseed*>(calibObject))) break;
1375 }
1376 for (Int_t l=0;(calibObject=friendTrack1->GetCalibObject(l));++l) {
1377 if ((seed1=dynamic_cast<AliTPCseed*>(calibObject))) break;
1378 }
1379 //
1380 if (pcstream){
1381 (*pcstream)<<"cosmicPairs"<<
1382 "run="<<fRun<< // run number
1383 "event="<<fEvent<< // event number
1384 "time="<<fTime<< // time stamp of event
1385 "trigger="<<fTrigger<< // trigger
1386 "triggerClass="<<&fTriggerClass<< // trigger
1387 "bz="<<fMagF<< // magnetic field
1388 //
1389 "nSPD="<<ntracksSPD<<
1390 "nTPC="<<ntracksTPC<<
1391 "vSPD.="<<vertexSPD<< //primary vertex -SPD
1392 "vTPC.="<<vertexTPC<< //primary vertex -TPC
1393 "t0.="<<track0<< //track0
1394 "t1.="<<track1<< //track1
1395 "ft0.="<<friendTrack0<< //track0
1396 "ft1.="<<friendTrack1<< //track1
1397 "s0.="<<seed0<< //track0
1398 "s1.="<<seed1<< //track1
1399 "\n";
1400 }
1401 if (!seed0) continue;
1402 if (!seed1) continue;
1403 Int_t nclA0=0, nclC0=0; // number of clusters
1404 Int_t nclA1=0, nclC1=0; // number of clusters
1405 //
1406 for (Int_t irow=0; irow<159; irow++){
1407 AliTPCclusterMI *cluster0=seed0->GetClusterPointer(irow);
1408 AliTPCclusterMI *cluster1=seed1->GetClusterPointer(irow);
1409 if (cluster0){
1410 if (cluster0->GetQ()>0 && cluster0->GetDetector()%36<18) nclA0++;
1411 if (cluster0->GetQ()>0 && cluster0->GetDetector()%36>=18) nclC0++;
1412 }
1413 if (cluster1){
1414 if (cluster1->GetQ()>0 && cluster1->GetDetector()%36<18) nclA1++;
1415 if (cluster1->GetQ()>0 && cluster1->GetDetector()%36>=18) nclC1++;
1416 }
1417 }
1418 Int_t cosmicType=0; // types of cosmic topology
1419 if ((nclA0>nclC0) && (nclA1>nclC1)) cosmicType=0; // AA side
1420 if ((nclA0<nclC0) && (nclA1<nclC1)) cosmicType=1; // CC side
1421 if ((nclA0>nclC0) && (nclA1<nclC1)) cosmicType=2; // AC side
1422 if ((nclA0<nclC0) && (nclA1>nclC1)) cosmicType=3; // CA side
1423 if (cosmicType<2) continue; // use only crossing tracks
1424 //
1425 Double_t deltaTimeCluster=0;
1426 deltaTimeCluster=0.5*(track1->GetZ()-track0->GetZ())/param->GetZWidth();
1427 if (nclA0>nclC0) deltaTimeCluster*=-1; // if A side track
1428 //
1429 for (Int_t irow=0; irow<159; irow++){
1430 AliTPCclusterMI *cluster0=seed0->GetClusterPointer(irow);
1431 if (cluster0 &&cluster0->GetX()>10){
1432 Double_t x0[3]={cluster0->GetRow(),cluster0->GetPad(),cluster0->GetTimeBin()+deltaTimeCluster};
1433 Int_t index0[1]={cluster0->GetDetector()};
1434 transform->Transform(x0,index0,0,1);
1435 cluster0->SetX(x0[0]);
1436 cluster0->SetY(x0[1]);
1437 cluster0->SetZ(x0[2]);
1438 //
1439 }
1440 AliTPCclusterMI *cluster1=seed1->GetClusterPointer(irow);
1441 if (cluster1&&cluster1->GetX()>10){
1442 Double_t x1[3]={cluster1->GetRow(),cluster1->GetPad(),cluster1->GetTimeBin()+deltaTimeCluster};
1443 Int_t index1[1]={cluster1->GetDetector()};
1444 transform->Transform(x1,index1,0,1);
1445 cluster1->SetX(x1[0]);
1446 cluster1->SetY(x1[1]);
1447 cluster1->SetZ(x1[2]);
1448 }
1449 }
1450 //
1451 //
1452 if (fPIDMatrix){
1453 (*fPIDMatrix)(itrack0,1)+=4; //
1454 (*fPIDMatrix)(itrack1,1)+=4; //
1455 }
1456 }
1457 }
1458}
1459
1460
1461
1462void AliTPCcalibGainMult::ProcessKinks(const AliESDEvent * event){
1463 //
1464 //
1465 //
1466 AliKFParticle::SetField(event->GetMagneticField());
1467 AliESDfriend *esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
1468 if (!esdFriend) {
1469 //Printf("ERROR: esdFriend not available");
1470 return;
1471 }
1472 // if (esdFriend->TestSkipBit()) return;
1473 //
1474 //
1475 const Double_t kChi2Cut=10;
1476 const Double_t kMinR=100;
1477 const Double_t kMaxR=230;
1478 const Int_t kMinNcl=110;
1479 //
1480 Int_t nkinks = event->GetNumberOfKinks();
1481 AliESDVertex *vertex= (AliESDVertex *)event->GetPrimaryVertex();
1482 AliKFVertex kfvertex=*vertex;
1483 TTreeSRedirector * pcstream = GetDebugStreamer();
1484 //
1485 for (Int_t ikink=0;ikink<nkinks;ikink++){
1486 AliESDkink *kink = event->GetKink(ikink);
1487 if (!kink) continue;
1488 if (kink->GetIndex(0)<0) continue;
1489 if (kink->GetIndex(1)<0) continue;
1490 if (TMath::Max(kink->GetIndex(1), kink->GetIndex(0))>event->GetNumberOfTracks()) continue;
1491 //
1492 //
1493 AliExternalTrackParam pd=kink->RefParamDaughter();
1494 AliExternalTrackParam pm=kink->RefParamMother();
1495 AliKFParticle kfpd( pd, 211 );
1496 AliKFParticle kfpm( pm, -13 );
1497 //
1498 AliKFParticle *v0KF = new AliKFParticle(kfpm,kfpd);
1499 v0KF->SetVtxGuess(kink->GetPosition()[0],kink->GetPosition()[1],kink->GetPosition()[2]);
1500 Double_t chi2 = v0KF->GetChi2();
1501 AliESDtrack * trackM = event->GetTrack(kink->GetIndex(0));
1502 AliESDtrack * trackD = event->GetTrack(kink->GetIndex(1));
1503 if (!trackM) continue;
1504 if (!trackD) continue;
1505 Int_t nclM= (Int_t)trackM->GetTPCClusterInfo(2,1);
1506 Int_t nclD= (Int_t)trackD->GetTPCClusterInfo(2,1);
1507 Double_t eta = TMath::Max(TMath::Abs(trackM->Eta()), TMath::Abs(trackD->Eta()));
1508 Double_t kx= v0KF->GetX();
1509 Double_t ky= v0KF->GetY();
1510 Double_t kz= v0KF->GetZ();
1511 Double_t ex= v0KF->GetErrX();
1512 Double_t ey= v0KF->GetErrY();
1513 Double_t ez= v0KF->GetErrZ();
1514 //
1515 Double_t radius=TMath::Sqrt(kx*kx+ky*ky);
1516 Double_t alpha=TMath::ATan2(ky,kx);
1517 if (!pd.Rotate(alpha)) continue;
1518 if (!pm.Rotate(alpha)) continue;
1519 if (!pd.PropagateTo(radius,event->GetMagneticField())) continue;
1520 if (!pm.PropagateTo(radius,event->GetMagneticField())) continue;
1521 Double_t pos[2]={0,kz};
1522 Double_t cov[3]={ex*ex+ey*ey,0,ez*ez};
1523 pd.Update(pos,cov);
1524 pm.Update(pos,cov);
1525 //
1526 if (pcstream){
1527 (*pcstream)<<"kinks"<<
1528 "eta="<<eta<<
1529 "nclM="<<nclM<<
1530 "nclD="<<nclD<<
1531 "kink.="<<kink<<
1532 "trackM.="<<trackM<<
1533 "trackD.="<<trackD<<
1534 "pm.="<<&pm<< //updated parameters
1535 "pd.="<<&pd<< // updated parameters
1536 "v0KF.="<<v0KF<<
1537 "chi2="<<chi2<<
1538 "\n";
1539 }
1540 /*
1541 TCut cutQ="chi2<10&&kink.fRr>90&&kink.fRr<220";
1542 TCut cutRD="20*sqrt(pd.fC[14])<abs(pd.fP[4])&&trackD.fTPCsignal>10&&trackD.fTPCsignalN>50";
1543
1544 */
1545 //
1546 if (chi2>kChi2Cut) continue;
1547 if (kink->GetR()<kMinR) continue;
1548 if (kink->GetR()>kMaxR) continue;
1549 if ((nclM+nclD)<kMinNcl) continue;
1550 if (TMath::Abs(eta)>1) continue;
1551 //
1552 //
1553 AliESDfriendTrack *friendTrackM = esdFriend->GetTrack(kink->GetIndex(0));
1554 AliESDfriendTrack *friendTrackD = esdFriend->GetTrack(kink->GetIndex(1));
1555 if (!friendTrackM) continue;
1556 if (!friendTrackD) continue;
1557 TObject *calibObject;
1558 AliTPCseed *seedM = 0;
1559 AliTPCseed *seedD = 0;
1560 for (Int_t l=0;(calibObject=friendTrackM->GetCalibObject(l));++l) {
1561 if ((seedM=dynamic_cast<AliTPCseed*>(calibObject))) break;
1562 }
1563 for (Int_t l=0;(calibObject=friendTrackD->GetCalibObject(l));++l) {
1564 if ((seedD=dynamic_cast<AliTPCseed*>(calibObject))) break;
1565 }
1566 }
1567}
1568
1569void AliTPCcalibGainMult::DumpHPT(const AliESDEvent * event){
1570 //
1571 // Function to select the HPT tracks and events
1572 // It is used to select event with HPT - list used later for the raw data downloading
1573 // - and reconstruction
1574 // Not actualy used for the calibration of the data
1575
1576 TTreeSRedirector * pcstream = GetDebugStreamer();
1577 AliKFParticle::SetField(event->GetMagneticField());
1578 AliESDfriend *esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
1579 if (!esdFriend) {
1580 //Printf("ERROR: esdFriend not available");
1581 return;
1582 }
1583 if (esdFriend->TestSkipBit()) return;
1584
1585 Int_t ntracks=event->GetNumberOfTracks();
1586 //
1587 for (Int_t i=0;i<ntracks;++i) {
1588 //
1589 AliESDtrack *track = event->GetTrack(i);
1590 if (!track) continue;
1591 if (track->Pt()<4) continue;
1592 UInt_t status = track->GetStatus();
1593 //
1594 AliExternalTrackParam * trackIn = (AliExternalTrackParam *)track->GetInnerParam();
1595 if (!trackIn) continue;
1596 if ((status&AliESDtrack::kTPCrefit)==0) continue;
1597 if ((status&AliESDtrack::kITSrefit)==0) continue;
1598 AliESDfriendTrack *friendTrack = esdFriend->GetTrack(i);
1599 if (!friendTrack) continue;
1600 AliExternalTrackParam * itsOut = (AliExternalTrackParam *)(friendTrack->GetITSOut());
1601 if (!itsOut) continue;
1602 AliExternalTrackParam * itsOut2 = (AliExternalTrackParam *)(friendTrack->GetITSOut()->Clone());
1603 AliExternalTrackParam * tpcIn2 = (AliExternalTrackParam *)(trackIn->Clone());
1604 if (!itsOut2->Rotate(trackIn->GetAlpha())) continue;
1605 //Double_t xmiddle=0.5*(itsOut2->GetX()+tpcIn2->GetX());
1606 Double_t xmiddle=(itsOut2->GetX());
1607 if (!itsOut2->PropagateTo(xmiddle,event->GetMagneticField())) continue;
1608 if (!tpcIn2->PropagateTo(xmiddle,event->GetMagneticField())) continue;
1609 //
1610 AliExternalTrackParam * tpcInner = (AliExternalTrackParam *)(track->GetTPCInnerParam());
1611 if (!tpcInner) continue;
1612 tpcInner->Rotate(track->GetAlpha());
1613 tpcInner->PropagateTo(track->GetX(),event->GetMagneticField());
1614 //
1615 // tpc constrained
1616 //
1617 AliExternalTrackParam * tpcInnerC = (AliExternalTrackParam *)(track->GetTPCInnerParam()->Clone());
1618 if (!tpcInnerC) continue;
1619 tpcInnerC->Rotate(track->GetAlpha());
1620 tpcInnerC->PropagateTo(track->GetX(),event->GetMagneticField());
1621 Double_t dz[2],cov[3];
1622 AliESDVertex *vtx= (AliESDVertex *)event->GetPrimaryVertex();
1623
1624 if (!tpcInnerC->PropagateToDCA(vtx, event->GetMagneticField(), 3, dz, cov)) continue;
1625 Double_t covar[6]; vtx->GetCovMatrix(covar);
1626 Double_t p[2]={tpcInnerC->GetParameter()[0]-dz[0],tpcInnerC->GetParameter()[1]-dz[1]};
1627 Double_t c[3]={covar[2],0.,covar[5]};
1628 //
1629 Double_t chi2C=tpcInnerC->GetPredictedChi2(p,c);
1630 tpcInnerC->Update(p,c);
1631
1632 if (pcstream){
1633 (*pcstream)<<"hpt"<<
1634 "run="<<fRun<<
1635 "time="<<fTime<<
1636 "vertex="<<vtx<<
1637 "bz="<<fMagF<<
1638 "track.="<<track<<
1639 "tpcInner.="<<tpcInner<<
1640 "tpcInnerC.="<<tpcInnerC<<
1641 "chi2C="<<chi2C<<
1642 //
1643 "its.="<<itsOut<<
1644 "its2.="<<itsOut2<<
1645 "tpc.="<<trackIn<<
1646 "tpc2.="<<tpcIn2<<
1647 "\n";
1648 }
1649 }
1650}
1651
1652
1653
1654void AliTPCcalibGainMult::ProcessTOF(const AliESDEvent * event){
1655 //
1656 // 1. Loop over tracks
1657 // 2. Fit T0
1658 // 3. Sign positivelly identified tracks
1659 //
1660 const Double_t kMaxDelta=1000;
1661 const Double_t kOrbit=50000; // distance in the time beween two orbits in the TOF time units - 50000=50 ns
1662 const Double_t kMaxD=20;
1663 const Double_t kRMS0=200;
1664 const Double_t kMaxDCAZ=10;
1665 AliESDVertex *vtx= (AliESDVertex *)event->GetPrimaryVertex();
1666 //
1667 TTreeSRedirector * pcstream = GetDebugStreamer();
1668 AliKFParticle::SetField(event->GetMagneticField());
1669 AliESDfriend *esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
1670 if (!esdFriend) {
1671 //Printf("ERROR: esdFriend not available");
1672 return;
1673 }
1674 //if (esdFriend->TestSkipBit()) return;
1675
1676 Int_t ntracks=event->GetNumberOfTracks();
1677 //
1678 Double_t deltaTPion[10000];
1679 Double_t medianT0=0;
1680 Double_t meanT0=0;
1681 Double_t rms=10000;
1682 Int_t counter=0;
1683 //
1684 // Get Median time for pion hypothesy
1685 //
1686 for (Int_t iter=0; iter<3; iter++){
1687 counter=0;
1688 for (Int_t i=0;i<ntracks;++i) {
1689 //
1690 AliESDtrack *track = event->GetTrack(i);
1691 if (!track) continue;
1692 if (!track->IsOn(AliESDtrack::kTIME)) continue;
1693 if (TMath::Abs(track->GetZ())>kMaxDCAZ) continue; // remove overlaped events
1694 if (TMath::Abs(track->GetTOFsignalDz())>kMaxD) continue;
1695 Double_t times[1000];
1696 track->GetIntegratedTimes(times);
1697 Int_t norbit=TMath::Nint((track->GetTOFsignal()-times[2])/kOrbit);
1698 Double_t torbit=norbit*kOrbit;
1699 if (iter==1 &&TMath::Abs(times[2]-times[3])<3*rms) continue; // skip umbigous points - kaon pion
1700 //
1701 Int_t indexBest=2;
1702 if (iter>1){
1703 for (Int_t j=3; j<5; j++)
1704 if (TMath::Abs(track->GetTOFsignal()-times[j]-torbit-medianT0)<TMath::Abs(track->GetTOFsignal()-times[j]-torbit-medianT0)) indexBest=j;
1705 }
1706 //
1707 if (iter>0) if (TMath::Abs(track->GetTOFsignal()-times[indexBest]-torbit-medianT0)>3*(kRMS0+rms)) continue;
1708 if (iter>0) if (TMath::Abs(track->GetTOFsignal()-times[indexBest]-torbit-medianT0)>kMaxDelta) continue;
1709 deltaTPion[counter]=track->GetTOFsignal()-times[indexBest]-torbit;
1710 counter++;
1711 }
1712 if (counter<2) return;
1713 medianT0=TMath::Median(counter,deltaTPion);
1714 meanT0=TMath::Median(counter,deltaTPion);
1715 rms=TMath::RMS(counter,deltaTPion);
1716 }
1717 if (counter<3) return;
1718 //
1719 // Dump
1720 //
1721 for (Int_t i=0;i<ntracks;++i) {
1722 //
1723 AliESDtrack *track = event->GetTrack(i);
1724 if (!track) continue;
1725 if (!track->IsOn(AliESDtrack::kTIME)) continue;
1726 if (TMath::Abs(track->GetZ())>kMaxDCAZ) continue; //remove overlapped events
1727 if (TMath::Abs(track->GetTOFsignalDz())>kMaxD) continue;
1728 Double_t times[1000];
1729 track->GetIntegratedTimes(times);
1730 Int_t norbit=TMath::Nint((track->GetTOFsignal()-times[2])/kOrbit);
1731 Double_t torbit=norbit*kOrbit;
1732 if (rms<=0) continue;
1733 //
1734 Double_t tPion = (track->GetTOFsignal()-times[2]-medianT0-torbit);
1735 Double_t tKaon = (track->GetTOFsignal()-times[3]-medianT0-torbit);
1736 Double_t tProton= (track->GetTOFsignal()-times[4]-medianT0-torbit);
1737 Double_t tElectron= (track->GetTOFsignal()-times[0]-medianT0-torbit);
1738 //
1739 Bool_t isPion = (TMath::Abs(tPion/rms)<6) && TMath::Abs(tPion)<(TMath::Min(TMath::Abs(tKaon), TMath::Abs(tProton))-rms);
1740 Bool_t isKaon = (TMath::Abs(tKaon/rms)<3) && TMath::Abs(tKaon)<0.2*(TMath::Min(TMath::Abs(tPion), TMath::Abs(tProton))-3*rms);
1741 Bool_t isProton = (TMath::Abs(tProton/rms)<6) && TMath::Abs(tProton)<0.5*(TMath::Min(TMath::Abs(tKaon), TMath::Abs(tPion))-rms);
1742 Bool_t isElectron = (TMath::Abs(tElectron/rms)<6) && TMath::Abs(tElectron)<0.2*(TMath::Min(TMath::Abs(tKaon), TMath::Abs(tPion))-rms) &&TMath::Abs(tElectron)<0.5*(TMath::Abs(tPion)-rms);
1743
1744 if (isPion) (*fPIDMatrix)(i,2)+=1;
1745 if (isKaon) (*fPIDMatrix)(i,3)+=1;
1746 if (isProton) (*fPIDMatrix)(i,4)+=1;
1747 // if (isElectron) (*fPIDMatrix)(i,0)+=1;
1748 //
1749 if (pcstream){
1750 // debug streamer to dump the information
1751 (*pcstream)<<"tof"<<
1752 "isPion="<<isPion<<
1753 "isKaon="<<isKaon<<
1754 "isProton="<<isProton<<
1755 "isElectron="<<isElectron<<
1756 //
1757 "counter="<<counter<<
1758 "torbit="<<torbit<<
1759 "norbit="<<norbit<<
1760 "medianT0="<<medianT0<<
1761 "meanT0="<<meanT0<<
1762 "rmsT0="<<rms<<
1763 "track.="<<track<<
1764 "vtx.="<<vtx<<
1765 "\n";
1766 }
1767
1768 }
1769 /*
1770 tof->SetAlias("isProton","(abs(track.fTOFsignal-track.fTrackTime[4]-medianT0-torbit)<(0.5*abs(track.fTOFsignal-track.fTrackTime[3]-medianT0-torbit)-rmsT0))");
1771 tof->SetAlias("isPion","(abs(track.fTOFsignal-track.fTrackTime[2]-medianT0-torbit)<(0.5*abs(track.fTOFsignal-track.fTrackTime[3]-medianT0-torbit)-rmsT0))");
1772 tof->SetAlias("isKaon","(abs(track.fTOFsignal-track.fTrackTime[3]-medianT0-torbit)<(0.5*abs(track.fTOFsignal-track.fTrackTime[2]-medianT0-torbit)-rmsT0))&&(abs(track.fTOFsignal-track.fTrackTime[3]-medianT0-torbit)<(0.5*abs(track.fTOFsignal-track.fTrackTime[4]-medianT0-torbit)-rmsT0))");
1773
1774 */
1775
1776}
1777
1778
1779// void AliTPCcalibGainMult::Terminate(){
1780// //
1781// // Terminate function
1782// // call base terminate + Eval of fitters
1783// //
1784// Info("AliTPCcalibGainMult","Terminate");
1785// TTreeSRedirector *pcstream = GetDebugStreamer();
1786// if (pcstream){
1787// TTreeStream &stream = (*pcstream)<<"dump";
1788// TTree* tree = stream.GetTree();
1789// if (tree) if ( tree->GetEntries()>0){
1790// TObjArray *array = tree->GetListOfBranches();
1791// for (Int_t i=0; i<array->GetEntries(); i++) {TBranch * br = (TBranch *)array->At(i); br->SetAddress(0);}
1792// gDirectory=gROOT;
1793// fdEdxTree=tree->CloneTree(10000);
1794// }
1795// }
1796// AliTPCcalibBase::Terminate();
1797// }
1798