From: martinez Date: Tue, 17 Jun 2008 16:39:51 +0000 (+0000) Subject: New class for applying cut in ESDMuonTracks. Example of how to use ESDMuonTrackCuts... X-Git-Url: http://git.uio.no/git/?a=commitdiff_plain;h=f70a1b5d1760bab2733a23076b618ee41b1b5f2f;hp=6240e6a4bdf51c8f3a0b84b37a29f30c497800b9;p=u%2Fmrichter%2FAliRoot.git New class for applying cut in ESDMuonTracks. Example of how to use ESDMuonTrackCuts class in the analysis tasks. Removing obsolete macro RunESDMuonFilter. (Roberta) --- diff --git a/PWG3/PWG3muonLinkDef.h b/PWG3/PWG3muonLinkDef.h index 4d05f326e97..47237195610 100644 --- a/PWG3/PWG3muonLinkDef.h +++ b/PWG3/PWG3muonLinkDef.h @@ -11,6 +11,7 @@ #pragma link C++ class AliAnalysisTaskLUT+; #pragma link C++ class AliAnalysisTaskESDMuonFilter+; #pragma link C++ class AliAnalysisTaskTrigChEff+; +#pragma link C++ class AliESDMuonTrackCuts+; #endif diff --git a/PWG3/libPWG3muon.pkg b/PWG3/libPWG3muon.pkg index 7c8db04ac51..ef4aeb0d525 100644 --- a/PWG3/libPWG3muon.pkg +++ b/PWG3/libPWG3muon.pkg @@ -6,7 +6,8 @@ SRCS:= muon/AliAnalysisTaskESDMuonFilter.cxx \ muon/AliAnalysisTaskLUT.cxx \ muon/AliAnalysisTaskTrigChEff.cxx \ muon/AliAODDimuon.cxx \ - muon/AliAODEventInfo.cxx + muon/AliAODEventInfo.cxx \ + muon/AliESDMuonTrackCuts.cxx HDRS:= $(SRCS:.cxx=.h) diff --git a/PWG3/muon/AliAnalysisTaskESDMuonFilter.cxx b/PWG3/muon/AliAnalysisTaskESDMuonFilter.cxx index c0b377b9ee3..0b580f86977 100644 --- a/PWG3/muon/AliAnalysisTaskESDMuonFilter.cxx +++ b/PWG3/muon/AliAnalysisTaskESDMuonFilter.cxx @@ -15,7 +15,7 @@ // Add the muon tracks to the generic AOD track branch during the // filtering of the ESD - R. Arnaldi 5/5/08 - + #include #include @@ -37,15 +37,15 @@ ClassImp(AliAnalysisTaskESDMuonFilter) //////////////////////////////////////////////////////////////////////// AliAnalysisTaskESDMuonFilter::AliAnalysisTaskESDMuonFilter(): - AliAnalysisTaskSE(), - fTrackFilter(0x0) + AliAnalysisTaskSE(), + fTrackFilter(0x0) { // Default constructor } AliAnalysisTaskESDMuonFilter::AliAnalysisTaskESDMuonFilter(const char* name): - AliAnalysisTaskSE(name), - fTrackFilter(0x0) + AliAnalysisTaskSE(name), + fTrackFilter(0x0) { // Constructor } @@ -53,13 +53,13 @@ AliAnalysisTaskESDMuonFilter::AliAnalysisTaskESDMuonFilter(const char* name): void AliAnalysisTaskESDMuonFilter::UserCreateOutputObjects() { // Create the output container - if (fTrackFilter) OutputTree()->GetUserInfo()->Add(fTrackFilter); + if (fTrackFilter) OutputTree()->GetUserInfo()->Add(fTrackFilter); } void AliAnalysisTaskESDMuonFilter::Init() { // Initialization - if (fDebug > 1) AliInfo("Init() \n"); + if (fDebug > 1) AliInfo("Init() \n"); } @@ -68,7 +68,7 @@ void AliAnalysisTaskESDMuonFilter::UserExec(Option_t */*option*/) // Execute analysis for current event Long64_t ientry = Entry(); printf("Muon Filter: Analysing event # %5d\n", (Int_t) ientry); - + ConvertESDtoAOD(); } @@ -76,87 +76,93 @@ void AliAnalysisTaskESDMuonFilter::ConvertESDtoAOD() { // ESD Muon Filter analysis task executed for each event AliESDEvent* esd = dynamic_cast(InputEvent()); - + // Define arrays for muons Double_t pos[3]; Double_t p[3]; Double_t pid[10]; - - // has to be changed once the muon pid is provided by the ESD - for (Int_t i = 0; i < 10; pid[i++] = 0.); - pid[AliAODTrack::kMuon]=1.; - + AliAODHeader* header = AODEvent()->GetHeader(); AliAODTrack *aodTrack = 0x0; - AliESDMuonTrack *esdMuTrack = 0x0; - + // Access to the AOD container of tracks TClonesArray &tracks = *(AODEvent()->GetTracks()); Int_t jTracks = header->GetRefMultiplicity(); // Read primary vertex from AOD event - AliAODVertex *primary = AODEvent()->GetPrimaryVertex(); + AliAODVertex * primary = AODEvent()->GetPrimaryVertex(); primary->Print(); - + // Loop on muon tracks to fill the AOD track branch Int_t nMuTracks = esd->GetNumberOfMuonTracks(); printf("Number of Muon Tracks=%d\n",nMuTracks); - - // Update number of positive and negative tracks from AOD event (M.G.) + + // Update number of positive and negative tracks from AOD event (M.G.) Int_t nPosTracks = header->GetRefMultiplicityPos(); Int_t nNegTracks = header->GetRefMultiplicityNeg(); - + for (Int_t nMuTrack = 0; nMuTrack < nMuTracks; ++nMuTrack) { - esdMuTrack = esd->GetMuonTrack(nMuTrack); - - if (!esdMuTrack->ContainTrackerData()) continue; - - p[0] = esdMuTrack->Px(); - p[1] = esdMuTrack->Py(); - p[2] = esdMuTrack->Pz(); - - pos[0] = esdMuTrack->GetNonBendingCoor(); - pos[1] = esdMuTrack->GetBendingCoor(); - pos[2] = esdMuTrack->GetZ(); - - aodTrack = new(tracks[jTracks++]) AliAODTrack(esdMuTrack->GetUniqueID(), // ID - 0, // label - p, // momentum - kTRUE, // cartesian coordinate system - pos, // position - kFALSE, // isDCA - 0x0, // covariance matrix - esdMuTrack->Charge(), // charge - 0, // ITSClusterMap - pid, // pid - primary, // primary vertex - kFALSE, // used for vertex fit? - kFALSE, // used for primary vertex fit? - AliAODTrack::kPrimary); // track type - - aodTrack->SetXYAtDCA(esdMuTrack->GetNonBendingCoorAtDCA(), esdMuTrack->GetBendingCoorAtDCA()); - aodTrack->SetPxPyPzAtDCA(esdMuTrack->PxAtDCA(), esdMuTrack->PyAtDCA(), esdMuTrack->PzAtDCA()); - aodTrack->ConvertAliPIDtoAODPID(); - aodTrack->SetChi2perNDF(esdMuTrack->GetChi2() / (2.*esdMuTrack->GetNHit() - 5.)); - aodTrack->SetChi2MatchTrigger(esdMuTrack->GetChi2MatchTrigger()); - aodTrack->SetHitsPatternInTrigCh(esdMuTrack->GetHitsPatternInTrigCh()); - aodTrack->SetMuonClusterMap(esdMuTrack->GetMuonClusterMap()); - aodTrack->SetMatchTrigger(esdMuTrack->GetMatchTrigger()); - - primary->AddDaughter(aodTrack); - - if (esdMuTrack->Charge() > 0) nPosTracks++; - else nNegTracks++; - } - + + AliESDMuonTrack *esdMuTrack = esd->GetMuonTrack(nMuTrack); + + UInt_t selectInfo = 0; + // Track selection + if (fTrackFilter) { + selectInfo = fTrackFilter->IsSelected(esdMuTrack); + if (!selectInfo) { + continue; + } + } + p[0] = esdMuTrack->Px(); + p[1] = esdMuTrack->Py(); + p[2] = esdMuTrack->Pz(); + pos[0] = primary->GetX(); + pos[1] = primary->GetY(); + pos[2] = primary->GetZ(); + + // has to be changed once the muon pid is provided by the ESD + for (Int_t i = 0; i < 10; pid[i++] = 0.); pid[AliAODTrack::kMuon]=1.; + primary->AddDaughter(aodTrack = + new(tracks[jTracks++]) AliAODTrack(0 , // no ID provided + 0, // no label provided + p, + kTRUE, + pos, + kFALSE, + NULL, // no covariance matrix provided + esdMuTrack->Charge(), + 0, // no ITSClusterMap + pid, + primary, + kTRUE, // check if this is right + kTRUE, // not used for vertex fit + AliAODTrack::kPrimary, + selectInfo) + ); + + if (esdMuTrack->Charge() > 0) nPosTracks++; + else nNegTracks++; + + aodTrack->ConvertAliPIDtoAODPID(); + aodTrack->SetHitsPatternInTrigCh(esdMuTrack->GetHitsPatternInTrigCh()); + Int_t track2Trigger = esdMuTrack->GetMatchTrigger(); + aodTrack->SetMatchTrigger(track2Trigger); + if (track2Trigger) + aodTrack->SetChi2MatchTrigger(esdMuTrack->GetChi2MatchTrigger()); + else + aodTrack->SetChi2MatchTrigger(0.); + } + header->SetRefMultiplicity(jTracks); header->SetRefMultiplicityPos(nPosTracks); header->SetRefMultiplicityNeg(nNegTracks); + + return; } void AliAnalysisTaskESDMuonFilter::Terminate(Option_t */*option*/) { - // Terminate analysis - // - if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n"); +// Terminate analysis +// + if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n"); } diff --git a/PWG3/muon/AliESDMuonTrackCuts.cxx b/PWG3/muon/AliESDMuonTrackCuts.cxx new file mode 100644 index 00000000000..1dfac25bce8 --- /dev/null +++ b/PWG3/muon/AliESDMuonTrackCuts.cxx @@ -0,0 +1,580 @@ +/************************************************************************** + * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * + * * + * Author: The ALICE Off-line Project. * + * Contributors are mentioned in the code where appropriate. * + * * + * Permission to use, copy, modify and distribute this software and its * + * documentation strictly for non-commercial purposes is hereby granted * + * without fee, provided that the above copyright notice appears in all * + * copies and that both the copyright notice and this permission notice * + * appear in the supporting documentation. The authors make no claims * + * about the suitability of this software for any purpose. It is * + * provided "as is" without express or implied warranty. * + **************************************************************************/ + +#include "AliESDMuonTrackCuts.h" + +#include +#include +#include +#include + +#include +#include +#include +#include + +//____________________________________________________________________ +ClassImp(AliESDMuonTrackCuts) + +// Cut names +const Char_t* AliESDMuonTrackCuts::fgkCutNames[kNCuts] = { + "p", + "p_{T}", + "p_{x}", + "p_{y}", + "p_{z}", + "y", + "eta" +}; + +//____________________________________________________________________ +AliESDMuonTrackCuts::AliESDMuonTrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title), + fPMin(0), + fPMax(0), + fPtMin(0), + fPtMax(0), + fPxMin(0), + fPxMax(0), + fPyMin(0), + fPyMax(0), + fPzMin(0), + fPzMax(0), + fEtaMin(0), + fEtaMax(0), + fRapMin(0), + fRapMax(0), + fHistogramsOn(0), + fhCutStatistics(0), + fhCutCorrelation(0) +{ + // + // constructor + // + Init(); + + //############################################################################## + // setting default cuts + SetPRange(); + SetPtRange(); + SetPxRange(); + SetPyRange(); + SetPzRange(); + SetEtaRange(); + SetRapRange(); + + SetHistogramsOn(); +} + +//_____________________________________________________________________________ +AliESDMuonTrackCuts::AliESDMuonTrackCuts(const AliESDMuonTrackCuts &c) : AliAnalysisCuts(c), + fPMin(0), + fPMax(0), + fPtMin(0), + fPtMax(0), + fPxMin(0), + fPxMax(0), + fPyMin(0), + fPyMax(0), + fPzMin(0), + fPzMax(0), + fEtaMin(0), + fEtaMax(0), + fRapMin(0), + fRapMax(0), + fHistogramsOn(0), + fhCutStatistics(0), + fhCutCorrelation(0) +{ + // + // copy constructor + // + ((AliESDMuonTrackCuts &) c).Copy(*this); +} + +AliESDMuonTrackCuts::~AliESDMuonTrackCuts() +{ + // + // destructor + // + for (Int_t i=0; i<2; i++) { + if (fhPt[i]) + delete fhPt[i]; + if (fhEta[i]) + delete fhEta[i]; + } + + if (fhCutStatistics) + delete fhCutStatistics; + if (fhCutCorrelation) + delete fhCutCorrelation; +} + +void AliESDMuonTrackCuts::Init() +{ + // + // sets everything to zero + // + fPMin = 0; + fPMax = 0; + fPtMin = 0; + fPtMax = 0; + fPxMin = 0; + fPxMax = 0; + fPyMin = 0; + fPyMax = 0; + fPzMin = 0; + fPzMax = 0; + fEtaMin = 0; + fEtaMax = 0; + fRapMin = 0; + fRapMax = 0; + + fHistogramsOn = kFALSE; + + for (Int_t i=0; i<2; ++i) + { + fhPt[i] = 0; + fhEta[i] = 0; + } + fhCutStatistics = 0; + fhCutCorrelation = 0; +} + +//_____________________________________________________________________________ +AliESDMuonTrackCuts &AliESDMuonTrackCuts::operator=(const AliESDMuonTrackCuts &c) +{ + // + // Assignment operator + // + if (this != &c) ((AliESDMuonTrackCuts &) c).Copy(*this); + return *this; +} + +//_____________________________________________________________________________ +void AliESDMuonTrackCuts::Copy(TObject &c) const +{ + // + // Copy function + // + + AliESDMuonTrackCuts& target = (AliESDMuonTrackCuts &) c; + + target.Init(); +// + target.fPMin = fPMin; + target.fPMax = fPMax; + target.fPtMin = fPtMin; + target.fPtMax = fPtMax; + target.fPxMin = fPxMin; + target.fPxMax = fPxMax; + target.fPyMin = fPyMin; + target.fPyMax = fPyMax; + target.fPzMin = fPzMin; + target.fPzMax = fPzMax; + target.fEtaMin = fEtaMin; + target.fEtaMax = fEtaMax; + target.fRapMin = fRapMin; + target.fRapMax = fRapMax; + + target.fHistogramsOn = fHistogramsOn; + + for (Int_t i=0; i<2; ++i) + { + if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone(); + if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone(); + } + + if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone(); + if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone(); + + TNamed::Copy(c); +} + +//_____________________________________________________________________________ +Long64_t AliESDMuonTrackCuts::Merge(TCollection* list) { + // Merge a list of AliESDMuonTrackCuts objects with this (needed for PROOF) + // Returns the number of merged objects (including this) + + if (!list) + return 0; + + if (list->IsEmpty()) + return 1; + + if (!fHistogramsOn) + return 0; + + TIterator* iter = list->MakeIterator(); + TObject* obj; + + // collection of measured and generated histograms + Int_t count = 0; + while ((obj = iter->Next())) { + + AliESDMuonTrackCuts* entry = dynamic_cast(obj); + if (entry == 0) + continue; + + if (!entry->fHistogramsOn) + continue; + + for (Int_t i=0; i<2; i++) { + fhPt[i]->Add(entry->fhPt[i]); + fhEta[i]->Add(entry->fhEta[i]); + } + + fhCutStatistics->Add(entry->fhCutStatistics); + fhCutCorrelation ->Add(entry->fhCutCorrelation); + + count++; + } + + return count+1; +} + +void AliESDMuonTrackCuts::EnableNeededBranches(TTree* tree) +{ + // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack + + tree->SetBranchStatus("fTracks.fFlags", 1); + tree->SetBranchStatus("fTracks.fP*", 1); + tree->SetBranchStatus("fTracks.fR*", 1); //detector response probability +} + +//____________________________________________________________________ +Bool_t +AliESDMuonTrackCuts::AcceptTrack(AliESDMuonTrack* esdMuTrack) { + // + // figure out if the tracks survives all the track cuts defined + // + // the different kinematic values are first + // retrieved from the track. then it is found out what cuts the + // track did not survive and finally the cuts are imposed. + + // getting the kinematic variables of the track + // (assuming the mass is known) + Double_t p[3]; + esdMuTrack->PxPyPz(p); + Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2)); + Float_t pt = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2)); + Float_t energy = TMath::Sqrt(TMath::Power(esdMuTrack->M(),2) + TMath::Power(momentum,2)); + + + //y-eta related calculations + Float_t eta = -100.; + Float_t y = -100.; + if((momentum != TMath::Abs(p[2]))&&(momentum != 0)) + eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2])); + if((energy != TMath::Abs(p[2]))&&(momentum != 0)) + y = 0.5*TMath::Log((energy + p[2])/(energy - p[2])); + + //######################################################################## + // cut the track? + + Bool_t cuts[kNCuts]; + for (Int_t i=0; i fPMax)) + cuts[0]=kTRUE; + if((pt < fPtMin) || (pt > fPtMax)) + cuts[1] = kTRUE; + if((p[0] < fPxMin) || (p[0] > fPxMax)) + cuts[2] = kTRUE; + if((p[1] < fPyMin) || (p[1] > fPyMax)) + cuts[3] = kTRUE; + if((p[2] < fPzMin) || (p[2] > fPzMax)) + cuts[4] = kTRUE; + if((eta < fEtaMin) || (eta > fEtaMax)) + cuts[5] = kTRUE; + if((y < fRapMin) || (y > fRapMax)) + cuts[6] = kTRUE; + + Bool_t cut=kFALSE; + for (Int_t i=0; iFill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks"))); + + if (cut) + fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks"))); + + for (Int_t i=0; iFill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i]))); + + for (Int_t j=i; jGetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i])); + Float_t yC = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j])); + fhCutCorrelation->Fill(xC, yC); + } + } + } + + fhPt[0]->Fill(pt); + fhEta[0]->Fill(eta); + + } + + //######################################################################## + // cut the track! + if (cut) return kFALSE; + + //######################################################################## + // filling histograms after cut + if (fHistogramsOn) { +// + fhPt[1]->Fill(pt); + fhEta[1]->Fill(eta); + + } + + return kTRUE; +} + +//____________________________________________________________________ +TObjArray* AliESDMuonTrackCuts::GetAcceptedTracks(AliESD* esd) +{ + // + // returns an array of all tracks that pass the cuts + // + + TObjArray* acceptedTracks = new TObjArray(); + + // loop over esd tracks + for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) { + AliESDMuonTrack* track = esd->GetMuonTrack(iTrack); + + if (AcceptTrack(track)) + acceptedTracks->Add(track); + } + + return acceptedTracks; +} + +//____________________________________________________________________ +Int_t AliESDMuonTrackCuts::CountAcceptedTracks(AliESD* esd) +{ + // + // returns an the number of tracks that pass the cuts + // + + Int_t count = 0; + + // loop over esd tracks + for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) { + AliESDMuonTrack* track = esd->GetMuonTrack(iTrack); + + if (AcceptTrack(track)) + count++; + } + + return count; +} + +//____________________________________________________________________ +TObjArray* AliESDMuonTrackCuts::GetAcceptedTracks(AliESDEvent* esd) +{ + // + // returns an array of all tracks that pass the cuts + // + + TObjArray* acceptedTracks = new TObjArray(); + + // loop over esd tracks + for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) { + AliESDMuonTrack* track = esd->GetMuonTrack(iTrack); + + if (AcceptTrack(track)) + acceptedTracks->Add(track); + } + + return acceptedTracks; +} + +//____________________________________________________________________ +Int_t AliESDMuonTrackCuts::CountAcceptedTracks(AliESDEvent* esd) +{ + // + // returns an the number of tracks that pass the cuts + // + + Int_t count = 0; + + // loop over esd tracks + for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) { + AliESDMuonTrack* track = esd->GetMuonTrack(iTrack); + + if (AcceptTrack(track)) + count++; + } + + return count; +} + +//____________________________________________________________________ + void AliESDMuonTrackCuts::DefineHistograms(Int_t color) { + // + // diagnostics histograms are defined + // + + fHistogramsOn=kTRUE; + + Bool_t oldStatus = TH1::AddDirectoryStatus(); + TH1::AddDirectory(kFALSE); + + //################################################################################### + // defining histograms + + fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5); + + fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks"); + fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks"); + + fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);; + + for (Int_t i=0; iGetXaxis()->SetBinLabel(i+4,fgkCutNames[i]); + fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]); + fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]); + } + + fhCutStatistics ->SetLineColor(color); + fhCutCorrelation ->SetLineColor(color); + fhCutStatistics ->SetLineWidth(2); + fhCutCorrelation ->SetLineWidth(2); + + Char_t str[256]; + for (Int_t i=0; i<2; i++) { + if (i==0) sprintf(str," "); + else sprintf(str,"_cut"); + + fhPt[i] = new TH1F(Form("pt%s",str) ,"p_{T} distribution;p_{T} (GeV/c)",500,0.0,100.0); + fhEta[i] = new TH1F(Form("eta%s",str) ,"#eta distribution;#eta",40,-2.0,2.0); + } + + TH1::AddDirectory(oldStatus); +} + +//____________________________________________________________________ +Bool_t AliESDMuonTrackCuts::LoadHistograms(const Char_t* dir) +{ + // + // loads the histograms from a file + // if dir is empty a directory with the name of this object is taken (like in SaveHistogram) + // + + if (!dir) + dir = GetName(); + + if (!gDirectory->cd(dir)) + return kFALSE; + + fhCutStatistics = dynamic_cast (gDirectory->Get("cut_statistics")); + fhCutCorrelation = dynamic_cast (gDirectory->Get("cut_correlation")); + + Char_t str[5]; + for (Int_t i=0; i<2; i++) { + if (i==0) + { + gDirectory->cd("before_cuts"); + str[0] = 0; + } + else + { + gDirectory->cd("after_cuts"); + sprintf(str,"_cut"); + } + + fhPt[i] = dynamic_cast (gDirectory->Get(Form("pt%s",str))); + fhEta[i] = dynamic_cast (gDirectory->Get(Form("eta%s",str))); + + gDirectory->cd("../"); + } + + gDirectory->cd(".."); + + return kTRUE; +} + +//____________________________________________________________________ +void AliESDMuonTrackCuts::SaveHistograms(const Char_t* dir) { + // + // saves the histograms in a directory (dir) + // + + if (!fHistogramsOn) { + AliDebug(0, "Histograms not on - cannot save histograms!!!"); + return; + } + + if (!dir) + dir = GetName(); + + gDirectory->mkdir(dir); + gDirectory->cd(dir); + + gDirectory->mkdir("before_cuts"); + gDirectory->mkdir("after_cuts"); + + fhCutStatistics->Write(); + fhCutCorrelation->Write(); + + for (Int_t i=0; i<2; i++) { + if (i==0) + gDirectory->cd("before_cuts"); + else + gDirectory->cd("after_cuts"); + + fhPt[i] ->Write(); + fhEta[i] ->Write(); + + gDirectory->cd("../"); + } + + gDirectory->cd("../"); +} + +//____________________________________________________________________ +void AliESDMuonTrackCuts::DrawHistograms() +{ + gStyle->SetPalette(1); + gStyle->SetFrameFillColor(10); + gStyle->SetCanvasColor(10); + + TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Cut Results", 800, 500); + canvas1->Divide(2, 1); + + canvas1->cd(1); + fhCutStatistics->SetStats(kFALSE); + fhCutStatistics->LabelsOption("v"); + gPad->SetBottomMargin(0.3); + fhCutStatistics->Draw(); + + canvas1->cd(2); + fhCutCorrelation->SetStats(kFALSE); + fhCutCorrelation->LabelsOption("v"); + gPad->SetBottomMargin(0.3); + gPad->SetLeftMargin(0.3); + fhCutCorrelation->Draw("COLZ"); + canvas1->Update(); + canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName())); + +} + diff --git a/PWG3/muon/AliESDMuonTrackCuts.h b/PWG3/muon/AliESDMuonTrackCuts.h new file mode 100644 index 00000000000..a19196aaac6 --- /dev/null +++ b/PWG3/muon/AliESDMuonTrackCuts.h @@ -0,0 +1,97 @@ +// +// Class for handling of ESD Muon track cuts +// (based on ANALYSIS/AliESDtrackCuts). +// +// The class manages some kinematic cuts. Two methods +// can be used to figure out if an ESD Muon track survives the cuts: +// AcceptTrack which takes a single AliESDMuonTrack as argument and +// returns kTRUE/kFALSE or GetAcceptedTracks which takes an AliESD +// object and returns an TObjArray (of AliESDMuonTracks) with the tracks +// in the ESD that survived the cuts. +// + +#ifndef ALIESDMUONTRACKCUTS_H +#define ALIESDMUONTRACKCUTS_H + +#include +#include +#include "AliAnalysisCuts.h" + +class AliESD; +class AliESDEvent; +class AliESDMuonTrack; +class AliLog; +class TTree; + +class AliESDMuonTrackCuts : public AliAnalysisCuts +{ +public: + AliESDMuonTrackCuts(const Char_t* name = "AliESDMuonTrackCuts", const Char_t* title = ""); + virtual ~AliESDMuonTrackCuts(); + Bool_t IsSelected(TObject* obj) + {return AcceptTrack((AliESDMuonTrack*)obj);} + Bool_t IsSelected(TList* /*list*/) {return kTRUE;} + Bool_t AcceptTrack(AliESDMuonTrack* esdMuTrack); + TObjArray* GetAcceptedTracks(AliESD* esd); + Int_t CountAcceptedTracks(AliESD* esd); + TObjArray* GetAcceptedTracks(AliESDEvent* esd); + Int_t CountAcceptedTracks(AliESDEvent* esd); + + virtual Long64_t Merge(TCollection* list); + virtual void Copy(TObject &c) const; + AliESDMuonTrackCuts(const AliESDMuonTrackCuts& pd); // Copy Constructor + AliESDMuonTrackCuts &operator=(const AliESDMuonTrackCuts &c); + + //###################################################### + + // track kinematic cut setters + void SetPRange(Float_t r1=0, Float_t r2=1e10) {fPMin=r1; fPMax=r2;} + void SetPtRange(Float_t r1=0, Float_t r2=1e10) {fPtMin=r1; fPtMax=r2;} + void SetPxRange(Float_t r1=-1e10, Float_t r2=1e10) {fPxMin=r1; fPxMax=r2;} + void SetPyRange(Float_t r1=-1e10, Float_t r2=1e10) {fPyMin=r1; fPyMax=r2;} + void SetPzRange(Float_t r1=-1e10, Float_t r2=1e10) {fPzMin=r1; fPzMax=r2;} + void SetEtaRange(Float_t r1=-1e10, Float_t r2=1e10) {fEtaMin=r1; fEtaMax=r2;} + void SetRapRange(Float_t r1=-1e10, Float_t r2=1e10) {fRapMin=r1; fRapMax=r2;} + + //###################################################### + + void SetHistogramsOn(Bool_t b=kFALSE) {fHistogramsOn = b;} + void DefineHistograms(Int_t color=1); + virtual Bool_t LoadHistograms(const Char_t* dir = 0); + void SaveHistograms(const Char_t* dir = 0); + void DrawHistograms(); + + static void EnableNeededBranches(TTree* tree); + +protected: + void Init(); // sets everything to 0 + + enum { kNCuts = 7 }; + + //###################################################### + static const Char_t* fgkCutNames[kNCuts]; //! names of cuts (for internal use) + + // esd kinematics cuts + Float_t fPMin, fPMax; // definition of the range of the P + Float_t fPtMin, fPtMax; // definition of the range of the Pt + Float_t fPxMin, fPxMax; // definition of the range of the Px + Float_t fPyMin, fPyMax; // definition of the range of the Py + Float_t fPzMin, fPzMax; // definition of the range of the Pz + Float_t fEtaMin, fEtaMax; // definition of the range of the eta + Float_t fRapMin, fRapMax; // definition of the range of the y + + //###################################################### + // diagnostics histograms + Bool_t fHistogramsOn; // histograms on/off + + TH1F* fhPt[2]; //-> pt of esd tracks + TH1F* fhEta[2]; //-> eta of esd tracks + + TH1F* fhCutStatistics; //-> statistics of what cuts the tracks did not survive + TH2F* fhCutCorrelation; //-> 2d statistics plot + + ClassDef(AliESDMuonTrackCuts, 1) +}; + + +#endif diff --git a/PWG3/muon/AnalysisTrainMuonCAF.C b/PWG3/muon/AnalysisTrainMuonCAF.C index 01db2fde2f6..8e27c339993 100644 --- a/PWG3/muon/AnalysisTrainMuonCAF.C +++ b/PWG3/muon/AnalysisTrainMuonCAF.C @@ -1,124 +1,136 @@ -void AnalysisTrainMuonCAF(char* fileout = "AliAOD.root", char *datasetname = "myDataSet", Int_t nev=1234567890) -{ -// Macro to produce a generic AOD starting from an ESD file. -// The AOD is filled with two tasks: -// 1- with the first one (AliAnalysisTaskESDfilter), -// all the branches of the AOD are filled apart from the muons. -// 2- with the second task (AliAnalysisTaskESDMuonFilter) -// muons tracks are added to the tracks branch -// This macro works on the CAF -// R. Arnaldi 4/5/08 - - gSystem->Load("libTree.so"); - gSystem->Load("libGeom.so"); - gSystem->Load("libVMC.so"); - gSystem->Load("libPhysics.so"); - - // Reset user processes if CAF if not responding anymore - // TProof::Reset("lxb6046"); - - // Connect to proof - TProof::Open("lxb6046"); // may be username@lxb6046 if user not the same as on local - - // Clear packages if changing ROOT version on CAF or local - // gProof->ClearPackages(); - // Enable proof debugging if needed - // gProof->SetLogLevel(5); - - // Common packages - gProof->UploadPackage("STEERBase.par"); - gProof->EnablePackage("STEERBase"); - gProof->UploadPackage("ESD.par"); - gProof->EnablePackage("ESD"); - gProof->UploadPackage("AOD.par"); - gProof->EnablePackage("AOD"); - gProof->UploadPackage("ANALYSIS.par"); - gProof->EnablePackage("ANALYSIS"); - gProof->UploadPackage("ANALYSISalice.par"); - gProof->EnablePackage("ANALYSISalice"); - // Analysis-specific - // --- Enable the PWG3base Package - gProof->UploadPackage("PWG3muon.par"); - gProof->EnablePackage("PWG3muon"); - - // Chain from files staged on CAF - // gROOT->LoadMacro("CreateESDChain.C"); - // TChain* chain = CreateESDChain("ESD1503X_v1.txt",3); - // TChain* chain = CreateESDChain("ESD82XX_30Kshort.txt", 10); - - // Chain from datasets - gROOT->LoadMacro("CreateChainFromDataSet.C"); - ds = gProof->GetDataSet(datasetname)->GetStagedSubset(); - chain = CreateChainFromDataSet(ds, "esdTree"); - - // Make the analysis manager - AliAnalysisManager *mgr = new AliAnalysisManager("Analysis Train", "Analysis train"); - - // ESD input handler - AliESDInputHandler *esdHandler = new AliESDInputHandler(); - esdHandler->SetInactiveBranches("FMD CaloCluster"); - // AOD output handler - AliAODHandler* aodHandler = new AliAODHandler(); - //aodHandler->SetOutputFileName(fileout); - aodHandler->SetOutputFileName("AOD.root"); - - mgr->SetInputEventHandler(esdHandler); - mgr->SetOutputEventHandler(aodHandler); - - // Set of cuts plugged into the ESD filter - // - // standard - AliESDtrackCuts* esdTrackCutsL = new AliESDtrackCuts("AliESDtrackCuts", "Loose"); - esdTrackCutsL->SetMinNClustersTPC(50); - esdTrackCutsL->SetMaxChi2PerClusterTPC(3.5); - esdTrackCutsL->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); - esdTrackCutsL->SetRequireTPCRefit(kTRUE); - esdTrackCutsL->SetMinNsigmaToVertex(3); - esdTrackCutsL->SetRequireSigmaToVertex(kTRUE); - esdTrackCutsL->SetAcceptKingDaughters(kFALSE); - // - // hard - AliESDtrackCuts* esdTrackCutsH = new AliESDtrackCuts("AliESDtrackCuts", "Hard"); - esdTrackCutsH->SetMinNClustersTPC(100); - esdTrackCutsH->SetMaxChi2PerClusterTPC(2.0); - esdTrackCutsH->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); - esdTrackCutsH->SetRequireTPCRefit(kTRUE); - esdTrackCutsH->SetMinNsigmaToVertex(2); - esdTrackCutsH->SetRequireSigmaToVertex(kTRUE); - esdTrackCutsH->SetAcceptKingDaughters(kFALSE); - // - AliAnalysisFilter* trackFilter = new AliAnalysisFilter("trackFilter"); - trackFilter->AddCuts(esdTrackCutsL); - trackFilter->AddCuts(esdTrackCutsH); - - // ESD filter task putting standard info to output generic AOD - AliAnalysisTaskESDfilter *esdfilter = new AliAnalysisTaskESDfilter("ESD Filter"); - //esdfilter->SetTrackFilter(trackFilter); - esdfilter->SetDebugLevel(10); - mgr->AddTask(esdfilter); - - // ESD filter task putting muon info to output generic AOD - AliAnalysisTaskESDMuonFilter *esdmuonfilter = new AliAnalysisTaskESDMuonFilter("ESD Muon Filter"); - mgr->AddTask(esdmuonfilter); - - // Containers for input/output - AliAnalysisDataContainer *cin_esd = mgr->CreateContainer("cESD",TChain::Class(), AliAnalysisManager::kInputContainer); - // Output AOD container. - AliAnalysisDataContainer *cout_aod = mgr->CreateContainer("cAOD", TTree::Class(), AliAnalysisManager::kOutputContainer, "default"); - - // Connect containers to tasks slots - mgr->ConnectInput (esdfilter, 0, cin_esd ); - mgr->ConnectOutput (esdfilter, 0, cout_aod ); - - mgr->ConnectInput (esdmuonfilter, 0, cin_esd); - mgr->ConnectOutput (esdmuonfilter, 0, cout_aod ); - - // - // Run the analysis - // - if (mgr->InitAnalysis()) { - mgr->PrintStatus(); - mgr->StartAnalysis("proof",chain,nev); - } -} - +void AnalysisTrainMuonCAF(char* fileout = "AliAOD.root", char *datasetname = "myDataSet", Int_t nev=1234567890) +{ +// Macro to produce a generic AOD starting from an ESD file. +// The AOD is filled with two tasks: +// 1- with the first one (AliAnalysisTaskESDfilter), +// all the branches of the AOD are filled apart from the muons. +// 2- with the second task (AliAnalysisTaskESDMuonFilter) +// muons tracks are added to the tracks branch +// This macro works on the CAF +// R. Arnaldi 4/5/08 + + gSystem->Load("libTree.so"); + gSystem->Load("libGeom.so"); + gSystem->Load("libVMC.so"); + gSystem->Load("libPhysics.so"); + + // Reset user processes if CAF if not responding anymore + // TProof::Reset("lxb6046"); + + // Connect to proof + TProof::Open("lxb6046"); // may be username@lxb6046 if user not the same as on local + + // Clear packages if changing ROOT version on CAF or local + // gProof->ClearPackages(); + // Enable proof debugging if needed + // gProof->SetLogLevel(5); + + // Common packages + gProof->UploadPackage("STEERBase.par"); + gProof->EnablePackage("STEERBase"); + gProof->UploadPackage("ESD.par"); + gProof->EnablePackage("ESD"); + gProof->UploadPackage("AOD.par"); + gProof->EnablePackage("AOD"); + gProof->UploadPackage("ANALYSIS.par"); + gProof->EnablePackage("ANALYSIS"); + gProof->UploadPackage("ANALYSISalice.par"); + gProof->EnablePackage("ANALYSISalice"); + // Analysis-specific + // --- Enable the PWG3base Package + gProof->UploadPackage("PWG3muon.par"); + gProof->EnablePackage("PWG3muon"); + + // Chain from files staged on CAF + // gROOT->LoadMacro("CreateESDChain.C"); + // TChain* chain = CreateESDChain("ESD1503X_v1.txt",3); + // TChain* chain = CreateESDChain("ESD82XX_30Kshort.txt", 10); + + // Chain from datasets + gROOT->LoadMacro("CreateChainFromDataSet.C"); + ds = gProof->GetDataSet(datasetname)->GetStagedSubset(); + chain = CreateChainFromDataSet(ds, "esdTree"); + + // Make the analysis manager + AliAnalysisManager *mgr = new AliAnalysisManager("Analysis Train", "Analysis train"); + + // ESD input handler + AliESDInputHandler *esdHandler = new AliESDInputHandler(); + esdHandler->SetInactiveBranches("FMD CaloCluster"); + + // AOD output handler + AliAODHandler* aodHandler = new AliAODHandler(); + aodHandler->SetOutputFileName(fileout); + //aodHandler->SetOutputFileName("AOD.root"); + + mgr->SetInputEventHandler(esdHandler); + mgr->SetOutputEventHandler(aodHandler); + + // Set of cuts plugged into the ESD filter + // + // standard + AliESDtrackCuts* esdTrackCutsL = new AliESDtrackCuts("AliESDtrackCuts", "Loose"); + esdTrackCutsL->SetMinNClustersTPC(50); + esdTrackCutsL->SetMaxChi2PerClusterTPC(3.5); + esdTrackCutsL->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); + esdTrackCutsL->SetRequireTPCRefit(kTRUE); + esdTrackCutsL->SetMinNsigmaToVertex(3); + esdTrackCutsL->SetRequireSigmaToVertex(kTRUE); + esdTrackCutsL->SetAcceptKingDaughters(kFALSE); + // + // hard cuts + AliESDtrackCuts* esdTrackCutsH = new AliESDtrackCuts("AliESDtrackCuts", "Hard"); + esdTrackCutsH->SetMinNClustersTPC(100); + esdTrackCutsH->SetMaxChi2PerClusterTPC(2.0); + esdTrackCutsH->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); + esdTrackCutsH->SetRequireTPCRefit(kTRUE); + esdTrackCutsH->SetMinNsigmaToVertex(2); + esdTrackCutsH->SetRequireSigmaToVertex(kTRUE); + esdTrackCutsH->SetAcceptKingDaughters(kFALSE); + esdTrackCutsH->SetPRange(0.,2.); + // + // muon cuts + AliESDMuonTrackCuts* esdMuonTrackCuts = new AliESDMuonTrackCuts("AliESDMuonTrackCuts", "test"); + esdMuonTrackCuts->SetPRange(0.,20.); + //esdMuonTrackCuts->SetPtRange(0.,0.5); // example of kinematic cuts that can be applied + + // track filter (to reject tracks not surviving the cuts - refers to all particles apart from muons) + AliAnalysisFilter* trackFilter = new AliAnalysisFilter("trackFilter"); + trackFilter->AddCuts(esdTrackCutsH); + + // muon track filter (to reject muon tracks not surviving the cuts) + AliAnalysisFilter* trackMuonFilter = new AliAnalysisFilter("trackMuonFilter"); + trackMuonFilter->AddCuts(esdMuonTrackCuts); + + // ESD filter task putting standard info to output generic AOD + AliAnalysisTaskESDfilter *esdfilter = new AliAnalysisTaskESDfilter("ESD Filter"); + //esdfilter->SetTrackFilter(trackFilter); + esdfilter->SetDebugLevel(10); + mgr->AddTask(esdfilter); + + // ESD filter task putting muon info to output generic AOD + AliAnalysisTaskESDMuonFilter *esdmuonfilter = new AliAnalysisTaskESDMuonFilter("ESD Muon Filter"); + esdmuonfilter->SetTrackFilter(trackMuonFilter); + mgr->AddTask(esdmuonfilter); + + // Containers for input/output + AliAnalysisDataContainer *cin_esd = mgr->CreateContainer("cESD",TChain::Class(), AliAnalysisManager::kInputContainer); + // Output AOD container. + AliAnalysisDataContainer *cout_aod = mgr->CreateContainer("cAOD", TTree::Class(), AliAnalysisManager::kOutputContainer, "default"); + + // Connect containers to tasks slots + mgr->ConnectInput (esdfilter, 0, cin_esd ); + mgr->ConnectOutput (esdfilter, 0, cout_aod ); + + mgr->ConnectInput (esdmuonfilter, 0, cin_esd); + mgr->ConnectOutput (esdmuonfilter, 0, cout_aod ); + + // + // Run the analysis + // + if (mgr->InitAnalysis()) { + mgr->PrintStatus(); + mgr->StartAnalysis("proof",chain,nev); + } +} + diff --git a/PWG3/muon/AnalysisTrainMuonLocal.C b/PWG3/muon/AnalysisTrainMuonLocal.C index 72b9f28a04e..402dae3c33b 100644 --- a/PWG3/muon/AnalysisTrainMuonLocal.C +++ b/PWG3/muon/AnalysisTrainMuonLocal.C @@ -1,138 +1,154 @@ -void AnalysisTrainMuonLocal(char* filein = "AliESDs.root", char* fileout = "AliAOD.root") - -// Macro to produce a generic AOD starting from an ESD file. -// The AOD is filled with two tasks: -// 1- with the first one (AliAnalysisTaskESDfilter), -// all the branches of the AOD are filled apart from the muons. -// 2- with the second task (AliAnalysisTaskESDMuonFilter) -// muons tracks are added to the tracks branch -// This macro works locally -// R. Arnaldi 5/5/08 - -{ - gSystem->Load("libTree.so"); - gSystem->Load("libGeom.so"); - gSystem->Load("libVMC.so"); - gSystem->Load("libPhysics.so"); - - // If analysis is .par based: - - // Common packages - SetupPar("STEERBase"); - SetupPar("ESD"); - SetupPar("AOD"); - SetupPar("ANALYSIS"); - SetupPar("ANALYSISalice"); - // Analysis-specific packages - SetupPar("PWG3muon"); - - // Input ESD file - TChain* chain = new TChain("esdTree"); - chain->Add(filein); - - // Define the analysis manager - AliAnalysisManager *mgr = new AliAnalysisManager("Analysis Train", "Analysis train"); - - // ESD input handler - AliESDInputHandler *esdHandler = new AliESDInputHandler(); - esdHandler->SetInactiveBranches("FMD CaloCluster"); - - // AOD output handler - AliAODHandler* aodHandler = new AliAODHandler(); - aodHandler->SetOutputFileName(fileout); - - mgr->SetInputEventHandler(esdHandler); - mgr->SetOutputEventHandler(aodHandler); - - // Set of cuts for the ESD filter - // - // standard cuts - AliESDtrackCuts* esdTrackCutsL = new AliESDtrackCuts("AliESDtrackCuts", "Loose"); - esdTrackCutsL->SetMinNClustersTPC(50); - esdTrackCutsL->SetMaxChi2PerClusterTPC(3.5); - esdTrackCutsL->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); - esdTrackCutsL->SetRequireTPCRefit(kTRUE); - esdTrackCutsL->SetMinNsigmaToVertex(3); - esdTrackCutsL->SetRequireSigmaToVertex(kTRUE); - esdTrackCutsL->SetAcceptKingDaughters(kFALSE); - // - // hard cuts - AliESDtrackCuts* esdTrackCutsH = new AliESDtrackCuts("AliESDtrackCuts", "Hard"); - esdTrackCutsH->SetMinNClustersTPC(100); - esdTrackCutsH->SetMaxChi2PerClusterTPC(2.0); - esdTrackCutsH->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); - esdTrackCutsH->SetRequireTPCRefit(kTRUE); - esdTrackCutsH->SetMinNsigmaToVertex(2); - esdTrackCutsH->SetRequireSigmaToVertex(kTRUE); - esdTrackCutsH->SetAcceptKingDaughters(kFALSE); - // - AliAnalysisFilter* trackFilter = new AliAnalysisFilter("trackFilter"); - trackFilter->AddCuts(esdTrackCutsL); - trackFilter->AddCuts(esdTrackCutsH); - - // ESD filter task putting standard info in the output generic AOD - AliAnalysisTaskESDfilter *esdfilter = new AliAnalysisTaskESDfilter("ESD Filter"); - //esdfilter->SetTrackFilter(trackFilter); - esdfilter->SetDebugLevel(10); - mgr->AddTask(esdfilter); - - // ESD filter task putting muon info in the output generic AOD - AliAnalysisTaskESDMuonFilter *esdmuonfilter = new AliAnalysisTaskESDMuonFilter("ESD Muon Filter"); - mgr->AddTask(esdmuonfilter); - - // Containers for input/output - AliAnalysisDataContainer *cin_esd = mgr->CreateContainer("cESD",TChain::Class(), - AliAnalysisManager::kInputContainer); - // Output AOD container. - AliAnalysisDataContainer *cout_aod = mgr->CreateContainer("cAOD", TTree::Class(), - AliAnalysisManager::kOutputContainer, "default"); - - // Connect containers to tasks slots - mgr->ConnectInput (esdfilter, 0, cin_esd ); - mgr->ConnectOutput (esdfilter, 0, cout_aod ); - - mgr->ConnectInput (esdmuonfilter, 0, cin_esd); - mgr->ConnectOutput (esdmuonfilter, 0, cout_aod ); - - // - // Run the analysis - // - if (mgr->InitAnalysis()) { - mgr->PrintStatus(); - mgr->StartAnalysis("local",chain); - } -} - -//______________________________________________________________________________ -void SetupPar(char* pararchivename) -{ - if (pararchivename) { - char processline[1024]; - sprintf(processline,".! tar xvzf %s.par",pararchivename); - gROOT->ProcessLine(processline); - TString ocwd = gSystem->WorkingDirectory(); - gSystem->ChangeDirectory(pararchivename); - - // check for BUILD.sh and execute - if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) { - printf("*******************************\n"); - printf("*** Building PAR archive ***\n"); - printf("*******************************\n"); - - if (gSystem->Exec("PROOF-INF/BUILD.sh")) { - Error("runProcess","Cannot Build the PAR Archive! - Abort!"); - return -1; - } - } - // check for SETUP.C and execute - if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) { - printf("*******************************\n"); - printf("*** Setup PAR archive ***\n"); - printf("*******************************\n"); - gROOT->Macro("PROOF-INF/SETUP.C"); - } - - gSystem->ChangeDirectory(ocwd.Data()); - printf("Current dir: %s\n", ocwd.Data()); - } -} +void AnalysisTrainMuonLocal(char* filein = "AliESDs.root", char* fileout = "AliAOD.root") + +// Macro to produce a generic AOD starting from an ESD file. +// The AOD is filled with two tasks: +// 1- with the first one (AliAnalysisTaskESDfilter), +// all the branches of the AOD are filled apart from the muons. +// 2- with the second task (AliAnalysisTaskESDMuonFilter) +// muons tracks are added to the tracks branch +// There is the possibility to apply cuts on the muon tracks in order +// to reject muons before filling the AOD +// This macro works locally +// R. Arnaldi 5/5/08 + +{ + gSystem->Load("libTree.so"); + gSystem->Load("libGeom.so"); + gSystem->Load("libVMC.so"); + gSystem->Load("libPhysics.so"); + + // If analysis is .par based: + + // Common packages + SetupPar("STEERBase"); + SetupPar("ESD"); + SetupPar("AOD"); + SetupPar("ANALYSIS"); + SetupPar("ANALYSISalice"); + // Analysis-specific packages + SetupPar("PWG3muon"); + + // Input ESD file + TChain* chain = new TChain("esdTree"); + chain->Add(filein); + + // Define the analysis manager + AliAnalysisManager *mgr = new AliAnalysisManager("Analysis Train", "Analysis train"); + + // ESD input handler + AliESDInputHandler *esdHandler = new AliESDInputHandler(); + esdHandler->SetInactiveBranches("FMD CaloCluster"); + + // AOD output handler + AliAODHandler* aodHandler = new AliAODHandler(); + aodHandler->SetOutputFileName(fileout); + + mgr->SetInputEventHandler(esdHandler); + mgr->SetOutputEventHandler(aodHandler); + + // Set of cuts for the ESD filter + // + // standard cut + AliESDtrackCuts* esdTrackCutsL = new AliESDtrackCuts("AliESDtrackCuts", "Loose"); + esdTrackCutsL->SetMinNClustersTPC(50); + esdTrackCutsL->SetMaxChi2PerClusterTPC(3.5); + esdTrackCutsL->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); + esdTrackCutsL->SetRequireTPCRefit(kTRUE); + esdTrackCutsL->SetMinNsigmaToVertex(3); + esdTrackCutsL->SetRequireSigmaToVertex(kTRUE); + esdTrackCutsL->SetAcceptKingDaughters(kFALSE); + // + // hard cuts + AliESDtrackCuts* esdTrackCutsH = new AliESDtrackCuts("AliESDtrackCuts", "Hard"); + esdTrackCutsH->SetMinNClustersTPC(100); + esdTrackCutsH->SetMaxChi2PerClusterTPC(2.0); + esdTrackCutsH->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); + esdTrackCutsH->SetRequireTPCRefit(kTRUE); + esdTrackCutsH->SetMinNsigmaToVertex(2); + esdTrackCutsH->SetRequireSigmaToVertex(kTRUE); + esdTrackCutsH->SetAcceptKingDaughters(kFALSE); + esdTrackCutsH->SetPRange(0.,2.); + // + // muon cuts + AliESDMuonTrackCuts* esdMuonTrackCuts = new AliESDMuonTrackCuts("AliESDMuonTrackCuts", "test"); + esdMuonTrackCuts->SetPRange(0.,20.); + //esdMuonTrackCuts->SetPtRange(0.,0.5); // examples of kinematic cuts that can be applied + esdMuonTrackCuts->SetHistogramsOn(kTRUE); // methods to draw control histos + esdMuonTrackCuts->DefineHistograms(); + esdMuonTrackCuts->DrawHistograms(); + + // track filter (to reject tracks not surviving the cuts - refers to all particles apart from muons) + AliAnalysisFilter* trackFilter = new AliAnalysisFilter("trackFilter"); + trackFilter->AddCuts(esdTrackCutsH); + + // muon track filter (to reject muon tracks not surviving the cuts) + AliAnalysisFilter* trackMuonFilter = new AliAnalysisFilter("trackMuonFilter"); + trackMuonFilter->AddCuts(esdMuonTrackCuts); + + // ESD filter task putting standard info in the output generic AOD + AliAnalysisTaskESDfilter *esdfilter = new AliAnalysisTaskESDfilter("ESD Filter"); + //esdfilter->SetTrackFilter(trackFilter); + esdfilter->SetDebugLevel(10); + mgr->AddTask(esdfilter); + + // ESD filter task putting muon info in the output generic AOD + AliAnalysisTaskESDMuonFilter *esdmuonfilter = new AliAnalysisTaskESDMuonFilter("ESD Muon Filter"); + esdmuonfilter->SetTrackFilter(trackMuonFilter); + mgr->AddTask(esdmuonfilter); + + // Containers for input/output + AliAnalysisDataContainer *cin_esd = mgr->CreateContainer("cESD",TChain::Class(), + AliAnalysisManager::kInputContainer); + // Output AOD container. + AliAnalysisDataContainer *cout_aod = mgr->CreateContainer("cAOD", TTree::Class(), + AliAnalysisManager::kOutputContainer, "default"); + + // Connect containers to tasks slots + mgr->ConnectInput (esdfilter, 0, cin_esd ); + mgr->ConnectOutput (esdfilter, 0, cout_aod ); + + mgr->ConnectInput (esdmuonfilter, 0, cin_esd); + mgr->ConnectOutput (esdmuonfilter, 0, cout_aod ); + + // + // Run the analysis + // + if (mgr->InitAnalysis()) { + mgr->PrintStatus(); + mgr->StartAnalysis("local",chain); + } +} + +//______________________________________________________________________________ +void SetupPar(char* pararchivename) +{ + if (pararchivename) { + char processline[1024]; + sprintf(processline,".! tar xvzf %s.par",pararchivename); + gROOT->ProcessLine(processline); + TString ocwd = gSystem->WorkingDirectory(); + gSystem->ChangeDirectory(pararchivename); + + // check for BUILD.sh and execute + if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) { + printf("*******************************\n"); + printf("*** Building PAR archive ***\n"); + printf("*******************************\n"); + + if (gSystem->Exec("PROOF-INF/BUILD.sh")) { + Error("runProcess","Cannot Build the PAR Archive! - Abort!"); + return -1; + } + } + // check for SETUP.C and execute + if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) { + printf("*******************************\n"); + printf("*** Setup PAR archive ***\n"); + printf("*******************************\n"); + gROOT->Macro("PROOF-INF/SETUP.C"); + } + + gSystem->ChangeDirectory(ocwd.Data()); + printf("Current dir: %s\n", ocwd.Data()); + } +} diff --git a/PWG3/muon/RunESDMuonFilter.C b/PWG3/muon/RunESDMuonFilter.C deleted file mode 100644 index 4a34ab3bcd8..00000000000 --- a/PWG3/muon/RunESDMuonFilter.C +++ /dev/null @@ -1,89 +0,0 @@ -// Macro to run AliAnalysisTaskESDMuonFilter -// - -void RunESDMuonFilter(char* filein = "AliESDs.root", char* fileout = "AliMuonAOD.root" ){ - - gSystem->Load("libTree.so"); - gSystem->Load("libGeom.so"); - gSystem->Load("libVMC.so"); - gSystem->Load("libPhysics"); - - // for analysis .par file based - - setupPar("STEERBase"); - setupPar("ESD"); - setupPar("AOD"); - setupPar("ANALYSIS"); - setupPar("ANALYSISalice"); - setupPar("PWG3base"); - - // Input ESD file - TChain* chain = new TChain("esdTree"); - chain->Add(filein); - AliESDInputHandler* esdHandler = new AliESDInputHandler(); - - // Make aod output handler - AliAODHandler* aodHandler = new AliAODHandler(); - aodHandler->SetOutputFileName(fileout); - - // Make the analysis manager - AliAnalysisManager *mgr = new AliAnalysisManager("Muon AOD Manager", "Muon AOD Manager"); - mgr->SetInputEventHandler(esdHandler); - mgr->SetOutputEventHandler(aodHandler); - mgr-> SetDebugLevel(10); - - // Task for MUON AOD generation - AliAnalysisTaskESDMuonFilter *esdfilter = new AliAnalysisTaskESDMuonFilter("ESD Muon Filter"); - mgr->AddTask(esdfilter); - - // Create containers for input/output - AliAnalysisDataContainer *cinput1 = mgr->CreateContainer("cESD",TChain::Class(), - AliAnalysisManager::kInputContainer); - AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("cAOD", TTree::Class(), - AliAnalysisManager::kOutputContainer, "default"); - - mgr->ConnectInput (esdfilter, 0, cinput1 ); - mgr->ConnectOutput (esdfilter, 0, coutput1 ); - - // Run the analysis - mgr->InitAnalysis(); - mgr->PrintStatus(); - mgr->StartAnalysis("local",chain); -} - - -Int_t setupPar(const char* pararchivename) { - /////////////////// - // Setup PAR File// - /////////////////// - if (pararchivename) { - char processline[1024]; - sprintf(processline,".! tar xvzf %s.par",pararchivename); - gROOT->ProcessLine(processline); - const char* ocwd = gSystem->WorkingDirectory(); - gSystem->ChangeDirectory(pararchivename); - - // check for BUILD.sh and execute - if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) { - printf("*******************************\n"); - printf("*** Building PAR archive ***\n"); - printf("*******************************\n"); - - if (gSystem->Exec("PROOF-INF/BUILD.sh")) { - Error("runAnalysis","Cannot Build the PAR Archive! - Abort!"); - return -1; - } - } - // check for SETUP.C and execute - if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) { - printf("*******************************\n"); - printf("*** Setup PAR archive ***\n"); - printf("*******************************\n"); - gROOT->Macro("PROOF-INF/SETUP.C"); - } - - gSystem->ChangeDirectory("../"); - } - - return 1; -}