1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
20 // 1. The laser track is associated with the mirror
21 // see function FindMirror
23 // 2. The laser track is accepted for the analysis under certain condition
24 // (see function Accpet laser)
26 // 3. The drift velocity and jitter is calculated event by event
27 // (see function drift velocity)
31 // To make laser scan the user interaction neccessary
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")
42 // laser track clasification;
44 TCut cutT("cutT","abs(Tr.fP[3])<0.06");
45 TCut cutPt("cutPt","abs(Tr.fP[4])<0.1");
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");
55 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
56 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
57 AliXRDPROOFtoolkit tool;
58 TChain * chain = tool.MakeChain("laserScan.txt","Mean",0,10200);
60 AliTPCcalibLaser::DumpScanInfo(chain)
61 TFile fscan("laserScan.root")
62 TTree * treeT = (TTree*)fscan.Get("Mean")
68 #include "TLinearFitter.h"
69 #include "AliTPCcalibLaser.h"
70 #include "AliExternalTrackParam.h"
71 #include "AliESDEvent.h"
72 #include "AliESDfriend.h"
73 #include "AliESDtrack.h"
74 #include "AliTPCTracklet.h"
77 #include "TTreeStream.h"
80 #include "TGraphErrors.h"
81 #include "AliTPCclusterMI.h"
82 #include "AliTPCseed.h"
83 #include "AliTracker.h"
84 #include "TClonesArray.h"
88 #include "TTreeStream.h"
91 #include "AliTPCLaserTrack.h"
95 ClassImp(AliTPCcalibLaser)
97 AliTPCcalibLaser::AliTPCcalibLaser():
103 fTracksEsdParam(336),
105 fDeltaZ(336), // array of histograms of delta z for each track
106 fDeltaPhi(336), // array of histograms of delta z for each track
107 fDeltaPhiP(336), // array of histograms of delta z for each track
108 fSignals(336), // array of dedx signals
109 fFitAside(new TVectorD(3)), // drift fit - A side
110 fFitCside(new TVectorD(3)), // drift fit - C- side
116 fTracksEsdParam.SetOwner(kTRUE);
119 AliTPCcalibLaser::AliTPCcalibLaser(const Text_t *name, const Text_t *title):
125 fTracksEsdParam(336),
127 fDeltaZ(336), // array of histograms of delta z for each track
128 fDeltaPhi(336), // array of histograms of delta z for each track
129 fDeltaPhiP(336), // array of histograms of delta z for each track
130 fSignals(336), // array of dedx signals
131 fFitAside(new TVectorD(3)), // drift fit - A side
132 fFitCside(new TVectorD(3)), // drift fit - C- side
140 fTracksEsdParam.SetOwner(kTRUE);
143 AliTPCcalibLaser::~AliTPCcalibLaser() {
151 void AliTPCcalibLaser::Process(AliESDEvent * event) {
154 // Loop over tracks and call Process function
160 fESDfriend=static_cast<AliESDfriend*>(fESD->FindListObject("AliESDfriend"));
166 fTracksEsdParam.Delete();
168 Int_t n=fESD->GetNumberOfTracks();
169 Int_t run = fESD->GetRunNumber();
171 for (Int_t i=0;i<n;++i) {
172 AliESDfriendTrack *friendTrack=fESDfriend->GetTrack(i);
173 AliESDtrack *track=fESD->GetTrack(i);
174 TObject *calibObject=0;
176 for (Int_t j=0;(calibObject=friendTrack->GetCalibObject(j));++j)
177 if ((seed=dynamic_cast<AliTPCseed*>(calibObject)))
179 if (track&&seed) FindMirror(track,seed);
186 for (Int_t id=0; id<336; id++){
189 if (!fTracksEsdParam.At(id)) continue;
196 void AliTPCcalibLaser::MakeDistHisto(){
200 for (Int_t id=0; id<336; id++){
203 if (!fTracksEsdParam.At(id)) continue;
204 if (!AcceptLaser(id)) continue;
207 TH1F * hisdz = (TH1F*)fDeltaZ.At(id);
208 TH1F * hisdphi = (TH1F*)fDeltaPhi.At(id);
209 TH1F * hisdphiP = (TH1F*)fDeltaPhiP.At(id);
210 TH1F * hisSignal = (TH1F*)fSignals.At(id);
213 hisdz = new TH1F(Form("hisdz%d",id),Form("hisdz%d",id),1000,-10,10);
214 hisdz->SetDirectory(0);
215 fDeltaZ.AddAt(hisdz,id);
217 hisdphi = new TH1F(Form("hisdphi%d",id),Form("hisdphi%d",id),1000,-1,1);
218 hisdphi->SetDirectory(0);
219 fDeltaPhi.AddAt(hisdphi,id);
221 hisdphiP = new TH1F(Form("hisdphiP%d",id),Form("hisdphiP%d",id),1000,-0.01,0.01);
222 hisdphiP->SetDirectory(0);
223 fDeltaPhiP.AddAt(hisdphiP,id);
224 hisSignal = new TH1F(Form("hisSignal%d",id),Form("hisSignal%d",id),1000,0,1000);
225 hisSignal->SetDirectory(0);
226 fSignals.AddAt(hisSignal,id);
229 AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
230 AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
231 AliESDtrack *track = (AliESDtrack*)fTracksEsd.At(id);
240 param->GetPxPyPz(pxyz);
242 ltrp->GetPxPyPz(lpxyz);
244 Float_t dz = param->GetZ()-ltrp->GetZ();
245 Float_t dphi = (TMath::ATan2(xyz[1],xyz[0])- TMath::ATan2(lxyz[1],lxyz[0]))*254.;
246 Float_t dphiP = param->GetParameter()[2]-ltrp->GetParameter()[2];
247 if (hisdz) hisdz->Fill(dz);
248 if (hisdphi) hisdphi->Fill(dphi);
249 if (hisdphiP) hisdphiP->Fill(dphiP);
250 if (hisSignal) hisSignal->Fill(track->GetTPCsignal());
254 void AliTPCcalibLaser::FitDriftV(){
256 // Fit drift velocity - linear approximation in the z and global y
258 static TLinearFitter fdriftA(3,"hyp2");
259 static TLinearFitter fdriftC(3,"hyp2");
260 fdriftA.ClearPoints();
261 fdriftC.ClearPoints();
263 for (Int_t id=0; id<336; id++){
264 if (!fTracksEsdParam.At(id)) continue;
265 if (!AcceptLaser(id)) continue;
266 AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
267 AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
273 param->GetPxPyPz(pxyz);
275 ltrp->GetPxPyPz(lpxyz);
276 Double_t xxx[2] = {lxyz[2],lxyz[1]};
277 if (ltrp->GetSide()==0){
278 fdriftA.AddPoint(xxx,xyz[2],1);
280 fdriftC.AddPoint(xxx,xyz[2],1);
288 if (fdriftA.GetNpoints()>10){
290 fdriftA.EvalRobust(0.8);
291 fdriftA.GetParameters(*fFitAside);
292 npointsA= fdriftA.GetNpoints();
293 chi2A = fdriftA.GetChisquare()/fdriftA.GetNpoints();
295 if (fdriftC.GetNpoints()>10){
297 fdriftC.EvalRobust(0.8);
298 fdriftC.GetParameters(*fFitCside);
299 npointsC= fdriftC.GetNpoints();
300 chi2C = fdriftC.GetChisquare()/fdriftC.GetNpoints();
304 TTreeSRedirector *cstream = GetDebugStreamer();
305 Int_t time = fESD->GetTimeStamp();
307 (*cstream)<<"driftv"<<
308 "driftA.="<<fFitAside<<
309 "driftC.="<<fFitCside<<
322 Bool_t AliTPCcalibLaser::AcceptLaser(Int_t id){
327 TCut cutT("cutT","abs(Tr.fP[3])<0.06");
328 TCut cutPt("cutPt","abs(Tr.fP[4])<0.1");
329 TCut cutN("cutN","fTPCncls>70");
330 TCut cutP("cutP","abs(atan2(x1,x0)-atan2(lx1,lx0))<0.03")
331 TCut cutA = cutT+cutPt+cutP;
333 AliExternalTrackParam *param =(AliExternalTrackParam*)fTracksEsdParam.At(id);
334 AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
335 AliESDtrack *track = (AliESDtrack*)fTracksEsd.At(id);
337 if (TMath::Abs(param->GetParameter()[4])>0.03) return kFALSE;
338 if (TMath::Abs(param->GetParameter()[3])>0.06) return kFALSE;
339 if (TMath::Abs(param->GetParameter()[2]-ltrp->GetParameter()[2])>0.06) return kFALSE;
340 if (TMath::Abs(param->GetParameter()[1]-ltrp->GetParameter()[1])>10) return kFALSE;
344 if (TMath::Abs(track->GetTPCsignal())<20) return kFALSE;
345 if (TMath::Abs(track->GetTPCsignal())>800) return kFALSE;
350 Int_t AliTPCcalibLaser::FindMirror(AliESDtrack *track, AliTPCseed *seed){
352 // Find corresponding mirror
353 // add the corresponding tracks
355 Float_t kRadius0 = 252;
356 Float_t kRadius = 253.4;
357 if (!track->GetOuterParam()) return -1;
358 AliExternalTrackParam param(*(track->GetOuterParam()));
359 AliTracker::PropagateTrackTo(¶m,kRadius0,0.10566,3,kTRUE);
360 AliTracker::PropagateTrackTo(¶m,kRadius,0.10566,0.1,kTRUE);
361 AliTPCLaserTrack ltr;
362 AliTPCLaserTrack *ltrp=0x0;
364 Int_t id = AliTPCLaserTrack::IdentifyTrack(¶m);
365 if (id!=-1 && (AliTPCLaserTrack::GetTracks()->UncheckedAt(id)))
366 ltrp=(AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(id);
373 Float_t radius=TMath::Abs(ltrp->GetX());
374 AliTracker::PropagateTrackTo(¶m,radius,0.10566,0.01,kTRUE);
376 if (!fTracksMirror.At(id)) fTracksMirror.AddAt(ltrp,id);
377 fTracksEsdParam.AddAt(param.Clone(),id);
378 fTracksEsd.AddAt(track,id);
379 fTracksTPC.AddAt(seed,id);
388 void AliTPCcalibLaser::DumpLaser(Int_t id) {
390 // Dump Laser info to the tree
392 AliESDtrack *track = (AliESDtrack*)fTracksEsd.At(id);
393 AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
394 AliTPCLaserTrack *ltrp = ( AliTPCLaserTrack*)fTracksMirror.At(id);
403 param->GetPxPyPz(pxyz);
405 ltrp->GetPxPyPz(lpxyz);
408 TTreeSRedirector *cstream = GetDebugStreamer();
409 Int_t time = fESD->GetTimeStamp();
410 Bool_t accept = AcceptLaser(id);
412 (*cstream)<<"Track"<<
416 "driftA.="<<fFitAside<<
417 "driftC.="<<fFitCside<<
442 void AliTPCcalibLaser::RefitLaser(Int_t id){
444 // Refit the track store residuals
447 AliTPCseed *track = (AliTPCseed*)fTracksTPC.At(id);
448 AliExternalTrackParam *param=(AliExternalTrackParam*)fTracksEsdParam.At(id);
449 AliTPCLaserTrack *ltrp = (AliTPCLaserTrack*)fTracksMirror.At(id);
452 static TLinearFitter fy2(3,"hyp2");
453 static TLinearFitter fz2(3,"hyp2");
454 static TLinearFitter fy1(2,"hyp1");
455 static TLinearFitter fz1(2,"hyp1");
456 static TVectorD vecy2,vecz2,vecy1,vecz1;
458 const Int_t kMinClusters=20;
461 for (Int_t i=0;i<72;++i) nclusters[i]=0;
463 for (Int_t i=0;i<160;++i) {
464 AliTPCclusterMI *c=track->GetClusterPointer(i);
465 if (c) nclusters[c->GetDetector()]++;
468 for (Int_t isec=0; isec<72;isec++){
469 if (nclusters[isec]<kMinClusters) continue;
475 for (Int_t irow=0;irow<160;++irow) {
476 AliTPCclusterMI *c=track->GetClusterPointer(irow);
477 //if (c && RejectCluster(c)) continue;
478 if (c&&c->GetDetector()==isec) {
479 Double_t xd = c->GetX()-120;;
480 Double_t x[2]={xd,xd*xd};
481 fy2.AddPoint(x,c->GetY());
482 fz2.AddPoint(x,c->GetZ());
484 fy1.AddPoint(x,c->GetY());
485 fz1.AddPoint(x,c->GetZ());
492 fy1.GetParameters(vecy1);
493 fy2.GetParameters(vecy2);
494 fz1.GetParameters(vecz1);
495 fz2.GetParameters(vecz2);
498 TTreeSRedirector *cstream = GetDebugStreamer();
500 Float_t dedx = track->GetdEdx();
501 (*cstream)<<"Tracklet"<<
505 "ncl="<<nclusters[isec]<<
519 // for (Int_t irow=0;irow<160;++irow) {
520 // AliTPCclusterMI *c=track->GetClusterPointer(irow);
521 // if (c && RejectCluster(c)) continue;
522 // if (c&&c->GetDetector()==isec) {
523 // Double_t x[2]={c->GetX(),c->GetX()*c->GetX()};
524 // fy2.AddPoint(&x,c->GetY());
525 // fz2.AddPoint(&x,c->GetZ());
527 // fy1.AddPoint(&x,c->GetY());
528 // fz1.AddPoint(&x,c->GetZ());
535 void AliTPCcalibLaser::DumpMeanInfo(Float_t bfield,Int_t minEntries){
537 // Dump information about laser beams
538 // isOK variable indicates usability of the beam
539 // Beam is not usable if:
540 // a. No entries in range (krmsCut0)
541 // b. Big sperad (krmscut1)
542 // c. RMSto fit sigma bigger then (kmultiCut)
543 // d. Too big angular spread
546 const Float_t krmsCut0=0.001;
547 const Float_t krmsCut1=0.16;
548 const Float_t kmultiCut=2;
549 const Float_t kcutP0=0.002;
551 AliTPCcalibLaser *laser = this;
552 TTreeSRedirector *pcstream = new TTreeSRedirector("laserMean.root");
556 for (Int_t id=0; id<336; id++){
558 TH1F * hisphi = (TH1F*)laser->fDeltaPhi.At(id);
559 TH1F * hisphiP = (TH1F*)laser->fDeltaPhiP.At(id);
560 TH1F * hisZ = (TH1F*)laser->fDeltaZ.At(id);
561 TH1F * hisS = (TH1F*)laser->fSignals.At(id);
562 if (!hisphi) continue;;
563 Double_t entries = hisphi->GetEntries();
564 if (entries<minEntries) continue;
566 AliTPCLaserTrack *ltrp = (AliTPCLaserTrack*)fTracksMirror.At(id);
568 AliTPCLaserTrack::LoadTracks();
569 ltrp =(AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(id);
571 Float_t meanphi = hisphi->GetMean();
572 Float_t rmsphi = hisphi->GetRMS();
574 Float_t meanphiP = hisphiP->GetMean();
575 Float_t rmsphiP = hisphiP->GetRMS();
576 Float_t meanZ = hisZ->GetMean();
577 Float_t rmsZ = hisZ->GetRMS();
578 hisphi->Fit(&fg,"","",hisphi->GetMean()-4*hisphi->GetRMS(),hisphi->GetMean()+4*hisphi->GetRMS());
579 Double_t gphi1 = fg.GetParameter(1);
580 Double_t gphi2 = fg.GetParameter(2);
581 hisphiP->Fit(&fg,"","",hisphiP->GetMean()-4*hisphiP->GetRMS(),hisphiP->GetMean()+4*hisphiP->GetRMS());
582 Double_t gphiP1 = fg.GetParameter(1);
583 Double_t gphiP2 = fg.GetParameter(2);
584 hisZ->Fit(&fg,"","",hisZ->GetMean()-4*hisZ->GetRMS(),hisZ->GetMean()+4*hisZ->GetRMS());
585 Double_t gz1 = fg.GetParameter(1);
586 Double_t gz2 = fg.GetParameter(2);
588 Float_t meanS=hisS->GetMean();
593 ltrp->GetPxPyPz(lpxyz);
595 if (rmsphi<krmsCut0) isOK=kFALSE; // empty in range - not entries inside
596 if (rmsphi>krmsCut1) isOK=kFALSE; // empty in range - not entries inside
597 if (rmsphi>krmsCut0) if (gphi2/rmsphi>kmultiCut) isOK=kFALSE; // multi peak structure
598 if (gphiP2>kcutP0) isOK=kFALSE;
600 (*pcstream)<<"Mean"<<
602 "entries="<<entries<< // number of entries
603 "bz="<<bfield<< // bfield
604 "LTr.="<<ltrp<< // refernece track
606 "lx0="<<lxyz[0]<< // reference x
607 "lx1="<<lxyz[1]<< // reference y
608 "lx2="<<lxyz[2]<< // refernece z
609 "lpx0="<<lpxyz[0]<< // reference x
610 "lpx1="<<lpxyz[1]<< // reference y
611 "lpx2="<<lpxyz[2]<< // refernece z
615 "mphi="<<meanphi<< //
616 "rmsphi="<<rmsphi<< //
620 "mphiP="<<meanphiP<< //
621 "rmsphiP="<<rmsphiP<< //
637 void AliTPCcalibLaser::DumpScanInfo(TTree * chain){
641 TTreeSRedirector *pcstream = new TTreeSRedirector("laserScan.root");
642 TFile * f = pcstream->GetFile();
651 char grnamefull[1000];
656 Double_t smphiP[100];
666 for (Int_t id=0; id<336; id++){
668 sprintf(cut,"isOK&&fId==%d",id);
669 Int_t entries = chain->Draw("bz",cut,"goff");
670 if (entries<3) continue;
671 AliTPCLaserTrack *ltrp = 0;;
672 if (!AliTPCLaserTrack::GetTracks()) AliTPCLaserTrack::LoadTracks();
673 ltrp =(AliTPCLaserTrack*)AliTPCLaserTrack::GetTracks()->UncheckedAt(id);
677 ltrp->GetPxPyPz(lpxyz);
679 chain->Draw("bz",cut,"goff");
680 memcpy(bz, chain->GetV1(), entries*sizeof(Double_t));
681 chain->Draw("0.01*abs(bz)+0.02",cut,"goff");
682 memcpy(sbz, chain->GetV1(), entries*sizeof(Double_t));
684 chain->Draw("gphi1",cut,"goff");
685 memcpy(mphi, chain->GetV1(), entries*sizeof(Double_t));
686 chain->Draw("0.05*abs(mphi)+gphi2",cut,"goff");
687 memcpy(smphi, chain->GetV1(), entries*sizeof(Double_t));
689 chain->Draw("gphiP1",cut,"goff");
690 memcpy(mphiP, chain->GetV1(), entries*sizeof(Double_t));
691 chain->Draw("0.05*abs(mphiP)+gphiP2",cut,"goff");
692 memcpy(smphiP, chain->GetV1(), entries*sizeof(Double_t));
694 chain->Draw("gz1",cut,"goff");
695 memcpy(mZ, chain->GetV1(), entries*sizeof(Double_t));
696 chain->Draw("0.01*abs(meanZ)+gz2",cut,"goff");
697 memcpy(smZ, chain->GetV1(), entries*sizeof(Double_t));
700 sprintf(grnamefull,"Side_%d_Bundle_%d_Rod_%d_Beam_%d",
701 ltrp->GetSide(), ltrp->GetBundle(), ltrp->GetRod(), ltrp->GetBeam());
705 TGraphErrors *grphi = new TGraphErrors(entries,bz,mphi,sbz,smphi);
708 pphi[0] = fp.GetParameter(0); // offset
709 pphi[1] = fp.GetParameter(1); // slope
710 pphi[2] = TMath::Sqrt(fp.GetChisquare()/(entries-2.)); // normalized chi2
711 sprintf(grname,"phi_id%d",id);
712 grphi->SetName(grname); grphi->SetTitle(grnamefull);
713 grphi->GetXaxis()->SetTitle("b_{z} (T)");
714 grphi->GetYaxis()->SetTitle("#Delta r#phi (cm)");
715 grphi->SetMaximum(1.2);
716 grphi->SetMinimum(-1.2);
720 gPad->SaveAs(Form("pic/phi/phi_%s.gif",grnamefull));
723 TGraphErrors *grphiP = new TGraphErrors(entries,bz,mphiP,sbz,smphiP);
726 pphiP[0] = fp.GetParameter(0); // offset
727 pphiP[1] = fp.GetParameter(1); // slope
728 pphiP[2] = TMath::Sqrt(fp.GetChisquare()/(entries-2.)); // normalized chi2
729 sprintf(grname,"phiP_id%d",id);
730 grphiP->SetName(grname); grphiP->SetTitle(grnamefull);
731 grphiP->GetXaxis()->SetTitle("b_{z} (T)");
732 grphiP->GetYaxis()->SetTitle("#Delta #phi (rad)");
733 grphiP->SetMaximum(pphiP[0]+0.005);
734 grphiP->SetMinimum(pphiP[0]-0.005);
736 gPad->SaveAs(Form("pic/phiP/phiP_%s.gif",grnamefull));
741 TGraphErrors *grmZ = new TGraphErrors(entries,bz,mZ,sbz,smZ);
744 pmZ[0] = fp.GetParameter(0); // offset
745 pmZ[1] = fp.GetParameter(1); // slope
746 pmZ[2] = TMath::Sqrt(fp.GetChisquare()/(entries-2.)); // normalized chi2
747 sprintf(grname,"mZ_id%d",id);
748 grmZ->SetName(grname); grmZ->SetTitle(grnamefull);
749 grmZ->GetXaxis()->SetTitle("b_{z} (T)");
750 grmZ->GetYaxis()->SetTitle("#Delta z (cm)");
752 gPad->SaveAs(Form("pic/z/z_%s.gif",grnamefull));
756 for (Int_t ientry=0; ientry<entries; ientry++){
757 (*pcstream)<<"Mean"<<
760 "entries="<<entries<<
762 "lx0="<<lxyz[0]<< // reference x
763 "lx1="<<lxyz[1]<< // reference y
764 "lx2="<<lxyz[2]<< // refernece z
765 "lpx0="<<lpxyz[0]<< // reference x
766 "lpx1="<<lpxyz[1]<< // reference y
767 "lpx2="<<lpxyz[2]<< // refernece z
769 "gphi1="<<mphi[ientry]<< // mean - from gaus fit
770 "pphi0="<<pphi[0]<< // offset
771 "pphi1="<<pphi[1]<< // mean
772 "pphi2="<<pphi[2]<< // norm chi2
774 "gphiP1="<<mphiP[ientry]<< // mean - from gaus fit
775 "pphiP0="<<pphiP[0]<< // offset
776 "pphiP1="<<pphiP[1]<< // mean
777 "pphiP2="<<pphiP[2]<< // norm chi2
780 "pmZ0="<<pmZ[0]<< // offset
781 "pmZ1="<<pmZ[1]<< // mean
782 "pmZ2="<<pmZ[2]<< // norm chi2
792 void AliTPCcalibLaser::Analyze(){
799 Long64_t AliTPCcalibLaser::Merge(TCollection *li) {
801 TIterator* iter = li->MakeIterator();
802 AliTPCcalibLaser* cal = 0;
804 while ((cal = (AliTPCcalibLaser*)iter->Next())) {
805 if (!cal->InheritsFrom(AliTPCcalibLaser::Class())) {
806 Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
810 // fHistNTracks->Add(cal->fHistNTracks);
811 // fClusters->Add(cal->fClusters);
812 // fModules->Add(cal->fModules);
813 // fHistPt->Add(cal->fHistPt);
814 // fPtResolution->Add(cal->fPtResolution);
815 // fDeDx->Add(cal->fDeDx);
821 for (Int_t id=0; id<336; id++){
822 // merge fDeltaZ histograms
823 hm = (TH1F*)cal->fDeltaZ.At(id);
824 h = (TH1F*)fDeltaZ.At(id);
826 h=new TH1F(Form("hisdz%d",id),Form("hisdz%d",id),1000,-10,10);
830 // merge fDeltaPhi histograms
831 hm = (TH1F*)cal->fDeltaPhi.At(id);
832 h = (TH1F*)fDeltaPhi.At(id);
834 h= new TH1F(Form("hisdphi%d",id),Form("hisdphi%d",id),1000,-1,1);
835 fDeltaPhi.AddAt(h,id);
838 // merge fDeltaPhiP histograms
839 hm = (TH1F*)cal->fDeltaPhiP.At(id);
840 h = (TH1F*)fDeltaPhiP.At(id);
842 h=new TH1F(Form("hisdphiP%d",id),Form("hisdphiP%d",id),1000,-0.01,0.01);
843 fDeltaPhiP.AddAt(h,id);
846 // merge fSignals histograms
847 hm = (TH1F*)cal->fSignals.At(id);
848 h = (TH1F*)fSignals.At(id);
850 h=new TH1F(Form("hisSignal%d",id),Form("hisSignal%d",id),1000,0,1000);
851 fSignals.AddAt(h,id);
862 gSystem->Load("libSTAT.so")
863 TStatToolkit toolkit;
868 TCut cutA("entries>2&&pphi2<3&&abs(gphiP1-pphiP0)<0.003&&abs(gz1)<6");
873 fstring+="(abs(LTr.fP[1]/250)^3-1)*bz++"; //1
874 fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*LTr.fP[2]++"; //2
875 fstring+="(abs(LTr.fP[1]/250)^1-1)*bz++"; //3
876 fstring+="(abs(LTr.fP[1]/250)-1)*bz*LTr.fP[2]++"; //4
878 fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*sin(atan2(lx1,lx0))++" //5
879 fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*sin(atan2(lx1,lx0))*LTr.fP[2]++" //6
880 fstring+="(abs(LTr.fP[1]/250)-1)*bz*sin(atan2(lx1,lx0))++" //7
881 fstring+="(abs(LTr.fP[1]/250)-1)*bz*sin(atan2(lx1,lx0))*LTr.fP[2]++" //8
883 fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*cos(atan2(lx1,lx0))++" //9
884 fstring+="(abs(LTr.fP[1]/250)^3-1)*bz*cos(atan2(lx1,lx0))*LTr.fP[2]++" //10
885 fstring+="(abs(LTr.fP[1]/250)-1)*bz*cos(atan2(lx1,lx0))++" //11
886 fstring+="(abs(LTr.fP[1]/250)-1)*bz*cos(atan2(lx1,lx0))*LTr.fP[2]++" //12
891 TString *strq0 = toolkit.FitPlane(treeT,"gphi1-pphi0",fstring->Data(), "fSide==1"+cutA, chi2,npoints,fitParam,covMatrix);
893 treeT->SetAlias("fit",strq0->Data());
896 TString *strqP = toolkit.FitPlane(treeT,"1000*(gphiP1-pphiP0)",fstring->Data(), "fSide==1"+cutA, chi2,npoints,fitParam,covMatrix);
898 treeT->SetAlias("fitP",strqP->Data());
901 TString *strqDrift = toolkit.FitPlane(treeT,"gz1","LTr.fP[1]++(1-2*(fSide==1))++lx1", cutA, chi2,npoints,fitParam,covMatrix);
902 treeT->SetAlias("fitD",strqDrift->Data());
905 treeT->Draw("fit:LTr.fP[1]","abs(bz+0.4)<0.05"+cutA,"");
907 for (Int_t i=0; i<6;i++){
908 treeT->SetLineColor(i+2);
909 treeT->SetMarkerSize(1);
910 treeT->SetMarkerStyle(22+i);
911 treeT->SetMarkerColor(i+2);
913 treeT->Draw("fit:LTr.fP[1]",Form("abs(bz+0.4)<0.05&fRod==%d",i)+cutA,"same");
921 TTree * tree = (TTree*)f.Get("FitModels");
923 TEventList listLFit0("listLFit0","listLFit0");
924 TEventList listLFit1("listLFit1","listLFit1");
926 tree->Draw(">>listLFit0","seed.fdEdx<200&&seed.fdEdx>40");
927 tree->SetEventList(&listLFit0);