]> 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
Remove self-assignment
[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 //
670 AliTPCROC *roc = AliTPCROC::Instance();
671 TDatabasePDG *pdg = TDatabasePDG::Instance();
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<<
1084 "medianMIP0="<<medianMIP0<< // median MIP position as estimated from the array of (not only) "MIPS"
1085 //dedx
1086 "truncUp.="<<&truncUp<<
1087 "truncDown.="<<&truncDown<<
1088 "vecAllMax.="<<&vecAllMax<<
1089 "vecIROCMax.="<<&vecIROCMax<<
1090 "vecOROCMax.="<<&vecOROCMax<<
1091 "vecOROC0Max.="<<&vecOROC0Max<<
1092 "vecOROC1Max.="<<&vecOROC1Max<<
1093 //
1094 "vecAllTot.="<<&vecAllTot<<
1095 "vecIROCTot.="<<&vecIROCTot<<
1096 "vecOROCTot.="<<&vecOROCTot<<
1097 "vecOROC0Tot.="<<&vecOROC0Tot<<
1098 "vecOROC1Tot.="<<&vecOROC1Tot<<
1099 //
1100 "vecAllTotLog.="<<&vecAllTotLog<<
1101 "vecIROCTotLog.="<<&vecIROCTotLog<<
1102 "vecOROCTotLog.="<<&vecOROCTotLog<<
1103 "vecOROC0TotLog.="<<&vecOROC0TotLog<<
1104 "vecOROC1TotLog.="<<&vecOROC1TotLog<<
1105 //
1106 "vecAllTotUp.="<<&vecAllTotUp<<
1107 "vecIROCTotUp.="<<&vecIROCTotUp<<
1108 "vecOROCTotUp.="<<&vecOROCTotUp<<
1109 "vecOROC0TotUp.="<<&vecOROC0TotUp<<
1110 "vecOROC1TotUp.="<<&vecOROC1TotUp<<
1111 //
1112 "vecAllTotDown.="<<&vecAllTotDown<<
1113 "vecIROCTotDown.="<<&vecIROCTotDown<<
1114 "vecOROCTotDown.="<<&vecOROCTotDown<<
1115 "vecOROC0TotDown.="<<&vecOROC0TotDown<<
1116 "vecOROC1TotDown.="<<&vecOROC1TotDown<<
1117 //
1118 "vecAllTotRMS.="<<&vecAllTotRMS<<
1119 "vecIROCTotRMS.="<<&vecIROCTotRMS<<
1120 "vecOROCTotRMS.="<<&vecOROCTotRMS<<
1121 "vecOROC0TotRMS.="<<&vecOROC0TotRMS<<
1122 "vecOROC1TotRMS.="<<&vecOROC1TotRMS<<
1123 //
1124 "vecAllTotM2.="<<&vecAllTotM2<<
1125 "vecIROCTotM2.="<<&vecIROCTotM2<<
1126 "vecOROCTotM2.="<<&vecOROCTotM2<<
1127 "vecOROC0TotM2.="<<&vecOROC0TotM2<<
1128 "vecOROC1TotM2.="<<&vecOROC1TotM2<<
1129 //
1130 "vecAllTotMS.="<<&vecAllTotMS<<
1131 "vecIROCTotMS.="<<&vecIROCTotMS<<
1132 "vecOROCTotMS.="<<&vecOROCTotMS<<
1133 "vecOROC0TotMS.="<<&vecOROC0TotMS<<
1134 "vecOROC1TotMS.="<<&vecOROC1TotMS<<
1135 "\n";
1136 }
1137}
1138
1139
1140
1141
1142void AliTPCcalibGainMult::ProcessV0s(AliESDEvent * event){
1143 //
1144 // Select the K0s and gamma - and sign daughter products
1145 //
1146 TTreeSRedirector * pcstream = GetDebugStreamer();
1147 AliKFParticle::SetField(event->GetMagneticField());
1148 AliESDfriend *esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
1149 if (!esdFriend) {
1150 //Printf("ERROR: esdFriend not available");
1151 return;
1152 }
1153 if (esdFriend->TestSkipBit()) return;
1154 //
1155 //
1156 TDatabasePDG *pdg = TDatabasePDG::Instance();
1157 const Double_t kChi2Cut=5;
1158 const Double_t kMinR=2;
1159 const Int_t kMinNcl=110;
1160 const Double_t kMinREl=5;
1161 const Double_t kMaxREl=70;
1162 //
1163 Int_t nv0 = event->GetNumberOfV0s();
1164 AliESDVertex *vertex= (AliESDVertex *)event->GetPrimaryVertex();
1165 AliKFVertex kfvertex=*vertex;
1166 //
1167 for (Int_t iv0=0;iv0<nv0;iv0++){
1168 AliESDv0 *v0 = event->GetV0(iv0);
1169 if (!v0) continue;
1170 if (v0->GetOnFlyStatus()<0.5) continue;
1171 if (v0->GetPindex()<0) continue;
1172 if (v0->GetNindex()<0) continue;
1173 if (TMath::Max(v0->GetPindex(), v0->GetNindex())>event->GetNumberOfTracks()) continue;
1174 //
1175 //
1176 AliExternalTrackParam pp=(v0->GetParamP()->GetSign()>0) ? (*(v0->GetParamP())):(*(v0->GetParamN()));
1177 AliExternalTrackParam pn=(v0->GetParamP()->GetSign()>0) ? (*(v0->GetParamN())):(*(v0->GetParamP()));
1178 AliKFParticle kfp1( pp, 211 );
1179 AliKFParticle kfp2( pn, -211 );
1180 //
1181 AliKFParticle *v0KFK0 = new AliKFParticle(kfp1,kfp2);
1182 AliKFParticle *v0KFK0CV = new AliKFParticle(*v0KFK0);
1183 v0KFK0CV->SetProductionVertex(kfvertex);
1184 v0KFK0CV->TransportToProductionVertex();
1185 Double_t chi2K0 = v0KFK0CV->GetChi2();
1186 if (chi2K0>kChi2Cut) continue;
1187 if (v0->GetRr()<kMinR) continue;
1188 Bool_t isOKC=TMath::Max(v0->GetCausalityP()[0],v0->GetCausalityP()[1])<0.7&&TMath::Min(v0->GetCausalityP()[2],v0->GetCausalityP()[3])>0.2;
1189 //
1190 Double_t effMass22=v0->GetEffMass(2,2);
1191 Double_t effMass42=v0->GetEffMass(4,2);
1192 Double_t effMass24=v0->GetEffMass(2,4);
1193 Double_t effMass00=v0->GetEffMass(0,0);
1194 AliKFParticle *v0KFK0CVM = new AliKFParticle(*v0KFK0CV);
1195 v0KFK0CVM->SetMassConstraint(pdg->GetParticle("K_S0")->Mass());
1196 Bool_t isV0= kFALSE;
1197 //
1198 Double_t d22 = TMath::Abs(effMass22-pdg->GetParticle("K_S0")->Mass());
1199 Double_t d42 = TMath::Abs(effMass42-pdg->GetParticle("Lambda0")->Mass());
1200 Double_t d24 = TMath::Abs(effMass24-pdg->GetParticle("Lambda0")->Mass());
1201 Double_t d00 = TMath::Abs(effMass00);
1202 //
1203 Bool_t isKaon = d22<0.01 && d22< 0.3 * TMath::Min(TMath::Min(d42,d24),d00);
1204 Bool_t isLambda = d42<0.01 && d42< 0.3 * TMath::Min(TMath::Min(d22,d24),d00);
1205 Bool_t isAntiLambda= d24<0.01 && d24< 0.3 * TMath::Min(TMath::Min(d22,d42),d00);
1206 Bool_t isGamma = d00<0.02 && d00< 0.3 * TMath::Min(TMath::Min(d42,d24),d22);
1207 //
1208 if (isGamma && (isKaon||isLambda||isAntiLambda)) continue;
1209 if (isLambda && (isKaon||isGamma||isAntiLambda)) continue;
1210 if (isKaon && (isLambda||isGamma||isAntiLambda)) continue;
1211 Double_t sign= v0->GetParamP()->GetSign()* v0->GetParamN()->GetSign();
1212 if (sign>0) continue;
1213 isV0=isKaon||isLambda||isAntiLambda||isGamma;
1214 if (!(isV0)) continue;
1215 if (isGamma&&v0->GetRr()<kMinREl) continue;
1216 if (isGamma&&v0->GetRr()>kMaxREl) continue;
1217 if (!isOKC) continue;
1218 //
1219 Int_t pindex = (v0->GetParamP()->GetSign()>0) ? v0->GetPindex() : v0->GetNindex();
1220 Int_t nindex = (v0->GetParamP()->GetSign()>0) ? v0->GetNindex() : v0->GetPindex();
1221 AliESDtrack * trackP = event->GetTrack(pindex);
1222 AliESDtrack * trackN = event->GetTrack(nindex);
1223 if (!trackN) continue;
1224 if (!trackP) continue;
1225 Int_t nclP= (Int_t)trackP->GetTPCClusterInfo(2,1);
1226 Int_t nclN= (Int_t)trackN->GetTPCClusterInfo(2,1);
1227 if (TMath::Min(nclP,nclN)<kMinNcl) continue;
1228 Double_t eta = TMath::Max(TMath::Abs(trackP->Eta()), TMath::Abs(trackN->Eta()));
1229 if (TMath::Abs(eta)>1) continue;
1230 //
1231 //
1232 AliESDfriendTrack *friendTrackP = esdFriend->GetTrack(pindex);
1233 AliESDfriendTrack *friendTrackN = esdFriend->GetTrack(nindex);
1234 if (!friendTrackP) continue;
1235 if (!friendTrackN) continue;
1236 TObject *calibObject;
1237 AliTPCseed *seedP = 0;
1238 AliTPCseed *seedN = 0;
1239 for (Int_t l=0;(calibObject=friendTrackP->GetCalibObject(l));++l) {
1240 if ((seedP=dynamic_cast<AliTPCseed*>(calibObject))) break;
1241 }
1242 for (Int_t l=0;(calibObject=friendTrackN->GetCalibObject(l));++l) {
1243 if ((seedN=dynamic_cast<AliTPCseed*>(calibObject))) break;
1244 }
1245 if (isGamma){
1246 if ( TMath::Abs((trackP->GetTPCsignal()/(trackN->GetTPCsignal()+0.0001)-1)>0.3)) continue;
1247 }
1248 if (isGamma) (*fPIDMatrix)(pindex, 0)+=2;
1249 if (isGamma) (*fPIDMatrix)(nindex, 0)+=2;
1250 //
1251 if (isKaon) (*fPIDMatrix)(pindex, 2)+=2;
1252 if (isKaon) (*fPIDMatrix)(nindex, 2)+=2;
1253 //
1254 //
1255 if (pcstream){
1256 (*pcstream)<<"v0s"<<
1257 "isGamma="<<isGamma<<
1258 "isKaon="<<isKaon<<
1259 "isLambda="<<isLambda<<
1260 "isAntiLambda="<<isAntiLambda<<
1261 "chi2="<<chi2K0<<
1262 "trackP.="<<trackP<<
1263 "trackN.="<<trackN<<
1264 "v0.="<<v0<<
1265 "\n";
1266 }
1267 }
1268}
1269
1270
1271
1272
1273void AliTPCcalibGainMult::ProcessCosmic(const AliESDEvent * event) {
1274 //
1275 // Find cosmic pairs trigger by random trigger
1276 //
1277 //
1278 AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;
1279 AliTPCParam *param = AliTPCcalibDB::Instance()->GetParameters();
1280
1281 AliESDVertex *vertexSPD = (AliESDVertex *)event->GetPrimaryVertexSPD();
1282 AliESDVertex *vertexTPC = (AliESDVertex *)event->GetPrimaryVertexTPC();
1283 AliESDfriend *esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
1284 const Double_t kMinPt=4;
1285 const Double_t kMinPtMax=0.8;
1286 const Double_t kMinNcl=159*0.5;
1287 const Double_t kMaxDelta[5]={2,600,0.02,0.02,0.1};
1288 Int_t ntracks=event->GetNumberOfTracks();
1289 const Double_t kMaxImpact=80;
1290 // Float_t dcaTPC[2]={0,0};
1291 // Float_t covTPC[3]={0,0,0};
1292
1293 UInt_t specie = event->GetEventSpecie(); // skip laser events
1294 if (specie==AliRecoParam::kCalib) return;
1295
1296
1297 for (Int_t itrack0=0;itrack0<ntracks;itrack0++) {
1298 AliESDtrack *track0 = event->GetTrack(itrack0);
1299 if (!track0) continue;
1300 if (!track0->IsOn(AliESDtrack::kTPCrefit)) continue;
1301
1302 if (TMath::Abs(AliTracker::GetBz())>1&&track0->Pt()<kMinPt) continue;
1303 if (track0->GetTPCncls()<kMinNcl) continue;
1304 if (TMath::Abs(track0->GetY())<2*kMaxDelta[0]) continue;
1305 if (TMath::Abs(track0->GetY())>kMaxImpact) continue;
1306 if (track0->GetKinkIndex(0)>0) continue;
1307 const Double_t * par0=track0->GetParameter(); //track param at rhe DCA
1308 //rm primaries
1309 //
1310 for (Int_t itrack1=itrack0+1;itrack1<ntracks;itrack1++) {
1311 AliESDtrack *track1 = event->GetTrack(itrack1);
1312 if (!track1) continue;
1313 if (!track1->IsOn(AliESDtrack::kTPCrefit)) continue;
1314 if (track1->GetKinkIndex(0)>0) continue;
1315 if (TMath::Abs(AliTracker::GetBz())>1&&track1->Pt()<kMinPt) continue;
1316 if (track1->GetTPCncls()<kMinNcl) continue;
1317 if (TMath::Abs(AliTracker::GetBz())>1&&TMath::Max(track1->Pt(), track0->Pt())<kMinPtMax) continue;
1318 if (TMath::Abs(track1->GetY())<2*kMaxDelta[0]) continue;
1319 if (TMath::Abs(track1->GetY())>kMaxImpact) continue;
1320 //
1321 const Double_t* par1=track1->GetParameter(); //track param at rhe DCA
1322 //
1323 Bool_t isPair=kTRUE;
1324 for (Int_t ipar=0; ipar<5; ipar++){
1325 if (ipar==4&&TMath::Abs(AliTracker::GetBz())<1) continue; // 1/pt not defined for B field off
1326 if (TMath::Abs(TMath::Abs(par0[ipar])-TMath::Abs(par1[ipar]))>kMaxDelta[ipar]) isPair=kFALSE;
1327 }
1328 if (!isPair) continue;
1329 if (TMath::Abs(TMath::Abs(track0->GetAlpha()-track1->GetAlpha())-TMath::Pi())>kMaxDelta[2]) isPair=kFALSE;
1330 //delta with correct sign
1331 if (TMath::Abs(par0[0]+par1[0])>kMaxDelta[0]) isPair=kFALSE; //delta y opposite sign
1332 if (TMath::Abs(par0[3]+par1[3])>kMaxDelta[3]) isPair=kFALSE; //delta tgl opposite sign
1333 if (TMath::Abs(AliTracker::GetBz())>1 && TMath::Abs(par0[4]+par1[4])>kMaxDelta[4]) isPair=kFALSE; //delta 1/pt opposite sign
1334 if (!isPair) continue;
1335 TString filename(AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile()->GetName());
1336 Int_t eventNumber = event->GetEventNumberInFile();
1337 Bool_t hasFriend=(esdFriend) ? (esdFriend->GetTrack(itrack0)!=0):0;
1338 Bool_t hasITS=(track0->GetNcls(0)+track1->GetNcls(0)>4);
1339 printf("DUMPHPTCosmic:%s|%f|%d|%d|%d\n",filename.Data(),(TMath::Min(track0->Pt(),track1->Pt())), eventNumber,hasFriend,hasITS);
1340 //
1341 //
1342 TTreeSRedirector * pcstream = GetDebugStreamer();
1343 Int_t ntracksSPD = vertexSPD->GetNContributors();
1344 Int_t ntracksTPC = vertexTPC->GetNContributors();
1345 //
1346 if (pcstream){
1347 (*pcstream)<<"cosmicPairsAll"<<
1348 "run="<<fRun<< // run number
1349 "event="<<fEvent<< // event number
1350 "time="<<fTime<< // time stamp of event
1351 "trigger="<<fTrigger<< // trigger
1352 "triggerClass="<<&fTriggerClass<< // trigger
1353 "bz="<<fMagF<< // magnetic field
1354 //
1355 "nSPD="<<ntracksSPD<<
1356 "nTPC="<<ntracksTPC<<
1357 "vSPD.="<<vertexSPD<< //primary vertex -SPD
1358 "vTPC.="<<vertexTPC<< //primary vertex -TPC
1359 "t0.="<<track0<< //track0
1360 "t1.="<<track1<< //track1
1361 "\n";
1362 }
1363 //
1364 AliESDfriendTrack *friendTrack0 = esdFriend->GetTrack(itrack0);
1365 if (!friendTrack0) continue;
1366 AliESDfriendTrack *friendTrack1 = esdFriend->GetTrack(itrack1);
1367 if (!friendTrack1) continue;
1368 TObject *calibObject;
1369 AliTPCseed *seed0 = 0;
1370 AliTPCseed *seed1 = 0;
1371 //
1372 for (Int_t l=0;(calibObject=friendTrack0->GetCalibObject(l));++l) {
1373 if ((seed0=dynamic_cast<AliTPCseed*>(calibObject))) break;
1374 }
1375 for (Int_t l=0;(calibObject=friendTrack1->GetCalibObject(l));++l) {
1376 if ((seed1=dynamic_cast<AliTPCseed*>(calibObject))) break;
1377 }
1378 //
1379 if (pcstream){
1380 (*pcstream)<<"cosmicPairs"<<
1381 "run="<<fRun<< // run number
1382 "event="<<fEvent<< // event number
1383 "time="<<fTime<< // time stamp of event
1384 "trigger="<<fTrigger<< // trigger
1385 "triggerClass="<<&fTriggerClass<< // trigger
1386 "bz="<<fMagF<< // magnetic field
1387 //
1388 "nSPD="<<ntracksSPD<<
1389 "nTPC="<<ntracksTPC<<
1390 "vSPD.="<<vertexSPD<< //primary vertex -SPD
1391 "vTPC.="<<vertexTPC<< //primary vertex -TPC
1392 "t0.="<<track0<< //track0
1393 "t1.="<<track1<< //track1
1394 "ft0.="<<friendTrack0<< //track0
1395 "ft1.="<<friendTrack1<< //track1
1396 "s0.="<<seed0<< //track0
1397 "s1.="<<seed1<< //track1
1398 "\n";
1399 }
1400 if (!seed0) continue;
1401 if (!seed1) continue;
1402 Int_t nclA0=0, nclC0=0; // number of clusters
1403 Int_t nclA1=0, nclC1=0; // number of clusters
1404 //
1405 for (Int_t irow=0; irow<159; irow++){
1406 AliTPCclusterMI *cluster0=seed0->GetClusterPointer(irow);
1407 AliTPCclusterMI *cluster1=seed1->GetClusterPointer(irow);
1408 if (cluster0){
1409 if (cluster0->GetQ()>0 && cluster0->GetDetector()%36<18) nclA0++;
1410 if (cluster0->GetQ()>0 && cluster0->GetDetector()%36>=18) nclC0++;
1411 }
1412 if (cluster1){
1413 if (cluster1->GetQ()>0 && cluster1->GetDetector()%36<18) nclA1++;
1414 if (cluster1->GetQ()>0 && cluster1->GetDetector()%36>=18) nclC1++;
1415 }
1416 }
1417 Int_t cosmicType=0; // types of cosmic topology
1418 if ((nclA0>nclC0) && (nclA1>nclC1)) cosmicType=0; // AA side
1419 if ((nclA0<nclC0) && (nclA1<nclC1)) cosmicType=1; // CC side
1420 if ((nclA0>nclC0) && (nclA1<nclC1)) cosmicType=2; // AC side
1421 if ((nclA0<nclC0) && (nclA1>nclC1)) cosmicType=3; // CA side
1422 if (cosmicType<2) continue; // use only crossing tracks
1423 //
1424 Double_t deltaTimeCluster=0;
1425 deltaTimeCluster=0.5*(track1->GetZ()-track0->GetZ())/param->GetZWidth();
1426 if (nclA0>nclC0) deltaTimeCluster*=-1; // if A side track
1427 //
1428 for (Int_t irow=0; irow<159; irow++){
1429 AliTPCclusterMI *cluster0=seed0->GetClusterPointer(irow);
1430 if (cluster0 &&cluster0->GetX()>10){
1431 Double_t x0[3]={cluster0->GetRow(),cluster0->GetPad(),cluster0->GetTimeBin()+deltaTimeCluster};
1432 Int_t index0[1]={cluster0->GetDetector()};
1433 transform->Transform(x0,index0,0,1);
1434 cluster0->SetX(x0[0]);
1435 cluster0->SetY(x0[1]);
1436 cluster0->SetZ(x0[2]);
1437 //
1438 }
1439 AliTPCclusterMI *cluster1=seed1->GetClusterPointer(irow);
1440 if (cluster1&&cluster1->GetX()>10){
1441 Double_t x1[3]={cluster1->GetRow(),cluster1->GetPad(),cluster1->GetTimeBin()+deltaTimeCluster};
1442 Int_t index1[1]={cluster1->GetDetector()};
1443 transform->Transform(x1,index1,0,1);
1444 cluster1->SetX(x1[0]);
1445 cluster1->SetY(x1[1]);
1446 cluster1->SetZ(x1[2]);
1447 }
1448 }
1449 //
1450 //
1451 if (fPIDMatrix){
1452 (*fPIDMatrix)(itrack0,1)+=4; //
1453 (*fPIDMatrix)(itrack1,1)+=4; //
1454 }
1455 }
1456 }
1457}
1458
1459
1460
1461void AliTPCcalibGainMult::ProcessKinks(const AliESDEvent * event){
1462 //
1463 //
1464 //
1465 AliKFParticle::SetField(event->GetMagneticField());
1466 AliESDfriend *esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
1467 if (!esdFriend) {
1468 //Printf("ERROR: esdFriend not available");
1469 return;
1470 }
1471 // if (esdFriend->TestSkipBit()) return;
1472 //
1473 //
1474 const Double_t kChi2Cut=10;
1475 const Double_t kMinR=100;
1476 const Double_t kMaxR=230;
1477 const Int_t kMinNcl=110;
1478 //
1479 Int_t nkinks = event->GetNumberOfKinks();
1480 AliESDVertex *vertex= (AliESDVertex *)event->GetPrimaryVertex();
1481 AliKFVertex kfvertex=*vertex;
1482 TTreeSRedirector * pcstream = GetDebugStreamer();
1483 //
1484 for (Int_t ikink=0;ikink<nkinks;ikink++){
1485 AliESDkink *kink = event->GetKink(ikink);
1486 if (!kink) continue;
1487 if (kink->GetIndex(0)<0) continue;
1488 if (kink->GetIndex(1)<0) continue;
1489 if (TMath::Max(kink->GetIndex(1), kink->GetIndex(0))>event->GetNumberOfTracks()) continue;
1490 //
1491 //
1492 AliExternalTrackParam pd=kink->RefParamDaughter();
1493 AliExternalTrackParam pm=kink->RefParamMother();
1494 AliKFParticle kfpd( pd, 211 );
1495 AliKFParticle kfpm( pm, -13 );
1496 //
1497 AliKFParticle *v0KF = new AliKFParticle(kfpm,kfpd);
1498 v0KF->SetVtxGuess(kink->GetPosition()[0],kink->GetPosition()[1],kink->GetPosition()[2]);
1499 Double_t chi2 = v0KF->GetChi2();
1500 AliESDtrack * trackM = event->GetTrack(kink->GetIndex(0));
1501 AliESDtrack * trackD = event->GetTrack(kink->GetIndex(1));
1502 if (!trackM) continue;
1503 if (!trackD) continue;
1504 Int_t nclM= (Int_t)trackM->GetTPCClusterInfo(2,1);
1505 Int_t nclD= (Int_t)trackD->GetTPCClusterInfo(2,1);
1506 Double_t eta = TMath::Max(TMath::Abs(trackM->Eta()), TMath::Abs(trackD->Eta()));
1507 Double_t kx= v0KF->GetX();
1508 Double_t ky= v0KF->GetY();
1509 Double_t kz= v0KF->GetZ();
1510 Double_t ex= v0KF->GetErrX();
1511 Double_t ey= v0KF->GetErrY();
1512 Double_t ez= v0KF->GetErrZ();
1513 //
1514 Double_t radius=TMath::Sqrt(kx*kx+ky*ky);
1515 Double_t alpha=TMath::ATan2(ky,kx);
1516 if (!pd.Rotate(alpha)) continue;
1517 if (!pm.Rotate(alpha)) continue;
1518 if (!pd.PropagateTo(radius,event->GetMagneticField())) continue;
1519 if (!pm.PropagateTo(radius,event->GetMagneticField())) continue;
1520 Double_t pos[2]={0,kz};
1521 Double_t cov[3]={ex*ex+ey*ey,0,ez*ez};
1522 pd.Update(pos,cov);
1523 pm.Update(pos,cov);
1524 //
1525 if (pcstream){
1526 (*pcstream)<<"kinks"<<
1527 "eta="<<eta<<
1528 "nclM="<<nclM<<
1529 "nclD="<<nclD<<
1530 "kink.="<<kink<<
1531 "trackM.="<<trackM<<
1532 "trackD.="<<trackD<<
1533 "pm.="<<&pm<< //updated parameters
1534 "pd.="<<&pd<< // updated parameters
1535 "v0KF.="<<v0KF<<
1536 "chi2="<<chi2<<
1537 "\n";
1538 }
1539 /*
1540 TCut cutQ="chi2<10&&kink.fRr>90&&kink.fRr<220";
1541 TCut cutRD="20*sqrt(pd.fC[14])<abs(pd.fP[4])&&trackD.fTPCsignal>10&&trackD.fTPCsignalN>50";
1542
1543 */
1544 //
1545 if (chi2>kChi2Cut) continue;
1546 if (kink->GetR()<kMinR) continue;
1547 if (kink->GetR()>kMaxR) continue;
1548 if ((nclM+nclD)<kMinNcl) continue;
1549 if (TMath::Abs(eta)>1) continue;
1550 //
1551 //
1552 AliESDfriendTrack *friendTrackM = esdFriend->GetTrack(kink->GetIndex(0));
1553 AliESDfriendTrack *friendTrackD = esdFriend->GetTrack(kink->GetIndex(1));
1554 if (!friendTrackM) continue;
1555 if (!friendTrackD) continue;
1556 TObject *calibObject;
1557 AliTPCseed *seedM = 0;
1558 AliTPCseed *seedD = 0;
1559 for (Int_t l=0;(calibObject=friendTrackM->GetCalibObject(l));++l) {
1560 if ((seedM=dynamic_cast<AliTPCseed*>(calibObject))) break;
1561 }
1562 for (Int_t l=0;(calibObject=friendTrackD->GetCalibObject(l));++l) {
1563 if ((seedD=dynamic_cast<AliTPCseed*>(calibObject))) break;
1564 }
1565 }
1566}
1567
1568void AliTPCcalibGainMult::DumpHPT(const AliESDEvent * event){
1569 //
1570 // Function to select the HPT tracks and events
1571 // It is used to select event with HPT - list used later for the raw data downloading
1572 // - and reconstruction
1573 // Not actualy used for the calibration of the data
1574
1575 TTreeSRedirector * pcstream = GetDebugStreamer();
1576 AliKFParticle::SetField(event->GetMagneticField());
1577 AliESDfriend *esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
1578 if (!esdFriend) {
1579 //Printf("ERROR: esdFriend not available");
1580 return;
1581 }
1582 if (esdFriend->TestSkipBit()) return;
1583
1584 Int_t ntracks=event->GetNumberOfTracks();
1585 //
1586 for (Int_t i=0;i<ntracks;++i) {
1587 //
1588 AliESDtrack *track = event->GetTrack(i);
1589 if (!track) continue;
1590 if (track->Pt()<4) continue;
1591 UInt_t status = track->GetStatus();
1592 //
1593 AliExternalTrackParam * trackIn = (AliExternalTrackParam *)track->GetInnerParam();
1594 if (!trackIn) continue;
1595 if ((status&AliESDtrack::kTPCrefit)==0) continue;
1596 if ((status&AliESDtrack::kITSrefit)==0) continue;
1597 AliESDfriendTrack *friendTrack = esdFriend->GetTrack(i);
1598 if (!friendTrack) continue;
1599 AliExternalTrackParam * itsOut = (AliExternalTrackParam *)(friendTrack->GetITSOut());
1600 if (!itsOut) continue;
1601 AliExternalTrackParam * itsOut2 = (AliExternalTrackParam *)(friendTrack->GetITSOut()->Clone());
1602 AliExternalTrackParam * tpcIn2 = (AliExternalTrackParam *)(trackIn->Clone());
1603 if (!itsOut2->Rotate(trackIn->GetAlpha())) continue;
1604 //Double_t xmiddle=0.5*(itsOut2->GetX()+tpcIn2->GetX());
1605 Double_t xmiddle=(itsOut2->GetX());
1606 if (!itsOut2->PropagateTo(xmiddle,event->GetMagneticField())) continue;
1607 if (!tpcIn2->PropagateTo(xmiddle,event->GetMagneticField())) continue;
1608 //
1609 AliExternalTrackParam * tpcInner = (AliExternalTrackParam *)(track->GetTPCInnerParam());
1610 if (!tpcInner) continue;
1611 tpcInner->Rotate(track->GetAlpha());
1612 tpcInner->PropagateTo(track->GetX(),event->GetMagneticField());
1613 //
1614 // tpc constrained
1615 //
1616 AliExternalTrackParam * tpcInnerC = (AliExternalTrackParam *)(track->GetTPCInnerParam()->Clone());
1617 if (!tpcInnerC) continue;
1618 tpcInnerC->Rotate(track->GetAlpha());
1619 tpcInnerC->PropagateTo(track->GetX(),event->GetMagneticField());
1620 Double_t dz[2],cov[3];
1621 AliESDVertex *vtx= (AliESDVertex *)event->GetPrimaryVertex();
1622
1623 if (!tpcInnerC->PropagateToDCA(vtx, event->GetMagneticField(), 3, dz, cov)) continue;
1624 Double_t covar[6]; vtx->GetCovMatrix(covar);
1625 Double_t p[2]={tpcInnerC->GetParameter()[0]-dz[0],tpcInnerC->GetParameter()[1]-dz[1]};
1626 Double_t c[3]={covar[2],0.,covar[5]};
1627 //
1628 Double_t chi2C=tpcInnerC->GetPredictedChi2(p,c);
1629 tpcInnerC->Update(p,c);
1630
1631 if (pcstream){
1632 (*pcstream)<<"hpt"<<
1633 "run="<<fRun<<
1634 "time="<<fTime<<
1635 "vertex="<<vtx<<
1636 "bz="<<fMagF<<
1637 "track.="<<track<<
1638 "tpcInner.="<<tpcInner<<
1639 "tpcInnerC.="<<tpcInnerC<<
1640 "chi2C="<<chi2C<<
1641 //
1642 "its.="<<itsOut<<
1643 "its2.="<<itsOut2<<
1644 "tpc.="<<trackIn<<
1645 "tpc2.="<<tpcIn2<<
1646 "\n";
1647 }
1648 }
1649}
1650
1651
1652
1653void AliTPCcalibGainMult::ProcessTOF(const AliESDEvent * event){
1654 //
1655 // 1. Loop over tracks
1656 // 2. Fit T0
1657 // 3. Sign positivelly identified tracks
1658 //
1659 const Double_t kMaxDelta=1000;
1660 const Double_t kOrbit=50000; // distance in the time beween two orbits in the TOF time units - 50000=50 ns
1661 const Double_t kMaxD=20;
1662 const Double_t kRMS0=200;
1663 const Double_t kMaxDCAZ=10;
1664 AliESDVertex *vtx= (AliESDVertex *)event->GetPrimaryVertex();
1665 //
1666 TTreeSRedirector * pcstream = GetDebugStreamer();
1667 AliKFParticle::SetField(event->GetMagneticField());
1668 AliESDfriend *esdFriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
1669 if (!esdFriend) {
1670 //Printf("ERROR: esdFriend not available");
1671 return;
1672 }
1673 //if (esdFriend->TestSkipBit()) return;
1674
1675 Int_t ntracks=event->GetNumberOfTracks();
1676 //
1677 Double_t deltaTPion[10000];
1678 Double_t medianT0=0;
1679 Double_t meanT0=0;
1680 Double_t rms=10000;
1681 Int_t counter=0;
1682 //
1683 // Get Median time for pion hypothesy
1684 //
1685 for (Int_t iter=0; iter<3; iter++){
1686 counter=0;
1687 for (Int_t i=0;i<ntracks;++i) {
1688 //
1689 AliESDtrack *track = event->GetTrack(i);
1690 if (!track) continue;
1691 if (!track->IsOn(AliESDtrack::kTIME)) continue;
1692 if (TMath::Abs(track->GetZ())>kMaxDCAZ) continue; // remove overlaped events
1693 if (TMath::Abs(track->GetTOFsignalDz())>kMaxD) continue;
1694 Double_t times[1000];
1695 track->GetIntegratedTimes(times);
1696 Int_t norbit=TMath::Nint((track->GetTOFsignal()-times[2])/kOrbit);
1697 Double_t torbit=norbit*kOrbit;
1698 if (iter==1 &&TMath::Abs(times[2]-times[3])<3*rms) continue; // skip umbigous points - kaon pion
1699 //
1700 Int_t indexBest=2;
1701 if (iter>1){
1702 for (Int_t j=3; j<5; j++)
1703 if (TMath::Abs(track->GetTOFsignal()-times[j]-torbit-medianT0)<TMath::Abs(track->GetTOFsignal()-times[j]-torbit-medianT0)) indexBest=j;
1704 }
1705 //
1706 if (iter>0) if (TMath::Abs(track->GetTOFsignal()-times[indexBest]-torbit-medianT0)>3*(kRMS0+rms)) continue;
1707 if (iter>0) if (TMath::Abs(track->GetTOFsignal()-times[indexBest]-torbit-medianT0)>kMaxDelta) continue;
1708 deltaTPion[counter]=track->GetTOFsignal()-times[indexBest]-torbit;
1709 counter++;
1710 }
1711 if (counter<2) return;
1712 medianT0=TMath::Median(counter,deltaTPion);
1713 meanT0=TMath::Median(counter,deltaTPion);
1714 rms=TMath::RMS(counter,deltaTPion);
1715 }
1716 if (counter<3) return;
1717 //
1718 // Dump
1719 //
1720 for (Int_t i=0;i<ntracks;++i) {
1721 //
1722 AliESDtrack *track = event->GetTrack(i);
1723 if (!track) continue;
1724 if (!track->IsOn(AliESDtrack::kTIME)) continue;
1725 if (TMath::Abs(track->GetZ())>kMaxDCAZ) continue; //remove overlapped events
1726 if (TMath::Abs(track->GetTOFsignalDz())>kMaxD) continue;
1727 Double_t times[1000];
1728 track->GetIntegratedTimes(times);
1729 Int_t norbit=TMath::Nint((track->GetTOFsignal()-times[2])/kOrbit);
1730 Double_t torbit=norbit*kOrbit;
1731 if (rms<=0) continue;
1732 //
1733 Double_t tPion = (track->GetTOFsignal()-times[2]-medianT0-torbit);
1734 Double_t tKaon = (track->GetTOFsignal()-times[3]-medianT0-torbit);
1735 Double_t tProton= (track->GetTOFsignal()-times[4]-medianT0-torbit);
1736 Double_t tElectron= (track->GetTOFsignal()-times[0]-medianT0-torbit);
1737 //
1738 Bool_t isPion = (TMath::Abs(tPion/rms)<6) && TMath::Abs(tPion)<(TMath::Min(TMath::Abs(tKaon), TMath::Abs(tProton))-rms);
1739 Bool_t isKaon = (TMath::Abs(tKaon/rms)<3) && TMath::Abs(tKaon)<0.2*(TMath::Min(TMath::Abs(tPion), TMath::Abs(tProton))-3*rms);
1740 Bool_t isProton = (TMath::Abs(tProton/rms)<6) && TMath::Abs(tProton)<0.5*(TMath::Min(TMath::Abs(tKaon), TMath::Abs(tPion))-rms);
1741 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);
1742
1743 if (isPion) (*fPIDMatrix)(i,2)+=1;
1744 if (isKaon) (*fPIDMatrix)(i,3)+=1;
1745 if (isProton) (*fPIDMatrix)(i,4)+=1;
1746 // if (isElectron) (*fPIDMatrix)(i,0)+=1;
1747 //
1748 if (pcstream){
1749 // debug streamer to dump the information
1750 (*pcstream)<<"tof"<<
1751 "isPion="<<isPion<<
1752 "isKaon="<<isKaon<<
1753 "isProton="<<isProton<<
1754 "isElectron="<<isElectron<<
1755 //
1756 "counter="<<counter<<
1757 "torbit="<<torbit<<
1758 "norbit="<<norbit<<
1759 "medianT0="<<medianT0<<
1760 "meanT0="<<meanT0<<
1761 "rmsT0="<<rms<<
1762 "track.="<<track<<
1763 "vtx.="<<vtx<<
1764 "\n";
1765 }
1766
1767 }
1768 /*
1769 tof->SetAlias("isProton","(abs(track.fTOFsignal-track.fTrackTime[4]-medianT0-torbit)<(0.5*abs(track.fTOFsignal-track.fTrackTime[3]-medianT0-torbit)-rmsT0))");
1770 tof->SetAlias("isPion","(abs(track.fTOFsignal-track.fTrackTime[2]-medianT0-torbit)<(0.5*abs(track.fTOFsignal-track.fTrackTime[3]-medianT0-torbit)-rmsT0))");
1771 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))");
1772
1773 */
1774
1775}
1776
1777
1778// void AliTPCcalibGainMult::Terminate(){
1779// //
1780// // Terminate function
1781// // call base terminate + Eval of fitters
1782// //
1783// Info("AliTPCcalibGainMult","Terminate");
1784// TTreeSRedirector *pcstream = GetDebugStreamer();
1785// if (pcstream){
1786// TTreeStream &stream = (*pcstream)<<"dump";
1787// TTree* tree = stream.GetTree();
1788// if (tree) if ( tree->GetEntries()>0){
1789// TObjArray *array = tree->GetListOfBranches();
1790// for (Int_t i=0; i<array->GetEntries(); i++) {TBranch * br = (TBranch *)array->At(i); br->SetAddress(0);}
1791// gDirectory=gROOT;
1792// fdEdxTree=tree->CloneTree(10000);
1793// }
1794// }
1795// AliTPCcalibBase::Terminate();
1796// }
1797