-/**************************************************************************
+/*************************************************************************
* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
* *
* Author: The ALICE Off-line Project. *
#include "AliMFTTrackerMU.h"
#include "TMath.h"
#include "AliRun.h"
+#include "AliRunLoader.h"
+#include "AliHeader.h"
+#include "AliGenEventHeader.h"
+#include "TArrayF.h"
#include "AliMFT.h"
#include "AliMUONTrackExtrap.h"
#include "AliMUONTrack.h"
#include "AliMUONESDInterface.h"
#include "AliMuonForwardTrack.h"
#include "AliMUONConstants.h"
+#include "TGrid.h"
ClassImp(AliMFTTrackerMU)
fMUONTrack(0),
fCurrentTrack(0),
fFinalBestCandidate(0),
- fVertexErrorX(0.015),
- fVertexErrorY(0.015),
- fVertexErrorZ(0.010),
+ fXExtrapVertex(0),
+ fYExtrapVertex(0),
+ fZExtrapVertex(0),
+ fXExtrapVertexError(0),
+ fYExtrapVertexError(0),
fBransonCorrection(kFALSE)
{
for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
AliInfo(Form("plane %02d: nClusters = %d\n", iPlane, fMFTClusterArray[iPlane]->GetEntries()));
}
+
+ AddClustersFromUnderlyingEvent();
+ AddClustersFromPileUpEvents();
+
SeparateFrontBackClusters();
return 0;
void AliMFTTrackerMU::UnloadClusters() {
-// //--------------------------------------------------------------------
-// // This function unloads MFT clusters
-// //--------------------------------------------------------------------
+ //--------------------------------------------------------------------
+ // This function unloads MFT clusters
+ //--------------------------------------------------------------------
for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
fMFTClusterArray[iPlane] -> Clear("C");
fESD = event;
- // get the ITS primary vertex
+ fXExtrapVertex = 0;
+ fYExtrapVertex = 0;
+ fZExtrapVertex = 0;
+ fXExtrapVertexError = AliMFTConstants::fXVertexTolerance;
+ fYExtrapVertexError = AliMFTConstants::fYVertexTolerance;
- Double_t vertex[3] = {0., 0., 0.};
+ // Imposing a fixed primary vertex, in order to keep memory of its position. Taking the primary vertex would imply the risk
+ // of loosing memory of its position when passing from ESD to AOD, due to possible refitting
const AliESDVertex* esdVert = fESD->GetVertex();
if (esdVert->GetNContributors() > 0 || !strcmp(esdVert->GetTitle(),"vertexer: smearMC")) {
- esdVert->GetXYZ(vertex);
- AliDebug(1,Form("found vertex (%e,%e,%e)",vertex[0],vertex[1],vertex[2]));
+ fXExtrapVertex = esdVert->GetX();
+ fYExtrapVertex = esdVert->GetY();
+ fZExtrapVertex = esdVert->GetZ();
+ fXExtrapVertexError = TMath::Max(AliMFTConstants::fXVertexTolerance, esdVert->GetXRes());
+ fYExtrapVertexError = TMath::Max(AliMFTConstants::fYVertexTolerance, esdVert->GetYRes());
+ AliInfo(Form("Found ESD vertex from %d contributors (%f +/- %f, %f +/- %f, %f)",
+ esdVert->GetNContributors(),fXExtrapVertex,fXExtrapVertexError,fYExtrapVertex,fYExtrapVertexError,fZExtrapVertex));
}
-
+ else GetVertexFromMC();
+
//----------- Read ESD MUON tracks -------------------
Int_t nTracksMUON = event->GetNumberOfMuonTracks();
if (fMUONTrack) delete fMUONTrack;
fMUONTrack = new AliMUONTrack();
- AliMUONESDInterface::ESDToMUON(*esdTrack, *fMUONTrack, kTRUE); // last arguments should be kTRUE or kFALSE??
+ AliMUONESDInterface::ESDToMUON(*esdTrack, *fMUONTrack, kFALSE);
+
+ if (!fMUONTrack->GetTrackParamAtCluster()->First()) {
+ AliInfo("Skipping track, no parameters available!!!");
+ iTrack++;
+ continue;
+ }
// the track we are going to build, starting from fMUONTrack and adding the MFT clusters
AliMuonForwardTrack *track = new ((*fCandidateTracks)[0]) AliMuonForwardTrack();
Int_t nCandidates = fCandidateTracks->GetEntriesFast();
for (Int_t iCandidate=0; iCandidate<nCandidates; iCandidate++) {
+
+ if (!(fCandidateTracks->UncheckedAt(iCandidate))) continue;
fCurrentTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(iCandidate);
// if the old track is compatible with the new cluster, the track is updated and inserted as new track in the array
for (Int_t iFinalCandidate=0; iFinalCandidate<nFinalTracks; iFinalCandidate++) {
+ if (!(fCandidateTracks->UncheckedAt(iFinalCandidate))) continue;
AliMuonForwardTrack *finalTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(iFinalCandidate);
Int_t nMFTClusters = finalTrack->GetNMFTClusters();
//----------------------- Save the information to the AliESDMuonGlobalTrack object
- newTrack -> EvalKinem(vertex[2]); // we evaluate the kinematics at the primary vertex
+ newTrack -> EvalKinem(AliMFTConstants::fZEvalKinem);
AliDebug(2,"Creating a new Muon Global Track");
AliESDMuonGlobalTrack *myESDTrack = event->NewMuonGlobalTrack();
myESDTrack -> SetPxPyPz(newTrack->Px(), newTrack->Py(), newTrack->Pz());
-// waiting for the commit of the new version of AliESDMuonGlobalTrack...
-
-// myESDTrack -> SetChi2OverNdf(newTrack->GetChi2OverNdf());
-// myESDTrack -> SetCharge(newTrack->GetCharge());
-// myESDTrack -> SetMatchTrigger(newTrack->GetMatchTrigger());
-// myESDTrack -> SetFirstTrackingPoint(newTrack->GetMFTCluster(0)->GetX(), newTrack->GetMFTCluster(0)->GetY(), newTrack->GetMFTCluster(0)->GetZ());
-// myESDTrack -> SetXYAtVertex(newTrack->GetOffsetX(vertex[0], vertex[2]), newTrack->GetOffsetX(vertex[1], vertex[2]));
-// myESDTrack -> SetRAtAbsorberEnd(newTrack->GetRAtAbsorberEnd());
-// myESDTrack -> SetChi2MatchTrigger(esdTrack->GetChi2MatchTrigger());
-// myESDTrack -> SetMuonClusterMap(esdTrack->GetMuonClusterMap());
-// myESDTrack -> SetHitsPatternInTrigCh(esdTrack->GetHitsPatternInTrigCh());
-// myESDTrack -> SetHitsPatternInTrigChTrk(esdTrack->GetHitsPatternInTrigChTrk());
-// myESDTrack -> Connected(esdTrack->IsConnected());
-
+ myESDTrack -> SetLabel(newTrack->GetMCLabel());
+ myESDTrack -> SetChi2OverNdf(newTrack->GetChi2OverNdf());
+ myESDTrack -> SetCharge(newTrack->GetCharge());
+ myESDTrack -> SetMatchTrigger(newTrack->GetMatchTrigger());
+ myESDTrack -> SetNMFTClusters(newTrack->GetNMFTClusters());
+ myESDTrack -> SetNWrongMFTClustersMC(newTrack->GetNWrongClustersMC());
+ myESDTrack -> SetFirstTrackingPoint(newTrack->GetMFTCluster(0)->GetX(), newTrack->GetMFTCluster(0)->GetY(), newTrack->GetMFTCluster(0)->GetZ());
+ myESDTrack -> SetXYAtVertex(newTrack->GetOffsetX(0., AliMFTConstants::fZEvalKinem), newTrack->GetOffsetY(0., AliMFTConstants::fZEvalKinem));
+ myESDTrack -> SetRAtAbsorberEnd(newTrack->GetRAtAbsorberEnd());
+ myESDTrack -> SetCovariances(newTrack->GetTrackParamAtMFTCluster(0)->GetCovariances());
+ myESDTrack -> SetChi2MatchTrigger(esdTrack->GetChi2MatchTrigger());
+ myESDTrack -> SetMuonClusterMap(esdTrack->GetMuonClusterMap());
+ myESDTrack -> SetHitsPatternInTrigCh(esdTrack->GetHitsPatternInTrigCh());
+ myESDTrack -> SetHitsPatternInTrigChTrk(esdTrack->GetHitsPatternInTrigChTrk());
+ myESDTrack -> Connected(esdTrack->IsConnected());
+
+ ULong_t mftClusterPattern = 0;
+ for (Int_t iCluster=0; iCluster<newTrack->GetNMFTClusters(); iCluster++) {
+ AliMFTCluster *localCluster = newTrack->GetMFTCluster(iCluster);
+ mftClusterPattern |= 1 << localCluster->GetPlane();
+ mftClusterPattern |= IsCorrectMatch(localCluster, newTrack->GetMCLabel()) << (fNMaxPlanes + localCluster->GetPlane());
+ }
+ myESDTrack -> SetMFTClusterPattern(mftClusterPattern);
+
//---------------------------------------------------------------------------------
}
currentParamBack = (*((AliMUONTrackParam*)(fMUONTrack->GetTrackParamAtCluster()->First())));
currentParamForResearchFront = currentParamFront;
currentParamForResearchBack = currentParamBack;
- Double_t xExtrap = gRandom->Gaus(0,fVertexErrorX);
- Double_t yExtrap = gRandom->Gaus(0,fVertexErrorY);
- Double_t zExtrap = gRandom->Gaus(0,fVertexErrorZ);
if (fBransonCorrection) {
- AliMUONTrackExtrap::ExtrapToVertex(¤tParamFront, xExtrap, yExtrap, zExtrap, fVertexErrorX, fVertexErrorY);
- AliMUONTrackExtrap::ExtrapToVertex(¤tParamBack, xExtrap, yExtrap, zExtrap, fVertexErrorX, fVertexErrorY);
+ AliMUONTrackExtrap::ExtrapToVertex(¤tParamFront, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
+ AliMUONTrackExtrap::ExtrapToVertex(¤tParamBack, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
}
else {
- AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(¤tParamFront, zExtrap);
- AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(¤tParamBack, zExtrap);
+ AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(¤tParamFront, fZExtrapVertex);
+ AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(¤tParamBack, fZExtrapVertex);
}
- AliMUONTrackExtrap::ExtrapToVertex(¤tParamForResearchFront, xExtrap, yExtrap, zExtrap, fVertexErrorX, fVertexErrorY);
- AliMUONTrackExtrap::ExtrapToVertex(¤tParamForResearchBack, xExtrap, yExtrap, zExtrap, fVertexErrorX, fVertexErrorY);
+ AliMUONTrackExtrap::ExtrapToVertex(¤tParamForResearchFront, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
+ AliMUONTrackExtrap::ExtrapToVertex(¤tParamForResearchBack, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
}
else { // MFT planes others than the last one: mult. scattering correction because of the upstream MFT planes is performed
currentParamFront = (*((AliMUONTrackParam*)(fCurrentTrack->GetTrackParamAtCluster()->First())));
//=========================================================================================================================================
+Bool_t AliMFTTrackerMU::IsCorrectMatch(AliMFTCluster *cluster, Int_t labelMC) {
+
+ Bool_t result = kFALSE;
+
+ // check if the cluster belongs to the correct MC track
+
+ for (Int_t iTrack=0; iTrack<cluster->GetNMCTracks(); iTrack++) {
+ if (cluster->GetMCLabel(iTrack)==labelMC) {
+ result = kTRUE;
+ break;
+ }
+ }
+
+ return result;
+
+}
+
+//======================================================================================================================================
+
+void AliMFTTrackerMU::GetVertexFromMC() {
+
+ AliRunLoader *runLoader = AliRunLoader::Open("galice.root");
+ if (!runLoader) {
+ AliError("no run loader found in file galice.root");
+ return;
+ }
+
+ runLoader->CdGAFile();
+ runLoader->LoadgAlice();
+ runLoader->LoadHeader();
+ runLoader->GetEvent(gAlice->GetEvNumber());
+
+ TArrayF vtx(3);
+ runLoader->GetHeader()->GenEventHeader()->PrimaryVertex(vtx);
+ AliInfo(Form("Primary vertex from MC found in (%f, %f, %f)\n",vtx[0], vtx[1], vtx[2]));
+
+ fXExtrapVertex = gRandom->Gaus(vtx[0], AliMFTConstants::fPrimaryVertexResX);
+ fYExtrapVertex = gRandom->Gaus(vtx[1], AliMFTConstants::fPrimaryVertexResY);
+ fZExtrapVertex = gRandom->Gaus(vtx[2], AliMFTConstants::fPrimaryVertexResZ);
+ fXExtrapVertexError = AliMFTConstants::fXVertexTolerance;
+ fYExtrapVertexError = AliMFTConstants::fYVertexTolerance;
+ AliInfo(Form("Set ESD vertex from MC (%f +/- %f, %f +/- %f, %f)",
+ fXExtrapVertex,fXExtrapVertexError,fYExtrapVertex,fYExtrapVertexError,fZExtrapVertex));
+
+}
+
+//======================================================================================================================================
+
+void AliMFTTrackerMU::AddClustersFromUnderlyingEvent() {
+
+ AliInfo("Adding clusters from underlying event");
+
+ if (!fMFT) return;
+
+ TGrid::Connect("alien://");
+
+ TFile* fileWithClustersToAdd = TFile::Open(fMFT->GetFileNameForUnderlyingEvent());
+ if (!fileWithClustersToAdd) return;
+ if (!(fileWithClustersToAdd->IsOpen())) return;
+ if (!(fileWithClustersToAdd->cd(Form("Event%d",fMFT->GetUnderlyingEventID())))) return;
+
+ TClonesArray *recPointsPerPlaneToAdd[AliMFTConstants::fNMaxPlanes] = {0};
+ TTree *treeIn = 0;
+
+ for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) recPointsPerPlaneToAdd[iPlane] = new TClonesArray("AliMFTCluster");
+
+ treeIn = (TTree*) gDirectory->Get("TreeR");
+
+ for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
+ if (!(treeIn->GetBranch(Form("Plane_%02d",iPlane)))) {
+ for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) delete recPointsPerPlaneToAdd[jPlane];
+ return;
+ }
+ else treeIn->SetBranchAddress(Form("Plane_%02d",iPlane), &(recPointsPerPlaneToAdd[iPlane]));
+ }
+
+ treeIn -> GetEntry(0);
+
+ for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
+ printf("plane %d -> before = %d ",iPlane,fMFTClusterArray[iPlane]->GetEntries());
+ Int_t nClusters = recPointsPerPlaneToAdd[iPlane]->GetEntries();
+ for (Int_t iCluster=0; iCluster<nClusters; iCluster++) {
+ AliMFTCluster *newCluster = (AliMFTCluster*) recPointsPerPlaneToAdd[iPlane]->At(iCluster);
+ for (Int_t iTrack=0; iTrack<newCluster->GetNMCTracks(); iTrack++) newCluster->SetMCLabel(iTrack, newCluster->GetMCLabel(iTrack)+AliMFTConstants::fLabelOffsetMC);
+ new ((*fMFTClusterArray[iPlane])[fMFTClusterArray[iPlane]->GetEntries()]) AliMFTCluster(*newCluster);
+ }
+ printf("after = %d\n",fMFTClusterArray[iPlane]->GetEntries());
+ }
+
+ for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) delete recPointsPerPlaneToAdd[jPlane];
+
+}
+
+//======================================================================================================================================
+
+void AliMFTTrackerMU::AddClustersFromPileUpEvents() {
+
+ AliInfo("Adding clusters from pile-up event(s)");
+
+ if (!fMFT) return;
+
+ TGrid::Connect("alien://");
+
+ TFile* fileWithClustersToAdd = TFile::Open(fMFT->GetFileNameForPileUpEvents());
+ if (!fileWithClustersToAdd) return;
+ if (!(fileWithClustersToAdd->IsOpen())) return;
+
+ TClonesArray *recPointsPerPlaneToAdd[AliMFTConstants::fNMaxPlanes] = {0};
+ TTree *treeIn = 0;
+
+ for (Int_t iPileUp=0; iPileUp<AliMFTConstants::fNMaxPileUpEvents; iPileUp++) {
+
+ if (!(fileWithClustersToAdd->cd(Form("Event%d",fMFT->GetPileUpEventID(iPileUp))))) continue;
+
+ for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) recPointsPerPlaneToAdd[iPlane] = new TClonesArray("AliMFTCluster");
+
+ treeIn = (TTree*) gDirectory->Get("TreeR");
+
+ for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
+ if (!(treeIn->GetBranch(Form("Plane_%02d",iPlane)))) {
+ for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) delete recPointsPerPlaneToAdd[jPlane];
+ return;
+ }
+ else treeIn->SetBranchAddress(Form("Plane_%02d",iPlane), &(recPointsPerPlaneToAdd[iPlane]));
+ }
+
+ treeIn -> GetEntry(0);
+
+ for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
+ printf("plane %d -> before = %d ",iPlane,fMFTClusterArray[iPlane]->GetEntries());
+ Int_t nClusters = recPointsPerPlaneToAdd[iPlane]->GetEntries();
+ for (Int_t iCluster=0; iCluster<nClusters; iCluster++) {
+ AliMFTCluster *newCluster = (AliMFTCluster*) recPointsPerPlaneToAdd[iPlane]->At(iCluster);
+ for (Int_t iTrack=0; iTrack<newCluster->GetNMCTracks(); iTrack++) newCluster->SetMCLabel(iTrack, newCluster->GetMCLabel(iTrack)+AliMFTConstants::fLabelOffsetMC);
+ new ((*fMFTClusterArray[iPlane])[fMFTClusterArray[iPlane]->GetEntries()]) AliMFTCluster(*newCluster);
+ }
+ printf("after = %d\n",fMFTClusterArray[iPlane]->GetEntries());
+ }
+
+ for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) delete recPointsPerPlaneToAdd[jPlane];
+
+ }
+
+}
+
+//======================================================================================================================================
+