]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TPC/AliTPCcalibLaser.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
M AliTPCcalibCalib.cxx - set after refitting also TPC out in friends
[u/mrichter/AliRoot.git] / TPC / AliTPCcalibLaser.cxx
CommitLineData
c6914c83 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/*
2b35e8e6 17 //
18 // FUNCTIONALITY:
19 //
20 // 1. The laser track is associated with the mirror
21 // see function FindMirror
22 //
23 // 2. The laser track is accepted for the analysis under certain condition
24 // (see function Accpet laser)
263d466a 25 //
2b35e8e6 26 // 3. The drift velocity and jitter is calculated event by event
263d466a 27 // (see function drift velocity)
2b35e8e6 28 //
be113f6e 29 // 4. The tracks are refitted at different sectors
30 // Fit model
31 // 4.a) line
32 // 4.b) parabola
33 // 4.c) parabola with common P2 for inner and outer
2b35e8e6 34 //
1fd56785 35 // To make laser scan the user interaction neccessary
36 //
880c3382 37 .x ~/NimStyle.C
2b35e8e6 38 gSystem->Load("libANALYSIS");
39 gSystem->Load("libTPCcalib");
82e5df1f 40 TFile fcalib("CalibObjectsTrain2.root");
41 AliTPCcalibLaser * laser = ( AliTPCcalibLaser *)fcalib->Get("laserTPC");
57dc06f2 42 laser->DumpMeanInfo(-0,0)
2b35e8e6 43 TFile fmean("laserMean.root")
1fd56785 44 //
2b35e8e6 45 // laser track clasification;
1fd56785 46 //
7b18d067 47 TCut cutT("cutT","abs(Tr.fP[3])<0.06");
c18f4385 48 TCut cutN("cutN","fTPCncls>70");
49 TCut cutP("cutP","abs(atan2(x1,x0)-atan2(lx1,lx0))<0.03")
50 TCut cutA = cutT+cutPt+cutP;
e5c59218 51 TChain * chainTrL = tool.MakeChain("laser.txt","Track",0,10200);
52
1fd56785 53 //
54 //
55 // Analyze LASER scan
56 //
e5c59218 57
2b35e8e6 58 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
59 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
0728a4f6 60 AliXRDPROOFtoolkit tool;
2b35e8e6 61 TChain * chain = tool.MakeChain("laserScan.txt","Mean",0,10200);
62 chain->Lookup();
1fd56785 63 AliTPCcalibLaser::DumpScanInfo(chain)
64 TFile fscan("laserScan.root")
65 TTree * treeT = (TTree*)fscan.Get("Mean")
95a0e09b 66 //
67 // Analyze laser
68 //
69 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
70 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
71 AliXRDPROOFtoolkit tool;
41a1f702 72 AliXRDPROOFtoolkit::FilterList("laser.txt","* driftvN",1)
73 TChain * chainDrift = tool.MakeChainRandom("laser.txt.Good","driftv",0,50);
57dc06f2 74 chainDrift->Lookup();
41a1f702 75 TChain * chainDriftN = tool.MakeChainRandom("laser.txt.Good","driftvN",0,300);
f6220766 76 chainDriftN->Lookup();
77
78
95a0e09b 79 TChain * chain = tool.MakeChain("laser.txt","Residuals",0,10200);
80 chain->Lookup();
e5c59218 81 TChain * chainFit = tool.MakeChain("laser.txt","FitModels",0,10200);
82 chainFit->Lookup();
82e5df1f 83 TChain * chainTrack = tool.MakeChainRandom("laser.txt","Track",0,30);
84 chainTrack->Lookup();
95a0e09b 85
c6914c83 86*/
87
88
89
90#include "TLinearFitter.h"
91#include "AliTPCcalibLaser.h"
92#include "AliExternalTrackParam.h"
e9f38a4b 93#include "AliESDEvent.h"
94#include "AliESDfriend.h"
c6914c83 95#include "AliESDtrack.h"
96#include "AliTPCTracklet.h"
97#include "TH1D.h"
ecb7e705 98#include "TH1F.h"
95a0e09b 99#include "TProfile.h"
c6914c83 100#include "TVectorD.h"
101#include "TTreeStream.h"
102#include "TFile.h"
103#include "TF1.h"
104#include "TGraphErrors.h"
105#include "AliTPCclusterMI.h"
106#include "AliTPCseed.h"
107#include "AliTracker.h"
263d466a 108#include "AliLog.h"
c6914c83 109#include "TClonesArray.h"
1fd56785 110#include "TPad.h"
e5c59218 111#include "TSystem.h"
112#include "TCut.h"
113#include "TChain.h"
114#include "TH2F.h"
115#include "TStatToolkit.h"
116#include "TROOT.h"
c6914c83 117
118
119#include "TTreeStream.h"
120#include <iostream>
121#include <sstream>
7b18d067 122#include "AliTPCLaserTrack.h"
108953e9 123#include "AliTPCcalibDB.h"
124#include "AliTPCParam.h"
f6220766 125#include "AliTPCParamSR.h"
3ac724f0 126#include "TTimeStamp.h"
31aa7c5c 127#include "AliDCSSensorArray.h"
128#include "AliDCSSensor.h"
7b18d067 129
c6914c83 130using namespace std;
131
132ClassImp(AliTPCcalibLaser)
133
134AliTPCcalibLaser::AliTPCcalibLaser():
e9f38a4b 135 AliTPCcalibBase(),
136 fESD(0),
592a0c8f 137 fESDfriend(0),
138 fNtracks(0),
c18f4385 139 fTracksMirror(336),
140 fTracksEsd(336),
141 fTracksEsdParam(336),
142 fTracksTPC(336),
3ac724f0 143 fFullCalib(kTRUE),
263d466a 144 fDeltaZ(336),
e5c59218 145 fDeltaP3(336),
146 fDeltaP4(336),
263d466a 147 fDeltaPhi(336),
148 fDeltaPhiP(336),
95a0e09b 149 fSignals(336),
be113f6e 150 //
d3ce44cb 151 fHisLaser(0), // N dim histogram of laser
be113f6e 152 fHisNclIn(0), //->Number of clusters inner
153 fHisNclOut(0), //->Number of clusters outer
154 fHisNclIO(0), //->Number of cluster inner outer
155 fHisLclIn(0), //->Level arm inner
156 fHisLclOut(0), //->Level arm outer
157 fHisLclIO(0), //->Number of cluster inner outer
158 fHisdEdx(0), //->dEdx histo
159 fHisdZfit(0), //->distance to the mirror after linear fit
160 //
161 //
162 fHisChi2YIn1(0), //->chi2 y inner - line
163 fHisChi2YOut1(0), //->chi2 y inner - line
164 fHisChi2YIn2(0), //->chi2 y inner - parabola
165 fHisChi2YOut2(0), //->chi2 y inner - parabola
166 fHisChi2YIO1(0), //->chi2 y IO - common
167 fHisChi2ZIn1(0), //->chi2 z inner - line
168 fHisChi2ZOut1(0), //->chi2 z inner - line
169 fHisChi2ZIn2(0), //->chi2 z inner - parabola
170 fHisChi2ZOut2(0), //->chi2 z inner - parabola
171 fHisChi2ZIO1(0), //->chi2 z IO - common
172 //
173 //
174 fHisPy1vP0(0), //-> delta y P0outer-P0inner - line
175 fHisPy2vP0(0), //-> delta y P0outer-P0inner - parabola
176 fHisPy3vP0(0), //-> delta y P0outer-P0inner - common parabola
177 fHisPy1vP1(0), //-> delta ky P1outer-P1inner - line
178 fHisPy2vP1(0), //-> delta ky P1outer-P1inner - parabola
179 fHisPy3vP1(0), //-> delta ky P1outer-P1inner - common parabola
180 fHisPy2vP2In(0), //-> Curv P2inner - parabola
181 fHisPy2vP2Out(0), //-> Curv P2outer - parabola
182 fHisPy3vP2IO(0), //-> Curv P2outerinner - common parabola
183 //
184 //
185 fHisPz1vP0(0), //-> delta z P0outer-P0inner - line
186 fHisPz2vP0(0), //-> delta z P0outer-P0inner - parabola
187 fHisPz3vP0(0), //-> delta z P0outer-P0inner - common parabola
188 fHisPz1vP1(0), //-> delta kz P1outer-P1inner - line
189 fHisPz2vP1(0), //-> delta kz P1outer-P1inner - parabola
190 fHisPz3vP1(0), //-> delta kz P1outer-P1inner - common parabola
191 fHisPz2vP2In(0), //-> Curv P2inner - parabola
192 fHisPz2vP2Out(0), //-> Curv P2outer - parabola
193 fHisPz3vP2IO(0), //-> Curv P2outerinner - common parabola
194 //
108953e9 195 fDeltaYres(336), //->2D histo of residuals
196 fDeltaZres(336), //->2D histo fo residuals
197 fDeltaYres2(336), //->2D histo of residuals
198 fDeltaZres2(336), //->2D histo fo residuals
41a1f702 199 fDeltaYresAbs(336), //->2D histo of residuals
200 fHisYAbsErrors(0), //-> total errors per beam in the abs res y analysis
201 fDeltaZresAbs(336), //->2D histo of residuals
202 fHisZAbsErrors(0), //-> total errors per beam in the abs res z analysis
108953e9 203 //fDeltaYres3(336), //->2D histo of residuals
204 //fDeltaZres3(336), //->2D histo fo residuals
592a0c8f 205 fFitAside(new TVectorD(5)),
206 fFitCside(new TVectorD(5)),
207 fFitACside(new TVectorD(6)),
be113f6e 208 fEdgeXcuts(3),
209 fEdgeYcuts(3),
263d466a 210 fNClCuts(5),
7af539c6 211 fNcuts(0),
212 fBeamSectorOuter(336),
213 fBeamSectorInner(336),
214 fBeamOffsetYOuter(336),
215 fBeamSlopeYOuter(336),
216 fBeamOffsetYInner(336),
217 fBeamSlopeYInner(336),
218 fBeamOffsetZOuter(336),
219 fBeamSlopeZOuter(336),
220 fBeamOffsetZInner(336),
221 fBeamSlopeZInner(336),
41a1f702 222 fInverseSlopeZ(kTRUE),
223 fUseFixedDriftV(0),
224 fFixedFitAside0(0.0),
225 fFixedFitAside1(1.0),
226 fFixedFitCside0(0.0),
227 fFixedFitCside1(1.0)
c6914c83 228{
229 //
230 // Constructor
231 //
e9f38a4b 232 fTracksEsdParam.SetOwner(kTRUE);
c6914c83 233}
234
3ac724f0 235AliTPCcalibLaser::AliTPCcalibLaser(const Text_t *name, const Text_t *title, Bool_t full):
e9f38a4b 236 AliTPCcalibBase(),
237 fESD(0),
238 fESDfriend(0),
592a0c8f 239 fNtracks(0),
c18f4385 240 fTracksMirror(336),
241 fTracksEsd(336),
242 fTracksEsdParam(336),
243 fTracksTPC(336),
3ac724f0 244 fFullCalib(full),
be113f6e 245 //
c18f4385 246 fDeltaZ(336), // array of histograms of delta z for each track
e5c59218 247 fDeltaP3(336), // array of histograms of delta z for each track
248 fDeltaP4(336), // array of histograms of P3 for each track
249 fDeltaPhi(336), // array of histograms of P4 for each track
c18f4385 250 fDeltaPhiP(336), // array of histograms of delta z for each track
2b35e8e6 251 fSignals(336), // array of dedx signals
be113f6e 252 //
253 //
d3ce44cb 254 fHisLaser(0), // N dim histogram of laser
be113f6e 255 fHisNclIn(0), //->Number of clusters inner
256 fHisNclOut(0), //->Number of clusters outer
257 fHisNclIO(0), //->Number of cluster inner outer
258 fHisLclIn(0), //->Level arm inner
259 fHisLclOut(0), //->Level arm outer
260 fHisLclIO(0), //->Number of cluster inner outer
261 fHisdEdx(0), //->dEdx histo
262 fHisdZfit(0), //->distance to the mirror after linear fit
263 //
264 //
265 fHisChi2YIn1(0), //->chi2 y inner - line
266 fHisChi2YOut1(0), //->chi2 y inner - line
267 fHisChi2YIn2(0), //->chi2 y inner - parabola
268 fHisChi2YOut2(0), //->chi2 y inner - parabola
269 fHisChi2YIO1(0), //->chi2 y IO - common
270 fHisChi2ZIn1(0), //->chi2 z inner - line
271 fHisChi2ZOut1(0), //->chi2 z inner - line
272 fHisChi2ZIn2(0), //->chi2 z inner - parabola
273 fHisChi2ZOut2(0), //->chi2 z inner - parabola
274 fHisChi2ZIO1(0), //->chi2 z IO - common
275 //
276 //
277 fHisPy1vP0(0), //-> delta y P0outer-P0inner - line
278 fHisPy2vP0(0), //-> delta y P0outer-P0inner - parabola
279 fHisPy3vP0(0), //-> delta y P0outer-P0inner - common parabola
280 fHisPy1vP1(0), //-> delta ky P1outer-P1inner - line
281 fHisPy2vP1(0), //-> delta ky P1outer-P1inner - parabola
282 fHisPy3vP1(0), //-> delta ky P1outer-P1inner - common parabola
283 fHisPy2vP2In(0), //-> Curv P2inner - parabola
284 fHisPy2vP2Out(0), //-> Curv P2outer - parabola
285 fHisPy3vP2IO(0), //-> Curv P2outerinner - common parabola
286 //
287 //
288 fHisPz1vP0(0), //-> delta z P0outer-P0inner - line
289 fHisPz2vP0(0), //-> delta z P0outer-P0inner - parabola
290 fHisPz3vP0(0), //-> delta z P0outer-P0inner - common parabola
291 fHisPz1vP1(0), //-> delta kz P1outer-P1inner - line
292 fHisPz2vP1(0), //-> delta kz P1outer-P1inner - parabola
293 fHisPz3vP1(0), //-> delta kz P1outer-P1inner - common parabola
294 fHisPz2vP2In(0), //-> Curv P2inner - parabola
295 fHisPz2vP2Out(0), //-> Curv P2outer - parabola
296 fHisPz3vP2IO(0), //-> Curv P2outerinner - common parabola
297 //
298 //
299 //
95a0e09b 300 fDeltaYres(336),
301 fDeltaZres(336),
108953e9 302 fDeltaYres2(336),
303 fDeltaZres2(336),
7af539c6 304 fDeltaYresAbs(336),
41a1f702 305 fHisYAbsErrors(0),
7af539c6 306 fDeltaZresAbs(336),
41a1f702 307 fHisZAbsErrors(0),
108953e9 308 // fDeltaYres3(336),
309 //fDeltaZres3(336),
592a0c8f 310 fFitAside(new TVectorD(5)), // drift fit - A side
311 fFitCside(new TVectorD(5)), // drift fit - C- side
312 fFitACside(new TVectorD(6)), // drift fit - AC- side
be113f6e 313 fEdgeXcuts(3), // cuts in local x direction; used in the refit of the laser tracks
314 fEdgeYcuts(3), // cuts in local y direction; used in the refit of the laser tracks
263d466a 315 fNClCuts(5), // cuts on the number of clusters per tracklet; used in the refit of the laser tracks
7af539c6 316 fNcuts(0), // number of cuts
317 fBeamSectorOuter(336),
318 fBeamSectorInner(336),
319 fBeamOffsetYOuter(336),
320 fBeamSlopeYOuter(336),
321 fBeamOffsetYInner(336),
322 fBeamSlopeYInner(336),
323 fBeamOffsetZOuter(336),
324 fBeamSlopeZOuter(336),
325 fBeamOffsetZInner(336),
326 fBeamSlopeZInner(336),
41a1f702 327 fInverseSlopeZ(kTRUE),
328 fUseFixedDriftV(0),
329 fFixedFitAside0(0.0),
330 fFixedFitAside1(1.0),
331 fFixedFitCside0(0.0),
332 fFixedFitCside1(1.0)
c6914c83 333{
334 SetName(name);
335 SetTitle(title);
336 //
337 // Constructor
338 //
263d466a 339 fTracksEsdParam.SetOwner(kTRUE);
c6914c83 340}
341
be113f6e 342AliTPCcalibLaser::AliTPCcalibLaser(const AliTPCcalibLaser& calibLaser):
343 AliTPCcalibBase(calibLaser),
344 fESD(0),
345 fESDfriend(0),
592a0c8f 346 fNtracks(0),
be113f6e 347 fTracksMirror(336),
348 fTracksEsd(336),
349 fTracksEsdParam(336),
350 fTracksTPC(336),
3ac724f0 351 fFullCalib(calibLaser.fFullCalib),
be113f6e 352 //
353 fDeltaZ(calibLaser.fDeltaZ), // array of histograms of delta z for each track
354 fDeltaP3(((calibLaser.fDeltaP3))), // array of histograms of delta z for each track
355 fDeltaP4(((calibLaser.fDeltaP4))), // array of histograms of P3 for each track
356 fDeltaPhi(((calibLaser.fDeltaPhi))), // array of histograms of P4 for each track
357 fDeltaPhiP(((calibLaser.fDeltaPhiP))), // array of histograms of delta z for each track
358 fSignals(((calibLaser.fSignals))), // array of dedx signals
359 //
360 //
d3ce44cb 361 fHisLaser(0), // N dim histogram of laser
be113f6e 362 fHisNclIn(new TH2F(*(calibLaser.fHisNclIn))), //->Number of clusters inner
363 fHisNclOut(new TH2F(*(calibLaser.fHisNclOut))), //->Number of clusters outer
364 fHisNclIO(new TH2F(*(calibLaser.fHisNclIO))), //->Number of cluster inner outer
365 fHisLclIn(new TH2F(*(calibLaser.fHisLclIn))), //->Level arm inner
366 fHisLclOut(new TH2F(*(calibLaser.fHisLclOut))), //->Level arm outer
367 fHisLclIO(new TH2F(*(calibLaser.fHisLclIO))), //->Number of cluster inner outer
368 fHisdEdx(new TH2F(*(calibLaser.fHisdEdx))), //
369 fHisdZfit(new TH2F(*(calibLaser.fHisdZfit))), //->distance to the mirror after linear fit
370 //
371 //
372 fHisChi2YIn1(new TH2F(*(calibLaser.fHisChi2YIn1))), //->chi2 y inner - line
373 fHisChi2YOut1(new TH2F(*(calibLaser.fHisChi2YOut1))), //->chi2 y inner - line
374 fHisChi2YIn2(new TH2F(*(calibLaser.fHisChi2YIn2))), //->chi2 y inner - parabola
375 fHisChi2YOut2(new TH2F(*(calibLaser.fHisChi2YOut2))), //->chi2 y inner - parabola
376 fHisChi2YIO1(new TH2F(*(calibLaser.fHisChi2YIO1))), //->chi2 y IO - common
377 fHisChi2ZIn1(new TH2F(*(calibLaser.fHisChi2ZIn1))), //->chi2 z inner - line
378 fHisChi2ZOut1(new TH2F(*(calibLaser.fHisChi2ZOut1))), //->chi2 z inner - line
379 fHisChi2ZIn2(new TH2F(*(calibLaser.fHisChi2ZIn2))), //->chi2 z inner - parabola
380 fHisChi2ZOut2(new TH2F(*(calibLaser.fHisChi2ZOut2))), //->chi2 z inner - parabola
381 fHisChi2ZIO1(new TH2F(*(calibLaser.fHisChi2ZIO1))), //->chi2 z IO - common
382 //
383 //
384 fHisPy1vP0(new TH2F(*(calibLaser.fHisPy1vP0))), //-> delta y P0outer-P0inner - line
385 fHisPy2vP0(new TH2F(*(calibLaser.fHisPy2vP0))), //-> delta y P0outer-P0inner - parabola
386 fHisPy3vP0(new TH2F(*(calibLaser.fHisPy3vP0))), //-> delta y P0outer-P0inner - common parabola
387 fHisPy1vP1(new TH2F(*(calibLaser.fHisPy1vP1))), //-> delta ky P1outer-P1inner - line
388 fHisPy2vP1(new TH2F(*(calibLaser.fHisPy2vP1))), //-> delta ky P1outer-P1inner - parabola
389 fHisPy3vP1(new TH2F(*(calibLaser.fHisPy3vP1))), //-> delta ky P1outer-P1inner - common parabola
390 fHisPy2vP2In(new TH2F(*(calibLaser.fHisPy2vP2In))), //-> Curv P2inner - parabola
391 fHisPy2vP2Out(new TH2F(*(calibLaser.fHisPy2vP2Out))), //-> Curv P2outer - parabola
392 fHisPy3vP2IO(new TH2F(*(calibLaser.fHisPy3vP2IO))), //-> Curv P2outerinner - common parabola
393 //
394 //
395 fHisPz1vP0(new TH2F(*(calibLaser.fHisPz1vP0))), //-> delta z P0outer-P0inner - line
396 fHisPz2vP0(new TH2F(*(calibLaser.fHisPz2vP0))), //-> delta z P0outer-P0inner - parabola
397 fHisPz3vP0(new TH2F(*(calibLaser.fHisPz3vP0))), //-> delta z P0outer-P0inner - common parabola
398 fHisPz1vP1(new TH2F(*(calibLaser. fHisPz1vP1))), //-> delta kz P1outer-P1inner - line
399 fHisPz2vP1(new TH2F(*(calibLaser.fHisPz2vP1))), //-> delta kz P1outer-P1inner - parabola
400 fHisPz3vP1(new TH2F(*(calibLaser.fHisPz3vP1))), //-> delta kz P1outer-P1inner - common parabola
401 fHisPz2vP2In(new TH2F(*(calibLaser.fHisPz2vP2In))), //-> Curv P2inner - parabola
402 fHisPz2vP2Out(new TH2F(*(calibLaser. fHisPz2vP2Out))), //-> Curv P2outer - parabola
403 fHisPz3vP2IO(new TH2F(*(calibLaser.fHisPz3vP2IO))), //-> Curv P2outerinner - common parabola
404 //
405 //
406 fDeltaYres(((calibLaser.fDeltaYres))),
407 fDeltaZres(((calibLaser.fDeltaZres))),
108953e9 408 fDeltaYres2(((calibLaser.fDeltaYres))),
409 fDeltaZres2(((calibLaser.fDeltaZres))),
7af539c6 410 fDeltaYresAbs(((calibLaser.fDeltaYresAbs))),
41a1f702 411 fHisYAbsErrors(new TH1F(*(calibLaser.fHisYAbsErrors))),
7af539c6 412 fDeltaZresAbs(((calibLaser.fDeltaZresAbs))),
41a1f702 413 fHisZAbsErrors(new TH1F(*(calibLaser.fHisZAbsErrors))),
108953e9 414 // fDeltaYres3(((calibLaser.fDeltaYres))),
415 //fDeltaZres3(((calibLaser.fDeltaZres))),
592a0c8f 416 fFitAside(new TVectorD(5)), // drift fit - A side
417 fFitCside(new TVectorD(5)), // drift fit - C- side
418 fFitACside(new TVectorD(6)), // drift fit - C- side
be113f6e 419 fEdgeXcuts(3), // cuts in local x direction; used in the refit of the laser tracks
420 fEdgeYcuts(3), // cuts in local y direction; used in the refit of the laser tracks
421 fNClCuts(5), // cuts on the number of clusters per tracklet; used in the refit of the laser tracks
7af539c6 422 fNcuts(0), // number of cuts
423 fBeamSectorOuter(336),
424 fBeamSectorInner(336),
425 fBeamOffsetYOuter(336),
426 fBeamSlopeYOuter(336),
427 fBeamOffsetYInner(336),
428 fBeamSlopeYInner(336),
429 fBeamOffsetZOuter(336),
430 fBeamSlopeZOuter(336),
431 fBeamOffsetZInner(336),
432 fBeamSlopeZInner(336),
41a1f702 433 fInverseSlopeZ(calibLaser.fInverseSlopeZ),
434 fUseFixedDriftV(calibLaser.fUseFixedDriftV),
435 fFixedFitAside0(calibLaser.fFixedFitAside0),
436 fFixedFitAside1(calibLaser.fFixedFitAside1),
437 fFixedFitCside0(calibLaser.fFixedFitCside0),
438 fFixedFitCside1(calibLaser.fFixedFitCside1)
be113f6e 439{
440 //
441 // copy constructor
442 //
443}
444
445
446
447AliTPCcalibLaser & AliTPCcalibLaser::operator=(const AliTPCcalibLaser& calibLaser){
448 //
449 // assgnment operator
450 //
451 if (this != &calibLaser) {
452 new (this) AliTPCcalibLaser(calibLaser);
453 }
454 return *this;
455
456}
457
458
459
460
c6914c83 461AliTPCcalibLaser::~AliTPCcalibLaser() {
462 //
463 // destructor
464 //
be113f6e 465 if ( fHisNclIn){
d3ce44cb 466 delete fHisLaser; //->
be113f6e 467 delete fHisNclIn; //->Number of clusters inner
468 delete fHisNclOut; //->Number of clusters outer
469 delete fHisNclIO; //->Number of cluster inner outer
470 delete fHisLclIn; //->Level arm inner
471 delete fHisLclOut; //->Level arm outer
472 delete fHisLclIO; //->Number of cluster inner outer
473 delete fHisdEdx;
474 delete fHisdZfit;
475 //
476 //
477 delete fHisChi2YIn1; //->chi2 y inner - line
478 delete fHisChi2YOut1; //->chi2 y inner - line
479 delete fHisChi2YIn2; //->chi2 y inner - parabola
480 delete fHisChi2YOut2; //->chi2 y inner - parabola
481 delete fHisChi2YIO1; //->chi2 y IO - common
482 delete fHisChi2ZIn1; //->chi2 z inner - line
483 delete fHisChi2ZOut1; //->chi2 z inner - line
484 delete fHisChi2ZIn2; //->chi2 z inner - parabola
485 delete fHisChi2ZOut2; //->chi2 z inner - parabola
486 delete fHisChi2ZIO1; //->chi2 z IO - common
487 //
488 //
489 delete fHisPy1vP0; //-> delta y P0outer-P0inner - line
490 delete fHisPy2vP0; //-> delta y P0outer-P0inner - parabola
491 delete fHisPy3vP0; //-> delta y P0outer-P0inner - common parabola
492 delete fHisPy1vP1; //-> delta ky P1outer-P1inner - line
493 delete fHisPy2vP1; //-> delta ky P1outer-P1inner - parabola
494 delete fHisPy3vP1; //-> delta ky P1outer-P1inner - common parabola
495 delete fHisPy2vP2In; //-> Curv P2inner - parabola
496 delete fHisPy2vP2Out; //-> Curv P2outer - parabola
497 delete fHisPy3vP2IO; //-> Curv P2outerinner - common parabola
498 //
499 //
500 delete fHisPz1vP0; //-> delta z P0outer-P0inner - line
501 delete fHisPz2vP0; //-> delta z P0outer-P0inner - parabola
502 delete fHisPz3vP0; //-> delta z P0outer-P0inner - common parabola
503 delete fHisPz1vP1; //-> delta kz P1outer-P1inner - line
504 delete fHisPz2vP1; //-> delta kz P1outer-P1inner - parabola
505 delete fHisPz3vP1; //-> delta kz P1outer-P1inner - common parabola
506 delete fHisPz2vP2In; //-> Curv P2inner - parabola
507 delete fHisPz2vP2Out; //-> Curv P2outer - parabola
508 delete fHisPz3vP2IO; //-> Curv P2outerinner - common parabola
509 }
5de4d13e 510 //
511 //
512 //
513 fDeltaZ.SetOwner(); //-> array of histograms of delta z for each track
514 fDeltaP3.SetOwner(); //-> array of histograms of P3 for each track
515 fDeltaP4.SetOwner(); //-> array of histograms of P4 for each track
516 fDeltaPhi.SetOwner(); //-> array of histograms of delta z for each track
517 fDeltaPhiP.SetOwner(); //-> array of histograms of delta z for each track
518 fSignals.SetOwner(); //->Array of dedx signals
519
520 fDeltaZ.Delete(); //-> array of histograms of delta z for each track
521 fDeltaP3.Delete(); //-> array of histograms of P3 for each track
522 fDeltaP4.Delete(); //-> array of histograms of P4 for each track
523 fDeltaPhi.Delete(); //-> array of histograms of delta z for each track
524 fDeltaPhiP.Delete(); //-> array of histograms of delta z for each track
525 fSignals.Delete(); //->Array of dedx signals
526
527 fDeltaYres.SetOwner();
528 fDeltaYres.Delete();
41a1f702 529 delete fHisYAbsErrors;
5de4d13e 530 fDeltaZres.SetOwner();
531 fDeltaZres.Delete();
41a1f702 532 delete fHisZAbsErrors;
5de4d13e 533 fDeltaYres2.SetOwner();
534 fDeltaYres2.Delete();
535 fDeltaZres2.SetOwner();
536 fDeltaZres2.Delete();
537
7af539c6 538 fDeltaYresAbs.SetOwner();
539 fDeltaYresAbs.Delete();
540 fDeltaZresAbs.SetOwner();
541 fDeltaZresAbs.Delete();
c6914c83 542}
543
e9f38a4b 544
545
546void AliTPCcalibLaser::Process(AliESDEvent * event) {
c6914c83 547 //
263d466a 548 //
e9f38a4b 549 // Loop over tracks and call Process function
550 //
3e55050f 551 const Int_t kMinTracks=20;
552 const Int_t kMinClusters=40;
38b1a1ca 553
e9f38a4b 554 fESD = event;
555 if (!fESD) {
556 return;
557 }
558 fESDfriend=static_cast<AliESDfriend*>(fESD->FindListObject("AliESDfriend"));
559 if (!fESDfriend) {
560 return;
561 }
108953e9 562 if (fESD->GetNumberOfTracks()<kMinTracks) return; //not enough tracks
263d466a 563 AliDebug(4,Form("Event number in current file: %d",event->GetEventNumberInFile()));
38b1a1ca 564 //
565 // find CE background if present
566 //
567 if (AliTPCLaserTrack::GetTracks()==0) AliTPCLaserTrack::LoadTracks();
82e5df1f 568 TH1D hisCE("hhisCE","hhisCE",100,-100,100);
569 for (Int_t i=0;i<fESD->GetNumberOfTracks();++i) {
570 AliESDtrack *track=fESD->GetTrack(i);
571 if (!track) continue;
572 hisCE.Fill(track->GetZ());
573 hisCE.Fill(track->GetZ()+2);
574 hisCE.Fill(track->GetZ()-2);
38b1a1ca 575 }
576 //
577 //
578
579
e9f38a4b 580 fTracksTPC.Clear();
581 fTracksEsd.Clear();
582 fTracksEsdParam.Delete();
be113f6e 583 for (Int_t id=0; id<336;id++) {
584 fCounter[id]=0;
585 fClusterCounter[id]=0;
586 fClusterSatur[id]=0;
587 }
e9f38a4b 588 //
589 Int_t n=fESD->GetNumberOfTracks();
108953e9 590 Int_t counter=0;
e9f38a4b 591 for (Int_t i=0;i<n;++i) {
592 AliESDfriendTrack *friendTrack=fESDfriend->GetTrack(i);
592a0c8f 593 if (!friendTrack) continue;
e9f38a4b 594 AliESDtrack *track=fESD->GetTrack(i);
82e5df1f 595 Double_t binC = hisCE.GetBinContent(hisCE.FindBin(track->GetZ()));
596 if (binC>336) continue; //remove CE background
e9f38a4b 597 TObject *calibObject=0;
598 AliTPCseed *seed=0;
599 for (Int_t j=0;(calibObject=friendTrack->GetCalibObject(j));++j)
600 if ((seed=dynamic_cast<AliTPCseed*>(calibObject)))
601 break;
3e55050f 602 if (track&&seed&&track->GetTPCNcls()>kMinClusters && seed->GetNumberOfClusters() >kMinClusters) {
5de4d13e 603 //filter CE tracks
108953e9 604 Int_t id = FindMirror(track,seed);
82e5df1f 605 if (id>=0) counter++;
be113f6e 606 }
e9f38a4b 607 //
592a0c8f 608 }
609 fNtracks=counter;
108953e9 610 if (counter<kMinTracks) return;
611
3e55050f 612 //FitDriftV();
82e5df1f 613 FitDriftV(0.3);
3ac724f0 614 if (!fFullCalib) return;
880c3382 615 static Bool_t init=kFALSE;
616 if (!init){
617 init = kTRUE; // way around for PROOF - to be investigated
618 MakeFitHistos();
619 }
e9f38a4b 620 //
be113f6e 621 for (Int_t id=0; id<336; id++){
e9f38a4b 622 //
623 if (!fTracksEsdParam.At(id)) continue;
be113f6e 624 if (fClusterSatur[id]>0.3) continue; // tracks in saturation
e9f38a4b 625 DumpLaser(id);
be113f6e 626 if ( AcceptLaser(id) && fFitZ[id]<0.5){
627 RefitLaserJW(id);
628 MakeDistHisto(id);
629 }
c18f4385 630 }
be113f6e 631
c18f4385 632}
633
be113f6e 634void AliTPCcalibLaser::MakeDistHisto(Int_t id){
c18f4385 635 //
636 //
637 //
be113f6e 638 // for (Int_t id=0; id<336; id++){
c18f4385 639 //
640 //
be113f6e 641 if (!fTracksEsdParam.At(id)) return;
642 if (!AcceptLaser(id)) return;
d3ce44cb 643 Double_t xhis[12]={0,0,0,0,0,0,0,0,0,0,0,0};
c18f4385 644 //
645 //
646 TH1F * hisdz = (TH1F*)fDeltaZ.At(id);
108953e9 647 if (!hisdz) MakeFitHistos();
648 hisdz = (TH1F*)fDeltaZ.At(id);
e5c59218 649 TH1F * hisP3 = (TH1F*)fDeltaP3.At(id);
650 TH1F * hisP4 = (TH1F*)fDeltaP4.At(id);
c18f4385 651 TH1F * hisdphi = (TH1F*)fDeltaPhi.At(id);
652 TH1F * hisdphiP = (TH1F*)fDeltaPhiP.At(id);
2b35e8e6 653 TH1F * hisSignal = (TH1F*)fSignals.At(id);
108953e9 654 //
c18f4385 655
656 AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
657 AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
2b35e8e6 658 AliESDtrack *track = (AliESDtrack*)fTracksEsd.At(id);
c03e3250 659 if (!param) return;
660 if (!ltrp) return;
661 if (!track) return;
c18f4385 662 Double_t xyz[3];
663 Double_t pxyz[3];
664 Double_t lxyz[3];
665 Double_t lpxyz[3];
666 param->GetXYZ(xyz);
667 param->GetPxPyPz(pxyz);
668 ltrp->GetXYZ(lxyz);
669 ltrp->GetPxPyPz(lpxyz);
670 //
671 Float_t dz = param->GetZ()-ltrp->GetZ();
672 Float_t dphi = (TMath::ATan2(xyz[1],xyz[0])- TMath::ATan2(lxyz[1],lxyz[0]))*254.;
673 Float_t dphiP = param->GetParameter()[2]-ltrp->GetParameter()[2];
c03e3250 674 if (hisdz) hisdz->Fill(dz);
e5c59218 675 if (hisP3) hisP3->Fill(param->GetParameter()[3]);
676 if (hisP4) hisP4->Fill(param->GetParameter()[4]);
c03e3250 677 if (hisdphi) hisdphi->Fill(dphi);
263d466a 678 if (hisdphiP) hisdphiP->Fill(dphiP);
e5c59218 679 if (hisSignal) hisSignal->Fill(TMath::Sqrt(TMath::Abs(track->GetTPCsignal())));
d3ce44cb 680 // fill HisLaser
681 xhis[0] = ltrp->GetId();
682 xhis[1] = ltrp->GetSide();
683 xhis[2] = ltrp->GetRod();
684 xhis[3] = ltrp->GetBundle();
685 xhis[4] = ltrp->GetBeam();
686 xhis[5] = dphi;
687 xhis[6] = fFitZ[id];
688 xhis[7] = param->GetParameter()[2]-ltrp->GetParameter()[2]; //dp2
689 xhis[8] = param->GetParameter()[3]-ltrp->GetParameter()[3]; //dp3
690 xhis[9] = param->GetParameter()[4];
691 xhis[10]= track->GetTPCNcls();
692 xhis[11]= TMath::Sqrt(TMath::Abs(track->GetTPCsignal()));
be113f6e 693 // }
d3ce44cb 694 fHisLaser->Fill(xhis);
c18f4385 695}
696
697void AliTPCcalibLaser::FitDriftV(){
698 //
108953e9 699 // Fit corrections to the drift velocity - linear approximation in the z and global y
700 //The transfromatiom from the drift time to the z position done in AliTPCTracnsform class
701 //
702 /*
703 Formulas:
704
705 z = s* (z0 - vd*(t-t0))
706
707 s - side -1 and +1
708 t0 - time 0
709 vd - nominal drift velocity
710 zs - miscalibrated position
711
712 zs = s*(z0 - vd*(1+vr)*(t-(t0+dt))
713 vr - relative change of the drift velocity
714 dzt - vd*dt
715 dr = zz0-s*z
716 ..
717 ==>
718 zs ~ z - s*vr*(z0-s*z)+s*dzt
719 --------------------------------
720 1. Correction function vr constant:
721
722
723 dz = zs-z = -s*vr *(z0-s*z)+s*dzt
724 dzs/dl = dz/dl +s*s*vr*dz/dl
725 d(dz/dl) = vr*dz/dl
726 */
727
728 const Float_t kZCut = 240; // remove the closest laser beam
be113f6e 729 const Float_t kSaturCut = 0.05; // remove saturated lasers - cut on fraction of saturated
108953e9 730 const Float_t kDistCut = 3; // distance sigma cut
731 const Float_t kDistCutAbs = 0.25;
be113f6e 732 const Float_t kMinClusters = 60; // minimal amount of the clusters
733 const Float_t kMinSignal = 16; // minimal mean height of the signal
734 const Float_t kChi2Cut = 0.1; // chi2 cut to accept drift fit
c18f4385 735 static TLinearFitter fdriftA(3,"hyp2");
736 static TLinearFitter fdriftC(3,"hyp2");
108953e9 737 static TLinearFitter fdriftAC(4,"hyp3");
738 TVectorD fitA(3),fitC(3),fitAC(4);
739
740 AliTPCcalibDB* calib=AliTPCcalibDB::Instance();
f6220766 741 AliTPCParam * tpcparam = calib->GetParameters();
108953e9 742 //
be113f6e 743 for (Int_t id=0; id<336; id++) fFitZ[id]=0;
744
108953e9 745 Float_t chi2A = 10;
746 Float_t chi2C = 10;
747 Float_t chi2AC = 10;
748 Int_t npointsA=0;
749 Int_t npointsC=0;
750 Int_t npointsAC=0;
be113f6e 751
108953e9 752
753 for (Int_t iter=0; iter<3; iter++){
be113f6e 754 fdriftA.ClearPoints();
755 fdriftC.ClearPoints();
108953e9 756 fdriftAC.ClearPoints();
be113f6e 757 //
758 for (Int_t id=0; id<336; id++){
759 if (!fTracksEsdParam.At(id)) continue;
760 if (!AcceptLaser(id)) continue;
be113f6e 761 if ( fClusterSatur[id]>kSaturCut) continue;
762 if ( fClusterCounter[id]<kMinClusters) continue;
108953e9 763 AliESDtrack *track = (AliESDtrack*)fTracksEsd.At(id);
764 if (track->GetTPCsignal()<kMinSignal) continue;
be113f6e 765 AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
766 AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
108953e9 767
be113f6e 768 Double_t xyz[3];
769 Double_t pxyz[3];
770 Double_t lxyz[3];
771 Double_t lpxyz[3];
772 param->GetXYZ(xyz);
773 param->GetPxPyPz(pxyz);
774 ltrp->GetXYZ(lxyz);
775 ltrp->GetPxPyPz(lpxyz);
776 if (TMath::Abs(lxyz[2])>kZCut) continue;
108953e9 777 Float_t sz = (ltrp->GetSide()==0) ? TMath::Sqrt(chi2A): TMath::Sqrt(chi2C);
778 if (npointsAC>0) sz =TMath::Sqrt(chi2AC);
779 if (TMath::Abs(fFitZ[id])>sz*kDistCut) continue;
780 if (iter>0 && TMath::Abs(fFitZ[id])>kDistCutAbs) continue;
781
782 // drift distance
783 Double_t zlength = tpcparam->GetZLength(0);
784 Double_t ldrift = zlength-TMath::Abs(lxyz[2]);
785 Double_t mdrift = zlength-TMath::Abs(xyz[2]);
786 Double_t xxx[2] = {ldrift,lxyz[1]*ldrift/(zlength*250.)};
be113f6e 787 if (ltrp->GetSide()==0){
108953e9 788 fdriftA.AddPoint(xxx,mdrift,1);
be113f6e 789 }else{
108953e9 790 fdriftC.AddPoint(xxx,mdrift,1);
be113f6e 791 }
108953e9 792 Double_t xxx2[3] = {ltrp->GetSide(), ldrift,lxyz[1]*ldrift/(zlength*250.)};
793 fdriftAC.AddPoint(xxx2,mdrift,1);
be113f6e 794 }
be113f6e 795 //
796 if (fdriftA.GetNpoints()>10){
108953e9 797 //
be113f6e 798 fdriftA.Eval();
108953e9 799 if (iter==0) fdriftA.EvalRobust(0.7);
800 else fdriftA.EvalRobust(0.8);
be113f6e 801 fdriftA.GetParameters(fitA);
802 npointsA= fdriftA.GetNpoints();
803 chi2A = fdriftA.GetChisquare()/fdriftA.GetNpoints();
880c3382 804 if (chi2A<kChi2Cut ||(*fFitAside)[0]==0 ) {
805 if (fFitAside->GetNoElements()<5) fFitAside->ResizeTo(5);
806 (*fFitAside)[0] = fitA[0];
807 (*fFitAside)[1] = fitA[1];
808 (*fFitAside)[2] = fitA[2];
809 (*fFitAside)[3] = fdriftA.GetNpoints();
810 (*fFitAside)[4] = chi2A;
811 }
be113f6e 812 }
813 if (fdriftC.GetNpoints()>10){
814 fdriftC.Eval();
108953e9 815 if (iter==0) fdriftC.EvalRobust(0.7);
816 else fdriftC.EvalRobust(0.8);
817 //
be113f6e 818 fdriftC.GetParameters(fitC);
819 npointsC= fdriftC.GetNpoints();
820 chi2C = fdriftC.GetChisquare()/fdriftC.GetNpoints();
880c3382 821 if (chi2C<kChi2Cut||(*fFitCside)[0]==0) {
822 if (fFitCside->GetNoElements()<5) fFitCside->ResizeTo(5);
823 (*fFitCside)[0] = fitC[0];
824 (*fFitCside)[1] = fitC[1];
825 (*fFitCside)[2] = fitC[2];
826 (*fFitCside)[3] = fdriftC.GetNpoints();
827 (*fFitCside)[4] = chi2C;
828 }
be113f6e 829 }
108953e9 830
831 if (fdriftAC.GetNpoints()>10&&fdriftC.GetNpoints()>10&&fdriftA.GetNpoints()>10){
832 fdriftAC.Eval();
833 if (iter==0) fdriftAC.EvalRobust(0.7);
834 else fdriftAC.EvalRobust(0.8);
835 //
836 fdriftAC.GetParameters(fitAC);
837 npointsAC= fdriftAC.GetNpoints();
838 chi2AC = fdriftAC.GetChisquare()/fdriftAC.GetNpoints();
839 if (chi2AC<kChi2Cut||(*fFitACside)[0]==0) (*fFitACside) = fitAC;
840 }
be113f6e 841
842 for (Int_t id=0; id<336; id++){
843 if (!fTracksEsdParam.At(id)) continue;
844 //
845 AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
846 AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
847 Double_t xyz[3];
848 Double_t pxyz[3];
849 Double_t lxyz[3];
850 Double_t lpxyz[3];
851 param->GetXYZ(xyz);
852 param->GetPxPyPz(pxyz);
853 ltrp->GetXYZ(lxyz);
108953e9 854 ltrp->GetPxPyPz(lpxyz);
855 Double_t zlength = tpcparam->GetZLength(0);
856 Double_t ldrift = zlength-TMath::Abs(lxyz[2]);
857 Double_t mdrift = zlength-TMath::Abs(xyz[2]);
858
be113f6e 859 Float_t fz =0;
860 if (ltrp->GetSide()==0){
108953e9 861 fz = (fitA)[0]+(fitA)[1]*ldrift+(fitA)[2]*lxyz[1]*ldrift/(zlength*250.);
be113f6e 862 }else{
108953e9 863 fz = (fitC)[0]+(fitC)[1]*ldrift+(fitC)[2]*lxyz[1]*ldrift/(zlength*250.);
be113f6e 864 }
108953e9 865 if (npointsAC>10){
866 fz = (fitAC)[0]+(fitAC)[1]*ltrp->GetSide()+(fitAC)[2]*ldrift+(fitAC)[3]*lxyz[1]*ldrift/(zlength*250.);
867 }
868 fFitZ[id]=mdrift-fz;
be113f6e 869 }
870 if (fStreamLevel>0){
871 TTreeSRedirector *cstream = GetDebugStreamer();
3ac724f0 872 TTimeStamp tstamp(fTime);
873 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
874 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
875 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
876 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
877 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
878 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
31aa7c5c 879 TVectorD vecGoofie(20);
880 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
881 if (goofieArray)
882 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
883 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
884 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
885 }
3ac724f0 886
be113f6e 887 if (cstream){
108953e9 888 (*cstream)<<"driftv"<<
889 "run="<<fRun<< // run number
890 "event="<<fEvent<< // event number
891 "time="<<fTime<< // time stamp of event
892 "trigger="<<fTrigger<< // trigger
893 "mag="<<fMagF<< // magnetic field
3ac724f0 894 // Environment values
895 "press0="<<valuePressure0<<
896 "press1="<<valuePressure1<<
897 "pt0="<<ptrelative0<<
898 "pt1="<<ptrelative1<<
899 "temp0="<<temp0<<
900 "temp1="<<temp1<<
31aa7c5c 901 "vecGoofie.="<<&vecGoofie<<
3ac724f0 902 //
903 //
be113f6e 904 "iter="<<iter<<
905 "driftA.="<<fFitAside<<
906 "driftC.="<<fFitCside<<
108953e9 907 "driftAC.="<<fFitACside<<
be113f6e 908 "chi2A="<<chi2A<<
909 "chi2C="<<chi2C<<
108953e9 910 "chi2AC="<<chi2AC<<
be113f6e 911 "nA="<<npointsA<<
912 "nC="<<npointsC<<
108953e9 913 "nAC="<<npointsAC<<
be113f6e 914 "\n";
915 }
c18f4385 916 }
917 }
be113f6e 918}
f6220766 919Bool_t AliTPCcalibLaser::FitDriftV(Float_t minFraction){
920 //
921 // Fit corrections to the drift velocity - linear approximation in the z and global y
922 //The transfromatiom from the drift time to the z position done in AliTPCTracnsform class
923 //
38b1a1ca 924 // Source of outlyers :
925 // 0. Track in the saturation - postpeak
926 // 1. gating grid close the part of the signal for first bundle
927
928 // The robust fit is performed in 2 itterations /robust fraction controlled by kFraction/
929 // 1. Robust fit is used in the itteration number 0
930 // only fraction of laser uted
931 // 2. Only the tracks close to the fit used in the second itteration
f6220766 932 /*
933 Formulas:
934
935 z = s* (z0 - vd*(t-t0))
936
937 s - side -1 and +1
938 t0 - time 0
939 vd - nominal drift velocity
940 zs - miscalibrated position
941
942 zs = s*(z0 - vd*(1+vr)*(t-(t0+dt))
943 vr - relative change of the drift velocity
944 dzt - vd*dt
945 dr = zz0-s*z
946 ..
947 ==>
948 zs ~ z - s*vr*(z0-s*z)+s*dzt
949 --------------------------------
950 1. Correction function vr constant:
951
952
953 dz = zs-z = -s*vr *(z0-s*z)+s*dzt
954 dzs/dl = dz/dl +s*s*vr*dz/dl
955 d(dz/dl) = vr*dz/dl
956 */
38b1a1ca 957 const Int_t knLaser = 336; //n laser tracks
958 const Float_t kFraction[2] = {0.70,0.95}; // robust fit fraction
959
960 const Float_t kSaturCut = 0.05; // remove saturated lasers - cut on fraction of saturated
961 const Float_t kDistCut = 3.; // distance sigma cut - 3 sigma
962 const Float_t kDistCutAbs = 1.; // absolute cut 1 cm
963 const Float_t kMinClusters = 60.; // minimal amount of the clusters
964 const Float_t kMinSignal = 10.; // minimal mean height of the signal
965 const Float_t kChi2Cut = 1.0; // chi2 cut to accept drift fit
f6220766 966 //
967 static TLinearFitter fdriftA(3,"hyp2");
968 static TLinearFitter fdriftC(3,"hyp2");
969 static TLinearFitter fdriftAC(4,"hyp3");
970 TVectorD fitA(3),fitC(3),fitAC(4);
971
972 AliTPCcalibDB* calib=AliTPCcalibDB::Instance();
973 AliTPCParam * tpcparam = calib->GetParameters();
974 //
975 // reset old data
976 //
977 for (Int_t id=0; id<336; id++) fFitZ[id]=0;
978 if (fFitAside->GetNoElements()<5) fFitAside->ResizeTo(5);
979 if (fFitCside->GetNoElements()<5) fFitCside->ResizeTo(5);
980 for (Int_t i=0;i<5; i++){
981 (*fFitCside)[i]=0;
982 (*fFitAside)[i]=0;
983 }
984 //
985 //
986 Float_t chi2A = 10;
987 Float_t chi2C = 10;
988 Float_t chi2AC = 10;
989 Int_t npointsA=0;
990 Int_t npointsC=0;
991 Int_t npointsAC=0;
992 Int_t nbA[4]={0,0,0,0};
993 Int_t nbC[4]={0,0,0,0};
82e5df1f 994 TVectorD vecZM(336); // measured z potion of laser
995 TVectorD vecA(336); // accepted laser
996 TVectorD vecZF(336); // fitted position
997 TVectorD vecDz(336); // deltaZ
998 TVectorD vecZS(336); // surveyed position of laser
999 // additional variable to cut
1000 TVectorD vecdEdx(336); // dEdx
1001 TVectorD vecSy(336); // shape y
1002 TVectorD vecSz(336); // shape z
1003 //
f6220766 1004 //
1005 for (Int_t id=0; id<336; id++){
82e5df1f 1006 Int_t reject=0;
1007 AliTPCLaserTrack *ltrp =
1008 (AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(id);
1009 AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
1010 AliTPCseed * seed = (AliTPCseed*)fTracksTPC.At(id);
1011 vecZM(id)= (param==0) ? 0:param->GetZ();
1012 vecZS(id)= ltrp->GetZ();
1013 vecDz(id)= 0;
1014 vecA(id)=1;
1015 vecdEdx(id)=(seed)?seed->GetdEdx():0;
1016 vecSy(id) =(seed)?seed->CookShape(1):0;
1017 vecSz(id) =(seed)?seed->CookShape(2):0;
1018 //
1019 fFitZ[id]=0;
1020 if (!fTracksEsdParam.At(id)) reject|=1;
1021 if (!param) continue;
1022 if (!AcceptLaser(id)) reject|=2;
1023 if ( fClusterSatur[id]>kSaturCut) reject|=4;
1024 if ( fClusterCounter[id]<kMinClusters) reject|=8;
1025 vecA(id)=reject;
1026 if (reject>0) continue;
f6220766 1027 if (ltrp->GetSide()==0){
1028 npointsA++;
1029 nbA[ltrp->GetBundle()]++;
1030 }
1031 if (ltrp->GetSide()==1){
82e5df1f 1032 npointsC++;
f6220766 1033 nbC[ltrp->GetBundle()]++;
1034 }
1035 }
1036 //
1037 // reject "bad events"
1038 //
1039 Bool_t isOK=kTRUE;
1040 Int_t naA=0;
1041 Int_t naC=0;
1042 if (npointsA<minFraction*0.5*knLaser && npointsC<minFraction*0.5*knLaser)
1043 isOK=kFALSE;
1044 for (Int_t i=0;i<4;i++){
1045 //count accepted for all layers
1046 if (nbA[i]>minFraction*0.5*0.25*knLaser) naA++;
1047 if (nbC[i]>minFraction*0.5*0.25*knLaser) naC++;
1048 }
1049 if (naA<3 &&naC<3) isOK=kFALSE;
1050 if (isOK==kFALSE) return kFALSE;
1051
1052 //
1053 //
1054 //
1055 for (Int_t iter=0; iter<2; iter++){
1056 fdriftA.ClearPoints();
1057 fdriftC.ClearPoints();
82e5df1f 1058 fdriftAC.ClearPoints();
1059 npointsA=0; npointsC=0; npointsAC=0;
f6220766 1060 //
1061 for (Int_t id=0; id<336; id++){
1062 if (!fTracksEsdParam.At(id)) continue;
1063 if (!AcceptLaser(id)) continue;
1064 if ( fClusterSatur[id]>kSaturCut) continue;
1065 if ( fClusterCounter[id]<kMinClusters) continue;
1066 AliESDtrack *track = (AliESDtrack*)fTracksEsd.At(id);
1067 if (track->GetTPCsignal()<kMinSignal) continue;
1068 AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
1069 AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
f6220766 1070 Double_t xyz[3];
1071 Double_t pxyz[3];
1072 Double_t lxyz[3];
1073 Double_t lpxyz[3];
1074 param->GetXYZ(xyz);
1075 param->GetPxPyPz(pxyz);
1076 ltrp->GetXYZ(lxyz);
1077 ltrp->GetPxPyPz(lpxyz);
1078 Float_t sz = (ltrp->GetSide()==0) ? TMath::Sqrt(chi2A): TMath::Sqrt(chi2C);
1079 if (npointsAC>0) sz =TMath::Sqrt(chi2AC);
82e5df1f 1080 if (iter>0 && TMath::Abs(fFitZ[id])>sz*kDistCut) continue;
f6220766 1081 if (iter>0 && TMath::Abs(fFitZ[id])>kDistCutAbs) continue;
1082
82e5df1f 1083// // drift distance
1084// Double_t zlength = tpcparam->GetZLength(0);
1085// Double_t ldrift = zlength-TMath::Abs(lxyz[2]);
1086// Double_t mdrift = zlength-TMath::Abs(xyz[2]);
1087 //
1088 Double_t zlength = (ltrp->GetSide()==0)? tpcparam->GetZLength(36): tpcparam->GetZLength(71);
1089 Double_t ldrift = (ltrp->GetSide()==0)?zlength-lxyz[2]:lxyz[2]+zlength;
1090 Double_t mdrift = (ltrp->GetSide()==0)?zlength-xyz[2]:xyz[2]+zlength;
1091
1092
f6220766 1093 Double_t xxx[2] = {ldrift,lxyz[1]*ldrift/(zlength*250.)};
82e5df1f 1094 if (iter==0 &&ltrp->GetBundle()==0) continue;
1095 // skip bundle 0 - can be wrong in the 0 iteration
f6220766 1096 if (ltrp->GetSide()==0){
1097 fdriftA.AddPoint(xxx,mdrift,1);
1098 }else{
1099 fdriftC.AddPoint(xxx,mdrift,1);
1100 }
1101 Double_t xxx2[3] = {ltrp->GetSide(), ldrift,lxyz[1]*ldrift/(zlength*250.)};
1102 fdriftAC.AddPoint(xxx2,mdrift,1);
1103 }
1104 //
1105 if (fdriftA.GetNpoints()>minFraction*0.5*knLaser){
1106 //
1107 fdriftA.Eval();
1108 npointsA= fdriftA.GetNpoints();
1109 chi2A = fdriftA.GetChisquare()/fdriftA.GetNpoints();
38b1a1ca 1110 fdriftA.EvalRobust(kFraction[iter]);
f6220766 1111 fdriftA.GetParameters(fitA);
1112 if (chi2A<kChi2Cut ||(*fFitAside)[0]==0 ) {
1113 if (fFitAside->GetNoElements()<5) fFitAside->ResizeTo(5);
1114 (*fFitAside)[0] = fitA[0];
1115 (*fFitAside)[1] = fitA[1];
1116 (*fFitAside)[2] = fitA[2];
1117 (*fFitAside)[3] = fdriftA.GetNpoints();
1118 (*fFitAside)[4] = chi2A;
1119 }
1120 }
1121 if (fdriftC.GetNpoints()>minFraction*0.5*knLaser){
1122 fdriftC.Eval();
1123 npointsC= fdriftC.GetNpoints();
1124 chi2C = fdriftC.GetChisquare()/fdriftC.GetNpoints();
38b1a1ca 1125 fdriftC.EvalRobust(kFraction[iter]);
f6220766 1126 fdriftC.GetParameters(fitC);
1127 if (chi2C<kChi2Cut||(*fFitCside)[0]==0) {
1128 if (fFitCside->GetNoElements()<5) fFitCside->ResizeTo(5);
1129 (*fFitCside)[0] = fitC[0];
1130 (*fFitCside)[1] = fitC[1];
1131 (*fFitCside)[2] = fitC[2];
1132 (*fFitCside)[3] = fdriftC.GetNpoints();
1133 (*fFitCside)[4] = chi2C;
1134 }
1135 }
1136
3e55050f 1137 if (fdriftAC.GetNpoints()>minFraction*knLaser &&npointsA>0.5*minFraction*knLaser&&npointsC>0.5*minFraction*knLaser){
f6220766 1138 fdriftAC.Eval();
1139 npointsAC= fdriftAC.GetNpoints();
1140 chi2AC = fdriftAC.GetChisquare()/fdriftAC.GetNpoints();
38b1a1ca 1141 fdriftAC.EvalRobust(kFraction[iter]);
f6220766 1142 fdriftAC.GetParameters(fitAC);
1143 if (chi2AC<kChi2Cut||(*fFitACside)[0]==0) (*fFitACside) = fitAC;
1144 (*fFitACside)[0] = fitAC[0];
1145 (*fFitACside)[1] = fitAC[1];
1146 (*fFitACside)[2] = fitAC[2];
1147 (*fFitACside)[3] = fdriftAC.GetNpoints();
1148 (*fFitACside)[4] = chi2AC;
1149 }
1150
1151 for (Int_t id=0; id<336; id++){
1152 if (!fTracksEsdParam.At(id)) continue;
1153 //
1154 AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
1155 AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
1156 Double_t xyz[3];
1157 Double_t pxyz[3];
1158 Double_t lxyz[3];
1159 Double_t lpxyz[3];
1160 param->GetXYZ(xyz);
1161 param->GetPxPyPz(pxyz);
1162 ltrp->GetXYZ(lxyz);
1163 ltrp->GetPxPyPz(lpxyz);
82e5df1f 1164 //Double_t zlength = tpcparam->GetZLength(0);
1165 //Double_t ldrift = zlength-TMath::Abs(lxyz[2]);
1166 //Double_t mdrift = zlength-TMath::Abs(xyz[2]);
1167 Double_t zlength = (ltrp->GetSide()==0)? tpcparam->GetZLength(36): tpcparam->GetZLength(71);
1168 Double_t ldrift = (ltrp->GetSide()==0)?zlength-lxyz[2]:lxyz[2]+zlength;
1169 Double_t mdrift = (ltrp->GetSide()==0)?zlength-xyz[2]:xyz[2]+zlength;
1170
f6220766 1171
1172 Float_t fz =0;
1173 if (ltrp->GetSide()==0){
1174 fz = (fitA)[0]+(fitA)[1]*ldrift+(fitA)[2]*lxyz[1]*ldrift/(zlength*250.);
1175 }else{
1176 fz = (fitC)[0]+(fitC)[1]*ldrift+(fitC)[2]*lxyz[1]*ldrift/(zlength*250.);
1177 }
1178 if (npointsAC>10){
82e5df1f 1179 //fz = (fitAC)[0]+(fitAC)[1]*ltrp->GetSide()+(fitAC)[2]*ldrift+(fitAC)[3]*lxyz[1]*ldrift/(zlength*250.);
f6220766 1180 }
1181 fFitZ[id]=mdrift-fz;
82e5df1f 1182 vecZF[id]=fz;
1183 vecDz[id]=mdrift-fz;
f6220766 1184 }
1185 if (fStreamLevel>0){
1186 TTreeSRedirector *cstream = GetDebugStreamer();
1187 TTimeStamp tstamp(fTime);
1188 Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
1189 Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
1190 Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
1191 Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
1192 Double_t temp0 = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
1193 Double_t temp1 = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
1194 TVectorD vecGoofie(20);
1195 AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
1196 if (goofieArray)
1197 for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
1198 AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
1199 if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
1200 }
1201
1202 if (cstream){
1203 (*cstream)<<"driftvN"<<
1204 "run="<<fRun<< // run number
1205 "event="<<fEvent<< // event number
1206 "time="<<fTime<< // time stamp of event
1207 "trigger="<<fTrigger<< // trigger
1208 "mag="<<fMagF<< // magnetic field
1209 // Environment values
1210 "press0="<<valuePressure0<<
1211 "press1="<<valuePressure1<<
1212 "pt0="<<ptrelative0<<
1213 "pt1="<<ptrelative1<<
1214 "temp0="<<temp0<<
1215 "temp1="<<temp1<<
1216 "vecGoofie.="<<&vecGoofie<<
1217 //
1218 //
82e5df1f 1219 "vecZM.="<<&vecZM<< // measured z position
1220 "vecZS.="<<&vecZS<< // surveyed z position
1221 "vecZF.="<<&vecZF<< // fitted z position
1222 "vecDz.="<<&vecDz<< // fitted z position
1223 "vecA.="<<&vecA<< // accept laser flag
1224 "vecdEdx.="<<&vecdEdx<< // dEdx - to cut on
1225 "vecSy.="<<&vecSy<< // shape y - to cut on
1226 "vecSz.="<<&vecSz<< // shape z - to cut on
1227 //
f6220766 1228 "iter="<<iter<<
1229 "driftA.="<<fFitAside<<
1230 "driftC.="<<fFitCside<<
1231 "driftAC.="<<fFitACside<<
1232 "chi2A="<<chi2A<<
1233 "chi2C="<<chi2C<<
1234 "chi2AC="<<chi2AC<<
1235 "nA="<<npointsA<<
1236 "nC="<<npointsC<<
1237 "nAC="<<npointsAC<<
1238 "\n";
82e5df1f 1239 /*
1240 //
1241 variables to check in debug mode:
1242 //
1243 chainDriftN->SetAlias("driftS","250-abs(vecZS.fElements)");
1244 chainDriftN->SetAlias("driftM","250-abs(vecZM.fElements)");
1245 chainDriftN->SetAlias("driftF","vecZF.fElements");
1246 chainDriftN->SetAlias("deltaZ","driftF-driftM"); //deltaZ
1247 TCut cutA="iter==1&&sqrt(chi2A)<0.1&&vecZS.fElements>0&&vecA.fElements==0";
1248 TCut cutC="iter==1&&sqrt(chi2C)<0.1&&vecZS.fElements<0&&vecA.fElements==0";
1249
1250 */
f6220766 1251 }
1252 }
1253 }
1254 return kTRUE;
1255}
263d466a 1256
be113f6e 1257Float_t AliTPCcalibLaser::GetDistance(AliExternalTrackParam *param, AliTPCLaserTrack *ltrp){
1258 //
1259 // get distance between mirror and laser track
1260 //
1261 //
1262 Double_t xyz[3];
1263 Double_t lxyz[3];
1264 param->GetXYZ(xyz);
1265 ltrp->GetXYZ(lxyz);
1266 //
1267 //
1268 Double_t dist = 0;
1269 //radial distance
1270 dist+=TMath::Abs((TMath::ATan2(xyz[1],xyz[0])-TMath::ATan2(lxyz[1],lxyz[0]))*param->GetX());
1271 //
1272 // z distance
1273 // apply drift correction if already exist
1274 //
1275 Float_t dz = 0;
1276 if (ltrp->GetSide()==0){
1277 if ((*fFitAside)[1]>0.) dz = ((*fFitAside)[0]+(*fFitAside)[1]*lxyz[2]+(*fFitAside)[2]*lxyz[1])-xyz[2];
1278 }else{
1279 if ((*fFitCside)[1]>0.) dz = ((*fFitCside)[0]+(*fFitCside)[1]*lxyz[2]+(*fFitCside)[2]*lxyz[1])-xyz[2];
c18f4385 1280 }
be113f6e 1281 if (TMath::Abs(dz)>3*(TMath::Abs(lxyz[2]-xyz[2])+1)) dz= TMath::Abs(lxyz[2]-xyz[2]);
1282 dist+=TMath::Abs(dz);
1283 //
1284 // phi dist - divergence on 50 cm
c18f4385 1285 //
be113f6e 1286 dist = TMath::Abs((param->GetParameter()[2]-ltrp->GetParameter()[2])*50);
1287 return dist;
c18f4385 1288}
1289
1290
1291Bool_t AliTPCcalibLaser::AcceptLaser(Int_t id){
1292 //
1293 //
1294 //
1295 /*
be113f6e 1296 TCut cutP0("cutP0","abs((atan2(x1,x0)-atan2(lx1,lx0))*254)<1.5");
e5c59218 1297 TCut cutP1("cutP1","abs(LTr.fP[1]-Tr.fP[1])<30");
1298 TCut cutP2("cutP2","abs(LTr.fP[2]-Tr.fP[2])<0.03");
1299 TCut cutP3("cutP3","abs(Tr.fP[3])<0.05");
e5c59218 1300
1301 TCut cutA = cutP0+cutP1+cutP2+cutP3+cutP4;
c18f4385 1302 */
1303 AliExternalTrackParam *param =(AliExternalTrackParam*)fTracksEsdParam.At(id);
1304 AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
43ad9d4c 1305 //AliESDtrack *track = (AliESDtrack*)fTracksEsd.At(id);
e5c59218 1306 Double_t xyz[3];
1307 Double_t lxyz[3];
1308 param->GetXYZ(xyz);
1309 ltrp->GetXYZ(lxyz);
be113f6e 1310 if (TMath::Abs((TMath::ATan2(xyz[1],xyz[0])-TMath::ATan2(lxyz[1],lxyz[0]))*254)>1.5) return kFALSE; //cut y- P0
e5c59218 1311 if (TMath::Abs(param->GetParameter()[1]-ltrp->GetParameter()[1])>30) return kFALSE; // cutZ -P1
1312 if (TMath::Abs(param->GetParameter()[2]-ltrp->GetParameter()[2])>0.03) return kFALSE; // cut -P2
1313 if (TMath::Abs(param->GetParameter()[3])>0.05) return kFALSE; // cut Tl -P3
c18f4385 1314 //
c18f4385 1315 //
e5c59218 1316
263d466a 1317 return kTRUE;
e9f38a4b 1318}
1319
1320Int_t AliTPCcalibLaser::FindMirror(AliESDtrack *track, AliTPCseed *seed){
1321 //
1322 // Find corresponding mirror
1323 // add the corresponding tracks
108953e9 1324
e9f38a4b 1325 if (!track->GetOuterParam()) return -1;
263d466a 1326
38b1a1ca 1327 Float_t kRadius0 = 252;
1328 Float_t kRadius = 254.2;
be113f6e 1329 Int_t countercl=0;
1330 Float_t counterSatur=0;
f6220766 1331 Int_t csideA =0;
1332 Int_t csideC =0;
be113f6e 1333 for (Int_t irow=158;irow>-1;--irow) {
1334 AliTPCclusterMI *c=seed->GetClusterPointer(irow);
1335 if (!c) continue;
1336 Double_t pedgeY = c->GetX()*TMath::DegToRad()*(10)-TMath::Abs(c->GetY());
1337 Double_t pedgeX = TMath::Min((irow)*0.75, (159.-irow)*1.5);
1338 if (pedgeY<3) continue;
1339 if (pedgeX<3) continue;
1340 countercl++;
f6220766 1341 if (c->GetDetector()%36<18) csideA++;
1342 if (c->GetDetector()%36>=18) csideC++;
be113f6e 1343 if (c->GetMax()>900) counterSatur++;
1344 }
1345 counterSatur/=(countercl+1);
1346 //
f6220766 1347 //
1348 //
38b1a1ca 1349 if (csideA<0.9*seed->GetNumberOfClusters() && csideC<0.9*seed->GetNumberOfClusters()) return 0; // cross laser track can not happen
1350
1351 Int_t side= 0;
1352 if (csideC>0.5*seed->GetNumberOfClusters()) side=1;
1353
1354
1355 AliExternalTrackParam param(*(track->GetOuterParam()));
1356 AliTracker::PropagateTrackTo(&param,kRadius0,0.10566,3,kTRUE);
1357 AliTracker::PropagateTrackTo(&param,kRadius,0.10566,0.1,kTRUE);
1358 AliTPCLaserTrack ltr;
1359 AliTPCLaserTrack *ltrp=0x0;
1360 // AliTPCLaserTrack *ltrpjw=0x0;
1361 //
1362 Int_t id = AliTPCLaserTrack::IdentifyTrack(&param,side);
1363 // Int_t idjw = AliTPCLaserTrack::IdentifyTrackJW(&param);
1364 //AliDebug(4,Form("Identified Track: %03d (%03d)",id,idjw));
1365
1366 if (id!=-1 && (AliTPCLaserTrack::GetTracks()->UncheckedAt(id)))
1367 ltrp=(AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(id);
1368 else
1369 ltrp=&ltr;
1370
1371 if (id<0) return -1;
1372 if (ltrp->GetSide()!=side) return -1;
1373 fCounter[id]++;
1374 //
1375 //
1376 //
f6220766 1377 //
be113f6e 1378 if (counterSatur>fClusterSatur[id]) fClusterSatur[id]=counterSatur;
1379 //
1380 //
1381 Float_t radius=TMath::Abs(ltrp->GetX());
880c3382 1382 param.Rotate(ltrp->GetAlpha());
38b1a1ca 1383 AliTracker::PropagateTrackTo(&param,radius,0.10566,0.01,kFALSE);
be113f6e 1384 //
1385 if (!fTracksMirror.At(id)) fTracksMirror.AddAt(ltrp,id);
1386 Bool_t accept=kTRUE;
1387 //
1388 // choose closer track
1389 //
1390 AliExternalTrackParam * param0 = (AliExternalTrackParam *)fTracksEsdParam.At(id);
1391 if (param0){
1392 Float_t dist0=GetDistance(param0,ltrp);
1393 Float_t dist1=GetDistance(&param,ltrp);
1394 if (dist0<dist1) accept=kFALSE;
1395 }
1396
1397 if (accept){
1398 fClusterCounter[id]=countercl;
e9f38a4b 1399 fTracksEsdParam.AddAt(param.Clone(),id);
1400 fTracksEsd.AddAt(track,id);
1401 fTracksTPC.AddAt(seed,id);
1402 }
1403 return id;
1404}
1405
1406
1407
1408void AliTPCcalibLaser::DumpLaser(Int_t id) {
1409 //
1410 // Dump Laser info to the tree
1411 //
1412 AliESDtrack *track = (AliESDtrack*)fTracksEsd.At(id);
1413 AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
1414 AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
1bd37cc0 1415 //
1416 // Fast laser ID
1417 //
1418 Double_t xyz[3];
1419 Double_t pxyz[3];
c18f4385 1420 Double_t lxyz[3];
1421 Double_t lpxyz[3];
e9f38a4b 1422 param->GetXYZ(xyz);
1423 param->GetPxPyPz(pxyz);
c18f4385 1424 ltrp->GetXYZ(lxyz);
1425 ltrp->GetPxPyPz(lpxyz);
be113f6e 1426 Float_t dist3D = GetDistance(param,ltrp);
1427 Float_t dist0 = (TMath::ATan2(xyz[1],xyz[0])-TMath::ATan2(lxyz[1],lxyz[0]))*param->GetX();
1428 Float_t distphi = (param->GetParameter()[2]-ltrp->GetParameter()[2])*50;
1429
1430
c6914c83 1431 if (fStreamLevel>0){
1432 TTreeSRedirector *cstream = GetDebugStreamer();
c18f4385 1433 Int_t time = fESD->GetTimeStamp();
1434 Bool_t accept = AcceptLaser(id);
c6914c83 1435 if (cstream){
1436 (*cstream)<<"Track"<<
108953e9 1437 //
1438 "run="<<fRun<< // run number
1439 "event="<<fEvent<< // event number
1440 "time="<<fTime<< // time stamp of event
1441 "trigger="<<fTrigger<< // trigger
1442 "mag="<<fMagF<< // magnetic field
1443
7b18d067 1444 "id="<<id<<
c18f4385 1445 "accept="<<accept<<
1446 "driftA.="<<fFitAside<<
1447 "driftC.="<<fFitCside<<
1448 "time="<<time<<
be113f6e 1449 "dist3D="<<dist3D<<
1450 "dist0="<<dist0<<
1451 "distphi="<<distphi<<
1452 //
1453 //
1454 "counter="<<fCounter[id]<<
1455 "clcounter="<<fClusterCounter[id]<<
1456 "satur="<<fClusterSatur[id]<<
1457 "fz="<<fFitZ[id]<<
7b18d067 1458 //
1459 "LTr.="<<ltrp<<
1460 "Esd.="<<track<<
e9f38a4b 1461 "Tr.="<<param<<
c6914c83 1462 "x0="<<xyz[0]<<
1463 "x1="<<xyz[1]<<
1464 "x2="<<xyz[2]<<
1465 "px0="<<pxyz[0]<<
1466 "px1="<<pxyz[1]<<
1467 "px2="<<pxyz[2]<<
c18f4385 1468 //
1469 "lx0="<<lxyz[0]<<
1470 "lx1="<<lxyz[1]<<
1471 "lx2="<<lxyz[2]<<
1472 "lpx0="<<lpxyz[0]<<
1473 "lpx1="<<lpxyz[1]<<
1474 "lpx2="<<lpxyz[2]<<
c6914c83 1475 "\n";
1476 }
1477 }
1478}
1479
263d466a 1480void AliTPCcalibLaser::RefitLaserJW(Int_t id){
1481 //
1482 // Refit the track with different tracklet models:
1483 // 1. Per ROC using the kalman filter, different edge cuts
1484 // 2. Per ROC linear in y and z
1485 // 3. Per ROC quadratic in y and z
1486 // 4. Per track offset for each sector, linear for each sector, common quadratic
1487 // store x, y, z information for all models and the cluster to calculate the residuals
1488 //
be113f6e 1489
108953e9 1490 // The clusters which do not fulfill given criteria are skipped
1491 //
1492 // Cluters removed from fit
1493 // 1. Extended shape > kShapeCut
1494 // 2. In saturation Max > kMax
1495 // 3. Distance to edge < cutEdge
1496 //
1497 // Clusters not used for the calibration:
1498 //
1499 // 1. Extended shape > kShapeCut
1500 // 2. In saturation Max > kMax
1501
1502
263d466a 1503 AliTPCseed *track = (AliTPCseed*)fTracksTPC.At(id);
1504 AliExternalTrackParam *extparam=(AliExternalTrackParam*)fTracksEsdParam.At(id);
1505 AliTPCLaserTrack *ltrp = (AliTPCLaserTrack*)fTracksMirror.At(id);
1506
1507 AliTPCclusterMI dummyCl;
1508
1509 //two tracklets
108953e9 1510 Int_t kMaxTracklets=2;
1511 Float_t kShapeCut = 1.3;
1512 Float_t kRatioCut = 2.;
1513
1514 Float_t kMax = 900.;
1515
1516
0728a4f6 1517 //=============================================//
1518 // Linear Fitters for the Different Approaches //
1519 //=============================================//
880c3382 1520 //linear fit model in y and z; inner - outer sector, combined with offset
0728a4f6 1521 static TLinearFitter fy1I(2,"hyp1");
1522 static TLinearFitter fy1O(2,"hyp1");
1523 static TLinearFitter fz1I(2,"hyp1");
1524 static TLinearFitter fz1O(2,"hyp1");
880c3382 1525 static TLinearFitter fy1IO(3,"hyp2");
1526 static TLinearFitter fz1IO(3,"hyp2");
0728a4f6 1527 //quadratic fit model in y and z; inner - sector
1528 static TLinearFitter fy2I(3,"hyp2");
1529 static TLinearFitter fy2O(3,"hyp2");
1530 static TLinearFitter fz2I(3,"hyp2");
1531 static TLinearFitter fz2O(3,"hyp2");
263d466a 1532 //common quadratic fit for IROC and OROC in y and z
1533 static TLinearFitter fy4(5,"hyp4");
1534 static TLinearFitter fz4(5,"hyp4");
1535
1536
1537 //set standard cuts
1538 if ( fNcuts==0 ){
1539 fNcuts=1;
1540 fEdgeXcuts[0]=4;
1541 fEdgeYcuts[0]=3;
1542 fNClCuts[0]=20;
1543 }
0728a4f6 1544 //=============================//
1545 // Loop over all Tracklet Cuts //
1546 //=============================//
263d466a 1547 for (Int_t icut=0; icut<fNcuts; icut++){
108953e9 1548 Float_t xinMin = 2500, xinMax=-90;
1549 Float_t xoutMin=2500, xoutMax=-90;
1550 Float_t msigmaYIn=0;
1551 Float_t msigmaYOut=0;
1552 Float_t mqratioIn=0;
1553 Float_t mqratioOut=0;
be113f6e 1554
1555
263d466a 1556 AliDebug(4,Form("Processing cut %d for track with ID %d",icut,id));
1557 //cut parameters
1558 Double_t edgeCutX = fEdgeXcuts[icut];
1559 Double_t edgeCutY = fEdgeYcuts[icut];
e5c59218 1560 Int_t nclCut = (Int_t)fNClCuts[icut];
0728a4f6 1561 //=========================//
1562 // Parameters to calculate //
1563 //=========================//
1564 //fit parameter inner, outer tracklet and combined fit
1565 TVectorD vecy1resInner(2),vecz1resInner(2); //pol1 fit parameters inner
1566 TVectorD vecy2resInner(3),vecz2resInner(3); //pol2 fit parameters inner
263d466a 1567 //
0728a4f6 1568 TVectorD vecy1resOuter(2),vecz1resOuter(2); //pol1 fit parameters outer
1569 TVectorD vecy2resOuter(3),vecz2resOuter(3); //pol2 fit parameters outer
263d466a 1570 TVectorD vecy4res(5),vecz4res(5);
880c3382 1571 TVectorD vecy1resIO(3),vecz1resIO(3);
263d466a 1572 // cluster and track positions for each row - used for residuals
108953e9 1573 TVectorD vecgX(159); // global X
1574 TVectorD vecgY(159); // global Y
1575 TVectorD vecgZ(159); // global Z
1576
263d466a 1577 TVectorD vecX(159); // x is the same for all (row center)
1578 TVectorD vecYkalman(159); // y from kalman fit
1579 TVectorD vecZkalman(159); // z from kalman fit
1580 TVectorD vecY1(159); // y from pol1 fit per ROC
1581 TVectorD vecZ1(159); // z from pol1 fit per ROC
880c3382 1582 TVectorD vecY1IO(159); // y from pol1 fit per ROC
1583 TVectorD vecZ1IO(159); // z from pol1 fit per ROC
263d466a 1584 TVectorD vecY2(159); // y from pol2 fit per ROC
1585 TVectorD vecZ2(159); // z from pol2 fit per ROC
1586 TVectorD vecY4(159); // y from sector fit
1587 TVectorD vecZ4(159); // z from sector fit
1588 TVectorD vecClY(159); // y cluster position
1589 TVectorD vecClZ(159); // z cluster position
1590 TVectorD vecSec(159); // sector for each row
108953e9 1591 TVectorD isReject(159); // flag - cluster to be rejected
0728a4f6 1592 //chi2 of fits
be113f6e 1593 Double_t chi2I1z=0; // chi2 of pol1 fit in z (inner)
1594 Double_t chi2I1y=0; // chi2 of pol1 fit in y (inner)
1595 Double_t chi2O1z=0; // chi2 of pol1 fit in z (outer)
1596 Double_t chi2O1y=0; // chi2 of pol1 fit in y (outer)
880c3382 1597 Double_t chi2IO1z=0; // chi2 of pol1 fit in z (outer)
1598 Double_t chi2IO1y=0; // chi2 of pol1 fit in y (outer)
be113f6e 1599 Double_t chi2I2z=0; // chi2 of pol2 fit in z (inner)
1600 Double_t chi2I2y=0; // chi2 of pol2 fit in y (inner)
1601 Double_t chi2O2z=0; // chi2 of pol2 fit in z (outer)
1602 Double_t chi2O2y=0; // chi2 of pol2 fit in y (outer)
1603 Double_t chi2IOz=0; // chi2 of hyp4 fit in z (inner+outer)
1604 Double_t chi2IOy=0; // chi2 of hyp4 fit in y (inner+outer)
0728a4f6 1605 //more
1606 Int_t innerSector = -1; // number of inner sector
1607 Int_t outerSector = -1; // number of outer sector
263d466a 1608 Int_t nclI=0; // number of clusters (inner)
1609 Int_t nclO=0; // number of clusters (outer)
0728a4f6 1610 //=================================================//
1611 // Perform the Kalman Fit using the Tracklet Class //
1612 //=================================================//
263d466a 1613 AliTPCTracklet::SetEdgeCut(edgeCutX,edgeCutY);
1614 TObjArray tracklets=
1615 AliTPCTracklet::CreateTracklets(track,AliTPCTracklet::kKalman,
1616 kFALSE,nclCut,kMaxTracklets);
263d466a 1617 // tracklet pointers
263d466a 1618 AliTPCTracklet *trInner = (AliTPCTracklet*)tracklets.At(0);
1619 AliTPCTracklet *trOuter = (AliTPCTracklet*)tracklets.At(1);
0728a4f6 1620 AliTPCTracklet *tr=0x0;
1621 AliTPCTracklet dummy;
1622 //continue if we didn't find a tracklet
263d466a 1623 if ( !trInner && !trOuter ) continue;
0728a4f6 1624 //================================================================================//
1625 // Swap Inner and Outer if necessary (inner sector is definde by smaller local x) //
1626 //================================================================================//
263d466a 1627 if ( trInner && trOuter ){
1628 if ( !trInner->GetInner() || !trOuter->GetInner() ) continue;
1629 if ( trInner->GetInner()->GetX() > trOuter->GetInner()->GetX() ){
1630 tr = trInner;
1631 trInner=trOuter;
1632 trOuter=tr;
1633 AliDebug(5,Form("Swapped Sectors: %02d (%f) <-> %02d (%f)", trOuter->GetSector(), trOuter->GetInner()->GetX(), trInner->GetSector(), trInner->GetInner()->GetX()));
1634 }
1635 } else {
1636 if ( trInner ){
1637 if ( !trInner->GetInner() ) continue;
1638 trOuter=&dummy;
1639 if ( trInner->GetSector()>35 ){
1640 trOuter=trInner;
1641 trInner=&dummy;
1642 }
1643 } else { //trOuter
1644 if ( !trOuter->GetInner() ) continue;
1645 trInner=&dummy;
1646 if ( trOuter->GetSector()<36 ){
1647 trInner=trOuter;
1648 trOuter=&dummy;
1649 }
1650 }
1651 }
1652 innerSector = trInner->GetSector();
1653 if ( innerSector>=0 ) AliDebug(5,Form("Found inner Sector %02d at X %.2f", innerSector, trInner->GetInner()->GetX()));
1654 outerSector = trOuter->GetSector();
1655 if ( outerSector>=0 ) AliDebug(5,Form("Found outer Sector %02d at X %.2f", outerSector, trOuter->GetInner()->GetX()));
1656
263d466a 1657 // array of clusters
1658 TClonesArray arrCl("AliTPCclusterMI",159);
1659 arrCl.ExpandCreateFast(159);
0728a4f6 1660 //=======================================//
1661 // fill fitters with cluster information //
1662 //=======================================//
1663 AliDebug(3,"Filing Cluster Information");
108953e9 1664
1665 //
1666 //
1667
263d466a 1668 for (Int_t irow=158;irow>-1;--irow) {
1669 AliTPCclusterMI *c=track->GetClusterPointer(irow);
0728a4f6 1670 AliTPCclusterMI &cl = (AliTPCclusterMI&) (*arrCl[irow]);
1671 cl=dummyCl;
be113f6e 1672 vecX[irow] = 0;
1673 vecClY[irow] = 0;
1674 vecClZ[irow] = 0;
1675 Float_t meanY=0, sumY=0;
1676 for (Int_t drow=-1;drow<=1;drow++) {
1677 if (irow+drow<0) continue;
1678 if (irow+drow>158) continue;
57dc06f2 1679 AliTPCclusterMI *ccurrent=track->GetClusterPointer(irow);
1680 if (!ccurrent) continue;
1681 Int_t roc = static_cast<Int_t>(ccurrent->GetDetector());
be113f6e 1682 if ( roc!=innerSector && roc!=outerSector ) continue;
57dc06f2 1683 if (ccurrent->GetX()<10) continue;
1684 if (ccurrent->GetY()==0) continue;
1685 meanY+=ccurrent->GetY();
be113f6e 1686 sumY++;
1687 }
1688 if (sumY>0) meanY/=sumY;
1689
0728a4f6 1690 //
263d466a 1691 vecSec[irow]=-1;
1692 if (!c) continue;
108953e9 1693 Double_t pedgeY = c->GetX()*TMath::Tan(TMath::DegToRad()*(10))-TMath::Abs(meanY);
e5c59218 1694 Double_t pedgeX = TMath::Min((irow)*0.75, (159.-irow)*1.5);
1695
0728a4f6 1696 //
1697 Int_t roc = static_cast<Int_t>(c->GetDetector());
1698 if ( roc!=innerSector && roc!=outerSector ) continue;
1699 vecSec[irow]=roc;
263d466a 1700 //store clusters in clones array
be113f6e 1701 cl=*c;
1702 //
1703 if (c->GetMax()<4) continue; // noise cluster?
1704 if (TMath::Abs(c->GetY())<0.0001) continue; // noise cluster?
263d466a 1705 //cluster position
1706 vecX[irow] = c->GetX();
1707 vecClY[irow] = c->GetY();
1708 vecClZ[irow] = c->GetZ();
108953e9 1709 //
31aa7c5c 1710// Float_t gxyz[3];
1711// c->GetGlobalXYZ(gxyz);
1712// vecgX[irow] = gxyz[0];
1713// vecgY[irow] = gxyz[1];
1714// vecgZ[irow] = gxyz[2];
263d466a 1715 //
0728a4f6 1716 Double_t x = vecX[irow]-133.4; //reference is between IROC and OROC
1717 Double_t y = vecClY[irow];
1718 Double_t z = vecClZ[irow];
1719 //
1720 Double_t x2[2]={x,x*x}; //linear and parabolic parameters
1721 Double_t x4[4]={0,0,0,0}; //hyp4 parameters
880c3382 1722 Double_t xIO[2]={0,x}; //common linear + offset IROC-OROC
263d466a 1723 if ( roc == innerSector ) {
0728a4f6 1724 x4[0]=1; //offset inner - outer sector
1725 x4[1]=x; //slope parameter inner sector
880c3382 1726 xIO[0]=1;
263d466a 1727 } else {
0728a4f6 1728 x4[2]=x; //slope parameter outer sector
263d466a 1729 }
0728a4f6 1730 x4[3]=x*x; //common parabolic shape
e5c59218 1731 if (pedgeX < fEdgeXcuts[icut]) continue;
1732 if (pedgeY < fEdgeYcuts[icut]) continue;
be113f6e 1733 if (c->GetMax()>900) continue; // cluster in saturation
263d466a 1734 //
0728a4f6 1735 if ( roc==innerSector ){
1736 fy1I.AddPoint(x2,y);
1737 fz1I.AddPoint(x2,z);
1738 fy2I.AddPoint(x2,y);
1739 fz2I.AddPoint(x2,z);
1740 ++nclI;
108953e9 1741 if (c->GetX()<xinMin) xinMin=c->GetX();
1742 if (c->GetX()>xinMax) xinMax=c->GetX();
1743
1744 msigmaYIn +=TMath::Sqrt(c->GetSigmaY2());
1745 mqratioIn +=c->GetMax()/c->GetQ();
1746
0728a4f6 1747 }
1748 if ( roc==outerSector ){
1749 fy1O.AddPoint(x2,y);
1750 fz1O.AddPoint(x2,z);
1751 fy2O.AddPoint(x2,y);
1752 fz2O.AddPoint(x2,z);
1753 ++nclO;
108953e9 1754 if (c->GetX()<xoutMin) xoutMin=c->GetX();
1755 if (c->GetX()>xoutMax) xoutMax=c->GetX();
1756 msigmaYOut +=TMath::Sqrt(c->GetSigmaY2());
1757 mqratioOut +=c->GetMax()/c->GetQ();
1758
0728a4f6 1759 }
263d466a 1760 fy4.AddPoint(x4,y);
1761 fz4.AddPoint(x4,z);
880c3382 1762 fy1IO.AddPoint(xIO,y);
1763 fz1IO.AddPoint(xIO,z);
263d466a 1764 }
108953e9 1765 if (nclI>0) {
1766 msigmaYIn/=nclI;
1767 mqratioIn/=nclI;
1768 }
1769 if (nclO>0) {
1770 msigmaYOut/=nclO;
1771 mqratioOut/=nclO;
1772 }
0728a4f6 1773 //======================================//
1774 // Evaluate and retrieve fit parameters //
1775 //======================================//
1776 AliDebug(5,Form("Evaluating fits with inner (outer) Sec: %02d (%02d)",innerSector,outerSector));
1777 //inner sector
1778 if ( (innerSector>-1) && (fy1I.GetNpoints()>0) ){
1779 fy1I.Eval();
1780 fz1I.Eval();
1781 fy2I.Eval();
1782 fz2I.Eval();
1783 fy1I.GetParameters(vecy1resInner);
1784 fz1I.GetParameters(vecz1resInner);
1785 fy2I.GetParameters(vecy2resInner);
1786 fz2I.GetParameters(vecz2resInner);
1787 chi2I1y=fy1I.GetChisquare()/(fy1I.GetNpoints()-2);
1788 chi2I1z=fz1I.GetChisquare()/(fz1I.GetNpoints()-2);
1789 chi2I2y=fy2I.GetChisquare()/(fy2I.GetNpoints()-3);
1790 chi2I2z=fz2I.GetChisquare()/(fz2I.GetNpoints()-3);
1791 }
1792 //outer sector
1793 if ( (outerSector>-1) && (fy1O.GetNpoints()>0) ){
1794 fy1O.Eval();
1795 fz1O.Eval();
1796 fy2O.Eval();
1797 fz2O.Eval();
1798 fy1O.GetParameters(vecy1resOuter);
1799 fz1O.GetParameters(vecz1resOuter);
1800 fy2O.GetParameters(vecy2resOuter);
1801 fz2O.GetParameters(vecz2resOuter);
1802 chi2O1y=fy1O.GetChisquare()/(fy1O.GetNpoints()-2);
1803 chi2O1z=fz1O.GetChisquare()/(fz1O.GetNpoints()-2);
1804 chi2O2y=fy2O.GetChisquare()/(fy2O.GetNpoints()-3);
1805 chi2O2z=fz2O.GetChisquare()/(fz2O.GetNpoints()-3);
1806 }
1807 //combined hyp4 fit
263d466a 1808 if ( innerSector>0 && outerSector>0 ){
0728a4f6 1809 if (fy4.GetNpoints()>0) {
1810 fy4.Eval();
1811 fy4.GetParameters(vecy4res);
1812 chi2IOy=fy4.GetChisquare()/(fy4.GetNpoints()-5);
1813 }
1814 if (fz4.GetNpoints()>0) {
1815 fz4.Eval();
1816 fz4.GetParameters(vecz4res);
1817 chi2IOz=fz4.GetChisquare()/(fz4.GetNpoints()-5);
1818 }
880c3382 1819 if (fy1IO.GetNpoints()>0) {
1820 fy1IO.Eval();
1821 fy1IO.GetParameters(vecy1resIO);
1822 chi2IO1y=fy1IO.GetChisquare()/(fy1IO.GetNpoints()-3);
1823 }
1824 if (fz1IO.GetNpoints()>0) {
1825 fz1IO.Eval();
1826 fz1IO.GetParameters(vecz1resIO);
1827 chi2IO1z=fz1IO.GetChisquare()/(fz1IO.GetNpoints()-3);
1828 }
263d466a 1829 }
0728a4f6 1830 //clear points
1831 fy4.ClearPoints(); fz4.ClearPoints();
1832 fy1I.ClearPoints(); fy1O.ClearPoints();
1833 fz1I.ClearPoints(); fz1O.ClearPoints();
1834 fy2I.ClearPoints(); fy2O.ClearPoints();
1835 fz2I.ClearPoints(); fz2O.ClearPoints();
880c3382 1836 fy1IO.ClearPoints(); fz1IO.ClearPoints();
0728a4f6 1837 //==============================//
1838 // calculate tracklet positions //
1839 //==============================//
263d466a 1840 AliDebug(4,"Calculate tracklet positions");
1841 for (Int_t irow=158;irow>-1;--irow) {
108953e9 1842 isReject[irow]=0;
1843 AliTPCclusterMI *c=track->GetClusterPointer(irow);
1844 if ( vecSec[irow]!=innerSector && vecSec[irow]!=outerSector ) { // no cluster in given sectors
1845 isReject[irow]+=1;
1846 }
1847
1848 if (!c) { //no cluster
1849 isReject[irow]+=2;
1850 }else{
1851 if (c->GetMax()>kMax){ //saturation
1852 isReject[irow]+=4;
1853 }
1854 if ( vecSec[irow] == outerSector ) { //extended shape
1855 if (c->GetMax()/c->GetQ()> mqratioOut*kRatioCut) isReject[irow]+=8;
1856 if (TMath::Sqrt(c->GetSigmaY2())>msigmaYOut*kShapeCut) isReject[irow]+=16;
1857 }else{
1858 if (c->GetMax()/c->GetQ()> mqratioIn*kRatioCut) isReject[irow]+=8;
1859 if (TMath::Sqrt(c->GetSigmaY2())>msigmaYIn*kShapeCut) isReject[irow]+=16;
1860 }
1861 }
1862
1863
1864
0728a4f6 1865 if ( vecSec[irow]==-1 ) continue; //no cluster info
263d466a 1866 if ( vecSec[irow]!=innerSector && vecSec[irow]!=outerSector ) continue;
1867 tr=&dummy;
1868 Double_t x=vecX[irow];
0728a4f6 1869 Double_t xref=x-133.4;
1870 //
263d466a 1871 Double_t yoffInner=0;
1872 Double_t zoffInner=0;
880c3382 1873 Double_t yoffInner1=0;
1874 Double_t zoffInner1=0;
263d466a 1875 Double_t yslopeInner=0;
1876 Double_t yslopeOuter=0;
1877 Double_t zslopeInner=0;
0728a4f6 1878 Double_t zslopeOuter=0;
1879 //positions of hyperplane fits
263d466a 1880 if ( vecSec[irow] == outerSector ) {
1881 tr=trOuter;
1882 vecY1[irow]=vecy1resOuter[0]+vecy1resOuter[1]*xref;
1883 vecZ1[irow]=vecz1resOuter[0]+vecz1resOuter[1]*xref;
1884 vecY2[irow]=vecy2resOuter[0]+vecy2resOuter[1]*xref+vecy2resOuter[2]*xref*xref;
1885 vecZ2[irow]=vecz2resOuter[0]+vecz2resOuter[1]*xref+vecz2resOuter[2]*xref*xref;
1886 yslopeOuter=vecy4res[3];
1887 zslopeOuter=vecz4res[3];
263d466a 1888 } else {
1889 tr=trInner;
1890 vecY1[irow]=vecy1resInner[0]+vecy1resInner[1]*xref;
1891 vecZ1[irow]=vecz1resInner[0]+vecz1resInner[1]*xref;
1892 vecY2[irow]=vecy2resInner[0]+vecy2resInner[1]*xref+vecy2resInner[2]*xref*xref;
1893 vecZ2[irow]=vecz2resInner[0]+vecz2resInner[1]*xref+vecz2resInner[2]*xref*xref;
1894 yoffInner=vecy4res[1];
880c3382 1895 zoffInner=vecz4res[1];
1896 yoffInner1=vecy1resIO[1];
1897 zoffInner1=vecz1resIO[1];
263d466a 1898 yslopeInner=vecy4res[2];
1899 zslopeInner=vecz4res[2];
263d466a 1900 }
880c3382 1901 vecY1IO[irow]=vecy1resIO[0]+yoffInner1+vecy1resIO[2]*xref;
1902 vecZ1IO[irow]=vecz1resIO[0]+zoffInner1+vecz1resIO[2]*xref;
263d466a 1903 vecY4[irow]=vecy4res[0]+yoffInner+yslopeInner*xref+yslopeOuter*xref+vecy4res[4]*xref*xref;
1904 vecZ4[irow]=vecz4res[0]+zoffInner+zslopeInner*xref+zslopeOuter*xref+vecz4res[4]*xref*xref;
0728a4f6 1905 //positions of kalman fits
263d466a 1906 Double_t gxyz[3],xyz[3];
1907 AliExternalTrackParam *param = 0x0;
0728a4f6 1908 //
263d466a 1909 param=tr->GetInner();
1910 if (param){
1911 param->GetXYZ(gxyz);
1912 Float_t bz = AliTracker::GetBz(gxyz);
1913 param->GetYAt(x, bz, xyz[1]);
1914 param->GetZAt(x, bz, xyz[2]);
1915 vecYkalman[irow]=xyz[1];
1916 vecZkalman[irow]=xyz[2];
1917 }
108953e9 1918 //
1919 //
1920 //
1921
263d466a 1922 }
0728a4f6 1923 //=====================================================================//
1924 // write results from the different tracklet fits with debug streamers //
1925 //=====================================================================//
263d466a 1926 if (fStreamLevel>4){
1927 TTreeSRedirector *cstream = GetDebugStreamer();
1928 if (cstream){
be113f6e 1929 Float_t dedx = track->GetdEdx();
1930 (*cstream)<<"FitModels"<<
108953e9 1931 "run="<<fRun<< // run number
1932 "event="<<fEvent<< // event number
1933 "time="<<fTime<< // time stamp of event
1934 "trigger="<<fTrigger<< // trigger
1935 "mag="<<fMagF<< // magnetic field
1936 //
be113f6e 1937 "cutNr=" << icut <<
1938 "edgeCutX=" << edgeCutX <<
1939 "edgeCutY=" << edgeCutY <<
1940 "nclCut=" << nclCut <<
1941 "innerSector="<< innerSector <<
1942 "outerSector="<< outerSector <<
1943 "dEdx=" << dedx <<
1944 "LTr.=" << ltrp <<
1945 "Tr.=" << extparam <<
880c3382 1946 "yPol1In.=" << &vecy1resInner <<
1947 "zPol1In.=" << &vecz1resInner <<
1948 "yPol1InOut.="<< &vecy1resIO <<
1949 "zPol1InOut.="<< &vecz1resIO <<
1950 "yPol2In.=" << &vecy2resInner <<
be113f6e 1951 "zPol2In.=" << &vecz2resInner <<
1952 "yPol1Out.=" << &vecy1resOuter <<
1953 "zPol1Out.=" << &vecz1resOuter <<
1954 "yPol2Out.=" << &vecy2resOuter <<
1955 "zPol2Out.=" << &vecz2resOuter <<
1956 "yInOut.=" << &vecy4res <<
1957 "zInOut.=" << &vecz4res <<
880c3382 1958 "chi2y1In=" << chi2I1y <<
1959 "chi2z1In=" << chi2I1z <<
1960 "chi2y1InOut="<< chi2IO1y <<
1961 "chi2z1InOut="<< chi2IO1z <<
1962 "chi2y1Out=" << chi2O1y <<
be113f6e 1963 "chi2z1Out=" << chi2O1z <<
1964 "chi2y2In=" << chi2I2y <<
1965 "chi2z2In=" << chi2I2z <<
1966 "chi2y2Out=" << chi2O2y <<
1967 "chi2z2Out=" << chi2O2z <<
1968 "chi2yInOut=" << chi2IOy <<
1969 "chi2zInOut=" << chi2IOz <<
1970 "trletIn.=" << trInner <<
1971 "trletOut.=" << trOuter <<
1972 "nclI=" << nclI <<
1973 "nclO=" << nclO <<
1974 "xinMin=" << xinMin<<
1975 "xinMax=" << xinMax<<
1976 "xoutMin=" << xoutMin<<
1977 "xoutMax=" << xoutMax<<
108953e9 1978 "msigmaYIn=" <<msigmaYIn<<
1979 "msigmaYOut=" <<msigmaYOut<<
1980 "mqratioIn=" <<mqratioIn<<
1981 "mqratioOut=" << mqratioOut <<
be113f6e 1982 "\n";
263d466a 1983 }
1984 }
1985
1986 // wirte residuals information
1987 if (fStreamLevel>5){
1988 TTreeSRedirector *cstream = GetDebugStreamer();
1989 if (cstream){
be113f6e 1990 Float_t dedx = track->GetdEdx();
1991 (*cstream)<<"Residuals"<<
108953e9 1992 "run="<<fRun<< // run number
1993 "event="<<fEvent<< // event number
1994 "time="<<fTime<< // time stamp of event
1995 "trigger="<<fTrigger<< // trigger
1996 "mag="<<fMagF<< // magnetic field
1997 //
be113f6e 1998 "cutNr=" << icut <<
1999 "edgeCutX=" << edgeCutX <<
2000 "edgeCutY=" << edgeCutY <<
2001 "nclCut=" << nclCut <<
2002 "LTr.=" << ltrp <<
2003 "Tr.=" << extparam<<
2004 "dEdx=" << dedx <<
2005 "Cl.=" << &arrCl <<
108953e9 2006 "vX.=" << &vecgX<< // global x
2007 "vY.=" << &vecgY<< // global y
2008 "vZ.=" << &vecgZ<< // global z
be113f6e 2009 "TrX.=" << &vecX <<
2010 "TrYpol1.=" << &vecY1 <<
2011 "TrZpol1.=" << &vecZ1 <<
2012 "TrYpol2.=" << &vecY2 <<
2013 "TrZpol2.=" << &vecZ2 <<
880c3382 2014 "TrYpol1InOut.="<< &vecY1IO <<
2015 "TrZpol1InOut.="<< &vecZ1IO <<
2016 "TrYInOut.=" << &vecY4 <<
2017 "TrZInOut.=" << &vecZ4 <<
2018 "ClY.=" << &vecClY <<
be113f6e 2019 "ClZ.=" << &vecClZ <<
108953e9 2020 "isReject.=" << &isReject<<
be113f6e 2021 "sec.=" << &vecSec <<
2022 "nclI=" << nclI <<
2023 "nclO=" << nclO <<
2024 "xinMin=" << xinMin<<
2025 "xinMax=" << xinMax<<
2026 "xoutMin=" << xoutMin<<
2027 "xoutMax=" << xoutMax<<
108953e9 2028 "msigmaYIn=" <<msigmaYIn<<
2029 "msigmaYOut=" <<msigmaYOut<<
2030 "mqratioIn=" <<mqratioIn<<
2031 "mqratioOut=" << mqratioOut <<
be113f6e 2032 "yInOut.=" << &vecy4res <<
2033 "zInOut.=" << &vecz4res <<
108953e9 2034 //chi2s
be113f6e 2035 "chi2y1In=" << chi2I1y << //
2036 "chi2z1In=" << chi2I1z <<
2037 "chi2y1Out=" << chi2O1y <<
2038 "chi2z1Out=" << chi2O1z <<
880c3382 2039 "chi2y1InOut="<< chi2IO1y <<
2040 "chi2z1InOut="<< chi2IO1z <<
2041 "chi2y2In=" << chi2I2y <<
be113f6e 2042 "chi2z2In=" << chi2I2z <<
2043 "chi2y2Out=" << chi2O2y <<
2044 "chi2z2Out=" << chi2O2z <<
2045 "chi2yInOut=" << chi2IOy <<
2046 "chi2zInOut=" << chi2IOz <<
108953e9 2047 // fit parameters
2048 "yPol1In.=" << &vecy1resInner <<
2049 "zPol1In.=" << &vecz1resInner <<
2050 "yPol2In.=" << &vecy2resInner <<
2051 "zPol2In.=" << &vecz2resInner <<
2052 "yPol1Out.=" << &vecy1resOuter <<
2053 "zPol1Out.=" << &vecz1resOuter <<
880c3382 2054 "yPol1InOut.="<< &vecy1resIO <<
2055 "zPol1InOut.="<< &vecz1resIO <<
2056 "yPol2Out.=" << &vecy2resOuter <<
108953e9 2057 "zPol2Out.=" << &vecz2resOuter <<
2058
be113f6e 2059 "\n";
263d466a 2060
2061 }
2062 }
95a0e09b 2063 //==========================//
2064 // Fill Residual Histograms //
2065 //==========================//
108953e9 2066 if (!fHisNclIn) MakeFitHistos();
2067
2068 TH2F *profy = (TH2F*)fDeltaYres.UncheckedAt(id);
2069 TH2F *profz = (TH2F*)fDeltaZres.UncheckedAt(id);
2070 TH2F *profy2 = (TH2F*)fDeltaYres2.UncheckedAt(id);
2071 TH2F *profz2 = (TH2F*)fDeltaZres2.UncheckedAt(id);
7af539c6 2072 TH2F *profyabs = (TH2F*)fDeltaYresAbs.UncheckedAt(id);
2073 TH2F *profzabs = (TH2F*)fDeltaZresAbs.UncheckedAt(id);
108953e9 2074 // TH2F *profy3 = (TH2F*)fDeltaYres3.UncheckedAt(id);
2075 //TH2F *profz3 = (TH2F*)fDeltaZres3.UncheckedAt(id);
2076 //
95a0e09b 2077 for (Int_t irow=158;irow>-1;--irow) {
108953e9 2078 if (vecSec[irow]==-1)continue; // no cluster info
2079 if (isReject[irow]>0.5) continue; //
95a0e09b 2080 Double_t ycl = vecClY[irow];
2081 Double_t yfit = vecY1[irow];
108953e9 2082 Double_t yfit2 = vecY2[irow];
7af539c6 2083 Double_t x = vecX[irow];
2084 Double_t yabsbeam = -1000;
f6220766 2085 if(vecSec[irow]==outerSector && outerSector==fBeamSectorOuter[id])
7af539c6 2086 yabsbeam = fBeamSlopeYOuter[id]*x + fBeamOffsetYOuter[id];
f6220766 2087 else if(innerSector==fBeamSectorInner[id])
7af539c6 2088 yabsbeam = fBeamSlopeYInner[id]*x + fBeamOffsetYInner[id];
2089
43ad9d4c 2090 // Double_t yfit3 = vecY2[irow];
95a0e09b 2091 Double_t zcl = vecClZ[irow];
2092 Double_t zfit = vecZ1[irow];
108953e9 2093 Double_t zfit2 = vecZ2[irow];
43ad9d4c 2094 //Double_t zfit3 = vecZ2[irow];
108953e9 2095
7af539c6 2096 // dz abs
2097 // The expressions for zcorrected has been obtained by
2098 // inverting the fits in the FitDriftV() method (ignoring the
2099 // global y dependence for now):
2100 // A side:
2101 // 250 - zmeasured = [0] + [1]*(250 - zreal) + .... (yglobal)
2102 // =>
2103 // zreal = (zmeasured + [0] - (1-[1])*250.0)/[1]
2104 //
2105 // C side:
2106 // 250 + zmeasured = [0] + [1]*(250+zreal) + .... (yglobal)
2107 // =>
2108 // zreal = (zmeasured - [0] + (1 - [1])*250.0)/[1]
2109
2110 Double_t dzabs = -1000;
41a1f702 2111 Double_t zcorrected = -1000;
7af539c6 2112 if (ltrp->GetSide()==0){
41a1f702 2113 if ((*fFitAside)[1]>0. || fUseFixedDriftV) {
7af539c6 2114 // ignore global y dependence for now
3e55050f 2115 zcorrected = 0;
41a1f702 2116 if(!fUseFixedDriftV)
2117 zcorrected = (zcl + (*fFitAside)[0] -
2118 (1.0-(*fFitAside)[1])*250.0)/(*fFitAside)[1];
2119 else
2120 zcorrected = (zcl + fFixedFitAside0 -
2121 (1.0-fFixedFitAside1)*250.0)/fFixedFitAside1;
7af539c6 2122 // zcorrected = zcl;
41a1f702 2123 if(vecSec[irow]==outerSector && outerSector==fBeamSectorOuter[id])
7af539c6 2124 dzabs = zcorrected -fBeamSlopeZOuter[id]*x -fBeamOffsetZOuter[id];
41a1f702 2125 else if(innerSector==fBeamSectorInner[id])
7af539c6 2126 dzabs = zcorrected -fBeamSlopeZInner[id]*x -fBeamOffsetZInner[id];
2127 }
2128 } else {
41a1f702 2129 if ((*fFitCside)[1]>0. || fUseFixedDriftV) {
2130
2131 if(!fUseFixedDriftV)
2132 zcorrected = (zcl - (*fFitCside)[0] +
2133 (1.0-(*fFitCside)[1])*250.0)/(*fFitCside)[1];
2134 else
2135 zcorrected = (zcl - fFixedFitCside0 +
2136 (1.0-fFixedFitCside1)*250.0)/fFixedFitCside1;
2137
7af539c6 2138 // zcorrected = zcl;
41a1f702 2139 if(vecSec[irow]==outerSector && outerSector==fBeamSectorOuter[id])
7af539c6 2140 dzabs = zcorrected -fBeamSlopeZOuter[id]*x -fBeamOffsetZOuter[id];
41a1f702 2141 else if(innerSector==fBeamSectorInner[id])
7af539c6 2142 dzabs = zcorrected -fBeamSlopeZInner[id]*x -fBeamOffsetZInner[id];
2143 }
2144 }
2145
be113f6e 2146 if (TMath::Abs(yfit-ycl)<2&&TMath::Abs(zfit-zcl)<2){
108953e9 2147 if (profy){
2148 profy->Fill(irow,ycl-yfit);
2149 profy2->Fill(irow,ycl-yfit2);
41a1f702 2150 if(yabsbeam<-100) {
2151 fHisYAbsErrors->Fill(id);
2152 // profyabs->Fill(irow,-0.99);
2153 } else
7af539c6 2154 profyabs->Fill(irow,ycl-yabsbeam);
2155
108953e9 2156 // profy3->Fill(irow,ycl-yfit3);
2157 }
2158 if (profz) {
2159 profz->Fill(irow,zcl-zfit);
2160 profz2->Fill(irow,zcl-zfit2);
2161 //profz3->Fill(irow,zcl-zfit3);
41a1f702 2162 if(dzabs<-100) {
2163
2164 fHisZAbsErrors->Fill(id);
2165 }else
7af539c6 2166 profzabs->Fill(irow,dzabs);
108953e9 2167 }
be113f6e 2168 }
95a0e09b 2169 }
be113f6e 2170 //
2171 //
2172 // Fill laser fit histograms
2173 //
be113f6e 2174 Float_t dedx = track->GetdEdx();
2175 if (nclI>20&&nclO>20){
2176 fHisNclIn->Fill(id,nclI); //->Number of clusters inner
2177 fHisNclOut->Fill(id,nclO); //->Number of clusters outer
2178 fHisNclIO->Fill(id,nclI+nclO); //->Number of cluster inner outer
2179 //
2180 fHisLclIn->Fill(id,xinMax-xinMin); //->Level arm inner
2181 fHisLclOut->Fill(id,xoutMax-xoutMin); //->Level arm outer
2182 fHisLclIO->Fill(id,xoutMax-xinMin); //->Number of cluster inner outer
2183 //
2184 fHisdEdx->Fill(id,TMath::Sqrt(TMath::Abs(dedx)));
2185 fHisdZfit->Fill(id,fFitZ[id]);
2186 //
2187 //
2188 fHisChi2YIn1->Fill(id,TMath::Sqrt(chi2I1y)); //->chi2 y inner - line
2189 fHisChi2YOut1->Fill(id,TMath::Sqrt(chi2O1y)); //->chi2 y inner - line
2190 fHisChi2YIn2->Fill(id,TMath::Sqrt(chi2I2y)); //->chi2 y inner - parabola
2191 fHisChi2YOut2->Fill(id,TMath::Sqrt(chi2O2y)); //->chi2 y inner - parabola
2192 fHisChi2YIO1->Fill(id,TMath::Sqrt(chi2IOy)); //->chi2 y IO - common
2193
2194
2195 fHisChi2ZIn1->Fill(id,TMath::Sqrt(chi2I1z)); //->chi2 z inner - line
2196 fHisChi2ZOut1->Fill(id,TMath::Sqrt(chi2O1z)); //->chi2 z inner - line
2197 fHisChi2ZIn2->Fill(id,TMath::Sqrt(chi2O2y)); //->chi2 z inner - parabola
2198 fHisChi2ZOut2->Fill(id,TMath::Sqrt(chi2O2z)); //->chi2 z inner - parabola
2199 fHisChi2ZIO1->Fill(id,TMath::Sqrt(chi2IOz)); //->chi2 z IO - common
2200 //
2201 //
2202 fHisPy1vP0->Fill(id,vecy1resOuter[0]-vecy1resInner[0]); //-> delta y P0outer-P0inner - line
2203 fHisPy2vP0->Fill(id,vecy2resOuter[0]-vecy2resInner[0]); //-> delta y P0outer-P0inner - parabola
2204 fHisPy3vP0->Fill(id,vecy4res[1]); //-> delta y P0outer-P0inner - common parabola
2205 //
2206 fHisPy1vP1->Fill(id,vecy1resOuter[1]-vecy1resInner[1]); //-> delta ky P1outer-P1inner - line
2207 fHisPy2vP1->Fill(id,vecy2resOuter[1]-vecy2resInner[1]); //-> delta ky P1outer-P1inner - parabola
2208 fHisPy3vP1->Fill(id,vecy4res[3]-vecy4res[2]); //-> delta ky P1outer-P1inner - common parabola
e9f38a4b 2209 //
be113f6e 2210 fHisPy3vP2IO->Fill(id,vecy4res[4]); //-> Curv P2outerinner - common parabola
2211 fHisPz1vP0->Fill(id,vecz1resOuter[0]-vecz1resInner[0]); //-> delta z P0outer-P0inner - line
2212 fHisPz2vP0->Fill(id,vecz2resOuter[0]-vecz2resInner[0]); //-> delta z P0outer-P0inner - parabola
2213 fHisPz3vP0->Fill(id,vecz4res[1]); //-> delta z P0outer-P0inner - common parabola
2214 //
2215 fHisPz1vP1->Fill(id,vecz1resOuter[1]-vecz1resInner[1]); //-> delta kz P1outer-P1inner - line
2216 fHisPz2vP1->Fill(id,vecz2resOuter[1]-vecz2resInner[1]); //-> delta kz P1outer-P1inner - parabola
2217 fHisPz3vP1->Fill(id,vecz4res[3]-vecz4res[2]); //-> delta kz P1outer-P1inner - common parabola
2218 fHisPz3vP2IO->Fill(id,vecz4res[4]); //-> Curv P2outerinner - common parabola
e9f38a4b 2219 }
be113f6e 2220 if(nclI>20){
2221 fHisPy2vP2In->Fill(id,vecy2resInner[2]); //-> Curv P2inner - parabola
2222 fHisPz2vP2In->Fill(id,vecz2resInner[2]); //-> Curv P2inner - parabola
e9f38a4b 2223 }
be113f6e 2224 //
2225 if (nclO>20){
2226 fHisPz2vP2Out->Fill(id,vecz2resOuter[2]); //-> Curv P2outer - parabola
2227 fHisPy2vP2Out->Fill(id,vecy2resOuter[2]); //-> Curv P2outer - parabola
2228 }
2229
e9f38a4b 2230 }
e9f38a4b 2231}
2232
2233
be113f6e 2234
57dc06f2 2235void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield, Int_t run){
2b35e8e6 2236 //
1fd56785 2237 // Dump information about laser beams
2238 // isOK variable indicates usability of the beam
2239 // Beam is not usable if:
2240 // a. No entries in range (krmsCut0)
2241 // b. Big sperad (krmscut1)
2242 // c. RMSto fit sigma bigger then (kmultiCut)
2243 // d. Too big angular spread
2244 //
2245
2246 const Float_t krmsCut0=0.001;
2247 const Float_t krmsCut1=0.16;
2248 const Float_t kmultiCut=2;
2249 const Float_t kcutP0=0.002;
2b35e8e6 2250 //
2251 AliTPCcalibLaser *laser = this;
2252 TTreeSRedirector *pcstream = new TTreeSRedirector("laserMean.root");
2253 TF1 fg("fg","gaus");
f6220766 2254
2255 // start set up for absolute residuals analysis
2256 //AliTPCcalibDB* calib=AliTPCcalibDB::Instance();
2257 // AliTPCParam * tpcparam = calib->GetParameters();
2258 AliTPCParam * tpcparam = 0;
2259 if (!tpcparam) tpcparam = new AliTPCParamSR;
2260 tpcparam->Update();
2261 SetBeamParameters(fBeamOffsetYOuter, fBeamSlopeYOuter, fBeamSectorOuter,0);
2262 SetBeamParameters(fBeamOffsetYInner, fBeamSlopeYInner, fBeamSectorInner,1);
2263 TLinearFitter lfabsyInner(2);
2264 lfabsyInner.SetFormula("1 ++ x");
2265 TLinearFitter lfabszInner(2);
2266 lfabszInner.SetFormula("1 ++ x");
2267
2268 TLinearFitter lfabsyOuter(2);
2269 lfabsyOuter.SetFormula("1 ++ x");
2270 TLinearFitter lfabszOuter(2);
2271 lfabszOuter.SetFormula("1 ++ x");
2272 // end set up for absolute residuals analysis
2273
2b35e8e6 2274 //
2275 //
2276 for (Int_t id=0; id<336; id++){
1fd56785 2277 Bool_t isOK=kTRUE;
2b35e8e6 2278 TH1F * hisphi = (TH1F*)laser->fDeltaPhi.At(id);
2279 TH1F * hisphiP = (TH1F*)laser->fDeltaPhiP.At(id);
2280 TH1F * hisZ = (TH1F*)laser->fDeltaZ.At(id);
e5c59218 2281 TH1F * hisP3 = (TH1F*)laser->fDeltaP3.At(id);
2282 TH1F * hisP4 = (TH1F*)laser->fDeltaP4.At(id);
1fd56785 2283 TH1F * hisS = (TH1F*)laser->fSignals.At(id);
be113f6e 2284 //if (!hisphi) continue;
108953e9 2285 Double_t entries = (hisphi==0)? 0: hisphi->GetEntries();
be113f6e 2286 //if (entries<minEntries) continue;
2b35e8e6 2287 //
2288 AliTPCLaserTrack *ltrp = (AliTPCLaserTrack*)fTracksMirror.At(id);
2289 if (!ltrp) {
2290 AliTPCLaserTrack::LoadTracks();
2291 ltrp =(AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(id);
2292 }
e5c59218 2293 pcstream->GetFile()->cd();
2294 if (hisphi) hisphi->Write();
2295 if (hisphiP) hisphiP->Write();
2296 if (hisZ) hisZ->Write();
2297 if (hisP3) hisP3->Write();
2298 if (hisP4) hisP4->Write();
2299
2b35e8e6 2300 Float_t meanphi = hisphi->GetMean();
2301 Float_t rmsphi = hisphi->GetRMS();
1fd56785 2302 //
2b35e8e6 2303 Float_t meanphiP = hisphiP->GetMean();
2304 Float_t rmsphiP = hisphiP->GetRMS();
2305 Float_t meanZ = hisZ->GetMean();
2306 Float_t rmsZ = hisZ->GetRMS();
d3ce44cb 2307 if (hisphi->GetRMS()>hisphi->GetBinWidth(1))
592a0c8f 2308 hisphi->Fit(&fg,"","",hisphi->GetMean()-4*hisphi->GetRMS(),hisphi->GetMean()+4*hisphi->GetRMS());
263d466a 2309 Double_t gphi1 = fg.GetParameter(1);
2310 Double_t gphi2 = fg.GetParameter(2);
592a0c8f 2311 if (hisphiP->GetRMS()>0)
2312 hisphiP->Fit(&fg,"","",hisphiP->GetMean()-4*hisphiP->GetRMS(),hisphiP->GetMean()+4*hisphiP->GetRMS());
263d466a 2313 Double_t gphiP1 = fg.GetParameter(1);
2314 Double_t gphiP2 = fg.GetParameter(2);
e5c59218 2315 //
d3ce44cb 2316 if (hisZ->GetRMS()>hisZ->GetBinWidth(1))
2317 hisZ->Fit(&fg,"","");
263d466a 2318 Double_t gz1 = fg.GetParameter(1);
2319 Double_t gz2 = fg.GetParameter(2);
1fd56785 2320 //
d3ce44cb 2321 if (hisP3->GetRMS()>hisP3->GetBinWidth(1))
592a0c8f 2322 hisP3->Fit(&fg,"","",hisP3->GetMean()-4*hisP3->GetRMS(),hisP3->GetMean()+4*hisP3->GetRMS());
e5c59218 2323 Double_t gp31 = fg.GetParameter(1);
2324 Double_t gp32 = fg.GetParameter(2);
2325 //
d3ce44cb 2326 if (hisP4->GetRMS()>hisP4->GetBinWidth(1))
592a0c8f 2327 hisP4->Fit(&fg,"","",hisP4->GetMean()-4*hisP4->GetRMS(),hisP4->GetMean()+4*hisP4->GetRMS());
e5c59218 2328 Double_t gp41 = fg.GetParameter(1);
2329 Double_t gp42 = fg.GetParameter(2);
2330 //
1fd56785 2331 Float_t meanS=hisS->GetMean();
2b35e8e6 2332 //
2333 Double_t lxyz[3];
2334 Double_t lpxyz[3];
2335 ltrp->GetXYZ(lxyz);
2336 ltrp->GetPxPyPz(lpxyz);
1fd56785 2337
2338 if (rmsphi<krmsCut0) isOK=kFALSE; // empty in range - not entries inside
2339 if (rmsphi>krmsCut1) isOK=kFALSE; // empty in range - not entries inside
2340 if (rmsphi>krmsCut0) if (gphi2/rmsphi>kmultiCut) isOK=kFALSE; // multi peak structure
2341 if (gphiP2>kcutP0) isOK=kFALSE;
be113f6e 2342 //
2343 //
2344 //
2345 TH1 *his =0;
2346 //
2347 his = fHisNclIn->ProjectionY("aaa",id+1,id+1);
2348 Float_t mnclIn = his->GetMean();
2349 delete his;
2350 his = fHisNclOut->ProjectionY("aaa",id+1,id+1);
2351 Float_t mnclOut = his->GetMean();
2352 delete his;
2353 his = fHisNclIO->ProjectionY("aaa",id+1,id+1);
2354 Float_t mnclIO = his->GetMean();
2355 delete his;
2356 his = fHisLclIn->ProjectionY("aaa",id+1,id+1);
2357 Float_t mLclIn = his->GetMean();
2358 delete his;
2359 his = fHisLclOut->ProjectionY("aaa",id+1,id+1);
2360 Float_t mLclOut = his->GetMean();
2361 delete his;
2362 his = fHisLclIO->ProjectionY("aaa",id+1,id+1);
2363 Float_t mLclIO = his->GetMean();
2364 delete his;
2365 //
2366 his = fHisdEdx->ProjectionY("aaa",id+1,id+1);
2367 Float_t mdEdx = his->GetMean();
2368 delete his;
2369 //
be113f6e 2370 //
2371 //
2372 //
2373 his = fHisChi2YIn1->ProjectionY("aaa",id+1,id+1); //->chi2 y inner - line
2374 Float_t mChi2YIn1= his->GetMean();
2375 delete his;
2376 his = fHisChi2YOut1->ProjectionY("aaa",id+1,id+1); //->chi2 y inner - line
2377 Float_t mChi2YOut1 = his->GetMean();
2378 delete his;
2379 his = fHisChi2YIn2->ProjectionY("aaa",id+1,id+1); //->chi2 y inner - parabola
2380 Float_t mChi2YIn2 = his->GetMean();
2381 delete his;
2382 his = fHisChi2YOut2->ProjectionY("aaa",id+1,id+1); //->chi2 y inner - parabola
2383 Float_t mChi2YOut2 = his->GetMean();
2384 delete his;
2385 his = fHisChi2YIO1->ProjectionY("aaa",id+1,id+1); //->chi2 y IO - common
2386 Float_t mChi2YIO1 = his->GetMean();
2387 delete his;
2388 his = fHisChi2ZIn1->ProjectionY("aaa",id+1,id+1); //->chi2 z inner - line
2389 Float_t mChi2ZIn1 = his->GetMean();
2390 delete his;
2391 his = fHisChi2ZOut1->ProjectionY("aaa",id+1,id+1); //->chi2 z inner - line
2392 Float_t mChi2ZOut1 = his->GetMean();
2393 delete his;
2394 his = fHisChi2ZIn2->ProjectionY("aaa",id+1,id+1); //->chi2 z inner - parabola
2395 Float_t mChi2ZIn2 = his->GetMean();
2396 delete his;
2397 his = fHisChi2ZOut2->ProjectionY("aaa",id+1,id+1); //->chi2 z inner - parabola
2398 Float_t mChi2ZOut2 = his->GetMean();
2399 delete his;
2400 his = fHisChi2ZIO1->ProjectionY("aaa",id+1,id+1); //->chi2 z IO - common
2401 Float_t mChi2ZIO1 = his->GetMean();
2402 delete his;
2403 //
2404 // fit res. histos
2405 //
108953e9 2406 his = fHisdZfit->ProjectionY("aaa",id+1,id+1);
2407 Float_t edZfit = his->GetEntries();
2408 Float_t mdZfit = his->GetMean();
2409 Float_t rdZfit = his->GetRMS();
2410 delete his;
2411
be113f6e 2412 his = fHisPy1vP0->ProjectionY("aaa",id+1,id+1); //-> delta y P0outer-P0inner - line
2413 Float_t ePy1vP0 = his->GetEntries();
2414 Float_t mPy1vP0 = his->GetMean();
2415 Float_t rPy1vP0 = his->GetRMS();
2416 delete his;
2417
2418 his = fHisPy2vP0->ProjectionY("aaa",id+1,id+1); //-> delta y P0outer-P0inner - parabola
2419 Float_t ePy2vP0 = his->GetEntries();
2420 Float_t mPy2vP0 = his->GetMean();
2421 Float_t rPy2vP0 = his->GetRMS();
2422 delete his;
2423
2424 his = fHisPy3vP0->ProjectionY("aaa",id+1,id+1); //-> delta y P0outer-P0inner - common parabola
2425 Float_t ePy3vP0 = his->GetEntries();
2426 Float_t mPy3vP0 = his->GetMean();
2427 Float_t rPy3vP0 = his->GetRMS();
2428 delete his;
2429
2430 his = fHisPy1vP1->ProjectionY("aaa",id+1,id+1); //-> delta ky P1outer-P1inner - line
2431 Float_t ePy1vP1 = his->GetEntries();
2432 Float_t mPy1vP1 = his->GetMean();
2433 Float_t rPy1vP1 = his->GetRMS();
2434 delete his;
2435
2436 his = fHisPy2vP1->ProjectionY("aaa",id+1,id+1); //-> delta ky P1outer-P1inner - parabola
2437 Float_t ePy2vP1 = his->GetEntries();
2438 Float_t mPy2vP1 = his->GetMean();
2439 Float_t rPy2vP1 = his->GetRMS();
2440 delete his;
2441
2442 his = fHisPy3vP1->ProjectionY("aaa",id+1,id+1); //-> delta ky P1outer-P1inner - common parabola
2443 Float_t ePy3vP1 = his->GetEntries();
2444 Float_t mPy3vP1 = his->GetMean();
2445 Float_t rPy3vP1 = his->GetRMS();
2446 delete his;
2447
2448 his = fHisPy2vP2In->ProjectionY("aaa",id+1,id+1); //-> Curv P2inner - parabola
2449 Float_t ePy2vP2In = his->GetEntries();
2450 Float_t mPy2vP2In = his->GetMean();
2451 Float_t rPy2vP2In = his->GetRMS();
2452 delete his;
2453
2454 his = fHisPy2vP2Out->ProjectionY("aaa",id+1,id+1); //-> Curv P2outer - parabola
2455 Float_t ePy2vP2Out = his->GetEntries();
2456 Float_t mPy2vP2Out = his->GetMean();
2457 Float_t rPy2vP2Out = his->GetRMS();
2458 delete his;
2459
2460 his = fHisPy3vP2IO->ProjectionY("aaa",id+1,id+1); //-> Curv P2outerinner - common parabola
2461 Float_t ePy3vP2IO = his->GetEntries();
2462 Float_t mPy3vP2IO = his->GetMean();
2463 Float_t rPy3vP2IO = his->GetRMS();
2464 delete his;
2465
2466 //
2467 //
2468 his = fHisPz1vP0->ProjectionY("aaa",id+1,id+1); //-> delta z P0outer-P0inner - line
2469 Float_t ePz1vP0 = his->GetEntries();
2470 Float_t mPz1vP0 = his->GetMean();
2471 Float_t rPz1vP0 = his->GetRMS();
2472 delete his;
2473
2474 his = fHisPz2vP0->ProjectionY("aaa",id+1,id+1); //-> delta z P0outer-P0inner - parabola
2475 Float_t ePz2vP0 = his->GetEntries();
2476 Float_t mPz2vP0 = his->GetMean();
2477 Float_t rPz2vP0 = his->GetRMS();
2478 delete his;
2479
2480 his = fHisPz3vP0->ProjectionY("aaa",id+1,id+1); //-> delta z P0outer-P0inner - common parabola
2481 Float_t ePz3vP0 = his->GetEntries();
2482 Float_t mPz3vP0 = his->GetMean();
2483 Float_t rPz3vP0 = his->GetRMS();
2484 delete his;
2485
2486 his = fHisPz1vP1->ProjectionY("aaa",id+1,id+1); //-> delta kz P1outer-P1inner - line
2487 Float_t ePz1vP1 = his->GetEntries();
2488 Float_t mPz1vP1 = his->GetMean();
2489 Float_t rPz1vP1 = his->GetRMS();
2490 delete his;
2491
2492 his = fHisPz2vP1->ProjectionY("aaa",id+1,id+1); //-> delta kz P1outer-P1inner - parabola
2493 Float_t ePz2vP1 = his->GetEntries();
2494 Float_t mPz2vP1 = his->GetMean();
2495 Float_t rPz2vP1 = his->GetRMS();
2496 delete his;
2497
2498 his = fHisPz3vP1->ProjectionY("aaa",id+1,id+1); //-> delta kz P1outer-P1inner - common parabola
2499 Float_t ePz3vP1 = his->GetEntries();
2500 Float_t mPz3vP1 = his->GetMean();
2501 Float_t rPz3vP1 = his->GetRMS();
2502 delete his;
2503
2504 his = fHisPz2vP2In->ProjectionY("aaa",id+1,id+1); //-> Curv P2inner - parabola
2505 Float_t ePz2vP2In = his->GetEntries();
2506 Float_t mPz2vP2In = his->GetMean();
2507 Float_t rPz2vP2In = his->GetRMS();
2508 delete his;
2509
2510 his = fHisPz2vP2Out->ProjectionY("aaa",id+1,id+1); //-> Curv P2outer - parabola
2511 Float_t ePz2vP2Out = his->GetEntries();
2512 Float_t mPz2vP2Out = his->GetMean();
2513 Float_t rPz2vP2Out = his->GetRMS();
2514 delete his;
2515
2516 his = fHisPz3vP2IO->ProjectionY("aaa",id+1,id+1); //-> Curv P2outerinner - common parabola
2517 Float_t ePz3vP2IO = his->GetEntries();
2518 Float_t mPz3vP2IO = his->GetMean();
2519 Float_t rPz3vP2IO = his->GetRMS();
2520 delete his;
f6220766 2521
2522 // Fit absolute laser residuals
2523 TH2F* histAbsY = (TH2F*)laser->fDeltaYresAbs[id];
2524 TH2F* histAbsZ = (TH2F*)laser->fDeltaZresAbs[id];
2525
2526 Int_t secInner = TMath::Nint(fBeamSectorInner[id]);
2527 Int_t secOuter = TMath::Nint(fBeamSectorOuter[id]);
2528
2529 TVectorD vecX(159); // X
2530 TVectorD vecY(159); // Y
2531 TVectorD vecR(159); // R
2532 TVectorD vecDY(159); // absolute residuals in Y
2533 TVectorD vecDZ(159); // absolute residuals in Z
2534 TVectorD vecN(159); // number of clusters
2535 TVectorD vecEy(159); //error y
2536 TVectorD vecEz(159); //error z
2537 TVectorD vecPhi(159); // local tangent
2538 TVectorD vecPhiR(159); // local tangent
2539
2540 lfabsyInner.ClearPoints();
2541 lfabszInner.ClearPoints();
2542 lfabsyOuter.ClearPoints();
2543 lfabszOuter.ClearPoints();
2544 // dummy fit values
2545 Int_t nClInY = 0;
2546 Float_t yAbsInOffset = -100;
2547 Float_t yAbsInSlope = -100;
2548 Float_t yAbsInDeltay = -100;
2549 Int_t nClInZ = 0;
2550 Float_t zAbsInOffset = -100;
2551 Float_t zAbsInSlope = -100;
2552 Float_t zAbsInDeltaz = -100;
2553 Int_t nClOutY = 0;
2554 Float_t yAbsOutOffset = -100;
2555 Float_t yAbsOutSlope = -100;
2556 Float_t yAbsOutDeltay = -100;
2557 Int_t nClOutZ = 0;
2558 Float_t zAbsOutOffset = -100;
2559 Float_t zAbsOutSlope = -100;
2560 Float_t zAbsOutDeltaz = -100;
2561
2562 Float_t lasTanPhiLocIn = -100;
2563 Float_t lasTanPhiLocOut = -100;
2564
2565 if(histAbsY->GetEntries()>0) {
2566
2567 Double_t rotAngOut = 10;
2568 Double_t rotAngIn = 10;
2569 if((secInner%36)!=(secOuter%36))
2570 rotAngIn += 20; // 30 degrees
2571
2572 // Calculate laser mirror X position in local frame
2573 Double_t laserposOut =
2574 TMath::Abs(ltrp->GetX()*TMath::Cos(rotAngOut*TMath::DegToRad()));
2575 Double_t laserposIn =
2576 TMath::Abs(ltrp->GetX()*TMath::Cos(rotAngIn*TMath::DegToRad()));
2577
2578 // Calculate laser tan phi in local frame
2579 lasTanPhiLocOut = TMath::ASin(ltrp->GetSnp());
2580 if(lasTanPhiLocOut<0) {
2581 lasTanPhiLocIn = lasTanPhiLocOut - rotAngIn*TMath::DegToRad();
2582 lasTanPhiLocOut -= rotAngOut*TMath::DegToRad();
2583 } else {
2584
2585 lasTanPhiLocIn = lasTanPhiLocOut + rotAngIn*TMath::DegToRad();
2586 lasTanPhiLocOut += rotAngOut*TMath::DegToRad();
2587 }
2588
2589 lasTanPhiLocIn = TMath::Tan(lasTanPhiLocIn);
2590 lasTanPhiLocOut = TMath::Tan(lasTanPhiLocOut);
2591
2592 TProfile* yprof = histAbsY->ProfileX("yprof");
2593 TProfile* zprof = histAbsZ->ProfileX("zprof");
2594
2595 for(Int_t bin = 1; bin<=159; bin++) {
2596
2597 if(yprof->GetBinEntries(bin)<10&&
2598 zprof->GetBinEntries(bin)<10) {
2599 continue;
2600 }
2601
2602 // There is a problem in defining inner and outer sectors for
2603 // the outer beams (0 and 6) where both sectors are OROCs. To
2604 // make sure there is no overlap row 94 to 99 are cutted.
2605 if(((ltrp->GetBeam()==0)||(ltrp->GetBeam()==6))&&bin>=95&&bin<=100)
2606 continue;
2607
2608 Int_t row = (bin-1);
2609 if(row>62)
2610 row -= 63;
2611
2612 Bool_t isOuter = kTRUE;
2613 Int_t sector = TMath::Nint(fBeamSectorOuter[id]);
2614
2615 if(bin<=62 ||
2616 (((ltrp->GetBeam()==0)||(ltrp->GetBeam()==6))&&bin<=95)) {
2617
2618 isOuter = kFALSE;
2619 sector = TMath::Nint(fBeamSectorInner[id]);
2620 }
2621
2622
2623 Double_t x = tpcparam->GetPadRowRadii(sector, row); // slope
2624 vecN[bin-1] =yprof->GetBinEntries(bin);
2625 vecEy[bin-1]=yprof->GetBinError(bin);
2626 vecEz[bin-1]=zprof->GetBinError(bin);
2627 vecX[bin-1] = x;
2628 vecDY[bin-1] = yprof->GetBinContent(bin);
2629 vecDZ[bin-1] = zprof->GetBinContent(bin);
2630 if(!isOuter) { // inner
2631 vecPhi[bin-1]=lasTanPhiLocIn;
2632 // Calculate local y from residual and database
2633 Double_t y = fBeamSlopeYInner[id]*x + fBeamOffsetYInner[id]
2634 + vecDY[bin-1];
2635 vecY[bin-1] = y;
2636 Double_t r = TMath::Sqrt(x*x + y*y);
2637 vecR[bin-1] = r;
2638 // Find angle between laser vector and R vector
2639 // cos angle = (x,y).(1,fBeamSlopeYInner)/R/sqrt(1+fBeamSlopeYInner**2)
2640 Double_t cosPhi = x + y*fBeamSlopeYInner[id];
2641 cosPhi /= r*TMath::Sqrt(1+fBeamSlopeYInner[id]*fBeamSlopeYInner[id]);
2642 vecPhiR[bin-1] = TMath::Tan(TMath::ACos(cosPhi));
2643 if(lasTanPhiLocIn<0)
2644 vecPhiR[bin-1]*=-1; // to have the same sign
2645
2646 if(yprof->GetBinEntries(bin)>=10) {
2647 lfabsyInner.AddPoint(&x, yprof->GetBinContent(bin),
2648 TMath::Max(yprof->GetBinError(bin), 0.001));
2649 }
2650 if(zprof->GetBinEntries(bin)>=10) {
2651 lfabszInner.AddPoint(&x, zprof->GetBinContent(bin),
2652 TMath::Max(zprof->GetBinError(bin), 0.001));
2653 }
2654 } else { // outer
2655 vecPhi[bin-1]=lasTanPhiLocOut;
2656 // Calculate local y from residual and database
2657 Double_t y = fBeamSlopeYOuter[id]*x + fBeamOffsetYOuter[id]
2658 + vecDY[bin-1];
2659 vecY[bin-1] = y;
2660 Double_t r = TMath::Sqrt(x*x + y*y);
2661 vecR[bin-1] = r;
2662
2663 Double_t cosPhi = x + y*fBeamSlopeYOuter[id];
2664 cosPhi /= r*TMath::Sqrt(1+fBeamSlopeYOuter[id]*fBeamSlopeYOuter[id]);
2665 vecPhiR[bin-1] = TMath::Tan(TMath::ACos(cosPhi));
2666 if(lasTanPhiLocOut<0)
2667 vecPhiR[bin-1]*=-1; // to have the same sign
2668
2669 if(yprof->GetBinEntries(bin)>=10) {
2670 lfabsyOuter.AddPoint(&x, yprof->GetBinContent(bin),
2671 TMath::Max(yprof->GetBinError(bin), 0.001));
2672 }
2673 if(zprof->GetBinEntries(bin)>=10) {
2674 lfabszOuter.AddPoint(&x, zprof->GetBinContent(bin),
2675 TMath::Max(zprof->GetBinError(bin), 0.001));
2676 }
2677 }
2678
2679 }
2680
2681 delete yprof; delete zprof;
2682
2683
2684 // Fit laser abs residuals with linear robust fit (exclude 5% outliers)
2685 nClInY = lfabsyInner.GetNpoints();
2686 if(lfabsyInner.GetNpoints()>10) {
2687 lfabsyInner.EvalRobust(0.95);
2688
2689 TVectorD result(2);
2690 lfabsyInner.GetParameters(result);
2691 yAbsInOffset = result[0];
2692 yAbsInSlope = result[1];
2693 yAbsInDeltay = yAbsInSlope*laserposIn + yAbsInOffset;
2694 }
2695 nClInZ = lfabszInner.GetNpoints();
2696 if(lfabszInner.GetNpoints()>10) {
2697 lfabszInner.EvalRobust(0.95);
2698
2699 TVectorD result(2);
2700 lfabszInner.GetParameters(result);
2701 zAbsInOffset = result[0];
2702 zAbsInSlope = result[1];
2703 zAbsInDeltaz = zAbsInSlope*laserposIn + zAbsInOffset;
2704 }
2705 nClOutY = lfabsyOuter.GetNpoints();
2706 if(lfabsyOuter.GetNpoints()>10) {
2707 lfabsyOuter.EvalRobust(0.95);
2708
2709 TVectorD result(2);
2710 lfabsyOuter.GetParameters(result);
2711 yAbsOutOffset = result[0];
2712 yAbsOutSlope = result[1];
2713 yAbsOutDeltay = yAbsOutSlope*laserposOut + yAbsOutOffset;
2714 }
2715 nClOutZ = lfabszOuter.GetNpoints();
2716 if(lfabszOuter.GetNpoints()>10) {
2717 lfabszOuter.EvalRobust(0.95);
2718
2719 TVectorD result(2);
2720 lfabszOuter.GetParameters(result);
2721 zAbsOutOffset = result[0];
2722 zAbsOutSlope = result[1];
2723 zAbsOutDeltaz = zAbsOutSlope*laserposOut + zAbsOutOffset;
2724 }
2725 }
2726
be113f6e 2727
f6220766 2728 Int_t itime=-1;
2729 Float_t coverIA=AliTPCcalibDB::GetCoverVoltage(run,0,itime);
2730 Float_t coverIC=AliTPCcalibDB::GetCoverVoltage(run,18,itime);
2731 Float_t coverOA=AliTPCcalibDB::GetCoverVoltage(run,36,itime);
2732 Float_t coverOC=AliTPCcalibDB::GetCoverVoltage(run,54,itime);
2733 //
2734 Float_t skirtA=AliTPCcalibDB::GetSkirtVoltage(run,0,itime);
2735 Float_t skirtC=AliTPCcalibDB::GetSkirtVoltage(run,18,itime);
2736 //
2737 Float_t ggOffA=AliTPCcalibDB::GetGGoffsetVoltage(run,0,itime);
2738 Float_t ggOffC=AliTPCcalibDB::GetGGoffsetVoltage(run,18,itime);
2739
2b35e8e6 2740 //
2741 (*pcstream)<<"Mean"<<
be113f6e 2742 "run="<<run<< //
f6220766 2743 //voltages
2744 "UcIA="<<coverIA<<
2745 "UcIC="<<coverIC<<
2746 "UcOA="<<coverOA<<
2747 "UcOC="<<coverOC<<
2748 "UsA="<<skirtA<<
2749 "UsC="<<skirtC<<
2750 "UggA="<<ggOffA<<
2751 "UggC="<<ggOffC<<
2752 //
be113f6e 2753 "isOK="<<isOK<< //
2754 "id="<<id<< // track id
2755 "entries="<<entries<< // number of entries
2756 "bz="<<bfield<< // bfield
2b35e8e6 2757 "LTr.="<<ltrp<< // refernece track
1fd56785 2758 //
be113f6e 2759 "mnclIn="<<mnclIn<< // mean number of clusters in inner
2760 "mnclOut="<<mnclOut<< // mean number of clusters in outer
2761 "mnclIO="<<mnclIO<< // mean number of clusters in inner+outer
2762 "mLclIn="<<mLclIn<< // mean number of clusters in inner
2763 "mLclOut="<<mLclOut<< // mean number of clusters in outer
2764 "mLclIO="<<mLclIO<< // mean number of clusters in inner+outer
2765 "mdEdx="<<mdEdx<< // mean dEdx
108953e9 2766 "edZfit="<<edZfit<< // entries z fit
be113f6e 2767 "mdZfit="<<mdZfit<< // mean z fit
108953e9 2768 "rdZfit="<<rdZfit<< // RMS z fit
be113f6e 2769 //
2770 //
2771 "mChi2YIn1="<<mChi2YIn1<< //->chi2 y inner - line
2772 "mChi2YOut1="<<mChi2YOut1<< //->chi2 y inner - line
2773 "mChi2YIn2="<<mChi2YIn2<< //->chi2 y inner - parabola
2774 "mChi2YOut2="<<mChi2YOut2<< //->chi2 y inner - parabola
2775 "mChi2YIO1="<<mChi2YIO1<< //->chi2 y IO - common
2776 "mChi2ZIn1="<<mChi2ZIn1<< //->chi2 z inner - line
2777 "mChi2ZOut1="<<mChi2ZOut1<< //->chi2 z inner - line
2778 "mChi2ZIn2="<<mChi2ZIn2<< //->chi2 z inner - parabola
2779 "mChi2ZOut2="<<mChi2ZOut2<< //->chi2 z inner - parabola
2780 "mChi2ZIO1="<<mChi2ZIO1<< //->chi2 z IO - common
2781 //
2782 //
2783 "ePy1vP0="<<ePy1vP0<<
2784 "mPy1vP0="<<mPy1vP0<<
2785 "rPy1vP0="<<rPy1vP0<<
2786 "ePy2vP0="<<ePy2vP0<<
2787 "mPy2vP0="<<mPy2vP0<<
2788 "rPy2vP0="<<rPy2vP0<<
2789 "ePy3vP0="<<ePy3vP0<<
2790 "mPy3vP0="<<mPy3vP0<<
2791 "rPy3vP0="<<rPy3vP0<<
2792 "ePy1vP1="<<ePy1vP1<<
2793 "mPy1vP1="<<mPy1vP1<<
2794 "rPy1vP1="<<rPy1vP1<<
2795 "ePy2vP1="<<ePy2vP1<<
2796 "mPy2vP1="<<mPy2vP1<<
2797 "rPy2vP1="<<rPy2vP1<<
2798 "ePy3vP1="<<ePy3vP1<<
2799 "mPy3vP1="<<mPy3vP1<<
2800 "rPy3vP1="<<rPy3vP1<<
2801 "ePy2vP2In="<<ePy2vP2In<<
2802 "mPy2vP2In="<<mPy2vP2In<<
2803 "rPy2vP2In="<<rPy2vP2In<<
2804 "ePy2vP2Out="<<ePy2vP2Out<<
2805 "mPy2vP2Out="<<mPy2vP2Out<<
2806 "rPy2vP2Out="<<rPy2vP2Out<<
2807 "ePy3vP2IO="<<ePy3vP2IO<<
2808 "mPy3vP2IO="<<mPy3vP2IO<<
2809 "rPy3vP2IO="<<rPy3vP2IO<<
2810 //
2811 //
2812 "ePz1vP0="<<ePz1vP0<<
2813 "mPz1vP0="<<mPz1vP0<<
2814 "rPz1vP0="<<rPz1vP0<<
2815 "ePz2vP0="<<ePz2vP0<<
2816 "mPz2vP0="<<mPz2vP0<<
2817 "rPz2vP0="<<rPz2vP0<<
2818 "ePz3vP0="<<ePz3vP0<<
2819 "mPz3vP0="<<mPz3vP0<<
2820 "rPz3vP0="<<rPz3vP0<<
2821 "ePz1vP1="<<ePz1vP1<<
2822 "mPz1vP1="<<mPz1vP1<<
2823 "rPz1vP1="<<rPz1vP1<<
2824 "ePz2vP1="<<ePz2vP1<<
2825 "mPz2vP1="<<mPz2vP1<<
2826 "rPz2vP1="<<rPz2vP1<<
2827 "ePz3vP1="<<ePz3vP1<<
2828 "mPz3vP1="<<mPz3vP1<<
2829 "rPz3vP1="<<rPz3vP1<<
2830 "ePz2vP2In="<<ePz2vP2In<<
2831 "mPz2vP2In="<<mPz2vP2In<<
2832 "rPz2vP2In="<<rPz2vP2In<<
2833 "ePz2vP2Out="<<ePz2vP2Out<<
2834 "mPz2vP2Out="<<mPz2vP2Out<<
2835 "rPz2vP2Out="<<rPz2vP2Out<<
2836 "ePz3vP2IO="<<ePz3vP2IO<<
2837 "mPz3vP2IO="<<mPz3vP2IO<<
2838 "rPz3vP2IO="<<rPz3vP2IO<<
2839 //
2840 //
2841 //
2842 "lx0="<<lxyz[0]<< // reference x
2843 "lx1="<<lxyz[1]<< // reference y
2844 "lx2="<<lxyz[2]<< // refernece z
2845 "lpx0="<<lpxyz[0]<< // reference x
2b35e8e6 2846 "lpx1="<<lpxyz[1]<< // reference y
263d466a 2847 "lpx2="<<lpxyz[2]<< // refernece z
1fd56785 2848 //
2849 "msig="<<meanS<<
2850 //
2b35e8e6 2851 "mphi="<<meanphi<< //
2852 "rmsphi="<<rmsphi<< //
2853 "gphi1="<<gphi1<<
2854 "gphi2="<<gphi2<<
1fd56785 2855 //
2b35e8e6 2856 "mphiP="<<meanphiP<< //
2857 "rmsphiP="<<rmsphiP<< //
2858 "gphiP1="<<gphiP1<<
2859 "gphiP2="<<gphiP2<<
1fd56785 2860 //
2b35e8e6 2861 "meanZ="<<meanZ<<
2862 "rmsZ="<<rmsZ<<
2863 "gz1="<<gz1<<
2864 "gz2="<<gz2<<
e5c59218 2865 //
2866 "gp31="<<gp31<< //gaus mean - tgl
2867 "gp32="<<gp32<< //gaus rms -tgl
2868 "gp41="<<gp41<< //gaus mean - P4
2869 "gp42="<<gp42<< //gaus rms - P4
f6220766 2870 // Parameters from abs res analysis
2871 "SecIn="<<secInner<< // inner sector
2872 "SecOut="<<secOuter<< // outer sector
2873 "lasTanPhiLocIn="<<lasTanPhiLocIn<< // laser tan phi in local frame (inner)
2874 "lasTanPhiLocOut="<<lasTanPhiLocOut<<// laser tan phi in local frame (outer)
2875 "X.="<<&vecX<< // local x
2876 "Y.="<<&vecY<< // local y
2877 "R.="<<&vecR<< // radius
2878 "dY.="<<&vecDY<< // abs y residuals
2879 "dZ.="<<&vecDZ<< // abs z residuals
2880 "eY.="<<&vecEy<< // error of y residuals
2881 "eZ.="<<&vecEz<< // error of z residuals
2882 "kY.="<<&vecPhi<< // local tangent y (fixed for sector)
2883 "kYR.="<<&vecPhiR<< // tangent between laser and R vector (varies inside sector)
2884 "nCl.="<<&vecN<< // number of clusters
2885 //
2886 "nClInY="<<nClInY<< // Number of clusters for inner y
2887 "yAbsInOffset="<<yAbsInOffset<< // fitted offset absres (inner y)
2888 "yAbsInSlope="<<yAbsInSlope << // fitted slope absres (inner y)
2889 "yAbsInDeltay="<<yAbsInDeltay<< // fitted laser offset absres (inner y)
2890 "nClInZ="<<nClInZ<< // Number of clusters for inner z
2891 "zAbsInOffset="<<zAbsInOffset<< // fitted offset absres (inner z)
2892 "zAbsInSlope="<<zAbsInSlope << // fitted slope absres (inner z)
2893 "zAbsInDeltaz="<<zAbsInDeltaz<< // fitted laser offset absres (inner z)
2894 //
2895 "nClOutY="<<nClOutY<< // Number of clusters for outer y
2896 "yAbsOutOffset="<<yAbsOutOffset<< // fitted offset absres (outer y)
2897 "yAbsOutSlope="<<yAbsOutSlope << // fitted slope absres (outer y)
2898 "yAbsOutDeltay="<<yAbsOutDeltay<< // fitted laser offset absres (outer y)
2899 "nClOutZ="<<nClOutZ<< // Number of clusters for outer z
2900 "zAbsOutOffset="<<zAbsOutOffset<< // fitted offset absres (outer z)
2901 "zAbsOutSlope="<<zAbsOutSlope << // fitted slope absres (outer z)
2902 "zAbsOutDeltaz="<<zAbsOutDeltaz<< // fitted laser offset absres (outer z)
2903 //
2b35e8e6 2904 "\n";
2905 }
2906 delete pcstream;
be113f6e 2907 /*
2908 Browse the content
2909 TFile fmean("laserMean.root");
2910
2911
2912 */
2913
2914
2b35e8e6 2915}
2916
1fd56785 2917
2918
592a0c8f 2919void AliTPCcalibLaser::DumpScanInfo(TTree * chain, const char * cutUser){
1fd56785 2920 //
2921 //
2922 //
2923 TTreeSRedirector *pcstream = new TTreeSRedirector("laserScan.root");
2924 TFile * f = pcstream->GetFile();
2925 f->mkdir("dirphi");
2926 f->mkdir("dirphiP");
2927 f->mkdir("dirZ");
2928 TF1 fp("p1","pol1");
2929 //
2930 //
2931 char cut[1000];
2932 char grname[1000];
2933 char grnamefull[1000];
2934
2935 Double_t mphi[100];
2936 Double_t mphiP[100];
2937 Double_t smphi[100];
2938 Double_t smphiP[100];
2939 Double_t mZ[100];
2940 Double_t smZ[100];
2941 Double_t bz[100];
2942 Double_t sbz[100];
2943 // fit parameters
2944 Double_t pphi[3];
2945 Double_t pphiP[3];
2946 Double_t pmZ[3];
0cd8bb8d 2947
1fd56785 2948 //
2949 for (Int_t id=0; id<336; id++){
2950 // id =205;
592a0c8f 2951 sprintf(cut,"fId==%d&&%s",id,cutUser);
1fd56785 2952 Int_t entries = chain->Draw("bz",cut,"goff");
2953 if (entries<3) continue;
e5c59218 2954 AliTPCLaserTrack *ltrp = 0;
1fd56785 2955 if (!AliTPCLaserTrack::GetTracks()) AliTPCLaserTrack::LoadTracks();
2956 ltrp =(AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(id);
2957 Double_t lxyz[3];
2958 Double_t lpxyz[3];
2959 ltrp->GetXYZ(lxyz);
2960 ltrp->GetPxPyPz(lpxyz);
2961
2962 chain->Draw("bz",cut,"goff");
2963 memcpy(bz, chain->GetV1(), entries*sizeof(Double_t));
2964 chain->Draw("0.01*abs(bz)+0.02",cut,"goff");
2965 memcpy(sbz, chain->GetV1(), entries*sizeof(Double_t));
2966 //
2967 chain->Draw("gphi1",cut,"goff");
2968 memcpy(mphi, chain->GetV1(), entries*sizeof(Double_t));
592a0c8f 2969 chain->Draw("0.05*abs(mphi)+abs(gphi2)*0.5+0.05",cut,"goff");
1fd56785 2970 memcpy(smphi, chain->GetV1(), entries*sizeof(Double_t));
2971 //
2972 chain->Draw("gphiP1",cut,"goff");
2973 memcpy(mphiP, chain->GetV1(), entries*sizeof(Double_t));
592a0c8f 2974 chain->Draw("0.05*abs(mphiP)+abs(gphiP2)*0.5+0.001",cut,"goff");
1fd56785 2975 memcpy(smphiP, chain->GetV1(), entries*sizeof(Double_t));
2976 //
2977 chain->Draw("gz1",cut,"goff");
2978 memcpy(mZ, chain->GetV1(), entries*sizeof(Double_t));
592a0c8f 2979 chain->Draw("0.01*abs(meanZ)+abs(gz2)*0.5+0.1",cut,"goff");
1fd56785 2980 memcpy(smZ, chain->GetV1(), entries*sizeof(Double_t));
2981 //
2982 //
2983 sprintf(grnamefull,"Side_%d_Bundle_%d_Rod_%d_Beam_%d",
2984 ltrp->GetSide(), ltrp->GetBundle(), ltrp->GetRod(), ltrp->GetBeam());
2985 // store data
2986 // phi
2987 f->cd("dirphi");
592a0c8f 2988 Float_t phiB0 =0;
2989 Float_t phiPB0=0;
2990 Float_t zB0=0;
2991 for (Int_t ientry=0; ientry<entries; ientry++){
2992 if (TMath::Abs(bz[ientry])<0.05){
2993 phiB0 = mphi[ientry];
2994 phiPB0 = mphiP[ientry];
2995 zB0 = mZ[ientry];
2996 }
2997 }
1fd56785 2998 TGraphErrors *grphi = new TGraphErrors(entries,bz,mphi,sbz,smphi);
2999 grphi->Draw("a*");
3000 grphi->Fit(&fp);
3001 pphi[0] = fp.GetParameter(0); // offset
3002 pphi[1] = fp.GetParameter(1); // slope
3003 pphi[2] = TMath::Sqrt(fp.GetChisquare()/(entries-2.)); // normalized chi2
3004 sprintf(grname,"phi_id%d",id);
3005 grphi->SetName(grname); grphi->SetTitle(grnamefull);
3006 grphi->GetXaxis()->SetTitle("b_{z} (T)");
3007 grphi->GetYaxis()->SetTitle("#Delta r#phi (cm)");
ae69f16f 3008 grphi->SetMaximum(1.2);
3009 grphi->SetMinimum(-1.2);
1fd56785 3010 grphi->Draw("a*");
3011
3012 grphi->Write();
3013 gPad->SaveAs(Form("pic/phi/phi_%s.gif",grnamefull));
3014 // phiP
3015 f->cd("dirphiP)");
3016 TGraphErrors *grphiP = new TGraphErrors(entries,bz,mphiP,sbz,smphiP);
3017 grphiP->Draw("a*");
3018 grphiP->Fit(&fp);
3019 pphiP[0] = fp.GetParameter(0); // offset
3020 pphiP[1] = fp.GetParameter(1); // slope
3021 pphiP[2] = TMath::Sqrt(fp.GetChisquare()/(entries-2.)); // normalized chi2
3022 sprintf(grname,"phiP_id%d",id);
3023 grphiP->SetName(grname); grphiP->SetTitle(grnamefull);
3024 grphiP->GetXaxis()->SetTitle("b_{z} (T)");
3025 grphiP->GetYaxis()->SetTitle("#Delta #phi (rad)");
ae69f16f 3026 grphiP->SetMaximum(pphiP[0]+0.005);
3027 grphiP->SetMinimum(pphiP[0]-0.005);
1fd56785 3028
3029 gPad->SaveAs(Form("pic/phiP/phiP_%s.gif",grnamefull));
3030 grphiP->Write();
3031 //
3032 //Z
3033 f->cd("dirZ");
3034 TGraphErrors *grmZ = new TGraphErrors(entries,bz,mZ,sbz,smZ);
3035 grmZ->Draw("a*");
3036 grmZ->Fit(&fp);
3037 pmZ[0] = fp.GetParameter(0); // offset
3038 pmZ[1] = fp.GetParameter(1); // slope
3039 pmZ[2] = TMath::Sqrt(fp.GetChisquare()/(entries-2.)); // normalized chi2
3040 sprintf(grname,"mZ_id%d",id);
3041 grmZ->SetName(grname); grmZ->SetTitle(grnamefull);
3042 grmZ->GetXaxis()->SetTitle("b_{z} (T)");
3043 grmZ->GetYaxis()->SetTitle("#Delta z (cm)");
3044
3045 gPad->SaveAs(Form("pic/z/z_%s.gif",grnamefull));
3046 grmZ->Write();
592a0c8f 3047 //
3048 // P4
3049 //
1fd56785 3050
3051 for (Int_t ientry=0; ientry<entries; ientry++){
3052 (*pcstream)<<"Mean"<<
3053 "id="<<id<<
3054 "LTr.="<<ltrp<<
3055 "entries="<<entries<<
3056 "bz="<<bz[ientry]<<
3057 "lx0="<<lxyz[0]<< // reference x
3058 "lx1="<<lxyz[1]<< // reference y
3059 "lx2="<<lxyz[2]<< // refernece z
3060 "lpx0="<<lpxyz[0]<< // reference x
3061 "lpx1="<<lpxyz[1]<< // reference y
3062 "lpx2="<<lpxyz[2]<< // refernece z
3063 //values
592a0c8f 3064 "phiB0="<<phiB0<< // position shift at 0 field
3065 "phiPB0="<<phiPB0<< // angular shift at 0 field
3066 "zB0="<<zB0<< // z shift for 0 field
3067 //
1fd56785 3068 "gphi1="<<mphi[ientry]<< // mean - from gaus fit
3069 "pphi0="<<pphi[0]<< // offset
592a0c8f 3070 "pphi1="<<pphi[1]<< // slope
1fd56785 3071 "pphi2="<<pphi[2]<< // norm chi2
3072 //
3073 "gphiP1="<<mphiP[ientry]<< // mean - from gaus fit
3074 "pphiP0="<<pphiP[0]<< // offset
592a0c8f 3075 "pphiP1="<<pphiP[1]<< // slope
1fd56785 3076 "pphiP2="<<pphiP[2]<< // norm chi2
3077 //
3078 "gz1="<<mZ[ientry]<<
3079 "pmZ0="<<pmZ[0]<< // offset
592a0c8f 3080 "pmZ1="<<pmZ[1]<< // slope
1fd56785 3081 "pmZ2="<<pmZ[2]<< // norm chi2
3082 "\n";
3083 }
3084 }
3085
3086 delete pcstream;
3087
3088}
3089
3090
c6914c83 3091void AliTPCcalibLaser::Analyze(){
3092 //
3093 //
3094 //
3095}
3096
3097
c03e3250 3098Long64_t AliTPCcalibLaser::Merge(TCollection *li) {
c6914c83 3099
c03e3250 3100 TIterator* iter = li->MakeIterator();
3101 AliTPCcalibLaser* cal = 0;
108953e9 3102 static Int_t counter0=0;
c03e3250 3103 while ((cal = (AliTPCcalibLaser*)iter->Next())) {
3104 if (!cal->InheritsFrom(AliTPCcalibLaser::Class())) {
5de4d13e 3105 // Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
c03e3250 3106 return -1;
3107 }
108953e9 3108 printf("Marging number %d\n", counter0);
3109 counter0++;
c03e3250 3110 //
be113f6e 3111 MergeFitHistos(cal);
c03e3250 3112 TH1F *h=0x0;
3113 TH1F *hm=0x0;
108953e9 3114 TH2F *h2=0x0;
3115 TH2F *h2m=0x0;
3116 // TProfile *hp=0x0;
3117 //TProfile *hpm=0x0;
c03e3250 3118
3119 for (Int_t id=0; id<336; id++){
3120 // merge fDeltaZ histograms
263d466a 3121 hm = (TH1F*)cal->fDeltaZ.At(id);
c03e3250 3122 h = (TH1F*)fDeltaZ.At(id);
3123 if (!h) {
3124 h=new TH1F(Form("hisdz%d",id),Form("hisdz%d",id),1000,-10,10);
5732d662 3125 h->SetDirectory(0);
c03e3250 3126 fDeltaZ.AddAt(h,id);
3127 }
3128 if (hm) h->Add(hm);
e5c59218 3129 // merge fP3 histograms
3130 hm = (TH1F*)cal->fDeltaP3.At(id);
3131 h = (TH1F*)fDeltaP3.At(id);
3132 if (!h) {
3133 h=new TH1F(Form("hisPar3v%d",id),Form("hisPar3v%d",id),400,-0.06,0.06);
3134 h->SetDirectory(0);
3135 fDeltaP3.AddAt(h,id);
3136 }
3137 if (hm) h->Add(hm);
3138 // merge fP4 histograms
3139 hm = (TH1F*)cal->fDeltaP4.At(id);
3140 h = (TH1F*)fDeltaP4.At(id);
3141 if (!h) {
3142 h=new TH1F(Form("hisPar4v%d",id),Form("hisPar4v%d",id),200,-0.06,0.06);
3143 h->SetDirectory(0);
3144 fDeltaP4.AddAt(h,id);
3145 }
3146 if (hm) h->Add(hm);
3147
3148 //
c03e3250 3149 // merge fDeltaPhi histograms
263d466a 3150 hm = (TH1F*)cal->fDeltaPhi.At(id);
c03e3250 3151 h = (TH1F*)fDeltaPhi.At(id);
3152 if (!h) {
3153 h= new TH1F(Form("hisdphi%d",id),Form("hisdphi%d",id),1000,-1,1);
5732d662 3154 h->SetDirectory(0);
c03e3250 3155 fDeltaPhi.AddAt(h,id);
3156 }
3157 if (hm) h->Add(hm);
3158 // merge fDeltaPhiP histograms
263d466a 3159 hm = (TH1F*)cal->fDeltaPhiP.At(id);
c03e3250 3160 h = (TH1F*)fDeltaPhiP.At(id);
3161 if (!h) {
3162 h=new TH1F(Form("hisdphiP%d",id),Form("hisdphiP%d",id),1000,-0.01,0.01);
5732d662 3163 h->SetDirectory(0);
c03e3250 3164 fDeltaPhiP.AddAt(h,id);
3165 }
3166 if (hm) h->Add(hm);
3167 // merge fSignals histograms
263d466a 3168 hm = (TH1F*)cal->fSignals.At(id);
c03e3250 3169 h = (TH1F*)fSignals.At(id);
3170 if (!h) {
e5c59218 3171 h=new TH1F(Form("hisSignal%d",id),Form("hisSignal%d",id),100,0,300);
5732d662 3172 h->SetDirectory(0);
c03e3250 3173 fSignals.AddAt(h,id);
3174 }
263d466a 3175 if (hm) h->Add(hm);
95a0e09b 3176 //
3177 //
108953e9 3178 // merge ProfileY histograms -0
3179 h2m = (TH2F*)cal->fDeltaYres.At(id);
3180 h2 = (TH2F*)fDeltaYres.At(id);
5de4d13e 3181 if (h2m&&h2) h2->Add(h2m);
95a0e09b 3182 //
108953e9 3183 h2m = (TH2F*)cal->fDeltaZres.At(id);
3184 h2 = (TH2F*)fDeltaZres.At(id);
5de4d13e 3185 if (h2m&&h2) h2->Add(h2m);
108953e9 3186 // merge ProfileY histograms - 2
3187 h2m = (TH2F*)cal->fDeltaYres2.At(id);
3188 h2 = (TH2F*)fDeltaYres2.At(id);
5de4d13e 3189 if (h2m&&h2) h2->Add(h2m);
108953e9 3190 //
3191 h2m = (TH2F*)cal->fDeltaZres2.At(id);
3192 h2 = (TH2F*)fDeltaZres2.At(id);
5de4d13e 3193 if (h2m&&h2) h2->Add(h2m);
7af539c6 3194
3195 // merge ProfileY histograms - abs
3196 h2m = (TH2F*)cal->fDeltaYresAbs.At(id);
3197 h2 = (TH2F*)fDeltaYresAbs.At(id);
3198 if (h2m&&h2) h2->Add(h2m);
3199
3200 h2m = (TH2F*)cal->fDeltaZresAbs.At(id);
3201 h2 = (TH2F*)fDeltaZresAbs.At(id);
3202 if (h2m&&h2) h2->Add(h2m);
108953e9 3203 // merge ProfileY histograms - 3
3204 //h2m = (TH2F*)cal->fDeltaYres3.At(id);
3205 //h2 = (TH2F*)fDeltaYres3.At(id);
3206 //if (h2m) h2->Add(h2m);
3207 //
3208 //h2m = (TH2F*)cal->fDeltaZres3.At(id);
3209 //h2 = (TH2F*)fDeltaZres3.At(id);
3210 //if (h2m) h->Add(h2m);
95a0e09b 3211 //
3212 //
c03e3250 3213 }
3214 }
3215 return 0;
3216}
3217
be113f6e 3218void AliTPCcalibLaser::MakeFitHistos(){
3219 //
3220 // Make a fit histograms
3221 //
3222 // Number of clusters
3223 //
3224 //TH2F *fHisNclIn; //->Number of clusters inner
3225 //TH2F *fHisNclOut; //->Number of clusters outer
3226 //TH2F *fHisNclIO; //->Number of cluster inner outer
3227 // TH2F *fHisdEdx; //->dEdx histo
3228 fHisNclIn = new TH2F("HisNclIn","HisNclIn",336,0,336,64,10,64);
3229 fHisNclOut = new TH2F("HisNclOut","HisNclOut",336,0,336,100,10,100);
3230 fHisNclIO = new TH2F("HisNclIO","HisNclIO",336,0,336,160,10,160);
3231 //
3232 fHisLclIn = new TH2F("HisLclIn","HisLclIn",336,0,336,64,10,64);
3233 fHisLclOut = new TH2F("HisLclOut","HisLclOut",336,0,336,100,10,150);
3234 fHisLclIO = new TH2F("HisLclIO","HisLclIO",336,0,336,160,10,160);
3235 //
3236 fHisdEdx = new TH2F("HisdEdx","HisdEdx",336,0,336,160,1,50);
108953e9 3237 fHisdZfit = new TH2F("HisdZfit","HisdZfit",336,0,336,300,-0.3,0.3);
be113f6e 3238
3239 //
3240 // Chi2
3241 //
3242 // TH2F *fHisChi2YIn1; //->chi2 y inner - line
3243 // TH2F *fHisChi2YOut1; //->chi2 y inner - line
3244 // TH2F *fHisChi2YIn2; //->chi2 y inner - parabola
3245 // TH2F *fHisChi2YOut2; //->chi2 y inner - parabola
3246 // TH2F *fHisChi2YIO1; //->chi2 y IO - common
3247 fHisChi2YIn1 = new TH2F("Chi2YIn1","Chi2YIn1",336,0,336,500,0.001,0.5);
3248 fHisChi2YOut1 = new TH2F("Chi2YOut1","Chi2YOut1",336,0,336,500,0.001,0.5);
3249 fHisChi2YIn2 = new TH2F("Chi2YIn2","Chi2YIn2",336,0,336,500,0.001,0.5);
3250 fHisChi2YOut2 = new TH2F("Chi2YOut2","Chi2YOut2",336,0,336,500,0.001,0.5);
3251 fHisChi2YIO1 = new TH2F("Chi2YIO1","Chi2YIO1",336,0,336,500,0.001,0.5);
3252 // TH2F *fHisChi2ZIn1; //->chi2 z inner - line
3253 // TH2F *fHisChi2ZOut1; //->chi2 z inner - line
3254 // TH2F *fHisChi2ZIn2; //->chi2 z inner - parabola
3255 // TH2F *fHisChi2ZOut2; //->chi2 z inner - parabola
3256 // TH2F *fHisChi2ZIO1; //->chi2 z IO - common
3257 fHisChi2ZIn1 = new TH2F("Chi2ZIn1","Chi2ZIn1",336,0,336,500,0.001,0.5);
3258 fHisChi2ZOut1 = new TH2F("Chi2ZOut1","Chi2ZOut1",336,0,336,500,0.001,0.5);
3259 fHisChi2ZIn2 = new TH2F("Chi2ZIn2","Chi2ZIn2",336,0,336,500,0.001,0.5);
3260 fHisChi2ZOut2 = new TH2F("Chi2ZOut2","Chi2ZOut2",336,0,336,500,0.001,0.5);
3261 fHisChi2ZIO1 = new TH2F("Chi2ZIO1","Chi2ZIO1",336,0,336,500,0.001,0.5);
3262 //
3263 // Fit
3264 //
3265 //
3266 // TH2F *fHisPy1vP0; //-> delta y P0outer-P0inner - line
3267 // TH2F *fHisPy2vP0; //-> delta y P0outer-P0inner - parabola
3268 // TH2F *fHisPy3vP0; //-> delta y P0outer-P0inner - common parabola
3269 // TH2F *fHisPy1vP1; //-> delta ky P1outer-P1inner - line
3270 // TH2F *fHisPy2vP1; //-> delta ky P1outer-P1inner - parabola
3271 // TH2F *fHisPy3vP1; //-> delta ky P1outer-P1inner - common parabola
3272 // TH2F *fHisPy2vP2In; //-> Curv P2inner - parabola
3273 // TH2F *fHisPy2vP2Out; //-> Curv P2outer - parabola
3274 // TH2F *fHisPy3vP2IO; //-> Curv P2outerinner - common parabola
3275 fHisPy1vP0 = new TH2F("HisPy1vP0", "HisPy1vP0",336,0,336,500,-0.3,0.3);
3276 fHisPy2vP0 = new TH2F("HisPy2vP0", "HisPy2vP0",336,0,336,500,-0.3,0.3);
3277 fHisPy3vP0 = new TH2F("HisPy3vP0", "HisPy3vP0",336,0,336,500,-0.3,0.3);
3278 fHisPy1vP1 = new TH2F("HisPy1vP1", "HisPy1vP1",336,0,336,500,-0.01,0.01);
3279 fHisPy2vP1 = new TH2F("HisPy2vP1", "HisPy2vP1",336,0,336,500,-0.01,0.01);
3280 fHisPy3vP1 = new TH2F("HisPy3vP1", "HisPy3vP1",336,0,336,500,-0.01,0.01);
3281 fHisPy2vP2In = new TH2F("HisPy2vP2In", "HisPy2vP2In",336,0,336,500,-0.0003,0.0003);
3282 fHisPy2vP2Out= new TH2F("HisPy2vP2Out","HisPy2vP2Out",336,0,336,500,-0.0003,0.0003);
3283 fHisPy3vP2IO = new TH2F("HisPy3vP2IO", "HisPy3vP2IO",336,0,336,500,-0.0003,0.0003);
3284 //
3285 //
3286 // TH2F *fHisPz1vP0; //-> delta z P0outer-P0inner - line
3287 // TH2F *fHisPz2vP0; //-> delta z P0outer-P0inner - parabola
3288 // TH2F *fHisPz3vP0; //-> delta z P0outer-P0inner - common parabola
3289 // TH2F *fHisPz1vP1; //-> delta kz P1outer-P1inner - line
3290 // TH2F *fHisPz2vP1; //-> delta kz P1outer-P1inner - parabola
3291 // TH2F *fHisPz3vP1; //-> delta kz P1outer-P1inner - common parabola
3292 // TH2F *fHisPz2vP2In; //-> Curv P2inner - parabola
3293 // TH2F *fHisPz2vP2Out; //-> Curv P2outer - parabola
3294 // TH2F *fHisPz3vP2IO; //-> Curv P2outerinner - common parabola
3295 fHisPz1vP0 = new TH2F("HisPz1vP0", "HisPz1vP0",336,0,336,500,-0.3,0.3);
3296 fHisPz2vP0 = new TH2F("HisPz2vP0", "HisPz2vP0",336,0,336,500,-0.3,0.3);
3297 fHisPz3vP0 = new TH2F("HisPz3vP0", "HisPz3vP0",336,0,336,500,-0.3,0.3);
3298 fHisPz1vP1 = new TH2F("HisPz1vP1", "HisPz1vP1",336,0,336,500,-0.01,0.01);
3299 fHisPz2vP1 = new TH2F("HisPz2vP1", "HisPz2vP1",336,0,336,500,-0.01,0.01);
3300 fHisPz3vP1 = new TH2F("HisPz3vP1", "HisPz3vP1",336,0,336,500,-0.01,0.01);
3301 fHisPz2vP2In = new TH2F("HisPz2vP2In", "HisPz2vP2In",336,0,336,500,-0.0003,0.0003);
3302 fHisPz2vP2Out= new TH2F("HisPz2vP2Out","HisPz2vP2Out",336,0,336,500,-0.0002,0.0002);
3303 fHisPz3vP2IO = new TH2F("HisPz3vP2IO", "HisPz3vP2IO",336,0,336,500,-0.0002,0.0002);
3304
41a1f702 3305 fHisYAbsErrors = new TH1F("HisYabsErrors", "Errors per beam (y)", 336, 0, 336);
3306 fHisZAbsErrors = new TH1F("HisZabsErrors", "Errors per beam (z)", 336, 0, 336);
3307
be113f6e 3308 fHisNclIn->SetDirectory(0); //->Number of clusters inner
3309 fHisNclOut->SetDirectory(0); //->Number of clusters outer
3310 fHisNclIO->SetDirectory(0); //->Number of cluster inner outer
3311 fHisLclIn->SetDirectory(0); //->Level arm inner
3312 fHisLclOut->SetDirectory(0); //->Level arm outer
3313 fHisLclIO->SetDirectory(0); //->Number of cluster inner outer
3314 fHisdEdx->SetDirectory(0); //->Number of cluster inner outer
3315 fHisdZfit->SetDirectory(0); //->Number of cluster inner outer
3316 //
3317 //
3318 fHisChi2YIn1->SetDirectory(0); //->chi2 y inner - line
3319 fHisChi2YOut1->SetDirectory(0); //->chi2 y inner - line
3320 fHisChi2YIn2->SetDirectory(0); //->chi2 y inner - parabola
3321 fHisChi2YOut2->SetDirectory(0); //->chi2 y inner - parabola
3322 fHisChi2YIO1->SetDirectory(0); //->chi2 y IO - common
3323 fHisChi2ZIn1->SetDirectory(0); //->chi2 z inner - line
3324 fHisChi2ZOut1->SetDirectory(0); //->chi2 z inner - line
3325 fHisChi2ZIn2->SetDirectory(0); //->chi2 z inner - parabola
3326 fHisChi2ZOut2->SetDirectory(0); //->chi2 z inner - parabola
3327 fHisChi2ZIO1->SetDirectory(0); //->chi2 z IO - common
3328 //
3329 //
3330 fHisPy1vP0->SetDirectory(0); //-> delta y P0outer-P0inner - line
3331 fHisPy2vP0->SetDirectory(0); //-> delta y P0outer-P0inner - parabola
3332 fHisPy3vP0->SetDirectory(0); //-> delta y P0outer-P0inner - common parabola
3333 fHisPy1vP1->SetDirectory(0); //-> delta ky P1outer-P1inner - line
3334 fHisPy2vP1->SetDirectory(0); //-> delta ky P1outer-P1inner - parabola
3335 fHisPy3vP1->SetDirectory(0); //-> delta ky P1outer-P1inner - common parabola
3336 fHisPy2vP2In->SetDirectory(0); //-> Curv P2inner - parabola
3337 fHisPy2vP2Out->SetDirectory(0); //-> Curv P2outer - parabola
3338 fHisPy3vP2IO->SetDirectory(0); //-> Curv P2outerinner - common parabola
3339 //
3340 //
3341 fHisPz1vP0->SetDirectory(0); //-> delta z P0outer-P0inner - line
3342 fHisPz2vP0->SetDirectory(0); //-> delta z P0outer-P0inner - parabola
3343 fHisPz3vP0->SetDirectory(0); //-> delta z P0outer-P0inner - common parabola
3344 fHisPz1vP1->SetDirectory(0); //-> delta kz P1outer-P1inner - line
3345 fHisPz2vP1->SetDirectory(0); //-> delta kz P1outer-P1inner - parabola
3346 fHisPz3vP1->SetDirectory(0); //-> delta kz P1outer-P1inner - common parabola
3347 fHisPz2vP2In->SetDirectory(0); //-> Curv P2inner - parabola
3348 fHisPz2vP2Out->SetDirectory(0); //-> Curv P2outer - parabola
3349 fHisPz3vP2IO->SetDirectory(0); //-> Curv P2outerinner - common parabola
3350
41a1f702 3351 fHisYAbsErrors->SetDirectory(0); //-> total errors per beam in the abs res y analysis
3352 fHisZAbsErrors->SetDirectory(0); //-> total errors per beam in the abs res z analysis
3353
3354
be113f6e 3355
108953e9 3356 //
3357 //
3358 //
3359 for (Int_t id=0; id<336;id++){
3360 TH2F *profy = (TH2F*)fDeltaYres.UncheckedAt(id);
3361 TH2F *profz = (TH2F*)fDeltaZres.UncheckedAt(id);
3362 TH2F *profy2 = (TH2F*)fDeltaYres2.UncheckedAt(id);
3363 TH2F *profz2 = (TH2F*)fDeltaZres2.UncheckedAt(id);
7af539c6 3364 TH2F *profyabs = (TH2F*)fDeltaYresAbs.UncheckedAt(id);
3365 TH2F *profzabs = (TH2F*)fDeltaYresAbs.UncheckedAt(id);
108953e9 3366 // TH2F *profy3 = (TH2F*)fDeltaYres3.UncheckedAt(id);
3367 //TH2F *profz3 = (TH2F*)fDeltaZres3.UncheckedAt(id);
3368 if (!profy){
3369 profy=new TH2F(Form("pry%03d",id),Form("Y Residuals for Laser Beam %03d -Linear",id),160,0,160,50,-0.5,0.5);
3370 profy->SetDirectory(0);
3371 fDeltaYres.AddAt(profy,id);
3372 profy2=new TH2F(Form("pry%03d",id),Form("Y Residuals for Laser Beam %03d -Parabolic",id),160,0,160,50,-0.5,0.5);
3373 profy2->SetDirectory(0);
3374 fDeltaYres2.AddAt(profy2,id);
41a1f702 3375 if(!fUseFixedDriftV)
3376 profyabs=new TH2F(Form("pryabs%03d",id),Form("Y Residuals for Laser Beam %03d -Absolute",id),160,0,160,100,-1.0,1.0); // has to be bigger based on earlier studies
3377 else
3378 profyabs=new TH2F(Form("pryabs%03d",id),Form("Y Residuals for Laser Beam %03d -Absolute",id),160,0,160,200,-2.0,2.0); // has to be bigger based on earlier studies
7af539c6 3379 profyabs->SetDirectory(0);
3380 fDeltaYresAbs.AddAt(profyabs,id);
108953e9 3381 //profy3=new TH2F(Form("pry%03d",id),Form("Y Residuals for Laser Beam %03d- Parabolic2",id),160,0,160,100,-0.5,0.5);
3382 //profy3->SetDirectory(0);
3383 //fDeltaYres3.AddAt(profy3,id);
3384 }
3385 if (!profz){
3386 profz=new TH2F(Form("prz%03d",id),Form("Z Residuals for Laser Beam %03d Linear",id),160,0,160,50,-0.5,0.5);
3387 profz->SetDirectory(0);
3388 fDeltaZres.AddAt(profz,id);
3389 profz2=new TH2F(Form("prz%03d",id),Form("Z Residuals for Laser Beam %03d - Parabolic",id),160,0,160,50,-0.5,0.5);
3390 profz2->SetDirectory(0);
3391 fDeltaZres2.AddAt(profz2,id);
41a1f702 3392 if(!fUseFixedDriftV)
3393 profzabs=new TH2F(Form("przabs%03d",id),Form("Z Residuals for Laser Beam %03d -Absolute",id),160,0,160,100,-1.0,1.0); // has to be bigger based on earlier studies
3394 else
3395 profzabs=new TH2F(Form("przabs%03d",id),Form("Z Residuals for Laser Beam %03d -Absolute",id),160,0,160,200,-2.0,2.0); // has to be bigger based on earlier studies
7af539c6 3396 profzabs->SetDirectory(0);
3397 fDeltaZresAbs.AddAt(profzabs,id);
108953e9 3398 //profz3=new TH2F(Form("prz%03d",id),Form("Z Residuals for Laser Beam %03d- Parabolic2",id),160,0,160,100,-0.5,0.5);
3399 //profz3->SetDirectory(0);
3400 //fDeltaZres3.AddAt(profz3,id);
3401 }
3402 }
3403 //
3404 //
3405 for (Int_t id=0; id<336;id++){
3406 TH1F * hisdz = (TH1F*)fDeltaZ.At(id);
3407 TH1F * hisP3 = (TH1F*)fDeltaP3.At(id);
3408 TH1F * hisP4 = (TH1F*)fDeltaP4.At(id);
3409
3410 TH1F * hisdphi = (TH1F*)fDeltaPhi.At(id);
3411 TH1F * hisdphiP = (TH1F*)fDeltaPhiP.At(id);
3412 TH1F * hisSignal = (TH1F*)fSignals.At(id);
3413
3414 if (!hisdz){
3415 hisdz = new TH1F(Form("hisdz%d",id),Form("hisdz%d",id),1000,-10,10);
3416 hisdz->SetDirectory(0);
3417 fDeltaZ.AddAt(hisdz,id);
3418
3419 hisP3 = new TH1F(Form("hisPar3v%d",id),Form("hisPar3v%d",id),400,-0.06,0.06);
3420 hisP3->SetDirectory(0);
3421 fDeltaP3.AddAt(hisP3,id);
3422 //
3423 hisP4 = new TH1F(Form("hisPar4v%d",id),Form("hisPar4v%d",id),200,-0.06,0.06);
3424 hisP4->SetDirectory(0);
3425 fDeltaP4.AddAt(hisP4,id);
be113f6e 3426
108953e9 3427 //
3428 hisdphi = new TH1F(Form("hisdphi%d",id),Form("hisdphi%d",id),1000,-1,1);
3429 hisdphi->SetDirectory(0);
3430 fDeltaPhi.AddAt(hisdphi,id);
3431 //
3432 hisdphiP = new TH1F(Form("hisdphiP%d",id),Form("hisdphiP%d",id),1000,-0.01,0.01);
3433 hisdphiP->SetDirectory(0);
3434 fDeltaPhiP.AddAt(hisdphiP,id);
3435 hisSignal = new TH1F(Form("hisSignal%d",id),Form("hisSignal%d",id),100,0,300);
3436 hisSignal->SetDirectory(0);
3437 fSignals.AddAt(hisSignal,id);
3438 }
3439 }
7af539c6 3440
3441 SetBeamParameters(fBeamOffsetZOuter, fBeamSlopeZOuter, fBeamSectorOuter,2);
3442 SetBeamParameters(fBeamOffsetZInner, fBeamSlopeZInner, fBeamSectorInner,3);
3443 SetBeamParameters(fBeamOffsetYOuter, fBeamSlopeYOuter, fBeamSectorOuter,0);
3444 SetBeamParameters(fBeamOffsetYInner, fBeamSlopeYInner, fBeamSectorInner,1);
d3ce44cb 3445 //
3446 // Make THnSparse
3447 //
3448 // id side rod bundle beam dP0 dP1 dP2 dP3 dP4 ncl dEdx
3449 Int_t binsLaser[12]= {336, //id
3450 2, //side
3451 6, //rod
3452 4, //bundle
3453 7, //beam
3454 300, //dP0
3455 300, //dP1
3456 300, //dP2
3457 300, //dP3
3458 300, //dP4
3459 80, //ncl
3460 50}; //dEdx
3461 Double_t xminLaser[12]= {0, //id
3462 0, //side
3463 0, //rod
3464 0, //bundle
3465 0, //beam
3466 -1, //dP0
3467 -1, //dP1
3468 -0.01, //dP2
3469 -0.01, //dP3
3470 -0.1, //dP4
3471 0, //ncl
3472 0}; //sqrt dEdx
3473 Double_t xmaxLaser[12]= {336, //id
3474 2, //side
3475 6, //rod
3476 4, //bundle
3477 7, //beam
3478 1, //dP0
3479 1, //dP1
3480 0.01, //dP2
3481 0.01, //dP3
3482 0.1, //dP4
3483 160, //ncl
3484 40}; //sqrt dEdx
3485
3486 TString nameLaser[12]= {"id",
3487 "side",
3488 "rod",
3489 "bundle",
3490 "beam",
3491 "dP0",
3492 "dP1",
3493 "dP2",
3494 "dP3",
3495 "dP4",
3496 "ncl",
3497 "sqrt dEdx"};
3498 TString titleLaser[12]= {"id",
3499 "side",
3500 "rod",
3501 "bundle",
3502 "beam",
3503 "#Delta_{P0}",
3504 "#Delta_{P1}",
3505 "#Delta_{P2}",
3506 "#Delta_{P3}",
3507 "#Delta_{P4}",
3508 "N_{cl}",
3509 "#sqrt{dEdx}"};
3510 fHisLaser = new THnSparseS("dLaser","#Delta_{Laser}", 12, binsLaser,xminLaser, xmaxLaser);
3511 for (Int_t iaxis=1; iaxis<12; iaxis++){
3512 fHisLaser->GetAxis(iaxis)->SetName(nameLaser[iaxis]);
3513 fHisLaser->GetAxis(iaxis)->SetTitle(titleLaser[iaxis]);
3514 }
be113f6e 3515}
3516
3517void AliTPCcalibLaser::MergeFitHistos(AliTPCcalibLaser * laser){
3518 //
3519 // Merge content of histograms
3520 //
3521 // Only first histogram is checked - all other should be the same
d3ce44cb 3522 if (fHisLaser &&laser->fHisLaser) fHisLaser->Add(laser->fHisLaser);
3523
be113f6e 3524 if (!laser->fHisNclIn) return; // empty histograms
3525 if (!fHisNclIn) MakeFitHistos();
3526 //
be113f6e 3527 fHisNclIn->Add(laser->fHisNclIn ); //->Number of clusters inner
3528 fHisNclOut->Add(laser->fHisNclOut ); //->Number of clusters outer
3529 fHisNclIO->Add(laser->fHisNclIO ); //->Number of cluster inner outer
3530 fHisLclIn->Add(laser->fHisLclIn ); //->Level arm inner
3531 fHisLclOut->Add(laser->fHisLclOut ); //->Level arm outer
3532 fHisLclIO->Add(laser->fHisLclIO ); //->Number of cluster inner outer
3533 fHisdEdx->Add(laser->fHisdEdx ); //->dedx
3534 fHisdZfit->Add(laser->fHisdZfit ); //->dedx
3535 //
3536 //
3537 fHisChi2YIn1->Add(laser->fHisChi2YIn1 ); //->chi2 y inner - line
3538 fHisChi2YOut1->Add(laser->fHisChi2YOut1 ); //->chi2 y inner - line
3539 fHisChi2YIn2->Add(laser->fHisChi2YIn2 ); //->chi2 y inner - parabola
3540 fHisChi2YOut2->Add(laser->fHisChi2YOut2 ); //->chi2 y inner - parabola
3541 fHisChi2YIO1->Add(laser->fHisChi2YIO1 ); //->chi2 y IO - common
3542 fHisChi2ZIn1->Add(laser->fHisChi2ZIn1 ); //->chi2 z inner - line
3543 fHisChi2ZOut1->Add(laser->fHisChi2ZOut1 ); //->chi2 z inner - line
3544 fHisChi2ZIn2->Add(laser->fHisChi2ZIn2 ); //->chi2 z inner - parabola
3545 fHisChi2ZOut2->Add(laser->fHisChi2ZOut2 ); //->chi2 z inner - parabola
3546 fHisChi2ZIO1->Add(laser->fHisChi2ZIO1 ); //->chi2 z IO - common
3547 //
3548 //
3549 fHisPy1vP0->Add(laser->fHisPy1vP0 ); //-> delta y P0outer-P0inner - line
3550 fHisPy2vP0->Add(laser->fHisPy2vP0 ); //-> delta y P0outer-P0inner - parabola
3551 fHisPy3vP0->Add(laser->fHisPy3vP0 ); //-> delta y P0outer-P0inner - common parabola
3552 fHisPy1vP1->Add(laser->fHisPy1vP1 ); //-> delta ky P1outer-P1inner - line
3553 fHisPy2vP1->Add(laser->fHisPy2vP1 ); //-> delta ky P1outer-P1inner - parabola
3554 fHisPy3vP1->Add(laser->fHisPy3vP1 ); //-> delta ky P1outer-P1inner - common parabola
3555 fHisPy2vP2In->Add(laser-> fHisPy2vP2In ); //-> Curv P2inner - parabola
3556 fHisPy2vP2Out->Add(laser->fHisPy2vP2Out ); //-> Curv P2outer - parabola
3557 fHisPy3vP2IO->Add(laser->fHisPy3vP2IO ); //-> Curv P2outerinner - common parabola
3558 //
3559 //
3560 fHisPz1vP0->Add(laser->fHisPz1vP0 ); //-> delta z P0outer-P0inner - line
3561 fHisPz2vP0->Add(laser->fHisPz2vP0 ); //-> delta z P0outer-P0inner - parabola
3562 fHisPz3vP0->Add(laser->fHisPz3vP0 ); //-> delta z P0outer-P0inner - common parabola
3563 fHisPz1vP1->Add(laser->fHisPz1vP1 ); //-> delta kz P1outer-P1inner - line
3564 fHisPz2vP1->Add(laser->fHisPz2vP1 ); //-> delta kz P1outer-P1inner - parabola
3565 fHisPz3vP1->Add(laser->fHisPz3vP1 ); //-> delta kz P1outer-P1inner - common parabola
3566 fHisPz2vP2In->Add(laser->fHisPz2vP2In ); //-> Curv P2inner - parabola
3567 fHisPz2vP2Out->Add(laser->fHisPz2vP2Out ); //-> Curv P2outer - parabola
3568 fHisPz3vP2IO->Add(laser->fHisPz3vP2IO ); //-> Curv P2outerinner - common parabola
41a1f702 3569 fHisYAbsErrors->Add(laser->fHisYAbsErrors); //-> total errors per beam in the abs res y analysis
3570 fHisZAbsErrors->Add(laser->fHisZAbsErrors); //-> total errors per beam in the abs res z analysis
be113f6e 3571 //
3572 //
3573 //
3574
3575
3576
3577
3578}
3579
3580
3581
3582
e5c59218 3583void AliTPCcalibLaser::DumpFitInfo(TTree * chainFit,Int_t id){
3584 //
3585 // Dump fit information - collect information from the streamers
3586 //
3587 /*
3588 TChain * chainFit=0;
3589 TChain * chainTrack=0;
3590 TChain * chain=0;
3591 //
3592 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
3593 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+");
3594 AliXRDPROOFtoolkit tool;
3595 chainTrack = tool.MakeChain("laser.txt","Track",0,10200);
3596 chainTrack->Lookup();
3597 chainTrack->SetProof(kTRUE);
be113f6e 3598 chainDrift = tool.MakeChain("laser.txt","driftv",0,10200);
3599 chainDrift->Lookup();
3600 chainDrift->SetProof(kTRUE);
3601
e5c59218 3602 chain = tool.MakeChain("laser.txt","Residuals",0,10200);
3603 chain->Lookup();
3604 chainFit = tool.MakeChain("laser.txt","FitModels",0,10200);
3605 chainFit->Lookup();
3606 chainFit->SetProof(kTRUE);
3607 chain->SetProof(kTRUE);
3608 AliTPCLaserTrack::LoadTracks();
3609 //AliTPCcalibLaser::DumpFitInfo(chainFit,0);
3610
3611 */
3612 //
3613 // Fit cuts
3614 //
3615 TCut cutP3("abs(Tr.fP[3])<0.1");
3616 TCut cutP4("abs(Tr.fP[4])<0.5");
3617 TCut cutPx = cutP3+cutP4;
3618 TCut cutChi2YOut("sqrt(chi2y2Out*dEdx)<5&&chi2y2Out>0");
3619 TCut cutChi2ZOut("sqrt(chi2z2Out*dEdx)<5&&chi2z2Out>0");
3620 TCut cutChi2YIn("sqrt(chi2y2In*dEdx)<5&&chi2y2In>0");
3621 TCut cutChi2ZIn("sqrt(chi2z2In*dEdx)<5&&chi2z2In>0");
3622 //
3623 TCut cutdEdx("sqrt(dEdx)>3");
3624 TCut cutDY("abs(yPol2In.fElements[2]*nclI*nclI/4.)<3");
3625 TCut cutN("nclO>20&&nclI>20");
be113f6e 3626 TCut cutA = cutChi2YOut+cutChi2ZOut+cutChi2YIn+cutChi2ZIn+cutN+cutdEdx+cutPx+"accept";
e5c59218 3627 //
3628 // Cluster cuts
3629 //
3630 TCut cutClY("abs(Cl[].fY-TrYpol2.fElements)<0.15");
3631 TCut cutClZ("abs(Cl[].fZ-TrZpol2.fElements)<0.15");
3632 TCut cutClX("abs(Cl[].fX)>10");
3633 TCut cutE("abs(Cl[].fY/Cl[].fX)<0.14");
3634 TCut cutSY("sqrt(Cl[].fSigmaY2)>0.05");
3635 TCut cutSZ("sqrt(Cl[].fSigmaZ2)>0.05");
3636 TCut cutQ("sqrt(Cl[].fMax)>4");
3637 TCut cutCl=cutClY+cutClZ+cutClX+cutE+cutSY+cutSZ+cutQ;
3638
3639
3640 TH1F * phisAl = 0;
3641 TH1F * phisAccept = 0;
3642 TH1F * phisOut = 0;
3643 TProfile * pdEdx = 0;
3644
3645 TProfile * pP0 = 0;
3646 TProfile * pP1 = 0;
3647 TProfile * pP2 = 0;
3648 TProfile * pP3 = 0;
3649 TProfile * pP4 = 0;
3650 //
3651 TProfile * pNclI = 0;
3652 TProfile * pNclO = 0;
3653 //
3654 TProfile * pchi2YIn =0;
3655 TProfile * pchi2ZIn =0;
3656 TProfile * pchi2YOut =0;
3657 TProfile * pchi2ZOut =0;
3658 TProfile * pchi2YInOut =0;
3659 TProfile * pchi2ZInOut =0;;
3660 // laser counters
3661 chainFit->Draw("LTr.fId>>hisAl(350,0,350)","LTr.fId<350");
3662 phisAl = (TH1F*)gROOT->FindObject("hisAl");
3663 chainFit->Draw("LTr.fId>>hisAccept(350,0,350)","LTr.fId<350"+cutA);
3664 phisAccept = (TH1F*)gROOT->FindObject("hisAccept");
3665 chainFit->Draw("LTr.fId>>hisOut(350,0,350)","LTr.fId<350"+!cutA);
3666 phisOut = (TH1F*)gROOT->FindObject("hisOut");
3667 //
3668 chainFit->Draw("sqrt(dEdx):LTr.fId>>hdEdx(350,0,350)","","prof");
3669 pdEdx = (TProfile*)gROOT->FindObject("hdEdx");
3670 // track param
3671 //
3672 chainFit->Draw("Tr.fP[0]:LTr.fId>>hP0(350,0,350)","Tr.fP[4]/sqrt(Tr.fC[14])<20"+cutA,"prof");
3673 pP0 = (TProfile*)gROOT->FindObject("hP0");
3674 chainFit->Draw("Tr.fP[1]:LTr.fId>>hP1(350,0,350)","Tr.fP[4]/sqrt(Tr.fC[14])<20"+cutA,"prof");
3675 pP1 = (TProfile*)gROOT->FindObject("hP1");
3676 chainFit->Draw("Tr.fP[2]:LTr.fId>>hP2(350,0,350)","Tr.fP[4]/sqrt(Tr.fC[14])<20"+cutA,"prof");
3677 pP2 = (TProfile*)gROOT->FindObject("hP2");
3678 chainFit->Draw("Tr.fP[3]:LTr.fId>>hP3(350,0,350)","Tr.fP[4]/sqrt(Tr.fC[14])<20"+cutA,"prof");
3679 pP3 = (TProfile*)gROOT->FindObject("hP3");
3680 chainFit->Draw("Tr.fP[4]:LTr.fId>>hP4(350,0,350)","Tr.fP[4]/sqrt(Tr.fC[14])<20"+cutA,"prof");
3681 pP4 = (TProfile*)gROOT->FindObject("hP4");
3682 //
3683 chainFit->Draw("nclI:LTr.fId>>hNclI(350,0,350)","","prof");
3684 pNclI = (TProfile*)gROOT->FindObject("hNclI");
3685 chainFit->Draw("nclO:LTr.fId>>hNclO(350,0,350)","","prof");
3686 pNclO = (TProfile*)gROOT->FindObject("hNclO");
3687 //
3688 //
3689 chainFit->Draw("sqrt(chi2y2In):LTr.fId>>hChi2YIn(350,0,350)",cutA+"","prof");
3690 pchi2YIn = (TProfile*)gROOT->FindObject("hChi2YIn");
3691 chainFit->Draw("sqrt(chi2y2Out):LTr.fId>>hChi2YOut(350,0,350)",cutA+"","prof");
3692 pchi2YOut = (TProfile*)gROOT->FindObject("hChi2YOut");
3693 chainFit->Draw("sqrt(chi2yInOut):LTr.fId>>hChi2YInOut(350,0,350)",cutA+"","prof");
3694 pchi2YInOut = (TProfile*)gROOT->FindObject("hChi2YInOut");
3695 chainFit->Draw("sqrt(chi2z2In):LTr.fId>>hChi2ZIn(350,0,350)",cutA+"","prof");
3696 pchi2ZIn = (TProfile*)gROOT->FindObject("hChi2ZIn");
3697 chainFit->Draw("sqrt(chi2z2Out):LTr.fId>>hChi2ZOut(350,0,350)",cutA+"","prof");
3698 pchi2ZOut = (TProfile*)gROOT->FindObject("hChi2ZOut");
3699 chainFit->Draw("sqrt(chi2zInOut):LTr.fId>>hChi2ZInOut(350,0,350)",cutA+"","prof");
3700 pchi2ZInOut = (TProfile*)gROOT->FindObject("hChi2ZInOut");
3701 //
3702 // second derivatives
3703 //
3704 TH2F * phisPy2In = new TH2F("Py2Inner","Py2Inner",350,0,350,100,-0.001,0.001);
3705 chainFit->Draw("yPol2In.fElements[2]:LTr.fId>>Py2Inner",cutA,"");
3706 TH2F * phisPy2Out = new TH2F("Py2Outer","Py2Outer",350,0,350,200,-0.0005,0.0005);
3707 chainFit->Draw("yPol2Out.fElements[2]:LTr.fId>>Py2Outer",cutA,"");
3708 TH2F * phisPy2InOut = new TH2F("Py2InOut","Py2InOut",350,0,350,200,-0.0005,0.0005);
3709 chainFit->Draw("yInOut.fElements[4]:LTr.fId>>Py2InOut",cutA,"");
3710 //
3711 phisPy2In->FitSlicesY();
3712 TH1D * phisPy2InEntries = (TH1D*)gROOT->FindObject("Py2Inner_0");
3713 TH1D * phisPy2InMean = (TH1D*)gROOT->FindObject("Py2Inner_1");
3714 TH1D * phisPy2InSigma = (TH1D*)gROOT->FindObject("Py2Inner_2");
3715 //
3716 phisPy2Out->FitSlicesY();
3717 TH1D * phisPy2OutEntries = (TH1D*)gROOT->FindObject("Py2Outer_0");
3718 TH1D * phisPy2OutMean = (TH1D*)gROOT->FindObject("Py2Outer_1");
3719 TH1D * phisPy2OutSigma = (TH1D*)gROOT->FindObject("Py2Outer_2");
3720 //
3721 phisPy2InOut->FitSlicesY();
3722 TH1D * phisPy2InOutEntries = (TH1D*)gROOT->FindObject("Py2InOut_0");
3723 TH1D * phisPy2InOutMean = (TH1D*)gROOT->FindObject("Py2InOut_1");
3724 TH1D * phisPy2InOutSigma = (TH1D*)gROOT->FindObject("Py2InOut_2");
3725 //
3726 TH2F * phisPz2In = new TH2F("Pz2Inner","Pz2Inner",350,0,350,100,-0.001,0.001);
3727 chainFit->Draw("zPol2In.fElements[2]:LTr.fId>>Pz2Inner",cutA,"");
3728 TH2F * phisPz2Out = new TH2F("Pz2Outer","Pz2Outer",350,0,350,200,-0.0005,0.0005);
3729 chainFit->Draw("zPol2Out.fElements[2]:LTr.fId>>Pz2Outer",cutA,"");
3730 TH2F * phisPz2InOut = new TH2F("Pz2InOut","Pz2InOut",350,0,350,200,-0.0005,0.0005);
3731 chainFit->Draw("zInOut.fElements[4]:LTr.fId>>Pz2InOut",cutA,"");
3732 //
3733 phisPz2In->FitSlicesY();
3734 TH1D * phisPz2InEntries = (TH1D*)gROOT->FindObject("Pz2Inner_0");
3735 TH1D * phisPz2InMean = (TH1D*)gROOT->FindObject("Pz2Inner_1");
3736 TH1D * phisPz2InSigma = (TH1D*)gROOT->FindObject("Pz2Inner_2");
3737 //
3738 phisPz2Out->FitSlicesY();
3739 TH1D * phisPz2OutEntries = (TH1D*)gROOT->FindObject("Pz2Outer_0");
3740 TH1D * phisPz2OutMean = (TH1D*)gROOT->FindObject("Pz2Outer_1");
3741 TH1D * phisPz2OutSigma = (TH1D*)gROOT->FindObject("Pz2Outer_2");
3742 //
3743 phisPz2InOut->FitSlicesY();
3744 TH1D * phisPz2InOutEntries = (TH1D*)gROOT->FindObject("Pz2InOut_0");
3745 TH1D * phisPz2InOutMean = (TH1D*)gROOT->FindObject("Pz2InOut_1");
3746 TH1D * phisPz2InOutSigma = (TH1D*)gROOT->FindObject("Pz2InOut_2");
3747 //
3748 //
3749 //
3750
3751
3752 {
3753 TTreeSRedirector *pcstream = new TTreeSRedirector("vscan.root");
3754 for (Int_t ilaser=0; ilaser<336; ilaser++){
3755 Float_t all=phisAl->GetBinContent(ilaser+1);
3756 Float_t accept=phisAccept->GetBinContent(ilaser+1);
3757 Float_t out=phisOut->GetBinContent(ilaser+1);
3758 Float_t sdedx = pdEdx->GetBinContent(ilaser+1);
3759 Float_t mP0 = pP0->GetBinContent(ilaser+1);
3760 Float_t mP1 = pP1->GetBinContent(ilaser+1);
3761 Float_t mP2 = pP2->GetBinContent(ilaser+1);
3762 Float_t mP3 = pP3->GetBinContent(ilaser+1);
3763 Float_t mP4 = pP4->GetBinContent(ilaser+1);
3764
3765
3766 Float_t nclI = pNclI->GetBinContent(ilaser+1);
3767 Float_t nclO = pNclO->GetBinContent(ilaser+1);
3768 //
3769 Float_t chi2YIn=pchi2YIn->GetBinContent(ilaser+1);
3770 Float_t chi2YOut=pchi2YOut->GetBinContent(ilaser+1);
3771 Float_t chi2YInOut=pchi2YInOut->GetBinContent(ilaser+1);
3772 Float_t chi2ZIn=pchi2ZIn->GetBinContent(ilaser+1);
3773 Float_t chi2ZOut=pchi2ZOut->GetBinContent(ilaser+1);
3774 Float_t chi2ZInOut=pchi2ZInOut->GetBinContent(ilaser+1);
3775 //
3776 Float_t entriesPy2In = phisPy2InEntries->GetBinContent(ilaser+1);
3777 Float_t meanPy2In = phisPy2InMean->GetBinContent(ilaser+1);
3778 Float_t sigmaPy2In = phisPy2InSigma->GetBinContent(ilaser+1);
3779 //
3780 Float_t entriesPy2Out = phisPy2OutEntries->GetBinContent(ilaser+1);
3781 Float_t meanPy2Out = phisPy2OutMean->GetBinContent(ilaser+1);
3782 Float_t sigmaPy2Out = phisPy2OutSigma->GetBinContent(ilaser+1);
3783 //
3784 Float_t entriesPy2InOut = phisPy2InOutEntries->GetBinContent(ilaser+1);
3785 Float_t meanPy2InOut = phisPy2InOutMean->GetBinContent(ilaser+1);
3786 Float_t sigmaPy2InOut = phisPy2InOutSigma->GetBinContent(ilaser+1);
3787 //
3788 Float_t entriesPz2In = phisPz2InEntries->GetBinContent(ilaser+1);
3789 Float_t meanPz2In = phisPz2InMean->GetBinContent(ilaser+1);
3790 Float_t sigmaPz2In = phisPz2InSigma->GetBinContent(ilaser+1);
3791 //
3792 Float_t entriesPz2Out = phisPz2OutEntries->GetBinContent(ilaser+1);
3793 Float_t meanPz2Out = phisPz2OutMean->GetBinContent(ilaser+1);
3794 Float_t sigmaPz2Out = phisPz2OutSigma->GetBinContent(ilaser+1);
3795 //
3796 Float_t entriesPz2InOut = phisPz2InOutEntries->GetBinContent(ilaser+1);
3797 Float_t meanPz2InOut = phisPz2InOutMean->GetBinContent(ilaser+1);
3798 Float_t sigmaPz2InOut = phisPz2InOutSigma->GetBinContent(ilaser+1);
3799
3800 AliTPCLaserTrack* ltrp =(AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(ilaser);
3801 (*pcstream)<<"Scan"<<
3802 "Run="<<id<<
3803 "LTr.="<<ltrp<<
3804 "all="<<all<<
3805 "accept="<<accept<<
3806 "out="<<out<<
3807 "sdedx="<<sdedx<<
3808 "mP0="<<mP0<<
3809 "mP1="<<mP1<<
3810 "mP2="<<mP2<<
3811 "mP3="<<mP3<<
3812 "mP4="<<mP4<<
3813 "nclI="<<nclI<<
3814 "nclO="<<nclO<<
3815 "chi2YIn="<<chi2YIn<<
3816 "chi2YOut="<<chi2YOut<<
3817 "chi2YInOut="<<chi2YInOut<<
3818 "chi2ZIn="<<chi2ZIn<<
3819 "chi2ZOut="<<chi2ZOut<<
3820 "chi2ZInOut="<<chi2ZInOut<<
3821 //
3822 "nPy2In="<<entriesPy2In<<
3823 "mPy2In="<<meanPy2In<<
3824 "sPy2In="<<sigmaPy2In<<
3825 //
3826 "nPy2Out="<<entriesPy2Out<<
3827 "mPy2Out="<<meanPy2Out<<
3828 "sPy2Out="<<sigmaPy2Out<<
3829 //
3830 "nPy2InOut="<<entriesPy2InOut<<
3831 "mPy2InOut="<<meanPy2InOut<<
3832 "sPy2InOut="<<sigmaPy2InOut<<
3833 //
3834 "nPz2In="<<entriesPz2In<<
3835 "mPz2In="<<meanPz2In<<
3836 "sPz2In="<<sigmaPz2In<<
3837 //
3838 "nPz2Out="<<entriesPz2Out<<
3839 "mPz2Out="<<meanPz2Out<<
3840 "sPz2Out="<<sigmaPz2Out<<
3841 //
3842 "nPz2InOut="<<entriesPz2InOut<<
3843 "mPz2InOut="<<meanPz2InOut<<
3844 "sPz2InOut="<<sigmaPz2InOut<<
3845 "\n";
3846 }
3847
3848 delete pcstream;
3849 }
7af539c6 3850}
3851
3852void AliTPCcalibLaser::SetBeamParameters(TVectorD& meanOffset,
3853 TVectorD& meanSlope,
3854 TVectorD& sectorArray,
3855 Int_t option)
3856{
3857 // This method should ideally go in AliTPCLaser
3858 // option == 0 (pads outer - closest to beam)
3859 // option == 1 (pads inner)
3860 // option == 2 (time outer)
3861 // option == 3 (time inner)
3862 Int_t nFailures = 0;
3863
3864 for(Int_t id = 0; id < 336; id++) {
3865
3866 if (!AliTPCLaserTrack::GetTracks())
3867 AliTPCLaserTrack::LoadTracks();
3868 AliTPCLaserTrack *ltrp =
3869 (AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(id);
3870
3871 AliExternalTrackParam trackParam(*ltrp);
3872
3873 Double_t deltaangle = 10.0; // sector is 1/2 a sector away from mirror
3874 if((option==1 || option==3)&& (ltrp->GetBeam()<=1 || ltrp->GetBeam()>=5))
3875 deltaangle = 30.0; // inner sector is 1 sector further away from mirror
3876
3877 Double_t angle = trackParam.GetAlpha();
3878 if(angle<0)
3879 angle += 2*TMath::Pi();
3880 if(trackParam.GetSnp()>0) // track points to sector "before"
3881 angle -= deltaangle*TMath::DegToRad();
3882 else // track points to sector "after"
3883 angle += deltaangle*TMath::DegToRad();
3884
3885 Bool_t success = trackParam.Rotate(angle);
3886
3887 if(!success) {
3888 // cout << "WARNING: Rotate failed for ID: " << id << endl;
3889 nFailures++;
3890 }
3891
3892 angle *= TMath::RadToDeg();
3893
3894 Int_t sector = TMath::Nint((angle-10.0)/20.0);
3895 if(sector<0)
3896 sector += 18;
3897 else if(sector>=18)
3898 sector -= 18;
3899 if(ltrp->GetSide()==1) // C side
3900 sector += 18;
f6220766 3901 if(option==0 || option==2)
3902 sector += 36;
3903 if((option==1||option==3)&&(ltrp->GetBeam()==0||ltrp->GetBeam()==6))
3904 sector += 36;
3905
7af539c6 3906 sectorArray[id] = sector;
3907
3908 const Double_t x0 = 0;
3909
3910 Double_t slopey = TMath::Tan(TMath::ASin(trackParam.GetSnp()));
3911 Double_t slopez = trackParam.GetTgl();
3912 // One needs a factor sqrt(1+slopey**2) to take into account the
3913 // longer path length
3914 slopez *= TMath::Sqrt(1.0 + slopey*slopey);
3915 if(fInverseSlopeZ) // wrong sign in database, should be fixed there
3916 slopez *= -1;
3917 // Double_t offsetz = trackParam.GetZ();
3918 Double_t offsety = trackParam.GetY() + slopey*(x0-trackParam.GetX());
3919 Double_t offsetz = trackParam.GetZ() + slopez*(x0-trackParam.GetX());
3920 if(option==2 || option==3) {
3921 meanOffset[id] = offsetz; meanSlope[id] = slopez;
3922 } else {
3923 meanOffset[id] = offsety; meanSlope[id] = slopey;
3924 }
3925 }
3926
3927 if(nFailures>0)
3928 AliWarning(Form("Rotate method failed %d times", nFailures));
3929}
3930
3931
e5c59218 3932 /*
3933 TFile f("vscan.root");
3934 */
3935
3936 /*
3937 pad binning effect
3938 chain->Draw("Cl[].fPad-int(Cl[].fPad)",cutA+cutCl+"Cl[].fZ>0&&Cl[].fPad>1","",10000);
3939 //
3940 chain->Draw("Cl[].fY-TrYpol1.fElements:Cl[].fPad-int(Cl[].fPad)",cutA+cutCl+"Cl[].fZ>0&&Cl[].fPad>1","prof",10000);
3941 //
3942
3943chain->Draw("Cl.fY-TrYpol1.fElements-AliTPCClusterParam::SPosCorrection(0,1,Cl[].fPad,Cl[].fTimeBin,Cl[].fZ,Cl[].fSigmaY2,Cl[].fSigmaZ2,Cl[].fMax):Cl[].fPad-int(Cl[].fPad)",cutA+cutCl+"Cl[].fZ>0&&Cl[].fPad>1","prof",10000);
3944
3945
3946chain->Draw("Cl[].fZ-TrZpol1.fElements-0*AliTPCClusterParam::SPosCorrection(1,1,Cl[].fPad,Cl[].fTimeBin,Cl[].fZ,Cl[].fSigmaY2,Cl[].fSigmaZ2,Cl[].fMax):Cl[].fTimeBin-int(Cl[].fTimeBin)",cutA+cutCl+"Cl[].fZ>0","prof",10000)
3947
3948 */
3949
3950
3951
3952
3953
3954 /*
3955 // check edge effects
3956 chain->Draw("Cl.fY-TrYpol1.fElements:Cl.fY/Cl.fX",""+cutA+cutCl,"prof",10000)
3957 //
3958 chain->Draw("Cl.fY-TrYpol2.fElements:Cl.fPad-int(Cl.fPad)","Cl.fZ>0"+cutA+cutCl+cutE,"prof",100000)
3959
3960 chain->Draw("Cl.fY-TrYpol2.fElements:Cl.fPad-int(Cl.fPad)","Cl.fX>80&&Cl.fZ>0&&Cl.fDetector>35"+cutA+cutCl+cutE,"prof",100000)
3961
3962
3963
3964 chainFit->Draw("yInOut.fElements[4]:LTr.fP[2]","LTr.fP[1]<0"+cutA,"prof",1000);
3965 chainFit->Draw("yPol2In.fElements[2]*90*90/4.:LTr.fP[2]","nclO>40&&LTr.fP[1]<0"+cutA+cutD,"prof")
3966
3967*/
3968
3969
3970
3971
3972
3973 /*
3974 Edge y effect
3975
3976 dedge = sign(Cl.fY)*(Cl.fX*tan(pi/18)-abs(Cl.fY))
3977
3978
3979 chain->Draw("sign(Cl.fY)*(Cl.fY-TrYpol1.fElements):pi/18-abs(Cl.fY/Cl.fX)>>hisYdphi(100,0,0.03)",""+cutA+cutCl,"prof",10000)
3980
3981 chain->Draw("sign(Cl.fY)*(Cl.fY-TrYpol1.fElements):Cl.fX*(pi/18-abs(Cl.fY/Cl.fX))>>hisYdy(100,0,5)",""+cutA+cutCl,"prof",10000)
3982
3983
3984
3985
3986
3987 chain->Draw("sign(Cl.fY)*(Cl.fY-TrYpol1.fElements):Cl.fX*(pi/18-abs(Cl.fY/Cl.fX))>>hisYdyIROC(100,0,5)","Cl.fDetector<36"+cutA+cutCl,"prof",100000)
3988
3989 chain->Draw("sign(Cl.fY)*(Cl.fY-TrYpol1.fElements):Cl.fX*(pi/18-abs(Cl.fY/Cl.fX))>>hisYdyOROC(100,0,5)","Cl.fDetector>36"+cutA+cutCl,"prof",100000)
3990
3991
3992
3993 chain->Draw("Cl.fY-TrYpol1.fElements:sign(Cl.fY)*(Cl.fX*tan(pi/18)-abs(Cl.fY))>>his(100,-5,5)",""+cutA+cutCl,"prof",100000)
3994
3995 chain->Draw("Cl.fY-TrYpol1.fElements:sign(Cl.fY)*(Cl.fX*tan(pi/18)-abs(Cl.fY))>>hisdyInner(100,-5,5)","Cl.fDetector<36"+cutA+cutCl,"prof",100000)
3996
3997
3998
3999*/
4000
4001
4002/*
4003
4004chainFit->Draw("yPol2Out.fElements[2]*90*90/4.:LTr.fP[2]","nclO>40&&LTr.fP[1]<0"+cutA+cutDY,"prof")
4005
4006chainFit->Draw("yPol2In.fElements[2]*64*64/4.:LTr.fP[2]","nclI>20&&LTr.fP[1]<0"+cutA+cutDY,"prof")
4007
4008
4009
4010chainFit->Draw("LTr.fId","nclI>10",100000)
4011
4012chainFit->Draw("yPol2In.fElements[2]:LTr.fId>>his(350,0,350,100,-0.002,0.002)","nclI>20","")
4013
4014chainFit->Draw("yPol2In.fElements[2]:LTr.fId>>hisPy2In0(350,0,350,100,-0.002,0.002)","nclI>20","");
4015
4016TH2 * phisPy2In = (TH2*) gROOT->FindObject("hisPy2In0")
4017
4018*/
4019
e5c59218 4020
4021
c03e3250 4022
4023
263d466a 4024
c03e3250 4025/*
4026 gSystem->Load("libSTAT.so")
4027 TStatToolkit toolkit;
4028 Double_t chi2;
4029 TVectorD fitParam;
4030 TMatrixD covMatrix;
4031 Int_t npoints;
e5c59218 4032
ae69f16f 4033 TCut cutA("entries>2&&pphi2<3&&abs(gphiP1-pphiP0)<0.003&&abs(gz1)<6");
c03e3250 4034
4035
4036TString fstring="";
1fd56785 4037//
4038fstring+="(abs(LTr.fP[1]/250)^3-1)*bz++"; //1
4039fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*LTr.fP[2]++"; //2
4040fstring+="(abs(LTr.fP[1]/250)^1-1)*bz++"; //3
4041fstring+="(abs(LTr.fP[1]/250)-1)*bz*LTr.fP[2]++"; //4
4042//
4043fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*sin(atan2(lx1,lx0))++" //5
4044fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*sin(atan2(lx1,lx0))*LTr.fP[2]++" //6
4045fstring+="(abs(LTr.fP[1]/250)-1)*bz*sin(atan2(lx1,lx0))++" //7
4046fstring+="(abs(LTr.fP[1]/250)-1)*bz*sin(atan2(lx1,lx0))*LTr.fP[2]++" //8
4047//
4048fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*cos(atan2(lx1,lx0))++" //9
4049fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*cos(atan2(lx1,lx0))*LTr.fP[2]++" //10
4050fstring+="(abs(LTr.fP[1]/250)-1)*bz*cos(atan2(lx1,lx0))++" //11
4051fstring+="(abs(LTr.fP[1]/250)-1)*bz*cos(atan2(lx1,lx0))*LTr.fP[2]++" //12
c5595838 4052
c5595838 4053
c5595838 4054
c5595838 4055
1fd56785 4056 TString *strq0 = toolkit.FitPlane(treeT,"gphi1-pphi0",fstring->Data(), "fSide==1"+cutA, chi2,npoints,fitParam,covMatrix);
c03e3250 4057
1fd56785 4058 treeT->SetAlias("fit",strq0->Data());
4059
c03e3250 4060
1fd56785 4061 TString *strqP = toolkit.FitPlane(treeT,"1000*(gphiP1-pphiP0)",fstring->Data(), "fSide==1"+cutA, chi2,npoints,fitParam,covMatrix);
c5595838 4062
1fd56785 4063 treeT->SetAlias("fitP",strqP->Data());
c5595838 4064
4065
ae69f16f 4066 TString *strqDrift = toolkit.FitPlane(treeT,"gz1","LTr.fP[1]++(1-2*(fSide==1))++lx1", cutA, chi2,npoints,fitParam,covMatrix);
4067 treeT->SetAlias("fitD",strqDrift->Data());
c5595838 4068
c03e3250 4069
ae69f16f 4070treeT->Draw("fit:LTr.fP[1]","abs(bz+0.4)<0.05"+cutA,"");
4071{
4072for (Int_t i=0; i<6;i++){
4073treeT->SetLineColor(i+2);
4074treeT->SetMarkerSize(1);
4075treeT->SetMarkerStyle(22+i);
4076treeT->SetMarkerColor(i+2);
c03e3250 4077
ae69f16f 4078treeT->Draw("fit:LTr.fP[1]",Form("abs(bz+0.4)<0.05&fRod==%d",i)+cutA,"same");
4079}
4080}
c03e3250 4081 */
68ff0583 4082
4083
4084
4085/*
4086 TTree * tree = (TTree*)f.Get("FitModels");
4087
4088 TEventList listLFit0("listLFit0","listLFit0");
4089 TEventList listLFit1("listLFit1","listLFit1");
68ff0583 4090 tree->Draw(">>listLFit0","seed.fdEdx<200&&seed.fdEdx>40");
4091 tree->SetEventList(&listLFit0);
4092
4093
95a0e09b 4094
4095
4096 gSystem->Load("libSTAT.so")
4097 TStatToolkit toolkit;
4098 Double_t chi2;
4099 TVectorD fitParam;
4100 TMatrixD covMatrix;
4101 Int_t npoints;
4102
4103 chain->SetAlias("dp","((Cl.fPad-int(Cl.fPad))*pi)");
4104 chain->SetAlias("dt","((Cl.fTimeBin-int(Cl.fTimeBin))*pi)");
4105
4106
4107 TString fstring="";
4108 fstring+="cos(dp)++";
4109 fstring+="sin(dp)++";
4110 fstring+="cos(dt)++";
4111 fstring+="sin(dt)++";
4112
4113 TString *str = toolkit.FitPlane(chain,"Cl.fZ-TrZInOut.fElements",fstring->Data(), "Cl.fDetector>35", chi2,npoints,fitParam,covMatrix,-1,0,200);
4114
4115
e5c59218 4116
4117*/
4118
4119
4120
4121/*
4122 Edge effects
4123 //
4124 //
4125 //
4126 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
4127 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
4128 AliXRDPROOFtoolkit tool;
4129 TChain * chainTrack = tool.MakeChain("laser.txt","Track",0,10200);
4130 chainTrack->Lookup();
4131 chainTrack->SetProof(kTRUE);
4132
4133 TChain * chain = tool.MakeChain("laser.txt","Residuals",0,10200);
4134 chain->Lookup();
4135 TChain * chainFit = tool.MakeChain("laser.txt","FitModels",0,10200);
4136 chainFit->Lookup();
4137 chainFit->SetProof(kTRUE);
4138 chain->SetProof(kTRUE);
4139 //
4140 // Fit cuts
4141 //
4142 TCut cutChi2YOut("sqrt(chi2y2Out*dEdx)<10");
4143 TCut cutChi2ZOut("sqrt(chi2z2Out*dEdx)<10");
4144 TCut cutChi2YIn("sqrt(chi2y2In*dEdx)<10");
4145 TCut cutChi2ZIn("sqrt(chi2z2In*dEdx)<10");
4146 //
4147 TCut cutdEdx("sqrt(dEdx)<30&&sqrt(dEdx)>3");
4148 TCut cutDY("abs(yPol2In.fElements[2]*nclO*nclO/4.)<3")
4149 TCut cutN("nclO>20&&nclI>20");
4150 TCut cutA = cutChi2YOut+cutChi2ZOut+cutChi2YIn+cutChi2ZIn+cutN+cutdEdx;
4151 //
4152 // Cluster cuts
4153 //
4154 TCut cutClY("abs(Cl.fY-TrYpol2.fElements)<0.2");
4155 TCut cutClZ("abs(Cl.fZ-TrZpol2.fElements)<0.4");
4156 TCut cutClX("abs(Cl.fX)>10");
4157 TCut cutE("abs(Cl.fY/Cl.fX)<0.14");
4158 TCut cutCl=cutClY+cutClZ+cutClX;
4159
4160
4161 // check edge effects
4162 chain->Draw("Cl.fY-TrYpol1.fElements:Cl.fY/Cl.fX",""+cutA+cutCl,"prof",10000)
4163 //
4164 chain->Draw("Cl.fY-TrYpol2.fElements:Cl.fPad-int(Cl.fPad)","Cl.fZ>0"+cutA+cutCl+cutE,"prof",100000)
4165
4166 chain->Draw("Cl.fY-TrYpol2.fElements:Cl.fPad-int(Cl.fPad)","Cl.fX>80&&Cl.fZ>0&&Cl.fDetector>35"+cutA+cutCl+cutE,"prof",100000)
4167
4168
4169
4170 chainFit->Draw("yInOut.fElements[4]:LTr.fP[2]","LTr.fP[1]<0"+cutA,"prof",1000);
4171 chainFit->Draw("yPol2In.fElements[2]*90*90/4.:LTr.fP[2]","nclO>40&&LTr.fP[1]<0"+cutA+cutD,"prof")
4172
68ff0583 4173*/
e5c59218 4174