]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - MUON/AliMUONTrackLight.cxx
Update to fix the fact that AliRawReader::Create does not work correctly in online...
[u/mrichter/AliRoot.git] / MUON / AliMUONTrackLight.cxx
index 29102fee522bbeca6e4139e369e48b1235dfc0b7..277cc2f178fbdc4c144feb3fda49646c2b651cca 100644 (file)
@@ -1,6 +1,6 @@
 /**************************************************************************
  * 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 "AliMUONTrackExtrap.h"
+#include "AliMUONTrackParam.h"
 
 #include "AliESDMuonTrack.h"
 #include "AliRunLoader.h"
@@ -40,7 +49,6 @@
 #include "AliHeader.h"
 
 #include "TDatabasePDG.h"
-#include "TClonesArray.h"
 #include "TParticle.h"
 #include "TString.h"
 
@@ -53,12 +61,14 @@ AliMUONTrackLight::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.); 
@@ -80,12 +90,14 @@ AliMUONTrackLight::AliMUONTrackLight(const AliMUONTrackLight &muonCopy)
     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]; 
@@ -105,20 +117,49 @@ AliMUONTrackLight::AliMUONTrackLight(AliESDMuonTrack* muonTrack)
     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); 
+  FillFromESD(muonTrack); 
+}
+
+//============================================
+AliMUONTrackLight::~AliMUONTrackLight()
+{
+/// Destructor
+} 
+
+//============================================
+
+void AliMUONTrackLight::FillFromAliMUONTrack(AliMUONTrack *trackReco,Double_t zvert){
+  /// this method sets the muon reconstructed momentum according to the value given by AliMUONTrack
+  AliMUONTrackParam trPar;
+  if (trackReco->GetTrackParamAtVertex()) trPar = *(trackReco->GetTrackParamAtVertex());
+  else {
+    trPar = *((AliMUONTrackParam*) trackReco->GetTrackParamAtCluster()->First());
+    AliMUONTrackExtrap::ExtrapToVertex(&trPar,0.,0.,0.,0.,0.);
+  }
+  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::ComputePRecAndChargeFromESD(AliESDMuonTrack* muonTrack){ 
+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();
@@ -132,6 +173,14 @@ void AliMUONTrackLight::ComputePRecAndChargeFromESD(AliESDMuonTrack* muonTrack){
   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();
 }
 
 //============================================
@@ -143,33 +192,58 @@ void AliMUONTrackLight::SetPxPyPz(Double_t px, Double_t py, Double_t pz){
 }
 
 //============================================
-TParticle* AliMUONTrackLight::FindRefTrack(AliMUONTrack* trackReco, TClonesArray* trackRefArray, AliRunLoader *runLoader){ 
-  /// find the MC particle that corresponds to a given rec track
+TParticle* AliMUONTrackLight::FindRefTrack(
+               AliMUONTrack* trackReco, AliMUONVTrackStore* trackRefArray,
+               AliStack* stack)
+{
+  /// Find the Monte Carlo (MC) particle that corresponds to a given reconstructed track.
+  /// @param trackReco  This is the reconstructed track for which we want to find a
+  ///           corresponding MC particle.
+  /// @param trackRefArray  The list of MC reference tracks as generated by
+  ///           AliMUONRecoCheck::ReconstructedTracks for example.
+  /// @param stack  The MC stack of simulated particles.
+  
   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);
+  const Double_t kSigmaCut = 4;  // 4 sigmas cut
+  Int_t compPart = 0;
+  Bool_t compTrack;
+  
+  TIter next(trackRefArray->CreateIterator());
+  AliMUONTrack* trackRef;
+  while ( (trackRef = static_cast<AliMUONTrack*>(next())) ) {
     // 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();
+    compTrack = kFALSE;
+    
+    if (trackReco->GetNClusters() > 1) {
+      
+      // check cluster by cluster if trackReco contain info at each cluster
+      Bool_t compTrackArray[10];
+      trackRef->CompatibleTrack(trackReco,kSigmaCut,compTrackArray);
+      if (TrackCheck(compTrackArray) == 4) compTrack = kTRUE;
+      
+    } else {
+      
+      // otherwise check only parameters at the z position of the first trackRef
+      AliMUONTrackParam *refParam = (AliMUONTrackParam*) trackRef->GetTrackParamAtCluster()->First();
+      AliMUONTrackParam recoParam(*((AliMUONTrackParam*) trackReco->GetTrackParamAtCluster()->First()));
+      AliMUONTrackExtrap::ExtrapToZCov(&recoParam, refParam->GetZ());
+      Double_t chi2;
+      if (refParam->CompatibleTrackParam(recoParam, kSigmaCut, chi2)) compTrack = kTRUE;
+      
+    }
+      
+    if (compTrack) {
+      compPart++;
+      Int_t trackID = trackRef->GetUniqueID();
       this->SetTrackPythiaLine(trackID);
-      part = ((AliStack *)(((AliHeader *) runLoader->GetHeader())->Stack()))->Particle(trackID);
+      part = 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); 
+  if (compPart>1) {
+    AliFatal("More than one particle compatible with the reconstructed track.");
   } 
-  return part; 
+  return part;
 }
 
 //============================================
@@ -177,6 +251,7 @@ 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;
   
@@ -192,20 +267,20 @@ Int_t AliMUONTrackLight::TrackCheck(Bool_t *compTrack){
 }
 
 //============================================
-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
@@ -217,9 +292,9 @@ void AliMUONTrackLight::FillMuonHistory(AliRunLoader *runLoader, TParticle *part
     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);
@@ -229,14 +304,25 @@ void AliMUONTrackLight::FillMuonHistory(AliRunLoader *runLoader, TParticle *part
     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
@@ -252,9 +338,11 @@ void AliMUONTrackLight::FillMuonHistory(AliRunLoader *runLoader, TParticle *part
        this->SetQuarkPDGCode(countP++, pdg);
        line = mother->GetFirstMother();
       }
-    }
+      this->PrintInfo("h");
+    }//mismatch
   }
 }
+
 //====================================
 void AliMUONTrackLight::ResetQuarkInfo(){
   /// resets parton information
@@ -346,5 +434,25 @@ void AliMUONTrackLight::PrintInfo(Option_t* opt){
           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){
+  /// 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;
+}