/**************************************************************************
* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * SigmaEffect_thetadegrees *
+ * SigmaEffect_thetadegrees *
* Author: The ALICE Off-line Project. *
* Contributors are mentioned in the code where appropriate. *
* *
/* $Id$ */
-//===================================================================
-//This class was prepared by INFN Cagliari, July 2006
-//(authors: H.Woehri, A.de Falco)
-//
-//
+//-----------------------------------------------------------------------------
// Compact information for the muon generated tracks in the MUON arm
// useful at the last stage of the analysis chain
// provides a link between the reconstructed track and the generated particle
// mother process
//
// To be used together with AliMUONPairLight
-//===================================================================
+//
+// This class was prepared by INFN Cagliari, July 2006
+// (authors: H.Woehri, A.de Falco)
+//-----------------------------------------------------------------------------
+
+// 13 Nov 2007:
+// Added a temporary fix to FindRefTrack to be able to handle reconstructed tracks
+// generated from ESD muon track information. The problem is that the ESD data at
+// the moment only contains the first hit on chamber 1. Hopefully in the near future
+// this will be fixed and all hit information will be available.
+// - Artur Szostak <artursz@iafrica.com>
#include "AliMUONTrackLight.h"
#include "AliMUONTrack.h"
#include "AliMUONConstants.h"
+#include "AliMUONVTrackStore.h"
+#include "AliMUONTrackParam.h"
#include "AliESDMuonTrack.h"
-#include "AliRunLoader.h"
#include "AliStack.h"
-#include "AliHeader.h"
+#include "AliLog.h"
#include "TDatabasePDG.h"
-#include "TClonesArray.h"
#include "TParticle.h"
#include "TString.h"
+#include <cstdio>
+
ClassImp(AliMUONTrackLight)
//===================================================================
fPrec(),
fIsTriggered(kFALSE),
fCharge(-999),
+ fChi2(-1),
fCentr(-1),
fPgen(),
fTrackPythiaLine(-999),
fTrackPDGCode(-999),
fOscillation(kFALSE),
- fNParents(0)
+ fNParents(0),
+ fWeight(1)
{
/// default constructor
fPgen.SetPxPyPzE(0.,0.,0.,0.);
fPrec(muonCopy.fPrec),
fIsTriggered(muonCopy.fIsTriggered),
fCharge(muonCopy.fCharge),
+ fChi2(muonCopy.fChi2),
fCentr(muonCopy.fCentr),
fPgen(muonCopy.fPgen),
fTrackPythiaLine(muonCopy.fTrackPythiaLine),
fTrackPDGCode(muonCopy.fTrackPDGCode),
fOscillation(muonCopy.fOscillation),
- fNParents(muonCopy.fNParents)
+ fNParents(muonCopy.fNParents),
+ fWeight(muonCopy.fWeight)
{
/// copy constructor
for (Int_t i=0; i<3; i++) fXYZ[i]=muonCopy.fXYZ[i];
fPrec(),
fIsTriggered(kFALSE),
fCharge(-999),
+ fChi2(-1),
fCentr(-1),
fPgen(),
fTrackPythiaLine(-999),
fTrackPDGCode(-999),
fOscillation(kFALSE),
- fNParents(0)
+ fNParents(0),
+ fWeight(1)
{
/// constructor
- //AliMUONTrackLight();
- ComputePRecAndChargeFromESD(muonTrack);
+ fPgen.SetPxPyPzE(0.,0.,0.,0.);
+ for (Int_t npar = 0; npar < fgkNParentsMax; npar++){
+ fParentPDGCode[npar] = -1;
+ fParentPythiaLine[npar] = -1;
+ }
+ for (Int_t i = 0; i < 4; i++){
+ fQuarkPDGCode[i] = -1;
+ fQuarkPythiaLine[i] = -1;
+ }
+ FillFromESD(muonTrack);
}
//============================================
-void AliMUONTrackLight::ComputePRecAndChargeFromESD(AliESDMuonTrack* muonTrack){
+AliMUONTrackLight::~AliMUONTrackLight()
+{
+/// Destructor
+}
+
+//============================================
+AliMUONTrackLight& AliMUONTrackLight::operator=(const AliMUONTrackLight& muonCopy)
+{
+ // check assignment to self
+ if (this == &muonCopy) return *this;
+
+ // base class assignment
+ TObject::operator=(muonCopy);
+
+ // assignment operator
+ fPrec = muonCopy.fPrec;
+ fIsTriggered = muonCopy.fIsTriggered;
+ fCharge = muonCopy.fCharge;
+ fChi2 = muonCopy.fChi2;
+ fCentr = muonCopy.fCentr;
+ fPgen = muonCopy.fPgen;
+ fTrackPythiaLine = muonCopy.fTrackPythiaLine;
+ fTrackPDGCode = muonCopy.fTrackPDGCode;
+ fOscillation = muonCopy.fOscillation;
+ fNParents = muonCopy.fNParents;
+ fWeight = muonCopy.fWeight;
+
+ for (Int_t i=0; i<3; i++) fXYZ[i]=muonCopy.fXYZ[i];
+ for (Int_t npar = 0; npar < fgkNParentsMax; npar++){
+ fParentPDGCode[npar] = muonCopy.fParentPDGCode[npar];
+ fParentPythiaLine[npar] = muonCopy.fParentPythiaLine[npar];
+ }
+ for (Int_t i = 0; i < 4; i++){
+ fQuarkPDGCode[i] = muonCopy.fQuarkPDGCode[i];
+ fQuarkPythiaLine[i] = muonCopy.fQuarkPythiaLine[i];
+ }
+
+ return *this;
+}
+
+//============================================
+
+void AliMUONTrackLight::FillFromAliMUONTrack(AliMUONTrack *trackReco,Double_t zvert){
+ /// this method sets the muon reconstructed momentum according to the value given by AliMUONTrack
+ AliMUONTrackParam* trPar = trackReco->GetTrackParamAtVertex();
+ if (!trPar) {
+ AliError("The track must contain the parameters at vertex");
+ return;
+ }
+ this->SetCharge(Int_t(TMath::Sign(1.,trPar->GetInverseBendingMomentum())));
+ this->SetPxPyPz(trPar->Px(),trPar->Py(), trPar->Pz());
+ this->SetTriggered(trackReco->GetMatchTrigger());
+
+ Double_t xyz[3] = { trPar->GetNonBendingCoor(),
+ trPar->GetBendingCoor(),
+ trPar->GetZ()};
+ if (zvert!=-9999) xyz[2] = zvert;
+ this->SetVertex(xyz);
+}
+
+//============================================
+void AliMUONTrackLight::FillFromESD(AliESDMuonTrack* muonTrack,Double_t zvert){
/// computes prec and charge from ESD track
Double_t mumass = TDatabasePDG::Instance()->GetParticle(13)->Mass();
Double_t thetaX = muonTrack->GetThetaX();
fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
Double_t energy = TMath::Sqrt(mumass * mumass + px*px + py*py + pz*pz);
fPrec.SetPxPyPzE(px,py,pz,energy);
+ // get the position
+ fXYZ[0] = muonTrack->GetNonBendingCoor();
+ fXYZ[1] = muonTrack->GetBendingCoor();
+ if (zvert==-9999) fXYZ[2] = muonTrack->GetZ();
+ else fXYZ[2] = zvert;
+ // get the chi2 per d.o.f.
+ fChi2 = muonTrack->GetChi2()/ (2.0 * muonTrack->GetNHit() - 5);
+ fIsTriggered = muonTrack->GetMatchTrigger();
}
//============================================
}
//============================================
-TParticle* AliMUONTrackLight::FindRefTrack(AliMUONTrack* trackReco, TClonesArray* trackRefArray, AliRunLoader *runLoader){
- /// find the MC particle that corresponds to a given rec track
- TParticle *part = 0;
- const Double_t kSigma2Cut = 16; // 4 sigmas cut, kSigma2Cut = 4*4
- Int_t nTrackRef = trackRefArray->GetEntriesFast();
- Int_t compPart = 0;
- for (Int_t iref = 0; iref < nTrackRef; iref++) {
- AliMUONTrack *trackRef = (AliMUONTrack *)trackRefArray->At(iref);
- // check if trackRef is compatible with trackReco:
- //routine returns for each chamber a yes/no information if the
- //hit of rec. track and hit of referenced track are compatible
- Bool_t *compTrack = trackRef->CompatibleTrack(trackReco,kSigma2Cut);
- Int_t iTrack = this->TrackCheck(compTrack); //returns number of validated conditions
- if (iTrack==4) {
- compPart++;
- Int_t trackID = trackRef->GetTrackID();
- this->SetTrackPythiaLine(trackID);
- part = ((AliStack *)(((AliHeader *) runLoader->GetHeader())->Stack()))->Particle(trackID);
- fTrackPDGCode = part->GetPdgCode();
- }
- }
- if (compPart>1) {
- printf ("<AliMUONTrackLight::FindRefTrack> ERROR: more than one particle compatible to the reconstructed track.\n");
- Int_t i=0, j=1/i;
- printf ("j=%d \n",j);
- }
- return part;
-}
-
-//============================================
-Int_t AliMUONTrackLight::TrackCheck(Bool_t *compTrack){
- /// Apply reconstruction requirements
- /// Return number of validated conditions
- /// If all the tests are verified then TrackCheck = 4 (good track)
- Int_t iTrack = 0;
- Int_t hitsInLastStations = 0;
-
- // apply reconstruction requirements
- if (compTrack[0] || compTrack[1]) iTrack++; // at least one hit in st. 0
- if (compTrack[2] || compTrack[3]) iTrack++; // at least one hit in st. 1
- if (compTrack[4] || compTrack[5]) iTrack++; // at least one hit in st. 2
- for (Int_t ch = 6; ch < AliMUONConstants::NTrackingCh(); ch++) {
- if (compTrack[ch]) hitsInLastStations++;
- }
- if (hitsInLastStations > 2) iTrack++; // at least 3 hits in st. 3 & 4
- return iTrack;
-}
-
-//============================================
-void AliMUONTrackLight::FillMuonHistory(AliRunLoader *runLoader, TParticle *part){
+void AliMUONTrackLight::FillMuonHistory(AliStack *stack, TParticle *part){
/// scans the muon history to determine parents pdg code and pythia line
Int_t countP = -1;
Int_t parents[10], parLine[10];
Int_t lineM = part->GetFirstMother();//line in the Pythia output of the particle's mother
- AliStack *stack = runLoader->GetHeader()->Stack();
TParticle *mother;
Int_t status=-1, pdg=-1;
while(lineM >= 0){
-
+
mother = stack->Particle(lineM); //direct mother of rec. track
pdg = mother->GetPdgCode();//store PDG code of first mother
- if(pdg == 92) break; // break if a string is found
+ // break if a string, gluon, quark or diquark is found
+ if(pdg == 92 || pdg == 21 || TMath::Abs(pdg) < 10 || IsDiquark(pdg)) break;
parents[++countP] = pdg;
parLine[countP] = lineM;
status = mother->GetStatusCode();//get its status code to check if oscillation occured
this->SetParentPDGCode(i,parents[countP-i]);
this->SetParentPythiaLine(i,parLine[countP-i]);
}
-
fNParents = countP+1;
countP = -1;
+
//and store the lines of the string and further quarks in another array:
while(lineM >= 0){
mother = stack->Particle(lineM);
this->SetQuarkPDGCode(countP, pdg);//store the pdg of the quarks in index 1,2
lineM = mother->GetFirstMother();
}
-
+
//check if in case of HF production, the string points to the correct end
//and correct it in case of need:
countP = 1;
+ for(int par = 0; par < 4; par++){
+ if(TMath::Abs(this->GetQuarkPDGCode(par)) < 6){
+ countP = par; //get the quark just before hadronisation
+ break;
+ }
+ }
if(this->GetQuarkPythiaLine(countP) > -1 && (this->GetParentFlavour(0)==4 || this->GetParentFlavour(0)==5)){
if(this->GetParentFlavour(0) != TMath::Abs(this->GetQuarkPDGCode(countP))){
- Int_t pdg = this->GetQuarkPDGCode(1), line = this->GetQuarkPythiaLine(1);
+ AliWarning(Form("quark flavour of parent and that of quark do not correspond: %d %d --> correcting\n",
+ this->GetParentFlavour(0), TMath::Abs(this->GetQuarkPDGCode(countP)))
+ );
+
+ pdg = this->GetQuarkPDGCode(countP);
+ Int_t line = this->GetQuarkPythiaLine(countP);
this->ResetQuarkInfo();
while(TMath::Abs(pdg) != this->GetParentFlavour(0)){//pdg of q,g in Pythia listing following the wrong string end
//must coincide with the flavour of the last fragmented mother
this->SetQuarkPDGCode(countP++, pdg);
line = mother->GetFirstMother();
}
- }
+ this->PrintInfo("h");
+ }//mismatch
}
}
+
//====================================
void AliMUONTrackLight::ResetQuarkInfo(){
/// resets parton information
}
//====================================
-void AliMUONTrackLight::PrintInfo(Option_t* opt){
+void AliMUONTrackLight::PrintInfo(const Option_t* opt){
/// prints information about the track:
/// - "H" muon's decay history
/// - "K" muon kinematics
TString pdg = "", line = "";
for(int i = 3; i >= 0; i--){
if(this->GetQuarkPythiaLine(i)>= 0){
- sprintf(name, "%4d --> ", this->GetQuarkPythiaLine(i));
+ snprintf(name, 100, "%4d --> ", this->GetQuarkPythiaLine(i));
line += name;
- sprintf(name, "%4d --> ", this->GetQuarkPDGCode(i));
+ snprintf(name, 100, "%4d --> ", this->GetQuarkPDGCode(i));
pdg += name;
}
}
for(int i = 0; i < fNParents; i++){
if(this->GetParentPythiaLine(i)>= 0){
- sprintf(name, "%7d --> ", this->GetParentPythiaLine(i));
+ snprintf(name, 100, "%7d --> ", this->GetParentPythiaLine(i));
line += name;
- sprintf(name, "%7d --> ", this->GetParentPDGCode(i));
+ snprintf(name, 100, "%7d --> ", this->GetParentPDGCode(i));
pdg += name;
}
}
- sprintf(name, "%4d", this->GetTrackPythiaLine()); line += name;
- sprintf(name, "%4d", this->GetTrackPDGCode()); pdg += name;
+ snprintf(name, 100, "%4d", this->GetTrackPythiaLine()); line += name;
+ snprintf(name, 100, "%4d", this->GetTrackPDGCode()); pdg += name;
printf("\nmuon's decay history:\n");
printf(" PDG: %s\n", pdg.Data());
momRec.Phi(), 180./TMath::Pi() * momRec.Phi());
}
}
-
-
+//====================================
+Bool_t AliMUONTrackLight::IsParentPionOrKaon(Int_t idparent){
+ /// checks if a muon comes from a pion or kaon or a particle that decays into one of these two
+ Int_t pdg = this->GetParentPDGCode(idparent);
+ if (TMath::Abs(pdg)==211 || //pi+
+ TMath::Abs(pdg)==321 || //K+
+ TMath::Abs(pdg)==213 || //rho+
+ TMath::Abs(pdg)==311 || //K0
+ TMath::Abs(pdg)==313 || //K*0
+ TMath::Abs(pdg)==323 //K*+
+ ) {
+ return kTRUE;
+ }
+ else return kFALSE;
+}
+//====================================
+Bool_t AliMUONTrackLight::IsDiquark(Int_t pdg) const{
+ /// check if the provided pdg code corresponds to a diquark
+ pdg = TMath::Abs(pdg);
+ if((pdg > 1000) && (pdg%100 < 10)) return kTRUE;
+ else return kFALSE;
+}