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