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 **************************************************************************/
16 //-------------------------------------------------------
17 // Implementation of the TPC combined trackfit
19 // Origin: Xianguo Lu lu@physi.uni-heidelberg.de Xianguo.Lu@cern.ch
21 //=========================================================================================
24 // In the default reconstruction in the ALICE the cosmic tracks are found as two independent particles.
26 // In general any of subtracks can be used for the physics studies. In order to avoid the double counting,
27 // the track from the upper hemisphere can be used.
29 // The momentum resolution is determined by the lever arm (1/L^2) and by the number of clusters
30 // used for the track fitting (1/sqrt(Ncl)).
31 // Combining/refitting the two segments together significantly better momentum resolution can be obtained.
32 // sigma_{1/pt} ~ 8x10^-3 - defaul tracking (e.g only upper track)
33 // sigma_{1/pt} ~ 8x10^-4 - combined tracking
34 //===========================================================================================
36 // Interface/Implementation:
37 // This class provide the functionality to combined two ESD tracks and perform the trackfit using all related track points.
38 // Input: AliESDtrack's (or AliTPCseed's, expert use)
39 // Output: the trackparams from the combined fit are stored as data member of this class, and they can be accessed via getter or copy.
41 //===========================================================================================
44 // 1. Read in two ESDtracks (see function CombineESDtracks())
45 // TPCseeds are got from ESDfriends
46 // Cuts are applied on the following quantities: #TPC cluster, relation between global-y coordinate of the two tracks, lever arm.
47 // 2. AliTPCCosmicUtils::CombinedFit is called and refit is performed using the TPCseeds from both tracks
48 // 3. The fitting quality is analyzed and stored and can be accessed via function GetStatus(). The detailed judgement of the quality is recorded as "enum CombineStatus"
51 //===========================================================================================
52 // Algorithm numerical debugging:
53 // The class can be used in the different debug/verbose level (see fDebugLevel)
54 // Several intermediate variables can be stored in the trees, printout, or draw.
55 // Given functionality (dumping of variables to the tree) was also used for the tuning of the pair
56 // selection criterias, and for validation of the fit functionality.
58 // ----- Debug output:
59 // for (debuglevel & 1)==1 && good fit, the following info saved:
61 (*fStreamer)<<"TrackProp"<<
62 "Tup.="<<fTrackparUp<< //AliExternalTrackParam at uppermost cluster obtained by upward propagation
63 "Tlow.="<<fTrackparLow<< //AliExternalTrackParam at lowermost cluster obtained by downward propagation
64 "icup.="<<&fInnerClusterUp<< //TVector3 position of the innermost cluster of the upper track
65 "iclow.="<<&fInnerClusterLow<<
66 "leverarm="<<fLeverArm<<
67 "ncl="<<fFitNcls<< //number of clusters used in successful propagation
68 "nmiss="<<fMissNcls<< //number of clusters failed in propagation, should always be 0 in this case.
69 "chi2="<<fPreChi2<< //chi2/nfit
70 "momup="<< momup << //momentum at uppermost cluster with upward propagation
71 "momlow="<< momlow << //momentum at lowermost cluster with downward propagation
76 // for (debuglevel & 2)==1, debug info in AliTPCCosmicUtils::FitKernel saved
78 //===========================================================================================
81 fCombinedTrackfit = new AliTPCCosmicTrackfit(debuglevel, "anystring");
83 //order not important; will be internally ordered (potinters modified due to &) such that track0 is the upper one
84 //kfit = kTRUE: good fit, kFALSE: bad fit
85 const Bool_t kfit = fCombinedTrackfit->CombineESDtracks(esdtrack0, esdtrack1);
87 //status = 0 for good fit (i.e. kfit=kTRUE), non-0 for bad fit (i.e. kfit=kFALSE), see "enum CombineStatus" definition in header file
88 const Int_t status = fCombinedTrackfit->GetStatus();
90 //===========================================================================================
93 // for 2011 Feb. cosmic data nch=2 events, the kfit and status look like:
95 kfit,status ( 0, 1): 16337 / 443901 = 3.680% //kFailGetTPCseeds
96 kfit,status ( 0, 2): 3692 / 443901 = 0.832% //not both tracks have ncl > AliTPCCosmicUtils::fgkNclsMin
97 kfit,status ( 0, 3): 6527 / 443901 = 1.470% //clusters in two tracks should be clearly separated in y, i.e. lowest cluster of upper track higher than highest cluster of lower track; otherwise fail
98 kfit,status ( 0, 4): 7033 / 443901 = 1.584% //fLeverArm<fgkCutLeverArm
99 kfit,status ( 0, 6): 4398 / 443901 = 0.991% //fail in propagation of at least one cluster
100 kfit,status ( 0, 7): 508 / 443901 = 0.114% //chi2/nfit > fgkMaxChi2
101 kfit,status ( 0, 8): 52 / 443901 = 0.012% //fail in geting impact parameters
102 kfit,status ( 1, 0): 405354 / 443901 = 91.316% //i.e. 91% of nch=2 events are successfully fitted.
105 //===========================================================================================
108 // for muon momentum small than 20 GeV, energy loss in muon filter is visable when compaing fTrackparUp and fTrackparLow; energy loss estimated as 5 MeV/cm.
109 // particle traversing muon filter can be rejected by requiring "fInnerClusterUp.fZ > -40 && fInnerClusterLow.fZ > -40"
110 // momentum resolution is estimated by comparing the trackfit result by upward propagation through odd pad rows and that by downward propagation through even pad rows. Number of clusters used in this case is only half of that in normal usage.
111 // RMS of log10 p = 0.01 at 10 GeV/c, 0.1 at 100 GeV/c, 0.5 at 1 TeV/c.
112 // muon filter deteriorates momentum resolution by about +0.01 (absolute value).
119 #include <TTreeStream.h>
120 #include <TVector3.h>
122 #include "AliESDtrack.h"
123 #include "AliESDfriendTrack.h"
124 #include "AliTPCseed.h"
125 #include "AliTrackerBase.h"
126 #include "AliTrackPointArray.h"
128 #include "AliTPCCosmicUtils.h"
129 #include "AliTPCCosmicTrackfit.h"
131 AliTPCCosmicTrackfit::AliTPCCosmicTrackfit(const Int_t dlev, const TString tag):
132 fStreamer(0x0), fDebugLevel(dlev)
133 , fSeedUp(0x0), fSeedLow(0x0), fTrackparUp(0x0), fTrackparLow(0x0), fIniTrackparUp(0x0), fIniTrackparLow(0x0)
136 , fFitNcls(-999), fMissNcls(-999), fPreChi2(-999), fFitLeverArm(-999), fImpactD(-999), fImpactZ(-999)
141 fInnerClusterUp.SetXYZ(-999,-999,-999);
142 fInnerClusterLow.SetXYZ(-999,-999,-999);
145 fStreamer = new TTreeSRedirector(Form("CombinedTrackfit_%s.root", tag.Data()));
151 AliTPCCosmicTrackfit::~AliTPCCosmicTrackfit()
161 delete fIniTrackparUp;
162 delete fIniTrackparLow;
165 Bool_t AliTPCCosmicTrackfit::CombineESDtracks(AliESDtrack * &trk0, AliESDtrack *&trk1)
168 //Get TPCseeds from the 2 ESDtracks, swap TPCseeds and ESDTracks (if necessary) according to y (0:upper 1:lower), perform trackfit using TPCseeds
169 //if fStatus==0, i.e. combine is successful, swap of the ESDtracks is kept since pointer *& is used
172 IniCombineESDtracks();
174 if(!GetTPCseeds(trk0, trk1)){
182 AliESDtrack * tmptrk = trk0;
192 Bool_t AliTPCCosmicTrackfit::CombineTPCseedsFast(AliTPCseed * tpcseeds[], const AliExternalTrackParam * trkpars[])
195 //do combined track fit for given TPC seeds and initial trackpar, [0]: upper, [1]: lower
198 IniCombineESDtracks();
200 fSeedUp = tpcseeds[0];
201 fSeedLow = tpcseeds[1];
203 fTrackparUp = new AliExternalTrackParam(*(trkpars[0]));
204 fTrackparLow = new AliExternalTrackParam(*(trkpars[1]));
206 AliExternalTrackParam * trackPars[]={fTrackparUp, fTrackparLow};
207 const AliTPCseed *seeds[]={fSeedUp, fSeedLow};
208 TTreeSRedirector * debugstreamer = 0x0;
210 debugstreamer = fStreamer;
213 AliTPCCosmicUtils::CombinedFit(trackPars, seeds, fRowStartShift, fRowStep, fXMin, fXMax, fFitNcls, fMissNcls, fPreChi2, fFitLeverArm, fImpactD, fImpactZ, debugstreamer);
224 Bool_t AliTPCCosmicTrackfit::CombineTPCseeds(AliTPCseed * &seed0, AliTPCseed *&seed1)
227 //same as AliTPCCosmicTrackfit::CombineESDtracks, except that the seeds are passed in from outside, which can be still unordered
228 //if fStatus==0, i.e. combine is successful, swap of the TPCseeds is kept since pointer *& is used
230 IniCombineESDtracks();
239 AliTPCseed * tmpseed = seed0;
249 void AliTPCCosmicTrackfit::Print() const
252 //print out variable values
254 printf("Status %2d NclsU %3d NclsD %3d ZinnerU %7.2f ZinnerD %7.2f LeverArm %7.2f\n", fStatus, fSeedUp->GetNumberOfClusters(), fSeedLow->GetNumberOfClusters(), fInnerClusterUp.Z(), fInnerClusterLow.Z(), fLeverArm);
258 Double_t AliTPCCosmicTrackfit::ImpactParameter2D() const
261 //calculate the 2D-impactparameter from (0,0)
264 const TVector3 p0(0,0,0);
265 const TVector3 l1(fInnerClusterUp.X(), fInnerClusterUp.Y(), 0);
266 const TVector3 l2(fInnerClusterLow.X(), fInnerClusterLow.Y(), 0);
268 return AliTPCCosmicUtils::Point2LineDist(p0, l1, l2);
271 Double_t AliTPCCosmicTrackfit::ImpactParameter3D() const
274 //calculate the 3D-impactparameter from (0,0,0)
277 const TVector3 p0(0,0,0);
278 const TVector3 l1(fInnerClusterUp);
279 const TVector3 l2(fInnerClusterLow);
281 return AliTPCCosmicUtils::Point2LineDist(p0, l1, l2);
285 Double_t AliTPCCosmicTrackfit::MinPhi() const
288 //the smaller phi of the two tracks w.r.t. horizon
290 Double_t fsp[] = {TMath::Abs(TMath::Sin(fTrackparUp->Phi())), TMath::Abs(TMath::Sin(fTrackparLow->Phi()))};;
291 return TMath::ASin(TMath::Min(fsp[0], fsp[1])) * TMath::RadToDeg();
293 //===================================================================================================
294 //===================================================================================================
296 void AliTPCCosmicTrackfit::IniCombineESDtracks()
299 //initialization, for reuse of the same AliTPCCosmicTrackfit instance
309 delete fIniTrackparUp;
310 delete fIniTrackparLow;
311 fIniTrackparUp = 0x0;
312 fIniTrackparLow = 0x0;
318 void AliTPCCosmicTrackfit::CombineTPCseeds()
321 //do combined trackfit using TPCseeds
326 || !AnaSeeds() // fSeedUp/Low swapped here! check by runTest.C
331 //AliExternalTrackParam object created
332 fTrackparUp = AliTPCCosmicUtils::MakeSeed(fSeedUp);
333 fTrackparLow = AliTPCCosmicUtils::MakeSeed(fSeedLow);
334 if(!fTrackparUp || !fTrackparLow){
335 fStatus = kFailMakeSeed;
339 fIniTrackparUp = new AliExternalTrackParam(*fTrackparUp);
340 fIniTrackparLow = new AliExternalTrackParam(*fTrackparLow);
342 AliExternalTrackParam * trackPars[]={fTrackparUp, fTrackparLow};
343 const AliTPCseed *seeds[]={fSeedUp, fSeedLow};
344 TTreeSRedirector * debugstreamer = 0x0;
346 debugstreamer = fStreamer;
349 AliTPCCosmicUtils::CombinedFit(trackPars, seeds, fRowStartShift, fRowStep, fXMin, fXMax, fFitNcls, fMissNcls, fPreChi2, fFitLeverArm, fImpactD, fImpactZ, debugstreamer);
356 void AliTPCCosmicTrackfit::Update()
359 //Update variables depending on the fit result
362 if(fMissNcls || fFitNcls==0){
363 fStatus = kFailPropagation;
367 fPreChi2 /= fFitNcls;
368 if(fPreChi2>fgkMaxChi2){
374 fStatus = kFailImpact;
378 if( fStatus == 0 && (fDebugLevel&1) ){
379 Double_t momup = fTrackparUp->P();
380 Double_t momlow = fTrackparLow->P();
381 Double_t ptup = fTrackparUp->Pt();
382 Double_t ptlow = fTrackparLow->Pt();
384 (*fStreamer)<<"TrackProp"<<
385 "Tup.="<<fTrackparUp<<
386 "Tlow.="<<fTrackparLow<<
387 "icup.="<<&fInnerClusterUp<<
388 "iclow.="<<&fInnerClusterLow<<
389 "leverarm="<<fLeverArm<<
391 "nmiss="<<fMissNcls<<
394 "momlow="<< momlow <<
401 Bool_t AliTPCCosmicTrackfit::CheckLeverArm()
404 //if lever arm is too short, no need to use combined track fit.
405 //On the other hand, short lever arm from two tracks mostly means they are fake pairs.
406 //lever arm extents over one quadrant, e.g. (0,250)-(250,0): 250*sqrt(2)~350
408 if(fLeverArm<fgkCutLeverArm){
409 fStatus = kFailLeverArm;
416 Bool_t AliTPCCosmicTrackfit::AnaSeeds()
419 //swap seeds (if necessary) so that (y of fSeedUp) > (y of fSeedLow)
422 //---------------------------------- navigate through all clusters ----------------------------------
423 AliTPCseed ** seeds[]={&fSeedUp, &fSeedLow};
425 //min, max according to y
426 TVector3 singlemin[2], singlemax[2];
427 for(Int_t ii=0; ii<2; ii++){
428 singlemin[ii].SetXYZ( 1e10, 1e10, 1e10);
429 singlemax[ii].SetXYZ(-1e10, -1e10, -1e10);
432 for(Int_t itrk=0; itrk<2; itrk++){
433 for(Int_t irow=0; irow<AliTPCCosmicUtils::fgkNRow; irow++){
434 const AliTPCclusterMI * cls = (*seeds[itrk])->GetClusterPointer(irow);
438 Float_t xyz[3]={-999,-999,-999};
439 cls->GetGlobalXYZ(xyz);
440 if(xyz[1]<singlemin[itrk].Y()){
441 singlemin[itrk].SetXYZ(xyz[0], xyz[1], xyz[2]);
443 if(xyz[1]>singlemax[itrk].Y()){
444 singlemax[itrk].SetXYZ(xyz[0], xyz[1], xyz[2]);
449 //--------------------------------
451 //kpass true if y of the two seeds clearly separate: min of one > max of the other
452 Bool_t kpass = kFALSE;
454 fInnerClusterUp.SetXYZ(-999,-999,-999);
455 fInnerClusterLow.SetXYZ(-999,-999,-999);
456 TVector3 combinedmin, combinedmax;
457 if(singlemin[0].Y()> singlemax[1].Y()){
458 fInnerClusterUp = singlemin[0];
459 fInnerClusterLow = singlemax[1];
466 combinedmax = singlemax[0];
467 combinedmin = singlemin[1];
469 else if(singlemin[1].Y()> singlemax[0].Y()){
470 fInnerClusterUp = singlemin[1];
471 fInnerClusterLow = singlemax[0];
475 AliTPCseed *tmp=*(seeds[0]);
476 *(seeds[0])=*(seeds[1]);
481 combinedmax = singlemax[1];
482 combinedmin = singlemin[0];
487 const TVector3 comdelta = combinedmax-combinedmin;
488 fLeverArm = comdelta.Pt();
491 fStatus = kFailSwapSeeds;
498 Bool_t AliTPCCosmicTrackfit::CheckNcls()
501 //check number of clusters in TPCseed, for too small number MakeSeed will fail
503 if( fSeedUp->GetNumberOfClusters()<AliTPCCosmicUtils::fgkNclsMin || fSeedLow->GetNumberOfClusters()<AliTPCCosmicUtils::fgkNclsMin ){
504 fStatus = kFailNclsMin;
511 Bool_t AliTPCCosmicTrackfit::GetTPCseeds(const AliESDtrack *trk0, const AliESDtrack *trk1)
514 //Get TPC seeds from ESDfriendTrack
516 fSeedUp = AliTPCCosmicUtils::GetTPCseed(trk0);
517 fSeedLow = AliTPCCosmicUtils::GetTPCseed(trk1);
519 if(!fSeedUp || !fSeedLow){
520 fStatus = kFailGetTPCseeds;