X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=MUON%2FMUONFakes.C;h=582761747c8c8beb11813e0222a177c2ccd4d166;hb=0249e9b86faf039c014b59bd206b610c9734aef8;hp=70ed78f1a8039fc1eba5db33f832473877f17c30;hpb=2303543417447e075cf208a9a42018f396b7ae74;p=u%2Fmrichter%2FAliRoot.git diff --git a/MUON/MUONFakes.C b/MUON/MUONFakes.C index 70ed78f1a80..582761747c8 100644 --- a/MUON/MUONFakes.C +++ b/MUON/MUONFakes.C @@ -1,30 +1,24 @@ #if !defined(__CINT__) || defined(__MAKECINT__) // ROOT includes -#include #include #include #include #include #include -#include -#include +#include #include // STEER includes #include "AliLog.h" -#include "AliMagF.h" #include "AliESDEvent.h" #include "AliESDMuonTrack.h" -#include "AliESDMuonCluster.h" -#include "AliCDBPath.h" -#include "AliCDBEntry.h" #include "AliCDBManager.h" // MUON includes +#include "AliMUONCDB.h" #include "AliMUONTrack.h" #include "AliMUONVTrackStore.h" #include "AliMUONTrackParam.h" -#include "AliMUONTrackExtrap.h" #include "AliMUONESDInterface.h" #include "AliMUONRecoCheck.h" #include "AliMUONVCluster.h" @@ -40,18 +34,17 @@ /// Results are saved in the root file Fakes.root /// Results are relevent provided that you use the same recoParams as for the reconstruction -void Prepare(AliMUONRecoParam *&recoParam, Double_t &sigmaCut); -TTree* GetESDTree(TFile *esdFile); -Bool_t TrackMatched(AliMUONTrack &track, AliMUONTrack &trackRef, Float_t &fractionOfMatchCluster, Double_t sigmaCut); -Bool_t IsRecontructible(AliMUONTrack &track, AliMUONRecoParam &recoParam); -AliMUONTrack* MatchWithTrackRef(AliESDMuonTrack &muonTrack, AliMUONVTrackStore &trackRefStore, - Float_t &fractionOfMatchCluster, Bool_t useLabel, Double_t sigmaCut); -Int_t RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStore &trackRefStore, AliMUONRecoParam &recoParam, - Bool_t useLabel, Double_t sigmaCut, TH1F &hFractionOfConnectedClusters); +UInt_t requestedStationMask = 0; +Bool_t request2ChInSameSt45 = kFALSE; +Double_t sigmaCut = -1.; + +Int_t RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStore &trackRefStore, + Bool_t useLabel, TH1F &hFractionOfConnectedClusters); //----------------------------------------------------------------------- void MUONFakes(Bool_t useLabel = kFALSE, Int_t FirstEvent = 0, Int_t LastEvent = -1, - const TString esdFileName = "AliESDs.root", const TString SimDir = "./generated/") + const TString esdFileName = "AliESDs.root", const TString SimDir = "./generated/", + const TString ocdbPath = "local://$ALICE_ROOT/OCDB") { //Reset ROOT and connect tree file @@ -86,19 +79,21 @@ void MUONFakes(Bool_t useLabel = kFALSE, Int_t FirstEvent = 0, Int_t LastEvent = TH1F *hPhiM = new TH1F("hPhiM"," matched track phi distribution",100,-1.,9.); TH1F *hPhiF = new TH1F("hPhiF"," fake track phi distribution",100,-1.,9.); - // prepare for analysis - AliMUONRecoParam *recoParam = 0x0; - Double_t sigmaCut = -1; - Prepare(recoParam, sigmaCut); + // link to reconstructed and simulated tracks + AliMUONRecoCheck rc(esdFileName, SimDir); - // link to reconstructed tracks - TFile* esdFile = TFile::Open(esdFileName); - TTree* esdTree = GetESDTree(esdFile); - AliESDEvent* esd = new AliESDEvent(); - esd->ReadFromTree(esdTree); + // load necessary data from OCDB + AliCDBManager::Instance()->SetDefaultStorage(ocdbPath); + AliCDBManager::Instance()->SetRun(rc.GetRunNumber()); + AliMUONRecoParam* recoParam = AliMUONCDB::LoadRecoParam(); + if (!recoParam) return; - // link to simulated tracks - AliMUONRecoCheck rc(esdFileName, SimDir); + // get sigma cut from recoParam to associate clusters with TrackRefs in case the label are not used + sigmaCut = (recoParam->ImproveTracks()) ? recoParam->GetSigmaCutForImprovement() : recoParam->GetSigmaCutForTracking(); + // compute the mask of requested stations from recoParam + for (Int_t i = 0; i < 5; i++) if (recoParam->RequestStation(i)) requestedStationMask |= ( 1 << i ); + // get from recoParam whether a track need 2 chambers hit in the same station (4 or 5) or not to be reconstructible + request2ChInSameSt45 = !recoParam->MakeMoreTrackCandidates(); // initialize global counters Int_t nReconstructibleTracks = 0; @@ -118,46 +113,45 @@ void MUONFakes(Bool_t useLabel = kFALSE, Int_t FirstEvent = 0, Int_t LastEvent = // Loop over ESD events FirstEvent = TMath::Max(0, FirstEvent); - LastEvent = (LastEvent>=0) ? TMath::Min((Int_t)esdTree->GetEntries() - 1, LastEvent) : (Int_t)esdTree->GetEntries() - 1; + LastEvent = (LastEvent>=0) ? TMath::Min(rc.NumberOfEvents() - 1, LastEvent) : rc.NumberOfEvents() - 1; for (Int_t iEvent = FirstEvent; iEvent <= LastEvent; iEvent++) { - // get the ESD of current event - esdTree->GetEvent(iEvent); - if (!esd) { - Error("CheckESD", "no ESD object found for event %d", iEvent); - return; - } - - // convert TrackRef to MUON tracks + // get reconstructed and simulated tracks + AliMUONVTrackStore* muonTrackStore = rc.ReconstructedTracks(iEvent, kFALSE); AliMUONVTrackStore* trackRefStore = rc.TrackRefs(iEvent); + if (!muonTrackStore || !trackRefStore) continue; // count the number of reconstructible tracks TIter next(trackRefStore->CreateIterator()); AliMUONTrack* trackRef; while ( ( trackRef = static_cast(next()) ) ) { - if (IsRecontructible(*trackRef,*recoParam)) nReconstructibleTracks++; + if (trackRef->IsValid(requestedStationMask, request2ChInSameSt45)) nReconstructibleTracks++; } // loop over ESD tracks Int_t nTrackerTracks = 0; trackReconstructedYet = kFALSE; AliMUONVTrackStore *fakeTrackStore = AliMUONESDInterface::NewTrackStore(); + const AliESDEvent* esd = rc.GetESDEvent(); Int_t nTracks = (Int_t)esd->GetNumberOfMuonTracks() ; for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) { - AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack); + AliESDMuonTrack* esdTrack = esd->GetMuonTrack(iTrack); // skip ghosts - if (!muonTrack->ContainTrackerData()) continue; + if (!esdTrack->ContainTrackerData()) continue; nTrackerTracks++; + // find the corresponding MUON track + AliMUONTrack* muonTrack = (AliMUONTrack*) muonTrackStore->FindObject(esdTrack->GetUniqueID()); + // get track info - Int_t nClusters = muonTrack->GetNClusters(); - Double_t normalizedChi2 = muonTrack->GetChi2() / (2. * muonTrack->GetNHit() - 5); - Double_t p = muonTrack->P(); - Double_t pT = muonTrack->Pt(); - Double_t eta = muonTrack->Eta(); - Double_t phi = muonTrack->Phi(); + Int_t nClusters = esdTrack->GetNClusters(); + Double_t normalizedChi2 = esdTrack->GetChi2() / (2. * esdTrack->GetNHit() - 5); + Double_t p = esdTrack->P(); + Double_t pT = esdTrack->Pt(); + Double_t eta = esdTrack->Eta(); + Double_t phi = esdTrack->Phi(); // fill global histograms hNumberOfClusters->Fill(nClusters); @@ -168,21 +162,21 @@ void MUONFakes(Bool_t useLabel = kFALSE, Int_t FirstEvent = 0, Int_t LastEvent = hPhi->Fill(phi); // try to match the reconstructed track with a simulated one - Float_t fractionOfMatchCluster = 0.; - AliMUONTrack* matchedTrackRef = MatchWithTrackRef(*muonTrack, *trackRefStore, fractionOfMatchCluster, useLabel, sigmaCut); + Int_t nMatchClusters = 0; + AliMUONTrack* matchedTrackRef = rc.FindCompatibleTrack(*muonTrack, *trackRefStore, nMatchClusters, useLabel, sigmaCut); // take actions according to matching result if (matchedTrackRef) { // global counter nTotMatchedTracks++; - if (!IsRecontructible(*matchedTrackRef,*recoParam)) { + if (!matchedTrackRef->IsValid(requestedStationMask, request2ChInSameSt45)) { trackReconstructedYet = kTRUE; nTotTracksReconstructedYet++; } // fill histograms - hFractionOfMatchedClusters->Fill(fractionOfMatchCluster); + hFractionOfMatchedClusters->Fill(((Float_t) nMatchClusters) / ((Float_t) nClusters)); hNumberOfClustersMC->Fill(matchedTrackRef->GetNClusters()); hNumberOfClustersM->Fill(nClusters); hChi2PerDofM->Fill(normalizedChi2); @@ -208,7 +202,7 @@ void MUONFakes(Bool_t useLabel = kFALSE, Int_t FirstEvent = 0, Int_t LastEvent = hPhiF->Fill(phi); // store fake tracks - AliMUONESDInterface::Add(*muonTrack, *fakeTrackStore); + fakeTrackStore->Add(*muonTrack); } @@ -223,8 +217,7 @@ void MUONFakes(Bool_t useLabel = kFALSE, Int_t FirstEvent = 0, Int_t LastEvent = if (fakeTrackStore->GetSize() > 0) { // remove the most connected fake tracks - Int_t nFreeMissingTracks = RemoveConnectedFakes(*fakeTrackStore, *trackRefStore, *recoParam, - useLabel, sigmaCut, *hFractionOfConnectedClusters); + Int_t nFreeMissingTracks = RemoveConnectedFakes(*fakeTrackStore, *trackRefStore, useLabel, *hFractionOfConnectedClusters); // remove the remaining free reconstructible tracks Int_t nAdditionalTracks = fakeTrackStore->GetSize() - nFreeMissingTracks; @@ -354,179 +347,11 @@ void MUONFakes(Bool_t useLabel = kFALSE, Int_t FirstEvent = 0, Int_t LastEvent = cout << "REMINDER: results are relevent provided that you use the same recoParams as for the reconstruction" << endl; cout << endl; - // clear memory - delete esd; - esdFile->Close(); - delete recoParam; - -} - -//----------------------------------------------------------------------- -void Prepare(AliMUONRecoParam *&recoParam, Double_t &sigmaCut) -{ - /// Set the magnetic field and return recoParam and sigmaCut to associate cluster and trackRef - - // prepare OCDB access - AliCDBManager* man = AliCDBManager::Instance(); - man->SetDefaultStorage("local://$ALICE_ROOT/OCDB"); - man->SetRun(0); - - // set mag field - // waiting for mag field in CDB - if (!TGeoGlobalMagField::Instance()->GetField()) { - printf("Loading field map...\n"); - AliMagF* field = new AliMagF("Maps","Maps",2,1.,1., 10.,AliMagF::k5kG); - TGeoGlobalMagField::Instance()->SetField(field); - } - // set the magnetic field for track extrapolations - AliMUONTrackExtrap::SetField(); - - // Load initial reconstruction parameters from OCDB - AliCDBPath path("MUON","Calib","RecoParam"); - AliCDBEntry *entry=man->Get(path.GetPath()); - if(entry) { - recoParam = dynamic_cast(entry->GetObject()); - entry->SetOwner(0); - AliCDBManager::Instance()->UnloadFromCache(path.GetPath()); - } - if (!recoParam) { - printf("Couldn't find RecoParam object in OCDB: create default one"); - recoParam = AliMUONRecoParam::GetLowFluxParam(); - } - - Info("MUONFakes", "\n recontruction parameters:"); - recoParam->Print("FULL"); - AliMUONESDInterface::ResetTracker(recoParam); - - // sigma cut to associate clusters with TrackRefs in case the label are not used - sigmaCut = (recoParam->ImproveTracks()) ? recoParam->GetSigmaCutForImprovement() : recoParam->GetSigmaCutForTracking(); - -} - -//----------------------------------------------------------------------- -TTree* GetESDTree(TFile *esdFile) -{ - /// Check that the file is properly open - /// Return pointer to the ESD Tree - - if (!esdFile || !esdFile->IsOpen()) { - Error("GetESDTree", "opening ESD file failed"); - exit(-1); - } - - TTree* tree = (TTree*) esdFile->Get("esdTree"); - if (!tree) { - Error("GetESDTree", "no ESD tree found"); - exit(-1); - } - - return tree; - -} - -//----------------------------------------------------------------------- -Bool_t TrackMatched(AliMUONTrack &track, AliMUONTrack &trackRef, Float_t &fractionOfMatchCluster, Double_t sigmaCut) -{ - /// Try to match 2 tracks - - Bool_t compTrack[10]; - Int_t nMatchClusters = track.CompatibleTrack(&trackRef, sigmaCut, compTrack); - fractionOfMatchCluster = ((Float_t)nMatchClusters) / ((Float_t)track.GetNClusters()); - - if ((compTrack[0] || compTrack[1] || compTrack[2] || compTrack[3]) && // at least 1 cluster matched in st 1 & 2 - (compTrack[6] || compTrack[7] || compTrack[8] || compTrack[9]) && // at least 1 cluster matched in st 4 & 5 - fractionOfMatchCluster > 0.5) return kTRUE; // more than 50% of clusters matched - else return kFALSE; - -} - -//----------------------------------------------------------------------- -Bool_t IsRecontructible(AliMUONTrack &track, AliMUONRecoParam &recoParam) -{ - /// Check il the track is reconstructible - Int_t nMinChHitInSt45 = (recoParam.MakeMoreTrackCandidates()) ? 2 : 3; - Int_t currentCh, previousCh = -1, nChHitInSt45 = 0; - Bool_t clusterInSt[5]; - for (Int_t iSt = 0; iSt < 5; iSt++) clusterInSt[iSt] = !recoParam.RequestStation(iSt); - - AliMUONTrackParam* trackParam = static_cast(track.GetTrackParamAtCluster()->First()); - while (trackParam) { - - currentCh = trackParam->GetClusterPtr()->GetChamberId(); - - clusterInSt[currentCh/2] = kTRUE; - - if (currentCh > 5 && currentCh != previousCh) { - nChHitInSt45++; - previousCh = currentCh; - } - - trackParam = static_cast(track.GetTrackParamAtCluster()->After(trackParam)); - } - - return (clusterInSt[0] && clusterInSt[1] && clusterInSt[2] && - clusterInSt[3] && clusterInSt[4] && nChHitInSt45 >= nMinChHitInSt45); - -} - -//----------------------------------------------------------------------- -AliMUONTrack* MatchWithTrackRef(AliESDMuonTrack &muonTrack, AliMUONVTrackStore &trackRefStore, - Float_t &fractionOfMatchCluster, Bool_t useLabel, Double_t sigmaCut) -{ - /// Return if the trackRef matched with the reconstructed track and the fraction of matched clusters - - AliMUONTrack *matchedTrackRef = 0x0; - fractionOfMatchCluster = 0.; - - if (useLabel) { // by using the MC label - - // get the corresponding simulated track if any - Int_t label = muonTrack.GetLabel(); - matchedTrackRef = (AliMUONTrack*) trackRefStore.FindObject(label); - - // get the fraction of matched clusters - if (matchedTrackRef) { - Int_t nMatchClusters = 0; - if (muonTrack.ClustersStored()) { - AliESDMuonCluster* cluster = (AliESDMuonCluster*) muonTrack.GetClusters().First(); - while (cluster) { - if (cluster->GetLabel() == label) nMatchClusters++; - cluster = (AliESDMuonCluster*) muonTrack.GetClusters().After(cluster); - } - } - fractionOfMatchCluster = ((Float_t)nMatchClusters) / ((Float_t)muonTrack.GetNClusters()); - } - - } else { // by comparing cluster/TrackRef positions - - // convert ESD track to MUON track - AliMUONTrack track; - AliMUONESDInterface::ESDToMUON(muonTrack,track); - - // look for the corresponding simulated track if any - TIter next(trackRefStore.CreateIterator()); - AliMUONTrack* trackRef; - while ( ( trackRef = static_cast(next()) ) ) { - - // check compatibility - Float_t f = 0.; - if (TrackMatched(track, *trackRef, f, sigmaCut)) { - matchedTrackRef = trackRef; - fractionOfMatchCluster = f; - break; - } - - } - - } - - return matchedTrackRef; - } //----------------------------------------------------------------------- -Int_t RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStore &trackRefStore, AliMUONRecoParam &recoParam, - Bool_t useLabel, Double_t sigmaCut, TH1F &hFractionOfConnectedClusters) +Int_t RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStore &trackRefStore, + Bool_t useLabel, TH1F &hFractionOfConnectedClusters) { /// loop over reconstructible TrackRef not associated with reconstructed track: /// for each of them, find and remove the most connected the fake track, if any, @@ -541,13 +366,13 @@ Int_t RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStor while ( ( trackRef = static_cast(next()) ) ) { // skip not reconstructible trackRefs - if (!IsRecontructible(*trackRef,recoParam)) continue; + if (!trackRef->IsValid(requestedStationMask, request2ChInSameSt45)) continue; Int_t label = trackRef->GetUniqueID(); // look for the most connected fake track AliMUONTrack *connectedFake = 0x0; - Float_t fractionOfConnectedClusters = 0.; + Double_t fractionOfConnectedClusters = 0.; TIter next2(fakeTrackStore.CreateIterator()); AliMUONTrack* fakeTrack; while ( ( fakeTrack = static_cast(next2()) ) ) { @@ -560,14 +385,14 @@ Int_t RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStor nConnectedClusters++; } else { // by comparing cluster/TrackRef positions Bool_t compTrack[10]; - nConnectedClusters = fakeTrack->CompatibleTrack(trackRef, sigmaCut, compTrack); + nConnectedClusters = fakeTrack->FindCompatibleClusters(*trackRef, sigmaCut, compTrack); } // skip non-connected fake tracks if (nConnectedClusters == 0) continue; // check if it is the most connected fake track - Float_t f = ((Float_t)nConnectedClusters) / ((Float_t)fakeTrack->GetNClusters()); + Double_t f = ((Double_t)nConnectedClusters) / ((Double_t)fakeTrack->GetNClusters()); if (f > fractionOfConnectedClusters) { connectedFake = fakeTrack; fractionOfConnectedClusters = f;