]>
Commit | Line | Data |
---|---|---|
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 | // |
29 | // | |
30 | // | |
1fd56785 | 31 | // To make laser scan the user interaction neccessary |
32 | // | |
33 | .x ~/UliStyle.C | |
2b35e8e6 | 34 | gSystem->Load("libANALYSIS"); |
35 | gSystem->Load("libTPCcalib"); | |
36 | TFile fcalib("CalibObjects.root"); | |
37 | TObjArray * array = (TObjArray*)fcalib.Get("TPCCalib"); | |
38 | AliTPCcalibLaser * laser = ( AliTPCcalibLaser *)array->FindObject("laserTPC"); | |
39 | laser->DumpMeanInfo(-0.4) | |
40 | TFile fmean("laserMean.root") | |
1fd56785 | 41 | // |
2b35e8e6 | 42 | // laser track clasification; |
1fd56785 | 43 | // |
7b18d067 | 44 | TCut cutT("cutT","abs(Tr.fP[3])<0.06"); |
45 | TCut cutPt("cutPt","abs(Tr.fP[4])<0.1"); | |
c18f4385 | 46 | TCut cutN("cutN","fTPCncls>70"); |
47 | TCut cutP("cutP","abs(atan2(x1,x0)-atan2(lx1,lx0))<0.03") | |
48 | TCut cutA = cutT+cutPt+cutP; | |
49 | TFile f("laserTPCDebug.root"); | |
50 | TTree * treeT = (TTree*)f.Get("Track"); | |
1fd56785 | 51 | // |
52 | // | |
53 | // Analyze LASER scan | |
54 | // | |
2b35e8e6 | 55 | gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros"); |
56 | gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+") | |
0728a4f6 | 57 | AliXRDPROOFtoolkit tool; |
2b35e8e6 | 58 | TChain * chain = tool.MakeChain("laserScan.txt","Mean",0,10200); |
59 | chain->Lookup(); | |
1fd56785 | 60 | AliTPCcalibLaser::DumpScanInfo(chain) |
61 | TFile fscan("laserScan.root") | |
62 | TTree * treeT = (TTree*)fscan.Get("Mean") | |
95a0e09b | 63 | // |
64 | // Analyze laser | |
65 | // | |
66 | gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros"); | |
67 | gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+") | |
68 | AliXRDPROOFtoolkit tool; | |
69 | TChain * chain = tool.MakeChain("laser.txt","Residuals",0,10200); | |
70 | chain->Lookup(); | |
71 | ||
c6914c83 | 72 | */ |
73 | ||
74 | ||
75 | ||
76 | #include "TLinearFitter.h" | |
77 | #include "AliTPCcalibLaser.h" | |
78 | #include "AliExternalTrackParam.h" | |
e9f38a4b | 79 | #include "AliESDEvent.h" |
80 | #include "AliESDfriend.h" | |
c6914c83 | 81 | #include "AliESDtrack.h" |
82 | #include "AliTPCTracklet.h" | |
83 | #include "TH1D.h" | |
ecb7e705 | 84 | #include "TH1F.h" |
95a0e09b | 85 | #include "TProfile.h" |
c6914c83 | 86 | #include "TVectorD.h" |
87 | #include "TTreeStream.h" | |
88 | #include "TFile.h" | |
89 | #include "TF1.h" | |
90 | #include "TGraphErrors.h" | |
91 | #include "AliTPCclusterMI.h" | |
92 | #include "AliTPCseed.h" | |
93 | #include "AliTracker.h" | |
263d466a | 94 | #include "AliLog.h" |
c6914c83 | 95 | #include "TClonesArray.h" |
1fd56785 | 96 | #include "TPad.h" |
c6914c83 | 97 | |
98 | ||
99 | #include "TTreeStream.h" | |
100 | #include <iostream> | |
101 | #include <sstream> | |
7b18d067 | 102 | #include "AliTPCLaserTrack.h" |
103 | ||
c6914c83 | 104 | using namespace std; |
105 | ||
106 | ClassImp(AliTPCcalibLaser) | |
107 | ||
108 | AliTPCcalibLaser::AliTPCcalibLaser(): | |
e9f38a4b | 109 | AliTPCcalibBase(), |
110 | fESD(0), | |
111 | fESDfriend(), | |
c18f4385 | 112 | fTracksMirror(336), |
113 | fTracksEsd(336), | |
114 | fTracksEsdParam(336), | |
115 | fTracksTPC(336), | |
263d466a | 116 | fDeltaZ(336), |
117 | fDeltaPhi(336), | |
118 | fDeltaPhiP(336), | |
95a0e09b | 119 | fSignals(336), |
120 | fDeltaYres(336), | |
121 | fDeltaZres(336), | |
ecb7e705 | 122 | fPol2Par2InY(336), |
123 | fDiffPar1InY(336), | |
124 | fPol2Par2OutY(336), | |
125 | fDiffPar1OutY(336), | |
126 | fPol2Par2InZ(336), | |
127 | fDiffPar1InZ(336), | |
128 | fPol2Par2OutZ(336), | |
129 | fDiffPar1OutZ(336), | |
130 | fFitAside(new TVectorD(3)), | |
263d466a | 131 | fFitCside(new TVectorD(3)), |
132 | fEdgeXcuts(5), | |
133 | fEdgeYcuts(5), | |
134 | fNClCuts(5), | |
135 | fNcuts(0), | |
136 | fRun(0), | |
137 | fEvent(0) | |
c6914c83 | 138 | { |
139 | // | |
140 | // Constructor | |
141 | // | |
e9f38a4b | 142 | fTracksEsdParam.SetOwner(kTRUE); |
c6914c83 | 143 | } |
144 | ||
145 | AliTPCcalibLaser::AliTPCcalibLaser(const Text_t *name, const Text_t *title): | |
e9f38a4b | 146 | AliTPCcalibBase(), |
147 | fESD(0), | |
148 | fESDfriend(0), | |
c18f4385 | 149 | fTracksMirror(336), |
150 | fTracksEsd(336), | |
151 | fTracksEsdParam(336), | |
152 | fTracksTPC(336), | |
153 | fDeltaZ(336), // array of histograms of delta z for each track | |
154 | fDeltaPhi(336), // array of histograms of delta z for each track | |
155 | fDeltaPhiP(336), // array of histograms of delta z for each track | |
2b35e8e6 | 156 | fSignals(336), // array of dedx signals |
95a0e09b | 157 | fDeltaYres(336), |
158 | fDeltaZres(336), | |
ecb7e705 | 159 | fPol2Par2InY(336), |
160 | fDiffPar1InY(336), | |
161 | fPol2Par2OutY(336), | |
162 | fDiffPar1OutY(336), | |
163 | fPol2Par2InZ(336), | |
164 | fDiffPar1InZ(336), | |
165 | fPol2Par2OutZ(336), | |
166 | fDiffPar1OutZ(336), | |
c18f4385 | 167 | fFitAside(new TVectorD(3)), // drift fit - A side |
168 | fFitCside(new TVectorD(3)), // drift fit - C- side | |
263d466a | 169 | fEdgeXcuts(5), // cuts in local x direction; used in the refit of the laser tracks |
170 | fEdgeYcuts(5), // cuts in local y direction; used in the refit of the laser tracks | |
171 | fNClCuts(5), // cuts on the number of clusters per tracklet; used in the refit of the laser tracks | |
172 | fNcuts(0), // number of cuts | |
173 | fRun(0), | |
174 | fEvent(0) | |
c6914c83 | 175 | { |
176 | SetName(name); | |
177 | SetTitle(title); | |
178 | // | |
179 | // Constructor | |
180 | // | |
263d466a | 181 | fTracksEsdParam.SetOwner(kTRUE); |
c6914c83 | 182 | } |
183 | ||
184 | AliTPCcalibLaser::~AliTPCcalibLaser() { | |
185 | // | |
186 | // destructor | |
187 | // | |
188 | } | |
189 | ||
e9f38a4b | 190 | |
191 | ||
192 | void AliTPCcalibLaser::Process(AliESDEvent * event) { | |
c6914c83 | 193 | // |
263d466a | 194 | // |
e9f38a4b | 195 | // Loop over tracks and call Process function |
196 | // | |
197 | fESD = event; | |
198 | if (!fESD) { | |
199 | return; | |
200 | } | |
201 | fESDfriend=static_cast<AliESDfriend*>(fESD->FindListObject("AliESDfriend")); | |
202 | if (!fESDfriend) { | |
203 | return; | |
204 | } | |
263d466a | 205 | AliDebug(4,Form("Event number in current file: %d",event->GetEventNumberInFile())); |
e9f38a4b | 206 | fTracksTPC.Clear(); |
207 | fTracksEsd.Clear(); | |
208 | fTracksEsdParam.Delete(); | |
209 | // | |
210 | Int_t n=fESD->GetNumberOfTracks(); | |
211 | Int_t run = fESD->GetRunNumber(); | |
212 | fRun = run; | |
213 | for (Int_t i=0;i<n;++i) { | |
214 | AliESDfriendTrack *friendTrack=fESDfriend->GetTrack(i); | |
215 | AliESDtrack *track=fESD->GetTrack(i); | |
216 | TObject *calibObject=0; | |
217 | AliTPCseed *seed=0; | |
218 | for (Int_t j=0;(calibObject=friendTrack->GetCalibObject(j));++j) | |
219 | if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) | |
220 | break; | |
221 | if (track&&seed) FindMirror(track,seed); | |
222 | // | |
223 | } | |
263d466a | 224 | |
e9f38a4b | 225 | FitDriftV(); |
c18f4385 | 226 | MakeDistHisto(); |
e9f38a4b | 227 | // |
c18f4385 | 228 | for (Int_t id=0; id<336; id++){ |
e9f38a4b | 229 | // |
230 | // | |
231 | if (!fTracksEsdParam.At(id)) continue; | |
232 | DumpLaser(id); | |
263d466a | 233 | // RefitLaser(id); |
234 | RefitLaserJW(id); | |
235 | ||
c18f4385 | 236 | } |
263d466a | 237 | // fEvent++; |
c18f4385 | 238 | } |
239 | ||
240 | void AliTPCcalibLaser::MakeDistHisto(){ | |
241 | // | |
242 | // | |
243 | // | |
244 | for (Int_t id=0; id<336; id++){ | |
245 | // | |
246 | // | |
263d466a | 247 | if (!fTracksEsdParam.At(id)) continue; |
c18f4385 | 248 | if (!AcceptLaser(id)) continue; |
249 | // | |
250 | // | |
251 | TH1F * hisdz = (TH1F*)fDeltaZ.At(id); | |
252 | TH1F * hisdphi = (TH1F*)fDeltaPhi.At(id); | |
253 | TH1F * hisdphiP = (TH1F*)fDeltaPhiP.At(id); | |
2b35e8e6 | 254 | TH1F * hisSignal = (TH1F*)fSignals.At(id); |
255 | ||
263d466a | 256 | if (!hisdz){ |
c03e3250 | 257 | hisdz = new TH1F(Form("hisdz%d",id),Form("hisdz%d",id),1000,-10,10); |
258 | hisdz->SetDirectory(0); | |
c18f4385 | 259 | fDeltaZ.AddAt(hisdz,id); |
260 | // | |
c03e3250 | 261 | hisdphi = new TH1F(Form("hisdphi%d",id),Form("hisdphi%d",id),1000,-1,1); |
262 | hisdphi->SetDirectory(0); | |
c18f4385 | 263 | fDeltaPhi.AddAt(hisdphi,id); |
264 | // | |
c03e3250 | 265 | hisdphiP = new TH1F(Form("hisdphiP%d",id),Form("hisdphiP%d",id),1000,-0.01,0.01); |
266 | hisdphiP->SetDirectory(0); | |
c18f4385 | 267 | fDeltaPhiP.AddAt(hisdphiP,id); |
c03e3250 | 268 | hisSignal = new TH1F(Form("hisSignal%d",id),Form("hisSignal%d",id),1000,0,1000); |
269 | hisSignal->SetDirectory(0); | |
2b35e8e6 | 270 | fSignals.AddAt(hisSignal,id); |
c18f4385 | 271 | } |
272 | ||
273 | AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id); | |
274 | AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id); | |
2b35e8e6 | 275 | AliESDtrack *track = (AliESDtrack*)fTracksEsd.At(id); |
c03e3250 | 276 | if (!param) return; |
277 | if (!ltrp) return; | |
278 | if (!track) return; | |
c18f4385 | 279 | Double_t xyz[3]; |
280 | Double_t pxyz[3]; | |
281 | Double_t lxyz[3]; | |
282 | Double_t lpxyz[3]; | |
283 | param->GetXYZ(xyz); | |
284 | param->GetPxPyPz(pxyz); | |
285 | ltrp->GetXYZ(lxyz); | |
286 | ltrp->GetPxPyPz(lpxyz); | |
287 | // | |
288 | Float_t dz = param->GetZ()-ltrp->GetZ(); | |
289 | Float_t dphi = (TMath::ATan2(xyz[1],xyz[0])- TMath::ATan2(lxyz[1],lxyz[0]))*254.; | |
290 | Float_t dphiP = param->GetParameter()[2]-ltrp->GetParameter()[2]; | |
c03e3250 | 291 | if (hisdz) hisdz->Fill(dz); |
292 | if (hisdphi) hisdphi->Fill(dphi); | |
263d466a | 293 | if (hisdphiP) hisdphiP->Fill(dphiP); |
c03e3250 | 294 | if (hisSignal) hisSignal->Fill(track->GetTPCsignal()); |
c18f4385 | 295 | } |
296 | } | |
297 | ||
298 | void AliTPCcalibLaser::FitDriftV(){ | |
299 | // | |
300 | // Fit drift velocity - linear approximation in the z and global y | |
301 | // | |
302 | static TLinearFitter fdriftA(3,"hyp2"); | |
303 | static TLinearFitter fdriftC(3,"hyp2"); | |
304 | fdriftA.ClearPoints(); | |
305 | fdriftC.ClearPoints(); | |
306 | // | |
307 | for (Int_t id=0; id<336; id++){ | |
263d466a | 308 | if (!fTracksEsdParam.At(id)) continue; |
c18f4385 | 309 | if (!AcceptLaser(id)) continue; |
310 | AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id); | |
311 | AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id); | |
312 | Double_t xyz[3]; | |
313 | Double_t pxyz[3]; | |
314 | Double_t lxyz[3]; | |
315 | Double_t lpxyz[3]; | |
316 | param->GetXYZ(xyz); | |
317 | param->GetPxPyPz(pxyz); | |
318 | ltrp->GetXYZ(lxyz); | |
319 | ltrp->GetPxPyPz(lpxyz); | |
320 | Double_t xxx[2] = {lxyz[2],lxyz[1]}; | |
321 | if (ltrp->GetSide()==0){ | |
322 | fdriftA.AddPoint(xxx,xyz[2],1); | |
323 | }else{ | |
324 | fdriftC.AddPoint(xxx,xyz[2],1); | |
325 | } | |
326 | } | |
327 | Float_t chi2A = 0; | |
328 | Float_t chi2C = 0; | |
329 | Int_t npointsA=0; | |
330 | Int_t npointsC=0; | |
331 | // | |
332 | if (fdriftA.GetNpoints()>10){ | |
333 | fdriftA.Eval(); | |
334 | fdriftA.EvalRobust(0.8); | |
335 | fdriftA.GetParameters(*fFitAside); | |
336 | npointsA= fdriftA.GetNpoints(); | |
337 | chi2A = fdriftA.GetChisquare()/fdriftA.GetNpoints(); | |
338 | } | |
339 | if (fdriftC.GetNpoints()>10){ | |
340 | fdriftC.Eval(); | |
341 | fdriftC.EvalRobust(0.8); | |
342 | fdriftC.GetParameters(*fFitCside); | |
343 | npointsC= fdriftC.GetNpoints(); | |
344 | chi2C = fdriftC.GetChisquare()/fdriftC.GetNpoints(); | |
e9f38a4b | 345 | } |
263d466a | 346 | |
c18f4385 | 347 | if (fStreamLevel>0){ |
348 | TTreeSRedirector *cstream = GetDebugStreamer(); | |
349 | Int_t time = fESD->GetTimeStamp(); | |
350 | if (cstream){ | |
351 | (*cstream)<<"driftv"<< | |
352 | "driftA.="<<fFitAside<< | |
353 | "driftC.="<<fFitCside<< | |
354 | "chi2A="<<chi2A<< | |
355 | "chi2C="<<chi2C<< | |
356 | "nA="<<npointsA<< | |
357 | "nC="<<npointsC<< | |
358 | "time="<<time<< | |
359 | "\n"; | |
360 | } | |
361 | } | |
362 | // | |
363 | } | |
364 | ||
365 | ||
366 | Bool_t AliTPCcalibLaser::AcceptLaser(Int_t id){ | |
367 | // | |
368 | // | |
369 | // | |
370 | /* | |
371 | TCut cutT("cutT","abs(Tr.fP[3])<0.06"); | |
372 | TCut cutPt("cutPt","abs(Tr.fP[4])<0.1"); | |
373 | TCut cutN("cutN","fTPCncls>70"); | |
374 | TCut cutP("cutP","abs(atan2(x1,x0)-atan2(lx1,lx0))<0.03") | |
375 | TCut cutA = cutT+cutPt+cutP; | |
376 | */ | |
377 | AliExternalTrackParam *param =(AliExternalTrackParam*)fTracksEsdParam.At(id); | |
378 | AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id); | |
379 | AliESDtrack *track = (AliESDtrack*)fTracksEsd.At(id); | |
380 | ||
381 | if (TMath::Abs(param->GetParameter()[4])>0.03) return kFALSE; | |
263d466a | 382 | if (TMath::Abs(param->GetParameter()[3])>0.06) return kFALSE; |
c18f4385 | 383 | if (TMath::Abs(param->GetParameter()[2]-ltrp->GetParameter()[2])>0.06) return kFALSE; |
384 | if (TMath::Abs(param->GetParameter()[1]-ltrp->GetParameter()[1])>10) return kFALSE; | |
385 | // | |
386 | // dedx cut | |
387 | // | |
388 | if (TMath::Abs(track->GetTPCsignal())<20) return kFALSE; | |
389 | if (TMath::Abs(track->GetTPCsignal())>800) return kFALSE; | |
390 | // | |
263d466a | 391 | return kTRUE; |
e9f38a4b | 392 | } |
393 | ||
394 | Int_t AliTPCcalibLaser::FindMirror(AliESDtrack *track, AliTPCseed *seed){ | |
395 | // | |
396 | // Find corresponding mirror | |
397 | // add the corresponding tracks | |
c6914c83 | 398 | // |
1bd37cc0 | 399 | Float_t kRadius0 = 252; |
c18f4385 | 400 | Float_t kRadius = 253.4; |
e9f38a4b | 401 | if (!track->GetOuterParam()) return -1; |
c6914c83 | 402 | AliExternalTrackParam param(*(track->GetOuterParam())); |
1bd37cc0 | 403 | AliTracker::PropagateTrackTo(¶m,kRadius0,0.10566,3,kTRUE); |
404 | AliTracker::PropagateTrackTo(¶m,kRadius,0.10566,0.1,kTRUE); | |
7b18d067 | 405 | AliTPCLaserTrack ltr; |
406 | AliTPCLaserTrack *ltrp=0x0; | |
5732d662 | 407 | // AliTPCLaserTrack *ltrpjw=0x0; |
1bd37cc0 | 408 | // |
263d466a | 409 | Int_t id = AliTPCLaserTrack::IdentifyTrack(¶m); |
410 | // Int_t idjw = AliTPCLaserTrack::IdentifyTrackJW(¶m); | |
411 | //AliDebug(4,Form("Identified Track: %03d (%03d)",id,idjw)); | |
412 | ||
413 | if (id!=-1 && (AliTPCLaserTrack::GetTracks()->UncheckedAt(id))) | |
1bd37cc0 | 414 | ltrp=(AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(id); |
263d466a | 415 | else |
1bd37cc0 | 416 | ltrp=<r; |
263d466a | 417 | /* |
418 | if (idjw!=-1 && (AliTPCLaserTrack::GetTracks()->UncheckedAt(idjw))) | |
419 | ltrpjw=(AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(idjw); | |
420 | else | |
421 | ltrpjw=<r; | |
422 | ||
423 | ||
424 | if (fStreamLevel>0){ | |
425 | TTreeSRedirector *cstream = GetDebugStreamer(); | |
426 | if (cstream){ | |
427 | (*cstream)<<"idcmp"<< | |
428 | "id=" << id << | |
429 | "idjw=" << idjw << | |
430 | "tr.=" << ltrp << | |
431 | "trjw.="<< ltrpjw << | |
432 | "seed.="<<seed<< | |
433 | "event="<<fEvent << | |
434 | "\n"; | |
435 | } | |
436 | } | |
437 | ||
438 | */ | |
439 | ||
440 | ||
e9f38a4b | 441 | if (id>=0){ |
c18f4385 | 442 | // |
443 | // | |
444 | Float_t radius=TMath::Abs(ltrp->GetX()); | |
445 | AliTracker::PropagateTrackTo(¶m,radius,0.10566,0.01,kTRUE); | |
446 | // | |
e9f38a4b | 447 | if (!fTracksMirror.At(id)) fTracksMirror.AddAt(ltrp,id); |
448 | fTracksEsdParam.AddAt(param.Clone(),id); | |
449 | fTracksEsd.AddAt(track,id); | |
450 | fTracksTPC.AddAt(seed,id); | |
c18f4385 | 451 | // |
e9f38a4b | 452 | } |
c18f4385 | 453 | |
e9f38a4b | 454 | return id; |
455 | } | |
456 | ||
457 | ||
458 | ||
459 | void AliTPCcalibLaser::DumpLaser(Int_t id) { | |
460 | // | |
461 | // Dump Laser info to the tree | |
462 | // | |
463 | AliESDtrack *track = (AliESDtrack*)fTracksEsd.At(id); | |
464 | AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id); | |
465 | AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id); | |
1bd37cc0 | 466 | // |
467 | // Fast laser ID | |
468 | // | |
469 | Double_t xyz[3]; | |
470 | Double_t pxyz[3]; | |
c18f4385 | 471 | Double_t lxyz[3]; |
472 | Double_t lpxyz[3]; | |
e9f38a4b | 473 | param->GetXYZ(xyz); |
474 | param->GetPxPyPz(pxyz); | |
c18f4385 | 475 | ltrp->GetXYZ(lxyz); |
476 | ltrp->GetPxPyPz(lpxyz); | |
477 | ||
c6914c83 | 478 | if (fStreamLevel>0){ |
479 | TTreeSRedirector *cstream = GetDebugStreamer(); | |
c18f4385 | 480 | Int_t time = fESD->GetTimeStamp(); |
481 | Bool_t accept = AcceptLaser(id); | |
c6914c83 | 482 | if (cstream){ |
483 | (*cstream)<<"Track"<< | |
e9f38a4b | 484 | "run="<<fRun<< |
7b18d067 | 485 | "id="<<id<< |
c18f4385 | 486 | "accept="<<accept<< |
487 | "driftA.="<<fFitAside<< | |
488 | "driftC.="<<fFitCside<< | |
489 | "time="<<time<< | |
7b18d067 | 490 | // |
491 | "LTr.="<<ltrp<< | |
492 | "Esd.="<<track<< | |
e9f38a4b | 493 | "Tr.="<<param<< |
c6914c83 | 494 | "x0="<<xyz[0]<< |
495 | "x1="<<xyz[1]<< | |
496 | "x2="<<xyz[2]<< | |
497 | "px0="<<pxyz[0]<< | |
498 | "px1="<<pxyz[1]<< | |
499 | "px2="<<pxyz[2]<< | |
c18f4385 | 500 | // |
501 | "lx0="<<lxyz[0]<< | |
502 | "lx1="<<lxyz[1]<< | |
503 | "lx2="<<lxyz[2]<< | |
504 | "lpx0="<<lpxyz[0]<< | |
505 | "lpx1="<<lpxyz[1]<< | |
506 | "lpx2="<<lpxyz[2]<< | |
c6914c83 | 507 | "\n"; |
508 | } | |
509 | } | |
510 | } | |
511 | ||
263d466a | 512 | void AliTPCcalibLaser::RefitLaserJW(Int_t id){ |
513 | // | |
514 | // Refit the track with different tracklet models: | |
515 | // 1. Per ROC using the kalman filter, different edge cuts | |
516 | // 2. Per ROC linear in y and z | |
517 | // 3. Per ROC quadratic in y and z | |
518 | // 4. Per track offset for each sector, linear for each sector, common quadratic | |
519 | // store x, y, z information for all models and the cluster to calculate the residuals | |
520 | // | |
521 | AliTPCseed *track = (AliTPCseed*)fTracksTPC.At(id); | |
522 | AliExternalTrackParam *extparam=(AliExternalTrackParam*)fTracksEsdParam.At(id); | |
523 | AliTPCLaserTrack *ltrp = (AliTPCLaserTrack*)fTracksMirror.At(id); | |
524 | ||
525 | AliTPCclusterMI dummyCl; | |
526 | ||
527 | //two tracklets | |
528 | Int_t kMaxTracklets=2; | |
0728a4f6 | 529 | //=============================================// |
530 | // Linear Fitters for the Different Approaches // | |
531 | //=============================================// | |
532 | //linear fit model in y and z; inner - outer sector | |
533 | static TLinearFitter fy1I(2,"hyp1"); | |
534 | static TLinearFitter fy1O(2,"hyp1"); | |
535 | static TLinearFitter fz1I(2,"hyp1"); | |
536 | static TLinearFitter fz1O(2,"hyp1"); | |
537 | //quadratic fit model in y and z; inner - sector | |
538 | static TLinearFitter fy2I(3,"hyp2"); | |
539 | static TLinearFitter fy2O(3,"hyp2"); | |
540 | static TLinearFitter fz2I(3,"hyp2"); | |
541 | static TLinearFitter fz2O(3,"hyp2"); | |
263d466a | 542 | //common quadratic fit for IROC and OROC in y and z |
543 | static TLinearFitter fy4(5,"hyp4"); | |
544 | static TLinearFitter fz4(5,"hyp4"); | |
545 | ||
546 | ||
547 | //set standard cuts | |
548 | if ( fNcuts==0 ){ | |
549 | fNcuts=1; | |
550 | fEdgeXcuts[0]=4; | |
551 | fEdgeYcuts[0]=3; | |
552 | fNClCuts[0]=20; | |
553 | } | |
0728a4f6 | 554 | //=============================// |
555 | // Loop over all Tracklet Cuts // | |
556 | //=============================// | |
263d466a | 557 | for (Int_t icut=0; icut<fNcuts; icut++){ |
558 | AliDebug(4,Form("Processing cut %d for track with ID %d",icut,id)); | |
559 | //cut parameters | |
560 | Double_t edgeCutX = fEdgeXcuts[icut]; | |
561 | Double_t edgeCutY = fEdgeYcuts[icut]; | |
562 | Int_t nclCut = fNClCuts[icut]; | |
0728a4f6 | 563 | //=========================// |
564 | // Parameters to calculate // | |
565 | //=========================// | |
566 | //fit parameter inner, outer tracklet and combined fit | |
567 | TVectorD vecy1resInner(2),vecz1resInner(2); //pol1 fit parameters inner | |
568 | TVectorD vecy2resInner(3),vecz2resInner(3); //pol2 fit parameters inner | |
263d466a | 569 | // |
0728a4f6 | 570 | TVectorD vecy1resOuter(2),vecz1resOuter(2); //pol1 fit parameters outer |
571 | TVectorD vecy2resOuter(3),vecz2resOuter(3); //pol2 fit parameters outer | |
263d466a | 572 | TVectorD vecy4res(5),vecz4res(5); |
573 | // cluster and track positions for each row - used for residuals | |
574 | TVectorD vecX(159); // x is the same for all (row center) | |
575 | TVectorD vecYkalman(159); // y from kalman fit | |
576 | TVectorD vecZkalman(159); // z from kalman fit | |
577 | TVectorD vecY1(159); // y from pol1 fit per ROC | |
578 | TVectorD vecZ1(159); // z from pol1 fit per ROC | |
579 | TVectorD vecY2(159); // y from pol2 fit per ROC | |
580 | TVectorD vecZ2(159); // z from pol2 fit per ROC | |
581 | TVectorD vecY4(159); // y from sector fit | |
582 | TVectorD vecZ4(159); // z from sector fit | |
583 | TVectorD vecClY(159); // y cluster position | |
584 | TVectorD vecClZ(159); // z cluster position | |
585 | TVectorD vecSec(159); // sector for each row | |
0728a4f6 | 586 | //chi2 of fits |
263d466a | 587 | Double_t chi2I1z=-1; // chi2 of pol1 fit in z (inner) |
588 | Double_t chi2I1y=-1; // chi2 of pol1 fit in y (inner) | |
589 | Double_t chi2O1z=-1; // chi2 of pol1 fit in z (outer) | |
590 | Double_t chi2O1y=-1; // chi2 of pol1 fit in y (outer) | |
591 | Double_t chi2I2z=-1; // chi2 of pol2 fit in z (inner) | |
592 | Double_t chi2I2y=-1; // chi2 of pol2 fit in y (inner) | |
593 | Double_t chi2O2z=-1; // chi2 of pol2 fit in z (outer) | |
594 | Double_t chi2O2y=-1; // chi2 of pol2 fit in y (outer) | |
595 | Double_t chi2IOz=-1; // chi2 of hyp4 fit in z (inner+outer) | |
596 | Double_t chi2IOy=-1; // chi2 of hyp4 fit in y (inner+outer) | |
0728a4f6 | 597 | //more |
598 | Int_t innerSector = -1; // number of inner sector | |
599 | Int_t outerSector = -1; // number of outer sector | |
263d466a | 600 | Int_t nclI=0; // number of clusters (inner) |
601 | Int_t nclO=0; // number of clusters (outer) | |
0728a4f6 | 602 | //=================================================// |
603 | // Perform the Kalman Fit using the Tracklet Class // | |
604 | //=================================================// | |
263d466a | 605 | AliTPCTracklet::SetEdgeCut(edgeCutX,edgeCutY); |
606 | TObjArray tracklets= | |
607 | AliTPCTracklet::CreateTracklets(track,AliTPCTracklet::kKalman, | |
608 | kFALSE,nclCut,kMaxTracklets); | |
263d466a | 609 | // tracklet pointers |
263d466a | 610 | AliTPCTracklet *trInner = (AliTPCTracklet*)tracklets.At(0); |
611 | AliTPCTracklet *trOuter = (AliTPCTracklet*)tracklets.At(1); | |
0728a4f6 | 612 | AliTPCTracklet *tr=0x0; |
613 | AliTPCTracklet dummy; | |
614 | //continue if we didn't find a tracklet | |
263d466a | 615 | if ( !trInner && !trOuter ) continue; |
0728a4f6 | 616 | //================================================================================// |
617 | // Swap Inner and Outer if necessary (inner sector is definde by smaller local x) // | |
618 | //================================================================================// | |
263d466a | 619 | if ( trInner && trOuter ){ |
620 | if ( !trInner->GetInner() || !trOuter->GetInner() ) continue; | |
621 | if ( trInner->GetInner()->GetX() > trOuter->GetInner()->GetX() ){ | |
622 | tr = trInner; | |
623 | trInner=trOuter; | |
624 | trOuter=tr; | |
625 | AliDebug(5,Form("Swapped Sectors: %02d (%f) <-> %02d (%f)", trOuter->GetSector(), trOuter->GetInner()->GetX(), trInner->GetSector(), trInner->GetInner()->GetX())); | |
626 | } | |
627 | } else { | |
628 | if ( trInner ){ | |
629 | if ( !trInner->GetInner() ) continue; | |
630 | trOuter=&dummy; | |
631 | if ( trInner->GetSector()>35 ){ | |
632 | trOuter=trInner; | |
633 | trInner=&dummy; | |
634 | } | |
635 | } else { //trOuter | |
636 | if ( !trOuter->GetInner() ) continue; | |
637 | trInner=&dummy; | |
638 | if ( trOuter->GetSector()<36 ){ | |
639 | trInner=trOuter; | |
640 | trOuter=&dummy; | |
641 | } | |
642 | } | |
643 | } | |
644 | innerSector = trInner->GetSector(); | |
645 | if ( innerSector>=0 ) AliDebug(5,Form("Found inner Sector %02d at X %.2f", innerSector, trInner->GetInner()->GetX())); | |
646 | outerSector = trOuter->GetSector(); | |
647 | if ( outerSector>=0 ) AliDebug(5,Form("Found outer Sector %02d at X %.2f", outerSector, trOuter->GetInner()->GetX())); | |
648 | ||
263d466a | 649 | // array of clusters |
650 | TClonesArray arrCl("AliTPCclusterMI",159); | |
651 | arrCl.ExpandCreateFast(159); | |
0728a4f6 | 652 | //=======================================// |
653 | // fill fitters with cluster information // | |
654 | //=======================================// | |
655 | AliDebug(3,"Filing Cluster Information"); | |
263d466a | 656 | for (Int_t irow=158;irow>-1;--irow) { |
657 | AliTPCclusterMI *c=track->GetClusterPointer(irow); | |
0728a4f6 | 658 | AliTPCclusterMI &cl = (AliTPCclusterMI&) (*arrCl[irow]); |
659 | cl=dummyCl; | |
660 | // | |
263d466a | 661 | vecSec[irow]=-1; |
662 | if (!c) continue; | |
0728a4f6 | 663 | // |
664 | Int_t roc = static_cast<Int_t>(c->GetDetector()); | |
665 | if ( roc!=innerSector && roc!=outerSector ) continue; | |
666 | vecSec[irow]=roc; | |
263d466a | 667 | //store clusters in clones array |
668 | cl=*c; | |
669 | //cluster position | |
670 | vecX[irow] = c->GetX(); | |
671 | vecClY[irow] = c->GetY(); | |
672 | vecClZ[irow] = c->GetZ(); | |
263d466a | 673 | // |
0728a4f6 | 674 | Double_t x = vecX[irow]-133.4; //reference is between IROC and OROC |
675 | Double_t y = vecClY[irow]; | |
676 | Double_t z = vecClZ[irow]; | |
677 | // | |
678 | Double_t x2[2]={x,x*x}; //linear and parabolic parameters | |
679 | Double_t x4[4]={0,0,0,0}; //hyp4 parameters | |
263d466a | 680 | if ( roc == innerSector ) { |
0728a4f6 | 681 | x4[0]=1; //offset inner - outer sector |
682 | x4[1]=x; //slope parameter inner sector | |
263d466a | 683 | } else { |
0728a4f6 | 684 | x4[2]=x; //slope parameter outer sector |
263d466a | 685 | } |
0728a4f6 | 686 | x4[3]=x*x; //common parabolic shape |
263d466a | 687 | // |
0728a4f6 | 688 | if ( roc==innerSector ){ |
689 | fy1I.AddPoint(x2,y); | |
690 | fz1I.AddPoint(x2,z); | |
691 | fy2I.AddPoint(x2,y); | |
692 | fz2I.AddPoint(x2,z); | |
693 | ++nclI; | |
694 | } | |
695 | if ( roc==outerSector ){ | |
696 | fy1O.AddPoint(x2,y); | |
697 | fz1O.AddPoint(x2,z); | |
698 | fy2O.AddPoint(x2,y); | |
699 | fz2O.AddPoint(x2,z); | |
700 | ++nclO; | |
701 | } | |
263d466a | 702 | fy4.AddPoint(x4,y); |
703 | fz4.AddPoint(x4,z); | |
704 | } | |
0728a4f6 | 705 | //======================================// |
706 | // Evaluate and retrieve fit parameters // | |
707 | //======================================// | |
708 | AliDebug(5,Form("Evaluating fits with inner (outer) Sec: %02d (%02d)",innerSector,outerSector)); | |
709 | //inner sector | |
710 | if ( (innerSector>-1) && (fy1I.GetNpoints()>0) ){ | |
711 | fy1I.Eval(); | |
712 | fz1I.Eval(); | |
713 | fy2I.Eval(); | |
714 | fz2I.Eval(); | |
715 | fy1I.GetParameters(vecy1resInner); | |
716 | fz1I.GetParameters(vecz1resInner); | |
717 | fy2I.GetParameters(vecy2resInner); | |
718 | fz2I.GetParameters(vecz2resInner); | |
719 | chi2I1y=fy1I.GetChisquare()/(fy1I.GetNpoints()-2); | |
720 | chi2I1z=fz1I.GetChisquare()/(fz1I.GetNpoints()-2); | |
721 | chi2I2y=fy2I.GetChisquare()/(fy2I.GetNpoints()-3); | |
722 | chi2I2z=fz2I.GetChisquare()/(fz2I.GetNpoints()-3); | |
723 | } | |
724 | //outer sector | |
725 | if ( (outerSector>-1) && (fy1O.GetNpoints()>0) ){ | |
726 | fy1O.Eval(); | |
727 | fz1O.Eval(); | |
728 | fy2O.Eval(); | |
729 | fz2O.Eval(); | |
730 | fy1O.GetParameters(vecy1resOuter); | |
731 | fz1O.GetParameters(vecz1resOuter); | |
732 | fy2O.GetParameters(vecy2resOuter); | |
733 | fz2O.GetParameters(vecz2resOuter); | |
734 | chi2O1y=fy1O.GetChisquare()/(fy1O.GetNpoints()-2); | |
735 | chi2O1z=fz1O.GetChisquare()/(fz1O.GetNpoints()-2); | |
736 | chi2O2y=fy2O.GetChisquare()/(fy2O.GetNpoints()-3); | |
737 | chi2O2z=fz2O.GetChisquare()/(fz2O.GetNpoints()-3); | |
738 | } | |
739 | //combined hyp4 fit | |
263d466a | 740 | if ( innerSector>0 && outerSector>0 ){ |
0728a4f6 | 741 | if (fy4.GetNpoints()>0) { |
742 | fy4.Eval(); | |
743 | fy4.GetParameters(vecy4res); | |
744 | chi2IOy=fy4.GetChisquare()/(fy4.GetNpoints()-5); | |
745 | } | |
746 | if (fz4.GetNpoints()>0) { | |
747 | fz4.Eval(); | |
748 | fz4.GetParameters(vecz4res); | |
749 | chi2IOz=fz4.GetChisquare()/(fz4.GetNpoints()-5); | |
750 | } | |
263d466a | 751 | } |
0728a4f6 | 752 | //clear points |
753 | fy4.ClearPoints(); fz4.ClearPoints(); | |
754 | fy1I.ClearPoints(); fy1O.ClearPoints(); | |
755 | fz1I.ClearPoints(); fz1O.ClearPoints(); | |
756 | fy2I.ClearPoints(); fy2O.ClearPoints(); | |
757 | fz2I.ClearPoints(); fz2O.ClearPoints(); | |
758 | //==============================// | |
759 | // calculate tracklet positions // | |
760 | //==============================// | |
263d466a | 761 | AliDebug(4,"Calculate tracklet positions"); |
762 | for (Int_t irow=158;irow>-1;--irow) { | |
0728a4f6 | 763 | if ( vecSec[irow]==-1 ) continue; //no cluster info |
263d466a | 764 | if ( vecSec[irow]!=innerSector && vecSec[irow]!=outerSector ) continue; |
765 | tr=&dummy; | |
766 | Double_t x=vecX[irow]; | |
0728a4f6 | 767 | Double_t xref=x-133.4; |
768 | // | |
263d466a | 769 | Double_t yoffInner=0; |
770 | Double_t zoffInner=0; | |
771 | Double_t yslopeInner=0; | |
772 | Double_t yslopeOuter=0; | |
773 | Double_t zslopeInner=0; | |
0728a4f6 | 774 | Double_t zslopeOuter=0; |
775 | //positions of hyperplane fits | |
263d466a | 776 | if ( vecSec[irow] == outerSector ) { |
777 | tr=trOuter; | |
778 | vecY1[irow]=vecy1resOuter[0]+vecy1resOuter[1]*xref; | |
779 | vecZ1[irow]=vecz1resOuter[0]+vecz1resOuter[1]*xref; | |
780 | vecY2[irow]=vecy2resOuter[0]+vecy2resOuter[1]*xref+vecy2resOuter[2]*xref*xref; | |
781 | vecZ2[irow]=vecz2resOuter[0]+vecz2resOuter[1]*xref+vecz2resOuter[2]*xref*xref; | |
782 | yslopeOuter=vecy4res[3]; | |
783 | zslopeOuter=vecz4res[3]; | |
263d466a | 784 | } else { |
785 | tr=trInner; | |
786 | vecY1[irow]=vecy1resInner[0]+vecy1resInner[1]*xref; | |
787 | vecZ1[irow]=vecz1resInner[0]+vecz1resInner[1]*xref; | |
788 | vecY2[irow]=vecy2resInner[0]+vecy2resInner[1]*xref+vecy2resInner[2]*xref*xref; | |
789 | vecZ2[irow]=vecz2resInner[0]+vecz2resInner[1]*xref+vecz2resInner[2]*xref*xref; | |
790 | yoffInner=vecy4res[1]; | |
791 | zoffInner=vecz4res[1]; | |
792 | yslopeInner=vecy4res[2]; | |
793 | zslopeInner=vecz4res[2]; | |
263d466a | 794 | } |
795 | vecY4[irow]=vecy4res[0]+yoffInner+yslopeInner*xref+yslopeOuter*xref+vecy4res[4]*xref*xref; | |
796 | vecZ4[irow]=vecz4res[0]+zoffInner+zslopeInner*xref+zslopeOuter*xref+vecz4res[4]*xref*xref; | |
0728a4f6 | 797 | //positions of kalman fits |
263d466a | 798 | Double_t gxyz[3],xyz[3]; |
799 | AliExternalTrackParam *param = 0x0; | |
0728a4f6 | 800 | // |
263d466a | 801 | param=tr->GetInner(); |
802 | if (param){ | |
803 | param->GetXYZ(gxyz); | |
804 | Float_t bz = AliTracker::GetBz(gxyz); | |
805 | param->GetYAt(x, bz, xyz[1]); | |
806 | param->GetZAt(x, bz, xyz[2]); | |
807 | vecYkalman[irow]=xyz[1]; | |
808 | vecZkalman[irow]=xyz[2]; | |
809 | } | |
263d466a | 810 | } |
0728a4f6 | 811 | //=====================================================================// |
812 | // write results from the different tracklet fits with debug streamers // | |
813 | //=====================================================================// | |
263d466a | 814 | if (fStreamLevel>4){ |
815 | TTreeSRedirector *cstream = GetDebugStreamer(); | |
816 | if (cstream){ | |
817 | Float_t dedx = track->GetdEdx(); | |
818 | (*cstream)<<"FitModels"<< | |
819 | "cutNr=" << icut << | |
820 | "edgeCutX=" << edgeCutX << | |
821 | "edgeCutY=" << edgeCutY << | |
822 | "nclCut=" << nclCut << | |
823 | "innerSector="<< innerSector << | |
824 | "outerSector="<< outerSector << | |
825 | "dEdx=" << dedx << | |
826 | "LTr.=" << ltrp << | |
827 | "Tr.=" << extparam << | |
828 | "yPol1In.=" << &vecy1resInner << | |
829 | "zPol1In.=" << &vecz1resInner << | |
830 | "yPol2In.=" << &vecy2resInner << | |
831 | "zPol2In.=" << &vecz2resInner << | |
832 | "yPol1Out.=" << &vecy1resOuter << | |
833 | "zPol1Out.=" << &vecz1resOuter << | |
834 | "yPol2Out.=" << &vecy2resOuter << | |
835 | "zPol2Out.=" << &vecz2resOuter << | |
836 | "yInOut.=" << &vecy4res << | |
837 | "zInOut.=" << &vecz4res << | |
838 | "chi2y1In=" << chi2I1y << | |
839 | "chi2z1In=" << chi2I1z << | |
840 | "chi2y1Out=" << chi2O1y << | |
841 | "chi2z1Out=" << chi2O1z << | |
842 | "chi2y2In=" << chi2I2y << | |
843 | "chi2z2In=" << chi2I2z << | |
844 | "chi2y2Out=" << chi2O2y << | |
845 | "chi2z2Out=" << chi2O2z << | |
846 | "chi2yInOut=" << chi2IOy << | |
847 | "chi2zInOut=" << chi2IOz << | |
848 | "trletIn.=" << trInner << | |
849 | "trletOut.=" << trOuter << | |
850 | "nclI=" << nclI << | |
851 | "nclO=" << nclO << | |
852 | "\n"; | |
853 | } | |
854 | } | |
855 | ||
856 | // wirte residuals information | |
857 | if (fStreamLevel>5){ | |
858 | TTreeSRedirector *cstream = GetDebugStreamer(); | |
859 | if (cstream){ | |
860 | Float_t dedx = track->GetdEdx(); | |
861 | (*cstream)<<"Residuals"<< | |
862 | "cutNr=" << icut << | |
863 | "edgeCutX=" << edgeCutX << | |
864 | "edgeCutY=" << edgeCutY << | |
865 | "nclCut=" << nclCut << | |
866 | "LTr.=" << ltrp << | |
867 | "Tr.=" << extparam<< | |
868 | "dEdx=" << dedx << | |
869 | "Cl.=" << &arrCl << | |
870 | "TrX.=" << &vecX << | |
871 | "TrYpol1.=" << &vecY1 << | |
872 | "TrZpol1.=" << &vecZ1 << | |
873 | "TrYpol2.=" << &vecY2 << | |
874 | "TrZpol2.=" << &vecZ2 << | |
875 | "TrYInOut.=" << &vecY4 << | |
876 | "TrZInOut.=" << &vecZ4 << | |
877 | "ClY.=" << &vecClY << | |
878 | "ClZ.=" << &vecClZ << | |
0728a4f6 | 879 | "sec.=" << &vecSec << |
263d466a | 880 | "nclI=" << nclI << |
881 | "nclO=" << nclO << | |
882 | "yInOut.=" << &vecy4res << | |
883 | "zInOut.=" << &vecz4res << | |
884 | "chi2y1In=" << chi2I1y << | |
885 | "chi2z1In=" << chi2I1z << | |
886 | "chi2y1Out=" << chi2O1y << | |
887 | "chi2z1Out=" << chi2O1z << | |
888 | "chi2y2In=" << chi2I2y << | |
889 | "chi2z2In=" << chi2I2z << | |
890 | "chi2y2Out=" << chi2O2y << | |
891 | "chi2z2Out=" << chi2O2z << | |
892 | "chi2yInOut=" << chi2IOy << | |
893 | "chi2zInOut=" << chi2IOz << | |
894 | "\n"; | |
895 | ||
896 | } | |
897 | } | |
95a0e09b | 898 | //==========================// |
899 | // Fill Residual Histograms // | |
900 | //==========================// | |
901 | TProfile *profy = (TProfile*)fDeltaYres.UncheckedAt(id); | |
902 | TProfile *profz = (TProfile*)fDeltaZres.UncheckedAt(id); | |
903 | if (!profy){ | |
904 | profy=new TProfile(Form("pry%03d",id),Form("Y Residuals for Laser Beam %03d",id),115,80,250); | |
5732d662 | 905 | profy->SetDirectory(0); |
95a0e09b | 906 | fDeltaYres.AddAt(profy,id); |
907 | } | |
908 | if (!profz){ | |
909 | profz=new TProfile(Form("prz%03d",id),Form("Z Residuals for Laser Beam %03d",id),115,80,250); | |
5732d662 | 910 | profz->SetDirectory(0); |
95a0e09b | 911 | fDeltaZres.AddAt(profz,id); |
912 | } | |
913 | for (Int_t irow=158;irow>-1;--irow) { | |
914 | if (vecSec[irow]==-1)continue; //no cluster info | |
915 | Double_t x = vecX[irow]; | |
916 | Double_t ycl = vecClY[irow]; | |
917 | Double_t yfit = vecY1[irow]; | |
918 | Double_t zcl = vecClZ[irow]; | |
919 | Double_t zfit = vecZ1[irow]; | |
5732d662 | 920 | if (profy) |
921 | if (profy->GetEntries()<1000000) | |
922 | profy->Fill(x,yfit-ycl); | |
923 | if (profz) | |
924 | if (profz->GetEntries()<1000000) | |
925 | profz->Fill(x,zfit-zcl); | |
95a0e09b | 926 | } |
ecb7e705 | 927 | //===============================// |
928 | // Fill Fit Parameter Histograms // | |
929 | //===============================// | |
930 | TH1F *pol2InnerY =(TH1F*)fPol2Par2InY.UncheckedAt(id); | |
931 | TH1F *diff1InnerY=(TH1F*)fDiffPar1InY.UncheckedAt(id); | |
932 | TH1F *pol2OuterY =(TH1F*)fPol2Par2OutY.UncheckedAt(id); | |
933 | TH1F *diff1OuterY=(TH1F*)fDiffPar1OutY.UncheckedAt(id); | |
934 | TH1F *pol2InnerZ =(TH1F*)fPol2Par2InZ.UncheckedAt(id); | |
935 | TH1F *diff1InnerZ=(TH1F*)fDiffPar1InZ.UncheckedAt(id); | |
936 | TH1F *pol2OuterZ =(TH1F*)fPol2Par2OutZ.UncheckedAt(id); | |
937 | TH1F *diff1OuterZ=(TH1F*)fDiffPar1OutZ.UncheckedAt(id); | |
938 | //create histograms if the do not already exist | |
939 | if (!pol2InnerY){ | |
940 | pol2InnerY =new TH1F(Form("pol2par2inY%03d",id), | |
941 | Form("2nd derivative from pol2 fit in Y for Laser beam %03d (inner sector)",id), | |
674ee45c | 942 | 500,-.005,.005); |
5732d662 | 943 | pol2InnerY->SetDirectory(0); |
ecb7e705 | 944 | pol2OuterY =new TH1F(Form("pol2par2outY%03d",id), |
945 | Form("2nd derivative from pol2 fit in Y for Laser beam %03d (outer sector)",id), | |
674ee45c | 946 | 500,0.01,.01); |
5732d662 | 947 | pol2OuterY->SetDirectory(0); |
948 | ||
ecb7e705 | 949 | diff1InnerY=new TH1F(Form("diff1inY%03d",id), |
950 | Form("diff of 1st derivative from pol1 and pol2 in Y fit for Laser beam %03d (inner sector)",id), | |
674ee45c | 951 | 500,-.5,.5); |
5732d662 | 952 | diff1InnerY->SetDirectory(0); |
953 | ||
ecb7e705 | 954 | diff1OuterY=new TH1F(Form("diff1outY%03d",id), |
955 | Form("diff of 1st derivative from pol1 and pol2 in Yfit for Laser beam %03d (outer sector)",id), | |
674ee45c | 956 | 500,-1,1); |
5732d662 | 957 | diff1OuterY->SetDirectory(0); |
958 | ||
ecb7e705 | 959 | pol2InnerZ =new TH1F(Form("pol2par2inZ%03d",id), |
960 | Form("2nd derivative from pol2 fit in Z for Laser beam %03d (inner sector)",id), | |
674ee45c | 961 | 500,-.002,.002); |
5732d662 | 962 | pol2InnerZ->SetDirectory(0); |
963 | ||
ecb7e705 | 964 | pol2OuterZ =new TH1F(Form("pol2par2outZ%03d",id), |
965 | Form("2nd derivative from pol2 fit in Z for Laser beam %03d (outer sector)",id), | |
674ee45c | 966 | 500,-.005,.005); |
5732d662 | 967 | pol2OuterZ->SetDirectory(0); |
ecb7e705 | 968 | diff1InnerZ=new TH1F(Form("diff1inZ%03d",id), |
969 | Form("diff of 1st derivative from pol1 and pol2 in Z fit for Laser beam %03d (inner sector)",id), | |
674ee45c | 970 | 500,-.02,.02); |
5732d662 | 971 | diff1InnerZ->SetDirectory(0); |
ecb7e705 | 972 | diff1OuterZ=new TH1F(Form("diff1outZ%03d",id), |
973 | Form("diff of 1st derivative from pol1 and pol2 in Zfit for Laser beam %03d (outer sector)",id), | |
674ee45c | 974 | 500,-.03,.03); |
5732d662 | 975 | diff1OuterZ->SetDirectory(0); |
ecb7e705 | 976 | //add |
977 | fPol2Par2InY.AddAt(pol2InnerY,id); | |
978 | fDiffPar1InY.AddAt(diff1InnerY,id); | |
979 | fPol2Par2OutY.AddAt(pol2OuterY,id); | |
980 | fDiffPar1OutY.AddAt(diff1OuterY,id); | |
981 | fPol2Par2InZ.AddAt(pol2InnerZ,id); | |
982 | fDiffPar1InZ.AddAt(diff1InnerZ,id); | |
983 | fPol2Par2OutZ.AddAt(pol2OuterZ,id); | |
984 | fDiffPar1OutZ.AddAt(diff1OuterZ,id); | |
985 | } | |
986 | //fill histograms | |
987 | pol2InnerY ->Fill(vecy2resInner[2]); | |
988 | pol2OuterY ->Fill(vecy2resOuter[2]); | |
989 | diff1InnerY->Fill(vecy2resInner[1]-vecy1resInner[1]); | |
990 | diff1OuterY->Fill(vecy2resOuter[1]-vecy1resOuter[1]); | |
991 | pol2InnerZ ->Fill(vecz2resInner[2]); | |
992 | pol2OuterZ ->Fill(vecz2resOuter[2]); | |
993 | diff1InnerZ->Fill(vecz2resInner[1]-vecz1resInner[1]); | |
994 | diff1OuterZ->Fill(vecz2resOuter[1]-vecz1resOuter[1]); | |
95a0e09b | 995 | |
263d466a | 996 | |
ecb7e705 | 997 | |
998 | } | |
263d466a | 999 | } |
1000 | ||
e9f38a4b | 1001 | |
1002 | void AliTPCcalibLaser::RefitLaser(Int_t id){ | |
1003 | // | |
1004 | // Refit the track store residuals | |
1005 | // | |
1006 | ||
1007 | AliTPCseed *track = (AliTPCseed*)fTracksTPC.At(id); | |
1008 | AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id); | |
1009 | AliTPCLaserTrack *ltrp = (AliTPCLaserTrack*)fTracksMirror.At(id); | |
263d466a | 1010 | |
1011 | //linear fit model in y and z per sector | |
e9f38a4b | 1012 | static TLinearFitter fy1(2,"hyp1"); |
1013 | static TLinearFitter fz1(2,"hyp1"); | |
263d466a | 1014 | //quadratic fit model in y and z per sector |
1015 | static TLinearFitter fy2(3,"hyp2"); | |
1016 | static TLinearFitter fz2(3,"hyp2"); | |
e9f38a4b | 1017 | static TVectorD vecy2,vecz2,vecy1,vecz1; |
1018 | ||
1019 | const Int_t kMinClusters=20; | |
263d466a | 1020 | Int_t nclusters[72]; |
e9f38a4b | 1021 | // |
1022 | for (Int_t i=0;i<72;++i) nclusters[i]=0; | |
1023 | ||
263d466a | 1024 | for (Int_t i=0;i<160;++i) { |
e9f38a4b | 1025 | AliTPCclusterMI *c=track->GetClusterPointer(i); |
1026 | if (c) nclusters[c->GetDetector()]++; | |
1027 | } | |
263d466a | 1028 | |
e9f38a4b | 1029 | for (Int_t isec=0; isec<72;isec++){ |
1030 | if (nclusters[isec]<kMinClusters) continue; | |
1031 | fy2.ClearPoints(); | |
1032 | fz2.ClearPoints(); | |
1033 | fy1.ClearPoints(); | |
1034 | fz1.ClearPoints(); | |
1035 | // | |
263d466a | 1036 | for (Int_t irow=0;irow<160;++irow) { |
e9f38a4b | 1037 | AliTPCclusterMI *c=track->GetClusterPointer(irow); |
1038 | //if (c && RejectCluster(c)) continue; | |
263d466a | 1039 | Double_t xd = c->GetX()-133.4; // reference x is beteen iroc and oroc |
e9f38a4b | 1040 | if (c&&c->GetDetector()==isec) { |
e9f38a4b | 1041 | Double_t x[2]={xd,xd*xd}; |
1042 | fy2.AddPoint(x,c->GetY()); | |
1043 | fz2.AddPoint(x,c->GetZ()); | |
1044 | // | |
1045 | fy1.AddPoint(x,c->GetY()); | |
1046 | fz1.AddPoint(x,c->GetZ()); | |
1047 | } | |
1048 | } | |
1049 | fy2.Eval(); | |
1050 | fz2.Eval(); | |
1051 | fy1.Eval(); | |
1052 | fz1.Eval(); | |
1053 | fy1.GetParameters(vecy1); | |
1054 | fy2.GetParameters(vecy2); | |
1055 | fz1.GetParameters(vecz1); | |
1056 | fz2.GetParameters(vecz2); | |
263d466a | 1057 | |
e9f38a4b | 1058 | if (fStreamLevel>0){ |
1059 | TTreeSRedirector *cstream = GetDebugStreamer(); | |
1060 | if (cstream){ | |
1061 | Float_t dedx = track->GetdEdx(); | |
1062 | (*cstream)<<"Tracklet"<< | |
1063 | "LTr.="<<ltrp<< | |
1064 | "Tr.="<<param<< | |
1065 | "sec="<<isec<< | |
1066 | "ncl="<<nclusters[isec]<< | |
1067 | "dedx="<<dedx<< | |
1068 | "dedx="<<dedx<< | |
1069 | "vy1.="<<&vecy1<< | |
1070 | "vy2.="<<&vecy2<< | |
1071 | "vz1.="<<&vecz1<< | |
1072 | "vz2.="<<&vecz2<< | |
1073 | "\n"; | |
1074 | } | |
1075 | } | |
1076 | } | |
1077 | // | |
1078 | // | |
1079 | // | |
263d466a | 1080 | // for (Int_t irow=0;irow<160;++irow) { |
e9f38a4b | 1081 | // AliTPCclusterMI *c=track->GetClusterPointer(irow); |
1082 | // if (c && RejectCluster(c)) continue; | |
1083 | // if (c&&c->GetDetector()==isec) { | |
1084 | // Double_t x[2]={c->GetX(),c->GetX()*c->GetX()}; | |
1085 | // fy2.AddPoint(&x,c->GetY()); | |
1086 | // fz2.AddPoint(&x,c->GetZ()); | |
1087 | // // | |
1088 | // fy1.AddPoint(&x,c->GetY()); | |
1089 | // fz1.AddPoint(&x,c->GetZ()); | |
1090 | // } | |
263d466a | 1091 | // } |
1092 | ||
e9f38a4b | 1093 | } |
1094 | ||
1095 | ||
1fd56785 | 1096 | void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield,Int_t minEntries){ |
2b35e8e6 | 1097 | // |
1fd56785 | 1098 | // Dump information about laser beams |
1099 | // isOK variable indicates usability of the beam | |
1100 | // Beam is not usable if: | |
1101 | // a. No entries in range (krmsCut0) | |
1102 | // b. Big sperad (krmscut1) | |
1103 | // c. RMSto fit sigma bigger then (kmultiCut) | |
1104 | // d. Too big angular spread | |
1105 | // | |
1106 | ||
1107 | const Float_t krmsCut0=0.001; | |
1108 | const Float_t krmsCut1=0.16; | |
1109 | const Float_t kmultiCut=2; | |
1110 | const Float_t kcutP0=0.002; | |
2b35e8e6 | 1111 | // |
1112 | AliTPCcalibLaser *laser = this; | |
1113 | TTreeSRedirector *pcstream = new TTreeSRedirector("laserMean.root"); | |
1114 | TF1 fg("fg","gaus"); | |
1115 | // | |
1116 | // | |
1117 | for (Int_t id=0; id<336; id++){ | |
1fd56785 | 1118 | Bool_t isOK=kTRUE; |
2b35e8e6 | 1119 | TH1F * hisphi = (TH1F*)laser->fDeltaPhi.At(id); |
1120 | TH1F * hisphiP = (TH1F*)laser->fDeltaPhiP.At(id); | |
1121 | TH1F * hisZ = (TH1F*)laser->fDeltaZ.At(id); | |
1fd56785 | 1122 | TH1F * hisS = (TH1F*)laser->fSignals.At(id); |
2b35e8e6 | 1123 | if (!hisphi) continue;; |
1124 | Double_t entries = hisphi->GetEntries(); | |
1fd56785 | 1125 | if (entries<minEntries) continue; |
2b35e8e6 | 1126 | // |
1127 | AliTPCLaserTrack *ltrp = (AliTPCLaserTrack*)fTracksMirror.At(id); | |
1128 | if (!ltrp) { | |
1129 | AliTPCLaserTrack::LoadTracks(); | |
1130 | ltrp =(AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(id); | |
1131 | } | |
1132 | Float_t meanphi = hisphi->GetMean(); | |
1133 | Float_t rmsphi = hisphi->GetRMS(); | |
1fd56785 | 1134 | // |
2b35e8e6 | 1135 | Float_t meanphiP = hisphiP->GetMean(); |
1136 | Float_t rmsphiP = hisphiP->GetRMS(); | |
1137 | Float_t meanZ = hisZ->GetMean(); | |
1138 | Float_t rmsZ = hisZ->GetRMS(); | |
1139 | hisphi->Fit(&fg,"","",hisphi->GetMean()-4*hisphi->GetRMS(),hisphi->GetMean()+4*hisphi->GetRMS()); | |
263d466a | 1140 | Double_t gphi1 = fg.GetParameter(1); |
1141 | Double_t gphi2 = fg.GetParameter(2); | |
2b35e8e6 | 1142 | hisphiP->Fit(&fg,"","",hisphiP->GetMean()-4*hisphiP->GetRMS(),hisphiP->GetMean()+4*hisphiP->GetRMS()); |
263d466a | 1143 | Double_t gphiP1 = fg.GetParameter(1); |
1144 | Double_t gphiP2 = fg.GetParameter(2); | |
2b35e8e6 | 1145 | hisZ->Fit(&fg,"","",hisZ->GetMean()-4*hisZ->GetRMS(),hisZ->GetMean()+4*hisZ->GetRMS()); |
263d466a | 1146 | Double_t gz1 = fg.GetParameter(1); |
1147 | Double_t gz2 = fg.GetParameter(2); | |
1fd56785 | 1148 | // |
1149 | Float_t meanS=hisS->GetMean(); | |
2b35e8e6 | 1150 | // |
1151 | Double_t lxyz[3]; | |
1152 | Double_t lpxyz[3]; | |
1153 | ltrp->GetXYZ(lxyz); | |
1154 | ltrp->GetPxPyPz(lpxyz); | |
1fd56785 | 1155 | |
1156 | if (rmsphi<krmsCut0) isOK=kFALSE; // empty in range - not entries inside | |
1157 | if (rmsphi>krmsCut1) isOK=kFALSE; // empty in range - not entries inside | |
1158 | if (rmsphi>krmsCut0) if (gphi2/rmsphi>kmultiCut) isOK=kFALSE; // multi peak structure | |
1159 | if (gphiP2>kcutP0) isOK=kFALSE; | |
2b35e8e6 | 1160 | // |
1161 | (*pcstream)<<"Mean"<< | |
1fd56785 | 1162 | "isOK="<<isOK<< |
2b35e8e6 | 1163 | "entries="<<entries<< // number of entries |
1164 | "bz="<<bfield<< // bfield | |
1165 | "LTr.="<<ltrp<< // refernece track | |
1fd56785 | 1166 | // |
2b35e8e6 | 1167 | "lx0="<<lxyz[0]<< // reference x |
1168 | "lx1="<<lxyz[1]<< // reference y | |
263d466a | 1169 | "lx2="<<lxyz[2]<< // refernece z |
2b35e8e6 | 1170 | "lpx0="<<lpxyz[0]<< // reference x |
1171 | "lpx1="<<lpxyz[1]<< // reference y | |
263d466a | 1172 | "lpx2="<<lpxyz[2]<< // refernece z |
1fd56785 | 1173 | // |
1174 | "msig="<<meanS<< | |
1175 | // | |
2b35e8e6 | 1176 | "mphi="<<meanphi<< // |
1177 | "rmsphi="<<rmsphi<< // | |
1178 | "gphi1="<<gphi1<< | |
1179 | "gphi2="<<gphi2<< | |
1fd56785 | 1180 | // |
2b35e8e6 | 1181 | "mphiP="<<meanphiP<< // |
1182 | "rmsphiP="<<rmsphiP<< // | |
1183 | "gphiP1="<<gphiP1<< | |
1184 | "gphiP2="<<gphiP2<< | |
1fd56785 | 1185 | // |
2b35e8e6 | 1186 | "meanZ="<<meanZ<< |
1187 | "rmsZ="<<rmsZ<< | |
1188 | "gz1="<<gz1<< | |
1189 | "gz2="<<gz2<< | |
1190 | ||
1191 | "\n"; | |
1192 | } | |
1193 | delete pcstream; | |
1194 | } | |
1195 | ||
1fd56785 | 1196 | |
1197 | ||
1198 | void AliTPCcalibLaser::DumpScanInfo(TTree * chain){ | |
1199 | // | |
1200 | // | |
1201 | // | |
1202 | TTreeSRedirector *pcstream = new TTreeSRedirector("laserScan.root"); | |
1203 | TFile * f = pcstream->GetFile(); | |
1204 | f->mkdir("dirphi"); | |
1205 | f->mkdir("dirphiP"); | |
1206 | f->mkdir("dirZ"); | |
1207 | TF1 fp("p1","pol1"); | |
1208 | // | |
1209 | // | |
1210 | char cut[1000]; | |
1211 | char grname[1000]; | |
1212 | char grnamefull[1000]; | |
1213 | ||
1214 | Double_t mphi[100]; | |
1215 | Double_t mphiP[100]; | |
1216 | Double_t smphi[100]; | |
1217 | Double_t smphiP[100]; | |
1218 | Double_t mZ[100]; | |
1219 | Double_t smZ[100]; | |
1220 | Double_t bz[100]; | |
1221 | Double_t sbz[100]; | |
1222 | // fit parameters | |
1223 | Double_t pphi[3]; | |
1224 | Double_t pphiP[3]; | |
1225 | Double_t pmZ[3]; | |
1226 | // | |
1227 | for (Int_t id=0; id<336; id++){ | |
1228 | // id =205; | |
1229 | sprintf(cut,"isOK&&fId==%d",id); | |
1230 | Int_t entries = chain->Draw("bz",cut,"goff"); | |
1231 | if (entries<3) continue; | |
1232 | AliTPCLaserTrack *ltrp = 0;; | |
1233 | if (!AliTPCLaserTrack::GetTracks()) AliTPCLaserTrack::LoadTracks(); | |
1234 | ltrp =(AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(id); | |
1235 | Double_t lxyz[3]; | |
1236 | Double_t lpxyz[3]; | |
1237 | ltrp->GetXYZ(lxyz); | |
1238 | ltrp->GetPxPyPz(lpxyz); | |
1239 | ||
1240 | chain->Draw("bz",cut,"goff"); | |
1241 | memcpy(bz, chain->GetV1(), entries*sizeof(Double_t)); | |
1242 | chain->Draw("0.01*abs(bz)+0.02",cut,"goff"); | |
1243 | memcpy(sbz, chain->GetV1(), entries*sizeof(Double_t)); | |
1244 | // | |
1245 | chain->Draw("gphi1",cut,"goff"); | |
1246 | memcpy(mphi, chain->GetV1(), entries*sizeof(Double_t)); | |
1247 | chain->Draw("0.05*abs(mphi)+gphi2",cut,"goff"); | |
1248 | memcpy(smphi, chain->GetV1(), entries*sizeof(Double_t)); | |
1249 | // | |
1250 | chain->Draw("gphiP1",cut,"goff"); | |
1251 | memcpy(mphiP, chain->GetV1(), entries*sizeof(Double_t)); | |
1252 | chain->Draw("0.05*abs(mphiP)+gphiP2",cut,"goff"); | |
1253 | memcpy(smphiP, chain->GetV1(), entries*sizeof(Double_t)); | |
1254 | // | |
1255 | chain->Draw("gz1",cut,"goff"); | |
1256 | memcpy(mZ, chain->GetV1(), entries*sizeof(Double_t)); | |
1257 | chain->Draw("0.01*abs(meanZ)+gz2",cut,"goff"); | |
1258 | memcpy(smZ, chain->GetV1(), entries*sizeof(Double_t)); | |
1259 | // | |
1260 | // | |
1261 | sprintf(grnamefull,"Side_%d_Bundle_%d_Rod_%d_Beam_%d", | |
1262 | ltrp->GetSide(), ltrp->GetBundle(), ltrp->GetRod(), ltrp->GetBeam()); | |
1263 | // store data | |
1264 | // phi | |
1265 | f->cd("dirphi"); | |
1266 | TGraphErrors *grphi = new TGraphErrors(entries,bz,mphi,sbz,smphi); | |
1267 | grphi->Draw("a*"); | |
1268 | grphi->Fit(&fp); | |
1269 | pphi[0] = fp.GetParameter(0); // offset | |
1270 | pphi[1] = fp.GetParameter(1); // slope | |
1271 | pphi[2] = TMath::Sqrt(fp.GetChisquare()/(entries-2.)); // normalized chi2 | |
1272 | sprintf(grname,"phi_id%d",id); | |
1273 | grphi->SetName(grname); grphi->SetTitle(grnamefull); | |
1274 | grphi->GetXaxis()->SetTitle("b_{z} (T)"); | |
1275 | grphi->GetYaxis()->SetTitle("#Delta r#phi (cm)"); | |
ae69f16f | 1276 | grphi->SetMaximum(1.2); |
1277 | grphi->SetMinimum(-1.2); | |
1fd56785 | 1278 | grphi->Draw("a*"); |
1279 | ||
1280 | grphi->Write(); | |
1281 | gPad->SaveAs(Form("pic/phi/phi_%s.gif",grnamefull)); | |
1282 | // phiP | |
1283 | f->cd("dirphiP)"); | |
1284 | TGraphErrors *grphiP = new TGraphErrors(entries,bz,mphiP,sbz,smphiP); | |
1285 | grphiP->Draw("a*"); | |
1286 | grphiP->Fit(&fp); | |
1287 | pphiP[0] = fp.GetParameter(0); // offset | |
1288 | pphiP[1] = fp.GetParameter(1); // slope | |
1289 | pphiP[2] = TMath::Sqrt(fp.GetChisquare()/(entries-2.)); // normalized chi2 | |
1290 | sprintf(grname,"phiP_id%d",id); | |
1291 | grphiP->SetName(grname); grphiP->SetTitle(grnamefull); | |
1292 | grphiP->GetXaxis()->SetTitle("b_{z} (T)"); | |
1293 | grphiP->GetYaxis()->SetTitle("#Delta #phi (rad)"); | |
ae69f16f | 1294 | grphiP->SetMaximum(pphiP[0]+0.005); |
1295 | grphiP->SetMinimum(pphiP[0]-0.005); | |
1fd56785 | 1296 | |
1297 | gPad->SaveAs(Form("pic/phiP/phiP_%s.gif",grnamefull)); | |
1298 | grphiP->Write(); | |
1299 | // | |
1300 | //Z | |
1301 | f->cd("dirZ"); | |
1302 | TGraphErrors *grmZ = new TGraphErrors(entries,bz,mZ,sbz,smZ); | |
1303 | grmZ->Draw("a*"); | |
1304 | grmZ->Fit(&fp); | |
1305 | pmZ[0] = fp.GetParameter(0); // offset | |
1306 | pmZ[1] = fp.GetParameter(1); // slope | |
1307 | pmZ[2] = TMath::Sqrt(fp.GetChisquare()/(entries-2.)); // normalized chi2 | |
1308 | sprintf(grname,"mZ_id%d",id); | |
1309 | grmZ->SetName(grname); grmZ->SetTitle(grnamefull); | |
1310 | grmZ->GetXaxis()->SetTitle("b_{z} (T)"); | |
1311 | grmZ->GetYaxis()->SetTitle("#Delta z (cm)"); | |
1312 | ||
1313 | gPad->SaveAs(Form("pic/z/z_%s.gif",grnamefull)); | |
1314 | grmZ->Write(); | |
1315 | ||
1316 | ||
1317 | for (Int_t ientry=0; ientry<entries; ientry++){ | |
1318 | (*pcstream)<<"Mean"<< | |
1319 | "id="<<id<< | |
1320 | "LTr.="<<ltrp<< | |
1321 | "entries="<<entries<< | |
1322 | "bz="<<bz[ientry]<< | |
1323 | "lx0="<<lxyz[0]<< // reference x | |
1324 | "lx1="<<lxyz[1]<< // reference y | |
1325 | "lx2="<<lxyz[2]<< // refernece z | |
1326 | "lpx0="<<lpxyz[0]<< // reference x | |
1327 | "lpx1="<<lpxyz[1]<< // reference y | |
1328 | "lpx2="<<lpxyz[2]<< // refernece z | |
1329 | //values | |
1330 | "gphi1="<<mphi[ientry]<< // mean - from gaus fit | |
1331 | "pphi0="<<pphi[0]<< // offset | |
1332 | "pphi1="<<pphi[1]<< // mean | |
1333 | "pphi2="<<pphi[2]<< // norm chi2 | |
1334 | // | |
1335 | "gphiP1="<<mphiP[ientry]<< // mean - from gaus fit | |
1336 | "pphiP0="<<pphiP[0]<< // offset | |
1337 | "pphiP1="<<pphiP[1]<< // mean | |
1338 | "pphiP2="<<pphiP[2]<< // norm chi2 | |
1339 | // | |
1340 | "gz1="<<mZ[ientry]<< | |
1341 | "pmZ0="<<pmZ[0]<< // offset | |
1342 | "pmZ1="<<pmZ[1]<< // mean | |
1343 | "pmZ2="<<pmZ[2]<< // norm chi2 | |
1344 | "\n"; | |
1345 | } | |
1346 | } | |
1347 | ||
1348 | delete pcstream; | |
1349 | ||
1350 | } | |
1351 | ||
1352 | ||
c6914c83 | 1353 | void AliTPCcalibLaser::Analyze(){ |
1354 | // | |
1355 | // | |
1356 | // | |
1357 | } | |
1358 | ||
1359 | ||
c03e3250 | 1360 | Long64_t AliTPCcalibLaser::Merge(TCollection *li) { |
c6914c83 | 1361 | |
c03e3250 | 1362 | TIterator* iter = li->MakeIterator(); |
1363 | AliTPCcalibLaser* cal = 0; | |
c6914c83 | 1364 | |
c03e3250 | 1365 | while ((cal = (AliTPCcalibLaser*)iter->Next())) { |
1366 | if (!cal->InheritsFrom(AliTPCcalibLaser::Class())) { | |
1367 | Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName()); | |
1368 | return -1; | |
1369 | } | |
1370 | // | |
1371 | // fHistNTracks->Add(cal->fHistNTracks); | |
1372 | // fClusters->Add(cal->fClusters); | |
1373 | // fModules->Add(cal->fModules); | |
1374 | // fHistPt->Add(cal->fHistPt); | |
1375 | // fPtResolution->Add(cal->fPtResolution); | |
1376 | // fDeDx->Add(cal->fDeDx); | |
263d466a | 1377 | |
c03e3250 | 1378 | |
1379 | TH1F *h=0x0; | |
1380 | TH1F *hm=0x0; | |
95a0e09b | 1381 | TProfile *hp=0x0; |
1382 | TProfile *hpm=0x0; | |
c03e3250 | 1383 | |
1384 | for (Int_t id=0; id<336; id++){ | |
1385 | // merge fDeltaZ histograms | |
263d466a | 1386 | hm = (TH1F*)cal->fDeltaZ.At(id); |
c03e3250 | 1387 | h = (TH1F*)fDeltaZ.At(id); |
1388 | if (!h) { | |
1389 | h=new TH1F(Form("hisdz%d",id),Form("hisdz%d",id),1000,-10,10); | |
5732d662 | 1390 | h->SetDirectory(0); |
c03e3250 | 1391 | fDeltaZ.AddAt(h,id); |
1392 | } | |
1393 | if (hm) h->Add(hm); | |
1394 | // merge fDeltaPhi histograms | |
263d466a | 1395 | hm = (TH1F*)cal->fDeltaPhi.At(id); |
c03e3250 | 1396 | h = (TH1F*)fDeltaPhi.At(id); |
1397 | if (!h) { | |
1398 | h= new TH1F(Form("hisdphi%d",id),Form("hisdphi%d",id),1000,-1,1); | |
5732d662 | 1399 | h->SetDirectory(0); |
c03e3250 | 1400 | fDeltaPhi.AddAt(h,id); |
1401 | } | |
1402 | if (hm) h->Add(hm); | |
1403 | // merge fDeltaPhiP histograms | |
263d466a | 1404 | hm = (TH1F*)cal->fDeltaPhiP.At(id); |
c03e3250 | 1405 | h = (TH1F*)fDeltaPhiP.At(id); |
1406 | if (!h) { | |
1407 | h=new TH1F(Form("hisdphiP%d",id),Form("hisdphiP%d",id),1000,-0.01,0.01); | |
5732d662 | 1408 | h->SetDirectory(0); |
c03e3250 | 1409 | fDeltaPhiP.AddAt(h,id); |
1410 | } | |
1411 | if (hm) h->Add(hm); | |
1412 | // merge fSignals histograms | |
263d466a | 1413 | hm = (TH1F*)cal->fSignals.At(id); |
c03e3250 | 1414 | h = (TH1F*)fSignals.At(id); |
1415 | if (!h) { | |
1416 | h=new TH1F(Form("hisSignal%d",id),Form("hisSignal%d",id),1000,0,1000); | |
5732d662 | 1417 | h->SetDirectory(0); |
c03e3250 | 1418 | fSignals.AddAt(h,id); |
1419 | } | |
263d466a | 1420 | if (hm) h->Add(hm); |
95a0e09b | 1421 | // |
1422 | // | |
1423 | // merge ProfileY histograms | |
1424 | hpm = (TProfile*)cal->fDeltaYres.At(id); | |
1425 | hp = (TProfile*)fDeltaYres.At(id); | |
1426 | if (!hp) { | |
5732d662 | 1427 | hp=new TProfile(Form("pry%03d",id),Form("Y Residuals for Laser Beam %03d",id),115,80,250); |
1428 | hp->SetDirectory(0); | |
95a0e09b | 1429 | fDeltaYres.AddAt(hp,id); |
1430 | } | |
1431 | if (hpm) hp->Add(hpm); | |
1432 | // | |
1433 | hpm = (TProfile*)cal->fDeltaZres.At(id); | |
1434 | hp = (TProfile*)fDeltaZres.At(id); | |
1435 | if (!hp) { | |
5732d662 | 1436 | hp=new TProfile(Form("prz%03d",id),Form("Z Residuals for Laser Beam %03d",id),115,80,250); |
1437 | hp->SetDirectory(0); | |
95a0e09b | 1438 | fDeltaZres.AddAt(hp,id); |
1439 | } | |
1440 | if (hpm) hp->Add(hpm); | |
1441 | // | |
1442 | // | |
ecb7e705 | 1443 | //merge fit param histograms |
1444 | //local hists | |
1445 | TH1F *pol2InnerY =(TH1F*)fPol2Par2InY.UncheckedAt(id); | |
1446 | TH1F *diff1InnerY=(TH1F*)fDiffPar1InY.UncheckedAt(id); | |
1447 | TH1F *pol2OuterY =(TH1F*)fPol2Par2OutY.UncheckedAt(id); | |
1448 | TH1F *diff1OuterY=(TH1F*)fDiffPar1OutY.UncheckedAt(id); | |
1449 | TH1F *pol2InnerZ =(TH1F*)fPol2Par2InZ.UncheckedAt(id); | |
1450 | TH1F *diff1InnerZ=(TH1F*)fDiffPar1InZ.UncheckedAt(id); | |
1451 | TH1F *pol2OuterZ =(TH1F*)fPol2Par2OutZ.UncheckedAt(id); | |
1452 | TH1F *diff1OuterZ=(TH1F*)fDiffPar1OutZ.UncheckedAt(id); | |
1453 | //hists to merge | |
1454 | TH1F *pol2InnerYm =(TH1F*)cal->fPol2Par2InY.UncheckedAt(id); | |
1455 | TH1F *diff1InnerYm=(TH1F*)cal->fDiffPar1InY.UncheckedAt(id); | |
1456 | TH1F *pol2OuterYm =(TH1F*)cal->fPol2Par2OutY.UncheckedAt(id); | |
1457 | TH1F *diff1OuterYm=(TH1F*)cal->fDiffPar1OutY.UncheckedAt(id); | |
1458 | TH1F *pol2InnerZm =(TH1F*)cal->fPol2Par2InZ.UncheckedAt(id); | |
1459 | TH1F *diff1InnerZm=(TH1F*)cal->fDiffPar1InZ.UncheckedAt(id); | |
1460 | TH1F *pol2OuterZm =(TH1F*)cal->fPol2Par2OutZ.UncheckedAt(id); | |
1461 | TH1F *diff1OuterZm=(TH1F*)cal->fDiffPar1OutZ.UncheckedAt(id); | |
1462 | //create histos if they do not exist | |
1463 | if (!pol2InnerY){ | |
1464 | pol2InnerY =new TH1F(Form("pol2par2inY%03d",id), | |
1465 | Form("2nd derivative from pol2 fit in Y for Laser beam %03d (inner sector)",id), | |
674ee45c | 1466 | 500,-.005,.005); |
5732d662 | 1467 | pol2InnerY->SetDirectory(0); |
1468 | ||
ecb7e705 | 1469 | pol2OuterY =new TH1F(Form("pol2par2outY%03d",id), |
1470 | Form("2nd derivative from pol2 fit in Y for Laser beam %03d (outer sector)",id), | |
674ee45c | 1471 | 500,0.01,.01); |
5732d662 | 1472 | pol2OuterY->SetDirectory(0); |
ecb7e705 | 1473 | diff1InnerY=new TH1F(Form("diff1inY%03d",id), |
1474 | Form("diff of 1st derivative from pol1 and pol2 in Y fit for Laser beam %03d (inner sector)",id), | |
674ee45c | 1475 | 500,-.5,.5); |
ecb7e705 | 1476 | diff1OuterY=new TH1F(Form("diff1outY%03d",id), |
1477 | Form("diff of 1st derivative from pol1 and pol2 in Yfit for Laser beam %03d (outer sector)",id), | |
674ee45c | 1478 | 500,-1,1); |
5732d662 | 1479 | diff1InnerY->SetDirectory(0); |
1480 | ||
1481 | ||
ecb7e705 | 1482 | pol2InnerZ =new TH1F(Form("pol2par2inZ%03d",id), |
1483 | Form("2nd derivative from pol2 fit in Z for Laser beam %03d (inner sector)",id), | |
674ee45c | 1484 | 500,-.002,.002); |
5732d662 | 1485 | pol2InnerZ->SetDirectory(0); |
1486 | ||
ecb7e705 | 1487 | pol2OuterZ =new TH1F(Form("pol2par2outZ%03d",id), |
1488 | Form("2nd derivative from pol2 fit in Z for Laser beam %03d (outer sector)",id), | |
674ee45c | 1489 | 500,-.005,.005); |
5732d662 | 1490 | pol2OuterZ->SetDirectory(0); |
ecb7e705 | 1491 | diff1InnerZ=new TH1F(Form("diff1inZ%03d",id), |
1492 | Form("diff of 1st derivative from pol1 and pol2 in Z fit for Laser beam %03d (inner sector)",id), | |
674ee45c | 1493 | 500,-.02,.02); |
5732d662 | 1494 | diff1InnerZ->SetDirectory(0); |
ecb7e705 | 1495 | diff1OuterZ=new TH1F(Form("diff1outZ%03d",id), |
1496 | Form("diff of 1st derivative from pol1 and pol2 in Zfit for Laser beam %03d (outer sector)",id), | |
674ee45c | 1497 | 500,-.03,.03); |
5732d662 | 1498 | diff1OuterZ->SetDirectory(0); |
ecb7e705 | 1499 | } |
1500 | pol2InnerYm =(TH1F*)cal->fPol2Par2InY.UncheckedAt(id); | |
1501 | diff1InnerYm=(TH1F*)cal->fDiffPar1InY.UncheckedAt(id); | |
1502 | pol2OuterYm =(TH1F*)cal->fPol2Par2OutY.UncheckedAt(id); | |
1503 | diff1OuterYm=(TH1F*)cal->fDiffPar1OutY.UncheckedAt(id); | |
1504 | pol2InnerZm =(TH1F*)cal->fPol2Par2InZ.UncheckedAt(id); | |
1505 | diff1InnerZm=(TH1F*)cal->fDiffPar1InZ.UncheckedAt(id); | |
1506 | pol2OuterZm =(TH1F*)cal->fPol2Par2OutZ.UncheckedAt(id); | |
1507 | diff1OuterZm=(TH1F*)cal->fDiffPar1OutZ.UncheckedAt(id); | |
c03e3250 | 1508 | } |
1509 | } | |
1510 | return 0; | |
1511 | } | |
1512 | ||
1513 | ||
1514 | ||
263d466a | 1515 | |
c03e3250 | 1516 | /* |
1517 | gSystem->Load("libSTAT.so") | |
1518 | TStatToolkit toolkit; | |
1519 | Double_t chi2; | |
1520 | TVectorD fitParam; | |
1521 | TMatrixD covMatrix; | |
1522 | Int_t npoints; | |
ae69f16f | 1523 | TCut cutA("entries>2&&pphi2<3&&abs(gphiP1-pphiP0)<0.003&&abs(gz1)<6"); |
c03e3250 | 1524 | |
1525 | ||
1526 | TString fstring=""; | |
1fd56785 | 1527 | // |
1528 | fstring+="(abs(LTr.fP[1]/250)^3-1)*bz++"; //1 | |
1529 | fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*LTr.fP[2]++"; //2 | |
1530 | fstring+="(abs(LTr.fP[1]/250)^1-1)*bz++"; //3 | |
1531 | fstring+="(abs(LTr.fP[1]/250)-1)*bz*LTr.fP[2]++"; //4 | |
1532 | // | |
1533 | fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*sin(atan2(lx1,lx0))++" //5 | |
1534 | fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*sin(atan2(lx1,lx0))*LTr.fP[2]++" //6 | |
1535 | fstring+="(abs(LTr.fP[1]/250)-1)*bz*sin(atan2(lx1,lx0))++" //7 | |
1536 | fstring+="(abs(LTr.fP[1]/250)-1)*bz*sin(atan2(lx1,lx0))*LTr.fP[2]++" //8 | |
1537 | // | |
1538 | fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*cos(atan2(lx1,lx0))++" //9 | |
1539 | fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*cos(atan2(lx1,lx0))*LTr.fP[2]++" //10 | |
1540 | fstring+="(abs(LTr.fP[1]/250)-1)*bz*cos(atan2(lx1,lx0))++" //11 | |
1541 | fstring+="(abs(LTr.fP[1]/250)-1)*bz*cos(atan2(lx1,lx0))*LTr.fP[2]++" //12 | |
c5595838 | 1542 | |
c5595838 | 1543 | |
c5595838 | 1544 | |
c5595838 | 1545 | |
1fd56785 | 1546 | TString *strq0 = toolkit.FitPlane(treeT,"gphi1-pphi0",fstring->Data(), "fSide==1"+cutA, chi2,npoints,fitParam,covMatrix); |
c03e3250 | 1547 | |
1fd56785 | 1548 | treeT->SetAlias("fit",strq0->Data()); |
1549 | ||
c03e3250 | 1550 | |
1fd56785 | 1551 | TString *strqP = toolkit.FitPlane(treeT,"1000*(gphiP1-pphiP0)",fstring->Data(), "fSide==1"+cutA, chi2,npoints,fitParam,covMatrix); |
c5595838 | 1552 | |
1fd56785 | 1553 | treeT->SetAlias("fitP",strqP->Data()); |
c5595838 | 1554 | |
1555 | ||
ae69f16f | 1556 | TString *strqDrift = toolkit.FitPlane(treeT,"gz1","LTr.fP[1]++(1-2*(fSide==1))++lx1", cutA, chi2,npoints,fitParam,covMatrix); |
1557 | treeT->SetAlias("fitD",strqDrift->Data()); | |
c5595838 | 1558 | |
c03e3250 | 1559 | |
ae69f16f | 1560 | treeT->Draw("fit:LTr.fP[1]","abs(bz+0.4)<0.05"+cutA,""); |
1561 | { | |
1562 | for (Int_t i=0; i<6;i++){ | |
1563 | treeT->SetLineColor(i+2); | |
1564 | treeT->SetMarkerSize(1); | |
1565 | treeT->SetMarkerStyle(22+i); | |
1566 | treeT->SetMarkerColor(i+2); | |
c03e3250 | 1567 | |
ae69f16f | 1568 | treeT->Draw("fit:LTr.fP[1]",Form("abs(bz+0.4)<0.05&fRod==%d",i)+cutA,"same"); |
1569 | } | |
1570 | } | |
c03e3250 | 1571 | */ |
68ff0583 | 1572 | |
1573 | ||
1574 | ||
1575 | /* | |
1576 | TTree * tree = (TTree*)f.Get("FitModels"); | |
1577 | ||
1578 | TEventList listLFit0("listLFit0","listLFit0"); | |
1579 | TEventList listLFit1("listLFit1","listLFit1"); | |
68ff0583 | 1580 | tree->Draw(">>listLFit0","seed.fdEdx<200&&seed.fdEdx>40"); |
1581 | tree->SetEventList(&listLFit0); | |
1582 | ||
1583 | ||
95a0e09b | 1584 | |
1585 | ||
1586 | gSystem->Load("libSTAT.so") | |
1587 | TStatToolkit toolkit; | |
1588 | Double_t chi2; | |
1589 | TVectorD fitParam; | |
1590 | TMatrixD covMatrix; | |
1591 | Int_t npoints; | |
1592 | ||
1593 | chain->SetAlias("dp","((Cl.fPad-int(Cl.fPad))*pi)"); | |
1594 | chain->SetAlias("dt","((Cl.fTimeBin-int(Cl.fTimeBin))*pi)"); | |
1595 | ||
1596 | ||
1597 | TString fstring=""; | |
1598 | fstring+="cos(dp)++"; | |
1599 | fstring+="sin(dp)++"; | |
1600 | fstring+="cos(dt)++"; | |
1601 | fstring+="sin(dt)++"; | |
1602 | ||
1603 | TString *str = toolkit.FitPlane(chain,"Cl.fZ-TrZInOut.fElements",fstring->Data(), "Cl.fDetector>35", chi2,npoints,fitParam,covMatrix,-1,0,200); | |
1604 | ||
1605 | ||
68ff0583 | 1606 | */ |