// - libCORRFW.so
// - libPWG3muondep.so
//
-// The macro reads ESDs and store outputs in standard output file (AnalysisResults.root)
-// It also needs to load magnetic field, mapping, geometry (+alignment), and recoonstruction parameters from the OCDB
+// It also needs to load magnetic field, mapping, geometry (+alignment), and reconstruction parameters from the OCDB
+//
+// The task reads ESDs
+// Intermediate results are stored in a file chamberResolution_step<#step>.root
+// Final results are stored in the file results.root
//
// Author: Philippe Pillot - SUBATECH Nantes
//--------------------------------------------------------------------------
-TString aliroot="VO_ALICE@AliRoot::v4-19-15-AN";
-
-enum {kLocal, kInteractif_xml, kInteractif_ESDList, kProof};
-
-// global declaration needed in Proof mode, certainly because of the interpretor
-TStopwatch* localTimer;
-TMultiGraph* mgClusterResXVsStep;
-TMultiGraph* mgClusterResYVsStep;
-Double_t clusterResNB[10];
-Double_t clusterResB[10];
-Double_t clusterResNBErr[10];
-Double_t clusterResBErr[10];
+void LoadAlirootLocally(TString& extraLibs);
-void RunMuonResolution(TString smode = "local", TString inputFileName = "AliESDs.root", Int_t nSteps = 10,
- Bool_t selectPhysics = kFALSE, Bool_t matchTrig = kFALSE, Double_t minMomentum = 0.,
- Bool_t correctForSystematics = kTRUE, Int_t extrapMode = 1)
+//______________________________________________________________________________
+void RunMuonResolution(TString smode = "local", TString inputFileName = "AliESDs.root",
+ TString alirootVersion = "VO_ALICE@AliRoot::v4-20-12-AN", Int_t nSteps = 5,
+ Bool_t selectPhysics = kFALSE, Bool_t selectTrigger = kFALSE, Bool_t matchTrig = kTRUE,
+ Bool_t applyAccCut = kTRUE, Double_t minMomentum = 0., Bool_t correctForSystematics = kTRUE,
+ Int_t extrapMode = 1, Int_t nevents = 1234567890)
{
/// Compute the cluster resolution by studying cluster-track residual, deconvoluting from track resolution
- /// if extrapMode == 0: extrapolate from the closest cluster
- /// if extrapMode == 1: extrapolate from the previous cluster except between stations 2-3-4
- /// if correctForSystematics == kTRUE: the systematic shifts of the residuals is included in the resolution
-
- // timer start...
- localTimer = new TStopwatch;
+ /// - smode = "local" or "proof"
+ /// - inputFileName = an ESD root file or a list of ESDs or a collection of ESDs or a dataset in proof mode
+ /// - alirootVersion = version of aliroot package to enable on AAF (only used in proof mode)
+ /// - nSteps = number of times to task is run (at each step it starts with the chamber resolution obtained in the previous one)
+ /// - selectPhysics : apply or not the physics selection
+ /// - selectTrigger : select only muon trigger events or not (the type of trigger can eventually be changed)
+ /// - matchTrigger : select only the tracks matching the trigger or not
+ /// - applyAccCut : select only the tracks passing the Rabs and the eta cut or not
+ /// - minMomentum : select only the tracks with a total momentum above this value
+ /// - if correctForSystematics == kTRUE: the systematic shifts of the residuals is included in the resolution
+ /// - if extrapMode == 0: extrapolate from the closest cluster
+ /// - if extrapMode == 1: extrapolate from the previous cluster except between stations 2-3-4
+ /// - nevents = maximum number of processed events
+
+ // Load libraries locally
+ TString extraLibs = "RAWDatabase:CDB:STEER:MUONcore:MUONmapping:MUONcalib:MUONgeometry:MUONtrigger:MUONraw:MUONbase:MUONrec:CORRFW:PWG3muondep";
+ LoadAlirootLocally(extraLibs);
+
+ // compile analysis macro locally
+ gROOT->LoadMacro("$ALICE_ROOT/PWG3/muondep/MuonResolution.C++g");
+ MuonResolution(smode, inputFileName, alirootVersion, nSteps, selectPhysics, selectTrigger, matchTrig,
+ applyAccCut, minMomentum, correctForSystematics, extrapMode, nevents, extraLibs);
- // check parameters
- nSteps = TMath::Max(nSteps,1);
- if (extrapMode != 0 && extrapMode != 1) {
- Error("RunMuonResolution","incorrect extrapolation mode!");
- return;
- }
-
- // Check runing mode
- Int_t mode = GetMode(smode, inputFileName);
- if(mode < 0){
- Error("RunMuonResolution","Please provide either an ESD root file a collection of ESDs or a dataset.");
- return;
- }
+}
+
+//______________________________________________________________________________
+void LoadAlirootLocally(TString& extraLibs)
+{
+ /// Load libraries locally
+ // Load common libraries
gSystem->Load("libVMC");
+ gSystem->Load("libTree.so");
+ gSystem->Load("libPhysics.so");
+ gSystem->Load("libMinuit.so");
+ gSystem->Load("libXMLParser.so");
+ gSystem->Load("libGui.so");
gSystem->Load("libSTEERBase");
gSystem->Load("libESD");
gSystem->Load("libAOD");
// Load additional libraries
gSystem->Load("libProofPlayer");
- TString extraLibs="Physics:Minuit:XMLParser:Gui:RAWDatabase:CDB:STEER:MUONcore:MUONmapping:MUONcalib:MUONgeometry:MUONtrigger:MUONraw:MUONbase:MUONrec:CORRFW:PWG3muondep";
TObjArray* libs = extraLibs.Tokenize(":");
for (Int_t i = 0; i < libs->GetEntriesFast(); i++)
gSystem->Load(Form("lib%s",static_cast<TObjString*>(libs->UncheckedAt(i))->GetName()));
+ delete libs;
- // check for old output file to removed
- char remove = '';
- if (!gSystem->Exec("ls chamberResolution_step*[0-9].root")) {
- cout<<"above files must be removed from the current directory. Delete? [y=yes, n=no] "<<flush;
- while (remove != 'y' && remove != 'n') cin>>remove;
- if (remove == 'y') gSystem->Exec("rm -f chamberResolution_step*[0-9].root");
- else {
- Error("RunMuonResolution","cannot proceed with these files there otherwise results will be mixed up!");
- return;
- }
- }
-
- // OCDB access
- if(mode == kLocal) AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
- else if (mode != kProof) {
- if (!TGrid::Connect("alien://")) return;
- if(mode == kInteractif_ESDList || !gSystem->AccessPathName("ConfigureCuts.C"))
- AliCDBManager::Instance()->SetDefaultStorage("raw://");
- }
-
- // Create input object
- TObject* inputObj = 0x0;
- if (mode == kProof) inputObj = new TObjString(inputFileName);
- else inputObj = CreateChain(mode, inputFileName);
- if (!inputObj) return;
-
- // set starting chamber resolution (if -1 they will be loaded from recoParam in the task)
- for (Int_t i = 0; i < AliMUONConstants::NTrackingCh(); i++) {
- clusterResNB[i] = -1.;
- clusterResB[i] = -1.;
- clusterResNBErr[i] = 0.;
- clusterResBErr[i] = 0.;
- }
-
- // output graphs
- mgClusterResXVsStep = new TMultiGraph("mgClusterResXVsStep","cluster X-resolution versus step;step;#sigma_{X} (cm)");
- mgClusterResYVsStep = new TMultiGraph("mgClusterResYVsStep","cluster Y-resolution versus step;step;#sigma_{Y} (cm)");
- TGraphErrors* gClusterResXVsStep[10];
- TGraphErrors* gClusterResYVsStep[10];
- for (Int_t i = 0; i < AliMUONConstants::NTrackingCh(); i++) {
- gClusterResXVsStep[i] = new TGraphErrors(nSteps+1);
- gClusterResXVsStep[i]->SetName(Form("gResX_ch%d",i+1));
- gClusterResXVsStep[i]->SetMarkerStyle(kFullDotMedium);
- gClusterResXVsStep[i]->SetMarkerColor(i+1+i/9);
- mgClusterResXVsStep->Add(gClusterResXVsStep[i],"lp");
-
- gClusterResYVsStep[i] = new TGraphErrors(nSteps+1);
- gClusterResYVsStep[i]->SetName(Form("gResY_ch%d",i+1));
- gClusterResYVsStep[i]->SetMarkerStyle(kFullDotMedium);
- gClusterResYVsStep[i]->SetMarkerColor(i+1+i/9);
- mgClusterResYVsStep->Add(gClusterResYVsStep[i],"lp");
- }
-
- // loop over step
- for (Int_t iStep = 0; iStep < nSteps; iStep++) {
- cout<<"step "<<iStep+1<<"/"<<nSteps<<endl;
-
- // Connect to proof if needed and prepare environment
- if (mode == kProof) {
-
- // set general environment and close previous session
- if (iStep == 0) {
- gEnv->SetValue("XSec.GSI.DelegProxy","2");
- TProof::AddEnvVar("ALIROOT_EXTRA_LIBS", extraLibs.Data());
- } else gProof->Close("s");
-
- // connect
- if (gSystem->Getenv("alien_API_USER") == NULL) TProof::Open("alice-caf.cern.ch","workers=20");
- else TProof::Open(Form("%s@alice-caf.cern.ch",gSystem->Getenv("alien_API_USER")),"workers=20");
- if (!gProof) return;
-
- // set environment and compile task on workers
- gProof->EnablePackage(aliroot.Data());
-
- // prepare OCDB access on workers
- gProof->Exec("AliCDBManager::Instance()->SetDefaultStorage(\"raw://\")");
-
- }
-
- // run one step
- run(mode, iStep, inputObj, extrapMode, correctForSystematics, minMomentum, matchTrig, selectPhysics, clusterResNB, clusterResB);
-
- // fill graph with starting resolutions from the task at first step
- if (iStep == 0) for (Int_t i = 0; i < AliMUONConstants::NTrackingCh(); i++) {
- gClusterResXVsStep[i]->SetPoint(0, 0, clusterResNB[i]);
- gClusterResXVsStep[i]->SetPointError(0, 0., clusterResNBErr[i]);
- gClusterResYVsStep[i]->SetPoint(0, 0, clusterResB[i]);
- gClusterResYVsStep[i]->SetPointError(0, 0., clusterResBErr[i]);
- }
-
- // read the chamber resolution from the output file
- if (!GetChamberResolution(iStep, clusterResNB, clusterResB, clusterResNBErr, clusterResBErr)) return;
-
- // fill graphs with computed resolutions
- for (Int_t i = 0; i < AliMUONConstants::NTrackingCh(); i++) {
- gClusterResXVsStep[i]->SetPoint(iStep+1, iStep+1, clusterResNB[i]);
- gClusterResXVsStep[i]->SetPointError(iStep+1, 0., clusterResNBErr[i]);
- gClusterResYVsStep[i]->SetPoint(iStep+1, iStep+1, clusterResB[i]);
- gClusterResYVsStep[i]->SetPointError(iStep+1, 0., clusterResBErr[i]);
- }
-
- }
-
- // copy final results in results.root file
- gSystem->Exec(Form("cp chamberResolution_step%d.root results.root", nSteps-1));
-
- // display convergence
- TCanvas* convergence = new TCanvas("convergence","convergence");
- convergence->Divide(1,2);
- convergence->cd(1);
- mgClusterResXVsStep->Draw("ap");
- convergence->cd(2);
- mgClusterResYVsStep->Draw("ap");
+ // Load lib for final mchview display
+ gSystem->Load("libMUONgraphics");
- // save convergence plots
- TFile* outFile = TFile::Open("results.root","UPDATE");
- if (!outFile || !outFile->IsOpen()) return;
- outFile->cd();
- mgClusterResXVsStep->Write();
- mgClusterResYVsStep->Write();
- convergence->Write();
- outFile->Close();
+ // Use AliRoot includes and compile our task
+ gROOT->ProcessLine(".include $ALICE_ROOT/include");
+ gROOT->ProcessLine(".include $ALICE_ROOT/MUON");
+ gROOT->ProcessLine(".include $ALICE_ROOT/MUON/mapping");
+ gROOT->ProcessLine(".include $ALICE_ROOT/ANALYSIS/macros");
+ gROOT->ProcessLine(".include $ALICE_ROOT/PWG3/muondep");
- // print results
- printf("\nchamber resolution:\n");
- printf(" - non-bending:");
- for (Int_t i = 0; i < AliMUONConstants::NTrackingCh(); i++) printf((i==0)?" %5.3f":", %5.3f",clusterResNB[i]);
- printf("\n - bending:");
- for (Int_t i = 0; i < AliMUONConstants::NTrackingCh(); i++) printf((i==0)?" %6.4f":", %6.4f",clusterResB[i]);
- printf("\n\n");
-
- // ...timer stop
- localTimer->Stop();
- localTimer->Print();
-
-}
-
-//______________________________________________________________________________
-void run(Int_t mode, Int_t iStep, TObject* input, Int_t extrapMode, Bool_t correctForSystematics,
- Double_t minMomentum, Bool_t matchTrig, Bool_t selectPhysics, Double_t clusterResNB[10], Double_t clusterResB[10])
-{
- /// launch the analysis with these parameters
-
- // Create the analysis manager
- AliAnalysisManager *mgr = new AliAnalysisManager("MuonResolutionAnalysis");
-
- // ESD input handler
- AliESDInputHandler* esdH = new AliESDInputHandler();
- esdH->SetReadFriends(kFALSE);
- mgr->SetInputEventHandler(esdH);
-
- // event selection
- if (selectPhysics) {
- gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C");
- AliPhysicsSelectionTask* physicsSelection = AddTaskPhysicsSelection();
- if(!physicsSelection) {
- Error("RunMuonResolution","AliPhysicsSelectionTask not created!");
- return;
- }
- physicsSelection->GetPhysicsSelection()->SetUseMuonTriggers();
- }
-
- // Muon Resolution analysis
- gROOT->LoadMacro("$ALICE_ROOT/PWG3/muondep/AddTaskMuonResolution.C");
- AliAnalysisManager::SetCommonFileName(Form("chamberResolution_step%d.root", iStep));
- AliAnalysisTaskMuonResolution* muonResolution = AddTaskMuonResolution(selectPhysics, matchTrig, minMomentum, correctForSystematics, extrapMode);
- if(!muonResolution) {
- Error("RunMuonResolution","AliAnalysisTaskMuonResolution not created!");
- return;
- }
- muonResolution->SetStartingResolution(clusterResNB, clusterResB);
-
- // Enable debug printouts
- //mgr->SetDebugLevel(2);
-
- // start local analysis
- if (mgr->InitAnalysis()) {
- mgr->PrintStatus();
- if (mode == kProof) mgr->StartAnalysis("proof", Form("%s",static_cast<TObjString*>(input)->GetName()));
- else mgr->StartAnalysis("local", static_cast<TChain*>(input));
- }
-
- // save the summary canvases
- if (muonResolution->GetCanvases()) {
- TFile* outFile = TFile::Open(Form("chamberResolution_step%d.root", iStep),"UPDATE");
- if (outFile && outFile->IsOpen()) {
- muonResolution->GetCanvases()->Write();
- outFile->Close();
- }
- }
-
- // return starting chamber resolution from the task
- muonResolution->GetStartingResolution(clusterResNB, clusterResB);
-
- // clean memory
- delete mgr;
- TObject::SetObjectStat(kFALSE);
-}
-
-//______________________________________________________________________________
-Bool_t GetChamberResolution(Int_t iStep, Double_t clusterResNB[10], Double_t clusterResB[10], Double_t clusterResNBErr[10], Double_t clusterResBErr[10])
-{
- /// read the chamber resolution from the output file
-
- TFile* outFile = TFile::Open(Form("chamberResolution_step%d.root", iStep),"READ");
-
- if (!outFile || !outFile->IsOpen()) {
- Error("GetChamberResolution","output file does not exist!");
- return kFALSE;
- }
-
- TObjArray* summary = static_cast<TObjArray*>(outFile->FindObjectAny("ChamberRes"));
- TGraphErrors* gCombinedResidualXPerChSigma = (summary) ? static_cast<TGraphErrors*>(summary->FindObject("gCombinedResidualXPerChSigma")) : 0x0;
- TGraphErrors* gCombinedResidualYPerChSigma = (summary) ? static_cast<TGraphErrors*>(summary->FindObject("gCombinedResidualYPerChSigma")) : 0x0;
-
- if (!gCombinedResidualXPerChSigma || !gCombinedResidualYPerChSigma) {
- Error("GetChamberResolution","resolution graphs do not exist!");
- return kFALSE;
- }
-
- Double_t dummy;
- for (Int_t i = 0; i < AliMUONConstants::NTrackingCh(); i++) {
- gCombinedResidualXPerChSigma->GetPoint(i, dummy, clusterResNB[i]);
- gCombinedResidualYPerChSigma->GetPoint(i, dummy, clusterResB[i]);
- clusterResNBErr[i] = gCombinedResidualXPerChSigma->GetErrorY(i);
- clusterResBErr[i] = gCombinedResidualYPerChSigma->GetErrorY(i);
- }
-
- outFile->Close();
-
- return kTRUE;
-}
-
-//______________________________________________________________________________
-Int_t GetMode(TString smode, TString input)
-{
- if (smode == "local") {
- if ( input.EndsWith(".xml") ) return kInteractif_xml;
- else if ( input.EndsWith(".txt") ) return kInteractif_ESDList;
- else if ( input.EndsWith(".root") ) return kLocal;
- } else if (smode == "proof") return kProof;
- return -1;
-}
-
-//______________________________________________________________________________
-TChain* CreateChainFromCollection(const char *xmlfile)
-{
- // Create a chain from the collection of tags.
- TAlienCollection* coll = TAlienCollection::Open(xmlfile);
- if (!coll) {
- Error("CreateChainFromCollection", Form("Cannot create an AliEn collection from %s!", xmlfile));
- return NULL;
- }
-
- TGridResult* tagResult = coll->GetGridResult("",kFALSE,kFALSE);
- AliTagAnalysis *tagAna = new AliTagAnalysis("ESD");
- tagAna->ChainGridTags(tagResult);
-
- AliRunTagCuts *runCuts = new AliRunTagCuts();
- AliLHCTagCuts *lhcCuts = new AliLHCTagCuts();
- AliDetectorTagCuts *detCuts = new AliDetectorTagCuts();
- AliEventTagCuts *evCuts = new AliEventTagCuts();
-
- // Check if the cuts configuration file was provided
- if (!gSystem->AccessPathName("ConfigureCuts.C")) {
- gROOT->LoadMacro("ConfigureCuts.C");
- ConfigureCuts(runCuts, lhcCuts, detCuts, evCuts);
- }
-
- TChain *chain = tagAna->QueryTags(runCuts, lhcCuts, detCuts, evCuts);
- if (!chain || !chain->GetNtrees()) return NULL;
- chain->ls();
- return chain;
-}
-
-//______________________________________________________________________________
-TChain* CreateChainFromFile(const char *rootfile)
-{
- // Create a chain using the root file.
- TChain* chain = new TChain("esdTree");
- chain->Add(rootfile);
- if (!chain->GetNtrees()) return NULL;
- chain->ls();
- return chain;
-}
-
-//______________________________________________________________________________
-TChain* CreateChainFromESDList(const char *esdList)
-{
- // Create a chain using tags from the run list.
- TChain* chain = new TChain("esdTree");
- ifstream inFile(esdList);
- TString inFileName;
- if (inFile.is_open()) {
- while (! inFile.eof() ) {
- inFileName.ReadLine(inFile,kFALSE);
- if(!inFileName.EndsWith(".root")) continue;
- chain->Add(inFileName.Data());
- }
- }
- inFile.close();
- if (!chain->GetNtrees()) return NULL;
- chain->ls();
- return chain;
-}
-
-//______________________________________________________________________________
-TChain* CreateChain(Int_t mode, TString input)
-{
- printf("*******************************\n");
- printf("*** Getting the Chain ***\n");
- printf("*******************************\n");
- if(mode == kInteractif_xml) return CreateChainFromCollection(input.Data());
- else if (mode == kInteractif_ESDList) return CreateChainFromESDList(input.Data());
- else if (mode == kLocal) return CreateChainFromFile(input.Data());
- else return NULL;
}