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