PWG/muon/AliVAnalysisMuon.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id: AliVAnalysisMuon.cxx 47782 2011-02-24 18:37:31Z martinez $ */
  17
  18 //-----------------------------------------------------------------------------
  19 /// \class AliVAnalysisMuon
  20 /// Base class with utilities for muon analysis
  21 ///
  22 /// \author Diego Stocco
  23 //-----------------------------------------------------------------------------
  24
  25 #include "AliVAnalysisMuon.h"
  26
  27 // ROOT includes
  28 #include "TROOT.h"
  29 #include "TH1.h"
  30 #include "TH2.h"
  31 #include "TAxis.h"
  32 #include "TCanvas.h"
  33 #include "TLegend.h"
  34 #include "TMath.h"
  35 #include "TObjString.h"
  36 #include "TObjArray.h"
  37 #include "THashList.h"
  38 #include "TStyle.h"
  39 //#include "TMCProcess.h"
  40 #include "TLorentzVector.h"
  41
  42 // STEER includes
  43 #include "AliInputEventHandler.h"
  44 #include "AliCentrality.h"
  45
  46 #include "AliAODEvent.h"
  47 #include "AliAODTrack.h"
  48 #include "AliAODMCParticle.h"
  49 #include "AliMCEvent.h"
  50 #include "AliMCParticle.h"
  51 //#include "AliStack.h"
  52 #include "AliESDEvent.h"
  53 #include "AliESDMuonTrack.h"
  54
  55 // ANALYSIS includes
  56 #include "AliAnalysisManager.h"
  57 #include "AliAnalysisTaskSE.h"
  58 #include "AliAnalysisDataSlot.h"
  59 #include "AliAnalysisDataContainer.h"
  60
  61 // CORRFW includes
  62 #include "AliCFGridSparse.h"
  63
  64 // PWG3 includes
  65 #include "AliCounterCollection.h"
  66 #include "AliMergeableCollection.h"
  67
  68 /// \cond CLASSIMP
  69 ClassImp(AliVAnalysisMuon) // Class implementation in ROOT context
  70 /// \endcond
  71
  72
  73 //________________________________________________________________________
  74 AliVAnalysisMuon::AliVAnalysisMuon() :
  75   AliAnalysisTaskSE(),
  76   fMuonTrackCuts(),
  77   fMuonPairCuts(),
  78   fESDEvent(0x0),
  79   fAODEvent(0x0),
  80   fTerminateOptions(0x0),
  81   fChargeKeys(0x0),
  82   fSrcKeys(0x0),
  83   fPhysSelKeys(0x0),
  84   fTriggerClasses(0x0),
  85   fCentralityClasses(0x0),
  86   fEventCounters(0x0),
  87   fMergeableCollection(0x0),
  88   fOutputList(0x0),
  89   fSelectedTrigPattern(0x0),
  90   fRejectedTrigPattern(0x0),
  91   fSelectedTrigLevel(0x0),
  92   fOutputPrototypeList(0x0)
  93 {
  94   /// Default ctor.
  95 }
  96
  97 //________________________________________________________________________
  98 AliVAnalysisMuon::AliVAnalysisMuon(const char *name, const AliMuonTrackCuts& trackCuts, const AliMuonPairCuts& pairCuts) :
  99   AliAnalysisTaskSE(name),
 100   fMuonTrackCuts(trackCuts),
 101   fMuonPairCuts(pairCuts),
 102   fESDEvent(0x0),
 103   fAODEvent(0x0),
 104   fTerminateOptions(0x0),
 105   fChargeKeys(0x0),
 106   fSrcKeys(0x0),
 107   fPhysSelKeys(0x0),
 108   fTriggerClasses(new THashList()),
 109   fCentralityClasses(0x0),
 110   fEventCounters(0x0),
 111   fMergeableCollection(0x0),
 112   fOutputList(0x0),
 113   fSelectedTrigPattern(new TObjArray()),
 114   fRejectedTrigPattern(new TObjArray()),
 115   fSelectedTrigLevel(new TObjArray()),
 116   fOutputPrototypeList(0x0)
 117 {
 118   //
 119   /// Constructor.
 120   //
 121
 122   fTriggerClasses->SetOwner();
 123   InitKeys();
 124   SetTrigClassPatterns();
 125   SetCentralityClasses();
 126
 127   DefineOutput(1, TObjArray::Class());
 128 }
 129
 130 //________________________________________________________________________
 131 AliVAnalysisMuon::AliVAnalysisMuon(const char *name, const AliMuonTrackCuts& trackCuts) :
 132   AliAnalysisTaskSE(name),
 133   fMuonTrackCuts(trackCuts),
 134   fMuonPairCuts(),
 135   fESDEvent(0x0),
 136   fAODEvent(0x0),
 137   fTerminateOptions(0x0),
 138   fChargeKeys(0x0),
 139   fSrcKeys(0x0),
 140   fPhysSelKeys(0x0),
 141   fTriggerClasses(new THashList()),
 142   fCentralityClasses(0x0),
 143   fEventCounters(0x0),
 144   fMergeableCollection(0x0),
 145   fOutputList(0x0),
 146   fSelectedTrigPattern(new TObjArray()),
 147   fRejectedTrigPattern(new TObjArray()),
 148   fSelectedTrigLevel(new TObjArray()),
 149   fOutputPrototypeList(0x0)
 150 {
 151   //
 152   /// Constructor.
 153   //
 154
 155   fTriggerClasses->SetOwner();
 156   InitKeys();
 157   SetTrigClassPatterns();
 158   SetCentralityClasses();
 159
 160   DefineOutput(1, TObjArray::Class());
 161 }
 162
 163
 164 //________________________________________________________________________
 165 AliVAnalysisMuon::AliVAnalysisMuon(const char *name, const AliMuonPairCuts& pairCuts) :
 166   AliAnalysisTaskSE(name),
 167   fMuonTrackCuts(),
 168   fMuonPairCuts(pairCuts),
 169   fESDEvent(0x0),
 170   fAODEvent(0x0),
 171   fTerminateOptions(0x0),
 172   fChargeKeys(0x0),
 173   fSrcKeys(0x0),
 174   fPhysSelKeys(0x0),
 175   fTriggerClasses(new THashList()),
 176   fCentralityClasses(0x0),
 177   fEventCounters(0x0),
 178   fMergeableCollection(0x0),
 179   fOutputList(0x0),
 180   fSelectedTrigPattern(new TObjArray()),
 181   fRejectedTrigPattern(new TObjArray()),
 182   fSelectedTrigLevel(new TObjArray()),
 183   fOutputPrototypeList(0x0)
 184 {
 185   //
 186   /// Constructor.
 187   //
 188   fTriggerClasses->SetOwner();
 189   InitKeys();
 190   SetTrigClassPatterns();
 191   SetCentralityClasses();
 192
 193   DefineOutput(1, TObjArray::Class());
 194 }
 195
 196
 197 //________________________________________________________________________
 198 AliVAnalysisMuon::~AliVAnalysisMuon()
 199 {
 200   //
 201   /// Destructor
 202   //
 203
 204   delete fTerminateOptions;
 205   delete fChargeKeys;
 206   delete fSrcKeys;
 207   delete fPhysSelKeys;
 208   delete fTriggerClasses;
 209   delete fCentralityClasses;
 210   delete fSelectedTrigPattern;
 211   delete fRejectedTrigPattern;
 212   delete fSelectedTrigLevel;
 213   delete fOutputPrototypeList;
 214
 215
 216   // For proof: do not delete output containers
 217   if ( ! AliAnalysisManager::GetAnalysisManager() || ! AliAnalysisManager::GetAnalysisManager()->IsProofMode() ) {
 218     delete fOutputList;
 219   }
 220 }
 221
 222 //___________________________________________________________________________
 223 void AliVAnalysisMuon::FinishTaskOutput()
 224 {
 225   //
 226   /// Remove empty histograms to reduce the number of histos to be merged
 227   //
 228
 229
 230   fMergeableCollection->PruneEmptyObjects();
 231
 232   TString objectName = "";
 233
 234   // Add stat. info from physics selection
 235   // (usefull when running on AODs)
 236   if ( fInputHandler ) {
 237     for ( Int_t istat=0; istat<2; istat++ ) {
 238       TString statType = ( istat == 0 ) ? "ALL" : "BIN0";
 239       TH2* hStat = dynamic_cast<TH2*>(fInputHandler->GetStatistics(statType.Data()));
 240       if ( hStat ) {
 241         objectName = Form("%s_%s", hStat->GetName(), GetName());
 242         TH2* cloneStat = static_cast<TH2*>(hStat->Clone(objectName.Data()));
 243         cloneStat->SetDirectory(0);
 244         fOutputList->Add(cloneStat);
 245       }
 246       else {
 247         AliWarning("Stat histogram not available");
 248         break;
 249       }
 250     } // loop on stat type
 251   }
 252 }
 253
 254
 255 //___________________________________________________________________________
 256 void AliVAnalysisMuon::NotifyRun()
 257 {
 258   /// Set run number for cuts
 259   if ( fMuonTrackCuts.GetFilterMask() ) fMuonTrackCuts.SetRun(fCurrentRunNumber);
 260   if ( fMuonPairCuts.GetFilterMask() ) fMuonPairCuts.SetRun(fCurrentRunNumber);
 261 }
 262
 263 //___________________________________________________________________________
 264 void AliVAnalysisMuon::UserCreateOutputObjects()
 265 {
 266   //
 267   /// Create output objects
 268   //
 269   AliInfo(Form("   CreateOutputObjects of task %s\n", GetName()));
 270
 271   fOutputList = new TObjArray();
 272   fOutputList->SetOwner();
 273
 274   fEventCounters = new AliCounterCollection("eventCounters");
 275
 276   if ( ! fCentralityClasses ) SetCentralityClasses();
 277   TString centralityClasses = "";
 278   for ( Int_t icent=1; icent<=fCentralityClasses->GetNbins(); ++icent ) {
 279     if ( ! centralityClasses.IsNull() ) centralityClasses += "/";
 280     centralityClasses += fCentralityClasses->GetBinLabel(icent);
 281   }
 282   fEventCounters->AddRubric("selected", "yes/no");
 283   fEventCounters->AddRubric("trigger", 100);
 284   fEventCounters->AddRubric("centrality", centralityClasses);
 285   fEventCounters->Init();
 286   fOutputList->Add(fEventCounters);
 287
 288   fMergeableCollection = new AliMergeableCollection("outputObjects");
 289   fOutputList->Add(fMergeableCollection);
 290
 291   PostData(1, fOutputList);
 292
 293   MyUserCreateOutputObjects();
 294 }
 295
 296
 297 //________________________________________________________________________
 298 void AliVAnalysisMuon::UserExec(Option_t * /*option*/)
 299 {
 300   //
 301   /// Main loop
 302   /// Called for each event
 303   //
 304
 305   fAODEvent = dynamic_cast<AliAODEvent*> (InputEvent());
 306   if ( ! fAODEvent )
 307     fESDEvent = dynamic_cast<AliESDEvent*> (InputEvent());
 308
 309   if ( ! fAODEvent && ! fESDEvent ) {
 310     AliError ("AOD or ESD event not found. Nothing done!");
 311     return;
 312   }
 313
 314   Int_t physSel = ( fInputHandler->IsEventSelected() & AliVEvent::kAny ) ? kPhysSelPass : kPhysSelReject;
 315
 316   //
 317   // Global event info
 318   //
 319
 320   TString firedTrigClasses = ( fAODEvent ) ? fAODEvent->GetFiredTriggerClasses() : fESDEvent->GetFiredTriggerClasses();
 321   firedTrigClasses.Prepend("ANY ");
 322   AliDebug(2, Form("Fired classes %s", firedTrigClasses.Data()));
 323   TObjArray* selectTrigClasses = BuildTriggerClasses(firedTrigClasses);
 324   if ( selectTrigClasses->GetEntries() == 0 ) {
 325     delete selectTrigClasses;
 326     return;
 327   }
 328
 329   Double_t centrality = InputEvent()->GetCentrality()->GetCentralityPercentile("V0M");
 330   Int_t centralityBin = fCentralityClasses->FindBin(centrality);
 331   TString centralityBinLabel = fCentralityClasses->GetBinLabel(centralityBin);
 332
 333   TString selKey = ( physSel == kPhysSelPass ) ? "yes" : "no";
 334   for ( Int_t itrig=0; itrig<selectTrigClasses->GetEntries(); ++itrig ) {
 335     TString trigName = selectTrigClasses->At(itrig)->GetName();
 336     fEventCounters->Count(Form("trigger:%s/selected:%s/centrality:%s", trigName.Data(), selKey.Data(), centralityBinLabel.Data()));
 337   }
 338
 339   ProcessEvent(fPhysSelKeys->At(physSel)->GetName(), *selectTrigClasses, fCentralityClasses->GetBinLabel(centralityBin));
 340
 341   delete selectTrigClasses;
 342
 343   // Post final data. It will be written to a file with option "RECREATE"
 344   PostData(1, fOutputList);
 345 }
 346
 347 //________________________________________________________________________
 348 void AliVAnalysisMuon::Terminate(Option_t *)
 349 {
 350   //
 351   /// Draw some histogram at the end.
 352   //
 353
 354   if ( gROOT->IsBatch() ) return;
 355
 356   fOutputList = dynamic_cast<TObjArray*>(GetOutputData(1));
 357   if ( ! fOutputList ) return;
 358   fEventCounters = static_cast<AliCounterCollection*>(fOutputList->FindObject("eventCounters"));
 359   fMergeableCollection = static_cast<AliMergeableCollection*>(fOutputList->FindObject("outputObjects"));
 360
 361   if ( ! fTerminateOptions ) SetTerminateOptions();
 362   if ( ! fMergeableCollection ) return;
 363   AliInfo(Form("Histogram collection size %g MB", fMergeableCollection->EstimateSize()/1024.0/1024.0));
 364   if ( fTerminateOptions->At(3) ) {
 365     TString sopt = fTerminateOptions->At(3)->GetName();
 366     if ( sopt.Contains("verbose") ) fMergeableCollection->Print("*");
 367   }
 368 }
 369
 370
 371
 372 //________________________________________________________________________
 373 Int_t AliVAnalysisMuon::GetNTracks()
 374 {
 375   //
 376   /// Return the number of tracks in event
 377   //
 378   return ( fAODEvent ) ? fAODEvent->GetNTracks() : fESDEvent->GetNumberOfMuonTracks();
 379 }
 380
 381
 382 //________________________________________________________________________
 383 AliVParticle* AliVAnalysisMuon::GetTrack(Int_t itrack)
 384 {
 385   //
 386   /// Get the current track
 387   //
 388   AliVParticle* track = 0x0;
 389   if ( fAODEvent ) track = fAODEvent->GetTrack(itrack);
 390   else track = fESDEvent->GetMuonTrack(itrack);
 391   return track;
 392 }
 393
 394 //________________________________________________________________________
 395 Double_t AliVAnalysisMuon::MuonMass2() const
 396 {
 397   /// A usefull constant
 398   static Double_t m2 = 1.11636129640000012e-02; // using a constant here as the line below is a problem for CINT...
 399   return m2;
 400 }
 401
 402 //________________________________________________________________________
 403 TLorentzVector AliVAnalysisMuon::GetTrackPair(AliVParticle* track1, AliVParticle* track2) const
 404 {
 405   //
 406   /// Get track pair
 407   //
 408
 409   AliVParticle* tracks[2] = {track1, track2};
 410
 411   TLorentzVector vec[2];
 412   for ( Int_t itrack=0; itrack<2; ++itrack ) {
 413     Double_t trackP = tracks[itrack]->P();
 414     Double_t energy = TMath::Sqrt(trackP*trackP + MuonMass2());
 415     vec[itrack].SetPxPyPzE(tracks[itrack]->Px(), tracks[itrack]->Py(), tracks[itrack]->Pz(), energy);
 416   }
 417
 418   TLorentzVector vecPair = vec[0] + vec[1];
 419   return vecPair;
 420 }
 421
 422
 423 //________________________________________________________________________
 424 Int_t AliVAnalysisMuon::GetNMCTracks()
 425 {
 426   //
 427   /// Return the number of MC tracks in event
 428   //
 429   Int_t nMCtracks = 0;
 430   if ( fMCEvent ) nMCtracks = fMCEvent->GetNumberOfTracks();
 431   else if ( fAODEvent ) {
 432     TClonesArray* mcArray = (TClonesArray*)fAODEvent->GetList()->FindObject(AliAODMCParticle::StdBranchName());
 433     if ( mcArray ) nMCtracks = mcArray->GetEntries();
 434   }
 435   return nMCtracks;
 436 }
 437
 438 //________________________________________________________________________
 439 AliVParticle* AliVAnalysisMuon::GetMCTrack(Int_t trackLabel)
 440 {
 441   //
 442   /// MC information can be provided by the MC input handler
 443   /// (mostly when analyising ESDs) or can be found inside AODs
 444   /// This method returns the correct one
 445   //
 446   AliVParticle* mcTrack = 0x0;
 447   if ( fMCEvent ) mcTrack = fMCEvent->GetTrack(trackLabel);
 448   else if ( fAODEvent ) {
 449     TClonesArray* mcArray = (TClonesArray*)fAODEvent->FindListObject(AliAODMCParticle::StdBranchName());
 450     if ( mcArray ) mcTrack =  (AliVParticle*)mcArray->At(trackLabel);
 451   }
 452   if ( ! mcTrack ) AliWarning(Form("No track with label %i!", trackLabel));
 453   return mcTrack;
 454 }
 455
 456 //________________________________________________________________________
 457 Bool_t AliVAnalysisMuon::IsMC()
 458 {
 459   //
 460   /// Contains MC info
 461   //
 462   return ( fMCEvent || ( fAODEvent && fAODEvent->FindListObject(AliAODMCParticle::StdBranchName()) ) );
 463 }
 464
 465
 466 //________________________________________________________________________
 467 Int_t AliVAnalysisMuon::GetParticleType(AliVParticle* track)
 468 {
 469   //
 470   /// Get particle type from mathced MC track
 471   //
 472
 473   Int_t trackSrc = kUnidentified;
 474   Int_t trackLabel = track->GetLabel();
 475   if ( trackLabel >= 0 ) {
 476     AliVParticle* matchedMCTrack = GetMCTrack(trackLabel);
 477     if ( matchedMCTrack ) trackSrc = RecoTrackMother(matchedMCTrack);
 478   } // track has MC label
 479   return trackSrc;
 480 }
 481
 482
 483 //________________________________________________________________________
 484 Int_t AliVAnalysisMuon::RecoTrackMother(AliVParticle* mcParticle)
 485 {
 486   //
 487   /// Find track mother from kinematics
 488   //
 489
 490   Int_t recoPdg = mcParticle->PdgCode();
 491
 492   // Track is not a muon
 493   if ( TMath::Abs(recoPdg) != 13 ) return kRecoHadron;
 494
 495   Int_t imother = ( fMCEvent ) ? ((AliMCParticle*)mcParticle)->GetMother() : ((AliAODMCParticle*)mcParticle)->GetMother();
 496
 497   Int_t den[3] = {100, 1000, 1};
 498
 499   Int_t motherType = kDecayMu;
 500   while ( imother >= 0 ) {
 501     AliVParticle* part = GetMCTrack(imother);
 502     //if ( ! part ) return motherType;
 503
 504     Int_t absPdg = TMath::Abs(part->PdgCode());
 505
 506     Bool_t isPrimary = ( fMCEvent ) ? ( imother < fMCEvent->GetNumberOfPrimaries() ) : ((AliAODMCParticle*)part)->IsPrimary();
 507
 508     if ( isPrimary ) {
 509       for ( Int_t idec=0; idec<3; idec++ ) {
 510         Int_t flv = (absPdg%100000)/den[idec];
 511         if ( flv > 0 && flv < 4 ) return kDecayMu;
 512         else if ( flv == 0 || flv > 5 ) continue;
 513         else {
 514           if ( den[idec] == 100 ) motherType = kQuarkoniumMu;
 515           else if ( flv == 4 ) motherType = kCharmMu;
 516           else motherType = kBeautyMu;
 517           break; // break loop on pdg code
 518           // but continue the global loop to find higher mass HF
 519         }
 520       } // loop on pdg code
 521       if ( absPdg < 10 ) break; // particle loop
 522     } // is primary
 523     else {
 524       if ( part->Zv() < -90. ) {
 525         // If hadronic process => secondary
 526         //if ( part->GetUniqueID() == kPHadronic ) {
 527         return kSecondaryMu;
 528         //}
 529       }
 530     } // is secondary
 531
 532     imother = ( fMCEvent ) ? ((AliMCParticle*)part)->GetMother() : ((AliAODMCParticle*)part)->GetMother();
 533
 534   } // loop on mothers
 535
 536   return motherType;
 537 }
 538
 539 //________________________________________________________________________
 540 Bool_t AliVAnalysisMuon::AddObjectToCollection(TObject* object, Int_t index)
 541 {
 542   //
 543   /// Add object to collection
 544   //
 545
 546   if ( ! fOutputPrototypeList ) {
 547     fOutputPrototypeList = new TObjArray();
 548     fOutputPrototypeList->SetOwner();
 549   }
 550   if ( fOutputPrototypeList->FindObject(object->GetName() ) ) {
 551     AliWarning(Form("Object with name %s already in the list", object->GetName()));
 552     return kFALSE;
 553   }
 554   if ( index < 0 ) fOutputPrototypeList->Add(object);
 555   else fOutputPrototypeList->AddAtAndExpand(object, index);
 556
 557   return kTRUE;
 558 }
 559
 560 //________________________________________________________________________
 561 TObject* AliVAnalysisMuon::GetMergeableObject(TString physSel, TString trigClassName, TString centrality, TString objectName)
 562 {
 563   //
 564   /// Get mergeable object
 565   /// (create collection if necessary)
 566   //
 567
 568   TString identifier = Form("/%s/%s/%s/", physSel.Data(), trigClassName.Data(), centrality.Data());
 569
 570   TObject* obj = fMergeableCollection->GetObject(identifier.Data(), objectName.Data());
 571   if ( ! obj ) {
 572     CreateMergeableObjects(physSel, trigClassName, centrality);
 573     obj = fMergeableCollection->GetObject(identifier.Data(), objectName.Data());
 574     AliInfo(Form("Mergeable object collection size %g MB", fMergeableCollection->EstimateSize()/1024.0/1024.0));
 575   }
 576   return obj;
 577 }
 578
 579 //________________________________________________________________________
 580 TObject* AliVAnalysisMuon::GetSum(TString physSel, TString trigClassNames, TString centrality, TString objectPattern)
 581 {
 582   //
 583   /// Sum objects
 584   /// - physSel, trigClassNames must be in the form: key1,key2
 585   /// - centrality must be in the form minValue_maxValue
 586   /// - objectPattern must be in the form match1&match2&match3,match4
 587   ///   meaning that the object name must contain match1 and match2 and either one of match3 and match4
 588
 589   if ( ! fMergeableCollection ) return 0x0;
 590
 591   // Get centrality range
 592   Int_t firstCentrality = 1;
 593   Int_t lastCentrality = fCentralityClasses->GetNbins();
 594
 595   TObjArray* centralityRange = centrality.Tokenize("_");
 596   Float_t range[2] = {0., 100.};
 597   if ( centralityRange->GetEntries() >= 2 ) {
 598     for ( Int_t irange=0; irange<2; ++irange ) {
 599       range[irange] = ((TObjString*)centralityRange->At(irange))->GetString().Atof();
 600     }
 601     firstCentrality = fCentralityClasses->FindBin(range[0]+0.0001);
 602     lastCentrality = fCentralityClasses->FindBin(range[1]-0.0001);
 603   }
 604   delete centralityRange;
 605
 606   TString sumCentralityString = "";
 607   for ( Int_t icent=firstCentrality; icent<=lastCentrality; ++icent ) {
 608     if ( ! sumCentralityString.IsNull() ) sumCentralityString += ",";
 609     sumCentralityString += fCentralityClasses->GetBinLabel(icent);
 610   }
 611
 612   objectPattern.ReplaceAll(" ","");
 613   TObjArray* objPatternList = objectPattern.Tokenize("&");
 614   TObjArray objPatternMatrix(objPatternList->GetEntries());
 615   objPatternMatrix.SetOwner();
 616   for ( Int_t ikey=0; ikey<objPatternList->GetEntries(); ikey++ ) {
 617     TObjArray* subKeyList = ((TObjString*)objPatternList->At(ikey))->GetString().Tokenize(",");
 618     objPatternMatrix.AddAt(subKeyList, ikey);
 619   }
 620   delete objPatternList;
 621
 622
 623   TObjArray objectNameInCollection;
 624   objectNameInCollection.SetOwner();
 625   TObjArray* physSelList = physSel.Tokenize(",");
 626   TObjArray* trigClassList = trigClassNames.Tokenize(",");
 627   TObjArray* centralityList = sumCentralityString.Tokenize(",");
 628   for ( Int_t isel=0; isel<physSelList->GetEntries(); isel++ ) {
 629     for ( Int_t itrig = 0; itrig<trigClassList->GetEntries(); itrig++ ) {
 630       for ( Int_t icent=0; icent<centralityList->GetEntries(); icent++ ) {
 631         TString currId = Form("/%s/%s/%s/", physSelList->At(isel)->GetName(), trigClassList->At(itrig)->GetName(),centralityList->At(icent)->GetName());
 632         TList* objNameList = fMergeableCollection->CreateListOfObjectNames(currId.Data());
 633         for ( Int_t iobj=0; iobj<objNameList->GetEntries(); iobj++ ) {
 634           TString objName = objNameList->At(iobj)->GetName();
 635           if ( ! objectNameInCollection.FindObject(objName.Data()) ) objectNameInCollection.Add(new TObjString(objName.Data()));
 636         }
 637         delete objNameList;
 638       }
 639     }
 640   }
 641   delete physSelList;
 642   delete trigClassList;
 643   delete centralityList;
 644
 645   TString matchingObjectNames = "";
 646   for ( Int_t iobj=0; iobj<objectNameInCollection.GetEntries(); iobj++ ) {
 647     TString objName = objectNameInCollection.At(iobj)->GetName();
 648     Bool_t matchAnd = kTRUE;
 649     for ( Int_t ikey=0; ikey<objPatternMatrix.GetEntries(); ikey++ ) {
 650       TObjArray*  subKeyList = (TObjArray*)objPatternMatrix.At(ikey);
 651       Bool_t matchOr = kFALSE;
 652       for ( Int_t isub=0; isub<subKeyList->GetEntries(); isub++ ) {
 653         TString subKeyString = ((TObjString*)subKeyList->At(isub))->GetString();
 654         if ( objName.Contains(subKeyString.Data()) ) {
 655           matchOr = kTRUE;
 656           break;
 657         }
 658       }
 659       if ( ! matchOr ) {
 660         matchAnd = kFALSE;
 661         break;
 662       }
 663     }
 664     if ( ! matchAnd ) continue;
 665     if ( ! matchingObjectNames.IsNull() ) matchingObjectNames.Append(",");
 666     matchingObjectNames += objName;
 667   }
 668
 669   TString idPattern = Form("/%s/%s/%s/%s", physSel.Data(), trigClassNames.Data(), sumCentralityString.Data(), matchingObjectNames.Data());
 670   idPattern.ReplaceAll(" ","");
 671
 672   AliDebug(1,Form("Sum pattern %s", idPattern.Data()));
 673
 674   return fMergeableCollection->GetSum(idPattern.Data());
 675 }
 676
 677 //___________________________________________________________________________
 678 void AliVAnalysisMuon::CreateMergeableObjects(TString physSel, TString trigClassName, TString centrality)
 679 {
 680   TObject* obj = 0x0;
 681   TString objectName = "";
 682   TString identifier = Form("/%s/%s/%s/", physSel.Data(), trigClassName.Data(), centrality.Data());
 683   for ( Int_t iobj=0; iobj<fOutputPrototypeList->GetEntries(); ++iobj ) {
 684     objectName = fOutputPrototypeList->At(iobj)->GetName();
 685     obj = fOutputPrototypeList->At(iobj)->Clone(objectName.Data());
 686     fMergeableCollection->Adopt(identifier, obj);
 687   } // loop on histos
 688 }
 689
 690
 691 //_______________________________________________________________________
 692 Bool_t AliVAnalysisMuon::SetSparseRange(AliCFGridSparse* gridSparse,
 693                                         Int_t ivar, TString labelName,
 694                                         Double_t varMin, Double_t varMax,
 695                                         TString option)
 696 {
 697   //
 698   /// Set range in a smart way.
 699   /// Allows to select a bin from the label.
 700   /// Check the bin limits.
 701   //
 702
 703   option.ToUpper();
 704   Int_t minVarBin = -1, maxVarBin = -1;
 705   TAxis* axis = gridSparse->GetAxis(ivar);
 706
 707   if ( ! axis ) {
 708     printf("Warning: Axis %i not found in %s", ivar, gridSparse->GetName());
 709     return kFALSE;
 710   }
 711
 712   if ( ! labelName.IsNull() ) {
 713     minVarBin = axis->FindBin(labelName.Data());
 714     maxVarBin = minVarBin;
 715     if ( minVarBin < 1 ) {
 716       printf("Warning: %s: label %s not found. Nothing done", gridSparse->GetName(), labelName.Data());
 717       return kFALSE;
 718     }
 719   }
 720   else if ( option.Contains( "USEBIN" ) ) {
 721     minVarBin = (Int_t)varMin;
 722     maxVarBin = (Int_t)varMax;
 723   }
 724   else {
 725     minVarBin = axis->FindBin(varMin);
 726     maxVarBin = axis->FindBin(varMax);
 727   }
 728
 729   if ( axis->GetFirst() == minVarBin && axis->GetLast() == maxVarBin ) return kFALSE;
 730
 731   gridSparse->SetRangeUser(ivar, axis->GetBinCenter(minVarBin), axis->GetBinCenter(maxVarBin));
 732
 733   return kTRUE;
 734 }
 735
 736 //________________________________________________________________________
 737 void AliVAnalysisMuon::SetTrigClassPatterns(TString pattern)
 738 {
 739   /// Set trigger classes
 740   ///
 741   /// Classes are filled dynamically according to the pattern
 742   /// - if name contains ! (without spaces): reject it
 743   /// - in the matching pattern it is also possible to specify the
 744   ///   pt cut level associated to the trigger
 745   /// example:
 746   /// SetTrigClassPatterns("CMBAC CPBI1MSL:Lpt CPBI1MSH:Hpt !ALLNOTRD")
 747   /// keeps classes containing CMBAC, CPBI1MSL and CPBI1MSH and not containing ALLNOTRD.
 748   /// In addition, it knows that the class matching CPBI1MSL requires an Lpt trigger
 749   /// and the one with CPBI1MSH requires a Hpt trigger.
 750   /// Hence, in the analysis, the function
 751   /// TrackPtCutMatchTrigClass(track, "CPBIMSL") returns true if track match Lpt
 752   /// TrackPtCutMatchTrigClass(track, "CPBIMSL") returns true if track match Hpt
 753   /// TrackPtCutMatchTrigClass(track, "CMBAC") always returns true
 754   ///
 755   /// CAVEAT: if you use an fCFContainer and you want an axis to contain the trigger classes,
 756   ///         please be sure that each pattern matches only 1 trigger class, or triggers will be mixed up
 757   ///         when merging different chuncks.
 758
 759   fSelectedTrigPattern->SetOwner();
 760   if ( fSelectedTrigPattern->GetEntries() > 0 ) fSelectedTrigPattern->Delete();
 761   fRejectedTrigPattern->SetOwner();
 762   if ( fRejectedTrigPattern->GetEntries() > 0 ) fRejectedTrigPattern->Delete();
 763   fSelectedTrigLevel->SetOwner();
 764   if ( fSelectedTrigLevel->GetEntries() > 0 ) fSelectedTrigLevel->Delete();
 765
 766   pattern.ReplaceAll("  "," ");
 767   pattern.ReplaceAll("! ","!");
 768   pattern.ReplaceAll(" :",":");
 769
 770   TObjArray* fullList = pattern.Tokenize(" ");
 771
 772   for ( Int_t ipat=0; ipat<fullList->GetEntries(); ++ipat ) {
 773     TString currPattern = fullList->At(ipat)->GetName();
 774     if ( currPattern.Contains("!") ) {
 775       currPattern.ReplaceAll("!","");
 776       fRejectedTrigPattern->AddLast(new TObjString(currPattern));
 777     }
 778     else {
 779       TObjArray* arr = currPattern.Tokenize(":");
 780       fSelectedTrigPattern->AddLast(new TObjString(arr->At(0)->GetName()));
 781       TString selTrigLevel = ( arr->At(1) ) ? arr->At(1)->GetName() : "none";
 782       selTrigLevel.ToUpper();
 783       fSelectedTrigLevel->AddLast(new TObjString(selTrigLevel));
 784       delete arr;
 785     }
 786   }
 787
 788   delete fullList;
 789 }
 790
 791 //________________________________________________________________________
 792 void AliVAnalysisMuon::SetCentralityClasses(Int_t nCentralityBins, Double_t* centralityBins)
 793 {
 794   //
 795   /// Set centrality classes
 796   //
 797   if ( ! centralityBins ) {
 798     Double_t defaultCentralityBins[] = {-5., 0., 5., 10., 20., 30., 40., 50., 60., 70., 80., 100., 105.};
 799     centralityBins = defaultCentralityBins;
 800     nCentralityBins = sizeof(defaultCentralityBins)/sizeof(defaultCentralityBins[0])-1;
 801   }
 802
 803   delete fCentralityClasses;
 804   fCentralityClasses = new TAxis(nCentralityBins, centralityBins);
 805   TString currClass = "";
 806   for ( Int_t ibin=1; ibin<=fCentralityClasses->GetNbins(); ++ibin ){
 807     currClass = Form("%.0f_%.0f",fCentralityClasses->GetBinLowEdge(ibin),fCentralityClasses->GetBinUpEdge(ibin));
 808     fCentralityClasses->SetBinLabel(ibin, currClass.Data());
 809   }
 810 }
 811
 812 //________________________________________________________________________
 813 void AliVAnalysisMuon::SetTerminateOptions(TString physSel, TString trigClass, TString centralityRange, TString furtherOpts)
 814 {
 815   //
 816   /// Set terminate options
 817   //
 818   if ( ! fTerminateOptions ) {
 819     fTerminateOptions = new TObjArray(4);
 820     fTerminateOptions->SetOwner();
 821   }
 822   fTerminateOptions->AddAt(new TObjString(physSel), 0);
 823   fTerminateOptions->AddAt(new TObjString(trigClass), 1);
 824   fTerminateOptions->AddAt(new TObjString(centralityRange),2);
 825   fTerminateOptions->AddLast(new TObjString(furtherOpts));
 826 }
 827
 828 //________________________________________________________________________
 829 void AliVAnalysisMuon::InitKeys()
 830 {
 831   TString chargeKeys = "MuMinus MuPlus";
 832   fChargeKeys = chargeKeys.Tokenize(" ");
 833
 834   TString srcKeys = "CharmMu BeautyMu QuarkoniumMu DecayMu SecondaryMu Hadron Unidentified";
 835   fSrcKeys = srcKeys.Tokenize(" ");
 836
 837   TString physSelKeys = "PhysSelPass PhysSelReject";
 838   fPhysSelKeys = physSelKeys.Tokenize(" ");
 839 }
 840
 841 //________________________________________________________________________
 842 TObjArray* AliVAnalysisMuon::BuildTriggerClasses(TString firedTrigClasses)
 843 {
 844   //
 845   /// Return the list of trigger classes to be considered
 846   /// for current event. Update the global list if necessary
 847   //
 848
 849   TObjArray* selectedTrigClasses = new TObjArray(0);
 850   selectedTrigClasses->SetOwner();
 851
 852   TObjArray* firedTrigClassesList = firedTrigClasses.Tokenize(" ");
 853
 854   for ( Int_t itrig=0; itrig<firedTrigClassesList->GetEntries(); ++itrig ) {
 855     TString trigName = ((TObjString*)firedTrigClassesList->At(itrig))->GetString();
 856     Bool_t rejectTrig = kFALSE;
 857     for ( Int_t ipat=0; ipat<fRejectedTrigPattern->GetEntries(); ++ipat ) {
 858       if ( trigName.Contains(fRejectedTrigPattern->At(ipat)->GetName() ) ) {
 859         rejectTrig = kTRUE;
 860         break;
 861       }
 862     } // loop on reject pattern
 863     if ( rejectTrig ) continue;
 864
 865     Int_t matchPatternIndex = -1;
 866     for ( Int_t ipat=0; ipat<fSelectedTrigPattern->GetEntries(); ++ipat ) {
 867       if ( trigName.Contains(fSelectedTrigPattern->At(ipat)->GetName() ) ) {
 868         matchPatternIndex = ipat;
 869         break;
 870       }
 871     } // loop on keep pattern
 872     if ( matchPatternIndex < 0 ) continue;
 873
 874     selectedTrigClasses->AddLast(new TObjString(trigName));
 875     if ( fTriggerClasses->FindObject(trigName.Data()) ) continue;
 876     Int_t trigLevel = 0;
 877     TString trigLevelString = fSelectedTrigLevel->At(matchPatternIndex)->GetName();
 878     if ( trigLevelString.Contains("APT") ) trigLevel = 1;
 879     else if ( trigLevelString.Contains("LPT") ) trigLevel = 2;
 880     else if ( trigLevelString.Contains("HPT") ) trigLevel = 3;
 881     AliInfo(Form("Adding %s to considered trigger classes",trigName.Data()));
 882     TObjString* addTrig = new TObjString(trigName);
 883     UInt_t uniqueId = trigLevel;
 884     addTrig->SetUniqueID(uniqueId);
 885     fTriggerClasses->Add(addTrig);
 886   } // loop on trigger classes
 887
 888   delete firedTrigClassesList;
 889
 890   return selectedTrigClasses;
 891 }
 892
 893
 894 //________________________________________________________________________
 895 Bool_t AliVAnalysisMuon::TrackPtCutMatchTrigClass(AliVParticle* track, TString trigClassName)
 896 {
 897   /// Check if track passes the trigger pt cut level used in the trigger class
 898   Int_t matchTrig = ( fAODEvent ) ? ((AliAODTrack*)track)->GetMatchTrigger() : ((AliESDMuonTrack*)track)->GetMatchTrigger();
 899   Int_t classMatchLevel = GetTrigClassPtCutLevel(trigClassName);
 900   return matchTrig >= classMatchLevel;
 901 }
 902
 903
 904 //________________________________________________________________________
 905 Int_t AliVAnalysisMuon::GetTrigClassPtCutLevel(TString trigClassName)
 906 {
 907   /// Get trigger class pt cut level for tracking/trigger matching
 908   TObject* obj = fTriggerClasses->FindObject(trigClassName.Data());
 909   if ( ! obj ) {
 910     AliWarning(Form("Class %s not in the list!", trigClassName.Data()));
 911     return -1;
 912   }
 913
 914   return obj->GetUniqueID();
 915 }