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