[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliNormalizationCounter.cxx

/**************************************************************************
 * Copyright(c) 1998-2008, 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.                  *
 **************************************************************************/

/* $Id$ */

//*************************************************************************
// Class AliNormalizationCounter
// Class to store the informations relevant for the normalization in the 
// barrel for each run
// Authors: G. Ortona, ortona@to.infn.it
// D. Caffarri, davide.caffarri@pd.to.infn.it
// with many thanks to P. Pillot
/////////////////////////////////////////////////////////////

#include "AliLog.h"
#include "AliNormalizationCounter.h"
#include <AliESDEvent.h>
#include <AliESDtrack.h>
#include <AliAODEvent.h>
#include <AliAODVZERO.h>
#include <AliVParticle.h>
#include <AliTriggerAnalysis.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TList.h>
#include <TString.h>
#include <TCanvas.h>
#include <AliPhysicsSelection.h>
#include <AliMultiplicity.h>

ClassImp(AliNormalizationCounter)

//____________________________________________
AliNormalizationCounter::AliNormalizationCounter(): 
TNamed(),
fCounters(),
fESD(kFALSE),
fMultiplicity(kFALSE),
fMultiplicityEtaRange(1.0),
fHistTrackFilterEvMult(0),
fHistTrackAnaEvMult(0),
fHistTrackFilterSpdMult(0),
fHistTrackAnaSpdMult(0)
{
  // empty constructor
}

//__________________________________________________				
AliNormalizationCounter::AliNormalizationCounter(const char *name): 
TNamed(name,name),
fCounters(name),
fESD(kFALSE),
fMultiplicity(kFALSE),
fMultiplicityEtaRange(1.0),
fHistTrackFilterEvMult(0),
fHistTrackAnaEvMult(0),
fHistTrackFilterSpdMult(0),
fHistTrackAnaSpdMult(0)
{
  ;
}

//______________________________________________
AliNormalizationCounter::~AliNormalizationCounter()
{
  //destructor
  if(fHistTrackFilterEvMult){
    delete fHistTrackFilterEvMult;
    fHistTrackFilterEvMult =0;
  }
  if(fHistTrackAnaEvMult){
    delete fHistTrackAnaEvMult;
    fHistTrackAnaEvMult=0;
  }
  if(fHistTrackFilterSpdMult){
    delete fHistTrackFilterSpdMult;
    fHistTrackFilterSpdMult=0;
  }
  if(fHistTrackAnaSpdMult){
    delete fHistTrackAnaSpdMult;
    fHistTrackAnaSpdMult=0;
  }
}

//______________________________________________
void AliNormalizationCounter::Init()
{
  //variables initialization
  fCounters.AddRubric("Event","triggered/V0AND/PileUp/PbPbC0SMH-B-NOPF-ALLNOTRD/Candles0.3/PrimaryV/countForNorm/noPrimaryV/zvtxGT10/!V0A&Candle03/!V0A&PrimaryV/Candid(Filter)/Candid(Analysis)/NCandid(Filter)/NCandid(Analysis)");
  if(fMultiplicity)  fCounters.AddRubric("Multiplicity", 5000);
  fCounters.AddRubric("Run", 1000000);
  fCounters.Init();
  fHistTrackFilterEvMult=new TH2F("FiltCandidvsTracksinEv","FiltCandidvsTracksinEv",10000,-0.5,9999.5,200,-0.5,199.5);
  fHistTrackFilterEvMult->GetYaxis()->SetTitle("NCandidates");
  fHistTrackFilterEvMult->GetXaxis()->SetTitle("NTracksinEvent");
  fHistTrackAnaEvMult=new TH2F("AnaCandidvsTracksinEv","AnaCandidvsTracksinEv",10000,-0.5,9999.5,100,-0.5,99.5);
  fHistTrackAnaEvMult->GetYaxis()->SetTitle("NCandidates");
  fHistTrackAnaEvMult->GetXaxis()->SetTitle("NTracksinEvent");
  fHistTrackFilterSpdMult=new TH2F("FilterCandidvsSpdMult","FilterCandidvsSpdMult",5000,-0.5,4999.5,200,-0.5,199.5);
  fHistTrackFilterSpdMult->GetYaxis()->SetTitle("NCandidates");
  fHistTrackFilterSpdMult->GetXaxis()->SetTitle("NSPDTracklets");
  fHistTrackAnaSpdMult=new TH2F("AnaCandidvsSpdMult","AnaCandidvsSpdMult",5000,-0.5,4999.5,100,-0.5,99.5);
  fHistTrackAnaSpdMult->GetYaxis()->SetTitle("NCandidates");
  fHistTrackAnaSpdMult->GetXaxis()->SetTitle("NSPDTracklets");
}

//______________________________________________
Long64_t AliNormalizationCounter::Merge(TCollection* list){
  if (!list) return 0;
  if (list->IsEmpty()) return 0;//(Long64_t)fCounters.Merge(list);

  TIter next(list);
  const TObject* obj = 0x0;
  while ((obj = next())) {
    
    // check that "obj" is an object of the class AliNormalizationCounter
    const AliNormalizationCounter* counter = dynamic_cast<const AliNormalizationCounter*>(obj);
    if (!counter) {
      AliError(Form("object named %s is not AliNormalizationCounter! Skipping it.", counter->GetName()));
      continue;
    }

    Add(counter);

  }
  
  return (Long64_t)1;//(Long64_t)fCounters->GetEntries();
}
//_______________________________________
void AliNormalizationCounter::Add(const AliNormalizationCounter *norm){
  fCounters.Add(&(norm->fCounters));
  fHistTrackFilterEvMult->Add(norm->fHistTrackFilterEvMult);
  fHistTrackAnaEvMult->Add(norm->fHistTrackAnaEvMult);
  fHistTrackFilterSpdMult->Add(norm->fHistTrackFilterSpdMult);
  fHistTrackAnaSpdMult->Add(norm->fHistTrackAnaSpdMult);
}
//_______________________________________
/*
Stores the variables used for normalization as function of run number
returns kTRUE if the event is to be counted for normalization
(pass event selection cuts OR has no primary vertex)
 */
void AliNormalizationCounter::StoreEvent(AliVEvent *event,AliRDHFCuts *rdCut,Bool_t mc, Int_t multiplicity){
  //

  Bool_t isEventSelected = rdCut->IsEventSelected(event);

  // events not passing physics selection. do nothing
  if(rdCut->IsEventRejectedDuePhysicsSelection()) return;

  Bool_t v0A=kFALSE; 
  Bool_t v0B=kFALSE;
  Bool_t flag03=kFALSE;
  Bool_t flagPV=kFALSE;

  //Run Number
  Int_t runNumber = event->GetRunNumber();
 
  // Evaluate the multiplicity
  if(multiplicity==-9999) Multiplicity(event);

  //Find CINT1B
  AliESDEvent *eventESD = (AliESDEvent*)event;
  if(!eventESD){AliError("ESD event not available");return;}
  if(mc&&event->GetEventType() != 0)return;
  //event must be either physics or MC
  if(!(event->GetEventType() == 7||event->GetEventType() == 0))return;
  
  if(fMultiplicity) 
    fCounters.Count(Form("Event:triggered/Run:%d/Multiplicity:%d",runNumber,multiplicity));
  else 
    fCounters.Count(Form("Event:triggered/Run:%d",runNumber));

  //Find V0AND
  AliTriggerAnalysis trAn; /// Trigger Analysis
  AliAODVZERO* aodV0 = (AliAODVZERO*)event->GetVZEROData();
  Bool_t isPP2012 = kFALSE;
  if(runNumber>=176326 && runNumber<=193766) isPP2012=kTRUE;
  if(aodV0 && !isPP2012){
    v0B = trAn.IsOfflineTriggerFired(eventESD , AliTriggerAnalysis::kV0C);
    v0A = trAn.IsOfflineTriggerFired(eventESD , AliTriggerAnalysis::kV0A);
  }
  if(v0A&&v0B){
    if(fMultiplicity) 
      fCounters.Count(Form("Event:V0AND/Run:%d/Multiplicity:%d",runNumber,multiplicity));
    else
      fCounters.Count(Form("Event:V0AND/Run:%d",runNumber));
  }
  
  //FindPrimary vertex  
  // AliVVertex *vtrc =  (AliVVertex*)event->GetPrimaryVertex();
  // if(vtrc && vtrc->GetNContributors()>0){
  //   fCounters.Count(Form("Event:PrimaryV/Run:%d",runNumber));
  //   flagPV=kTRUE;
  // }

  //trigger
  AliAODEvent *eventAOD = (AliAODEvent*)event;
  TString trigclass=eventAOD->GetFiredTriggerClasses();
  if(trigclass.Contains("C0SMH-B-NOPF-ALLNOTRD")||trigclass.Contains("C0SMH-B-NOPF-ALL")){
    if(fMultiplicity) 
      fCounters.Count(Form("Event:PbPbC0SMH-B-NOPF-ALLNOTRD/Run:%d/Multiplicity:%d",runNumber,multiplicity));
    else 
      fCounters.Count(Form("Event:PbPbC0SMH-B-NOPF-ALLNOTRD/Run:%d",runNumber));
  }

  //FindPrimary vertex  
  if(isEventSelected){
    if(fMultiplicity) 
      fCounters.Count(Form("Event:PrimaryV/Run:%d/Multiplicity:%d",runNumber,multiplicity));
    else
      fCounters.Count(Form("Event:PrimaryV/Run:%d",runNumber));
    flagPV=kTRUE;
  }else{
    if(rdCut->GetWhyRejection()==0){
      if(fMultiplicity) 
	fCounters.Count(Form("Event:noPrimaryV/Run:%d/Multiplicity:%d",runNumber,multiplicity));
      else
	fCounters.Count(Form("Event:noPrimaryV/Run:%d",runNumber));
    }
    //find good vtx outside range
    if(rdCut->GetWhyRejection()==6){
      if(fMultiplicity) {
	fCounters.Count(Form("Event:zvtxGT10/Run:%d/Multiplicity:%d",runNumber,multiplicity));
	fCounters.Count(Form("Event:PrimaryV/Run:%d/Multiplicity:%d",runNumber,multiplicity));
      } else {
	fCounters.Count(Form("Event:zvtxGT10/Run:%d",runNumber));
	fCounters.Count(Form("Event:PrimaryV/Run:%d",runNumber));
      }
      flagPV=kTRUE;
    }
    if(rdCut->GetWhyRejection()==1){
      if(fMultiplicity) 
	fCounters.Count(Form("Event:PileUp/Run:%d/Multiplicity:%d",runNumber,multiplicity));
      else
	fCounters.Count(Form("Event:PileUp/Run:%d",runNumber));
    }
  }
  //to be counted for normalization
  if(rdCut->CountEventForNormalization()){
    if(fMultiplicity) 
      fCounters.Count(Form("Event:countForNorm/Run:%d/Multiplicity:%d",runNumber,multiplicity));
    else 
      fCounters.Count(Form("Event:countForNorm/Run:%d",runNumber));
  }


  //Find Candle
  Int_t trkEntries = (Int_t)event->GetNumberOfTracks();
  for(Int_t i=0;i<trkEntries&&!flag03;i++){
    AliAODTrack *track=(AliAODTrack*)event->GetTrack(i);
    if((track->Pt()>0.3)&&(!flag03)){
      if(fMultiplicity) 
	fCounters.Count(Form("Event:Candles0.3/Run:%d/Multiplicity:%d",runNumber,multiplicity));
      else
	fCounters.Count(Form("Event:Candles0.3/Run:%d",runNumber));
      flag03=kTRUE;
      break;
    }
  }
  
  if(!(v0A&&v0B)&&(flag03)){ 
    if(fMultiplicity) 
      fCounters.Count(Form("Event:!V0A&Candle03/Run:%d/Multiplicity:%d",runNumber,multiplicity));
    else 
      fCounters.Count(Form("Event:!V0A&Candle03/Run:%d",runNumber));
  }
  if(!(v0A&&v0B)&&flagPV){
    if(fMultiplicity) 
      fCounters.Count(Form("Event:!V0A&PrimaryV/Run:%d/Multiplicity:%d",runNumber,multiplicity));
    else
      fCounters.Count(Form("Event:!V0A&PrimaryV/Run:%d",runNumber));
  }
  
  return;
}
//_____________________________________________________________________
void AliNormalizationCounter::StoreCandidates(AliVEvent *event,Int_t nCand,Bool_t flagFilter){
  
  Int_t ntracks=event->GetNumberOfTracks();
  if(flagFilter)fHistTrackFilterEvMult->Fill(ntracks,nCand);
  else fHistTrackAnaEvMult->Fill(ntracks,nCand);
  Int_t nSPD=0;
  if(fESD){
    AliESDEvent *ESDevent=(AliESDEvent*)event;
    const AliMultiplicity *alimult = ESDevent->GetMultiplicity();
    nSPD = alimult->GetNumberOfTracklets();

  }else{
    AliAODEvent *aodEvent =(AliAODEvent*)event;
    AliAODTracklets *trklets=aodEvent->GetTracklets();
    nSPD = trklets->GetNumberOfTracklets();
  }
  if(flagFilter)fHistTrackFilterSpdMult->Fill(nSPD,nCand);
  else fHistTrackAnaSpdMult->Fill(nSPD,nCand);
  
  Int_t runNumber = event->GetRunNumber();
  Int_t multiplicity = Multiplicity(event);
  if(nCand==0)return;
  if(flagFilter){
    if(fMultiplicity) 
      fCounters.Count(Form("Event:Candid(Filter)/Run:%d/Multiplicity:%d",runNumber,multiplicity));
    else 
      fCounters.Count(Form("Event:Candid(Filter)/Run:%d",runNumber));
    for(Int_t i=0;i<nCand;i++){ 
      if(fMultiplicity) 
	fCounters.Count(Form("Event:NCandid(Filter)/Run:%d/Multiplicity:%d",runNumber,multiplicity));
      else 
	fCounters.Count(Form("Event:NCandid(Filter)/Run:%d",runNumber));
    }
  }else{
    if(fMultiplicity) 
      fCounters.Count(Form("Event:Candid(Analysis)/Run:%d/Multiplicity:%d",runNumber,multiplicity));
    else
      fCounters.Count(Form("Event:Candid(Analysis)/Run:%d",runNumber));
    for(Int_t i=0;i<nCand;i++){ 
      if(fMultiplicity) 
	fCounters.Count(Form("Event:NCandid(Analysis)/Run:%d/Multiplicity:%d",runNumber,multiplicity));
      else
	fCounters.Count(Form("Event:NCandid(Analysis)/Run:%d",runNumber));
    }
  }
  return;
}
//_______________________________________________________________________
TH1D* AliNormalizationCounter::DrawAgainstRuns(TString candle,Bool_t drawHist){
  //
  fCounters.SortRubric("Run");
  TString selection;
  selection.Form("event:%s",candle.Data());
  TH1D* histoneD = fCounters.Get("run",selection.Data());

  histoneD->Sumw2();
  if(drawHist)histoneD->DrawClone();
  return histoneD;
}
//___________________________________________________________________________
TH1D* AliNormalizationCounter::DrawRatio(TString candle1,TString candle2){
  //
  fCounters.SortRubric("Run");
  TString name;

  name.Form("%s/%s",candle1.Data(),candle2.Data());
  TH1D* num=DrawAgainstRuns(candle1.Data(),kFALSE);
  TH1D* den=DrawAgainstRuns(candle2.Data(),kFALSE);

  den->SetTitle(candle2.Data());
  den->SetName(candle2.Data());
  num->Divide(num,den,1,1,"B");
  num->SetTitle(name.Data());
  num->SetName(name.Data());
  num->DrawClone();
  return num;
}
//___________________________________________________________________________
void AliNormalizationCounter::PrintRubrics(){
  fCounters.PrintKeyWords();
}
//___________________________________________________________________________
Double_t AliNormalizationCounter::GetSum(TString candle){
  TString selection="event:";
  selection.Append(candle);
  return fCounters.GetSum(selection.Data());
}
//___________________________________________________________________________
TH2F* AliNormalizationCounter::GetHist(Bool_t filtercuts,Bool_t spdtracklets,Bool_t drawHist){
  if(filtercuts){
    if(spdtracklets){
      if(drawHist)fHistTrackFilterSpdMult->DrawCopy("LEGO2Z 0");
      return fHistTrackFilterSpdMult;
    }else{
      if(drawHist)fHistTrackFilterEvMult->DrawCopy("LEGO2Z 0");
      return fHistTrackFilterEvMult;
    }
  }else{
    if(spdtracklets){
      if(drawHist)fHistTrackAnaSpdMult->DrawCopy("LEGO2Z 0");
      return fHistTrackAnaSpdMult;
    }else{
      if(drawHist)fHistTrackAnaEvMult->DrawCopy("LEGO2Z 0");
      return fHistTrackAnaEvMult;
    }
  }
}
//___________________________________________________________________________
Double_t AliNormalizationCounter::GetNEventsForNorm(){
  Double_t noVtxzGT10=GetSum("noPrimaryV")*GetSum("zvtxGT10")/GetSum("PrimaryV");
  return GetSum("countForNorm")-noVtxzGT10;
}
//___________________________________________________________________________
Double_t AliNormalizationCounter::GetNEventsForNorm(Int_t runnumber){
  TString listofruns = fCounters.GetKeyWords("RUN");
  if(!listofruns.Contains(Form("%d",runnumber))){
    printf("WARNING: %d is not a valid run number\n",runnumber);
    fCounters.Print("Run","",kTRUE);
    return 0.;
  }
  TString suffix;suffix.Form("/RUN:%d",runnumber);
  TString zvtx;zvtx.Form("zvtxGT10%s",suffix.Data());
  TString noPV;noPV.Form("noPrimaryV%s",suffix.Data());
  TString pV;pV.Form("PrimaryV%s",suffix.Data());
  TString tbc;tbc.Form("countForNorm%s",suffix.Data());
  Double_t noVtxzGT10=GetSum(noPV.Data())*GetSum(zvtx.Data())/GetSum(pV.Data());
  return GetSum(tbc.Data())-noVtxzGT10;
}

//___________________________________________________________________________
Double_t AliNormalizationCounter::GetNEventsForNorm(Int_t minmultiplicity, Int_t maxmultiplicity){

  if(!fMultiplicity) {
    AliInfo("Sorry, you didn't activate the multiplicity in the counter!");
    return 0.;
  }

  TString listofruns = fCounters.GetKeyWords("Multiplicity");

  Int_t nmultbins = maxmultiplicity - minmultiplicity;
  Double_t sumnoPV=0., sumZvtx=0., sumPv=0., sumEvtNorm=0.;
  for (Int_t ibin=0; ibin<=nmultbins; ibin++) {
    //    cout << " Looking at bin "<< ibin+minmultiplicity<<endl;
    if(!listofruns.Contains(Form("%d",ibin+minmultiplicity))){
      AliInfo(Form("WARNING: %d is not a valid multiplicity number. \n",ibin+minmultiplicity));
      continue;
    }
    TString suffix;suffix.Form("/Multiplicity:%d",ibin+minmultiplicity);
    TString zvtx;zvtx.Form("zvtxGT10%s",suffix.Data());
    TString noPV;noPV.Form("noPrimaryV%s",suffix.Data());
    TString pV;pV.Form("PrimaryV%s",suffix.Data());
    TString tbc;tbc.Form("countForNorm%s",suffix.Data());
    sumnoPV += GetSum(noPV.Data());
    sumZvtx += GetSum(zvtx.Data());
    sumPv += GetSum(pV.Data());
    sumEvtNorm += GetSum(tbc.Data());
  }
  Double_t noVtxzGT10 = sumPv>0. ? sumnoPV * sumZvtx / sumPv : 0.;
  return sumEvtNorm - noVtxzGT10;
}

//___________________________________________________________________________
TH1D* AliNormalizationCounter::DrawNEventsForNorm(Bool_t drawRatio){
  //usare algebra histos
  fCounters.SortRubric("Run");
  TString selection;

  selection.Form("event:noPrimaryV");
  TH1D* hnoPrimV = fCounters.Get("run",selection.Data());
  hnoPrimV->Sumw2();

  selection.Form("event:zvtxGT10");
  TH1D*  hzvtx= fCounters.Get("run",selection.Data());
  hzvtx->Sumw2();

  selection.Form("event:PrimaryV");
  TH1D* hPrimV = fCounters.Get("run",selection.Data());
  hPrimV->Sumw2();

  hzvtx->Multiply(hnoPrimV);
  hzvtx->Divide(hPrimV);

  selection.Form("event:countForNorm");
  TH1D* hCountForNorm = fCounters.Get("run",selection.Data());
  hCountForNorm->Sumw2();

  hCountForNorm->Add(hzvtx,-1.);

  if(drawRatio){
    selection.Form("event:triggered");
    TH1D* htriggered = fCounters.Get("run",selection.Data());
    htriggered->Sumw2();
    hCountForNorm->Divide(htriggered);
  }

  hCountForNorm->DrawClone();
  return hCountForNorm;
}

//___________________________________________________________________________
Int_t AliNormalizationCounter::Multiplicity(AliVEvent* event){

  Int_t multiplicity = 0;
  AliAODEvent *eventAOD = (AliAODEvent*)event;
  AliAODTracklets * aodTracklets = (AliAODTracklets*)eventAOD->GetTracklets();
  Int_t ntracklets = (Int_t)aodTracklets->GetNumberOfTracklets();
  for(Int_t i=0;i<ntracklets; i++){
    Double_t theta = aodTracklets->GetTheta(i);
    Double_t eta = -TMath::Log( TMath::Tan(theta/2.) ); // check the formula
    if(TMath::Abs(eta)<fMultiplicityEtaRange){ // set the proper cut on eta
      multiplicity++;
    }
  }

  return multiplicity;
}
Commit	Line	Data
8eea8388	1	/**************************************************************************
	2	* Copyright(c) 1998-2008, 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$ */
	17
	18	//*************************************************************************
	19	// Class AliNormalizationCounter
	20	// Class to store the informations relevant for the normalization in the
	21	// barrel for each run
	22	// Authors: G. Ortona, ortona@to.infn.it
	23	// D. Caffarri, davide.caffarri@pd.to.infn.it
	24	// with many thanks to P. Pillot
	25	/////////////////////////////////////////////////////////////
	26
	27	#include "AliLog.h"
	28	#include "AliNormalizationCounter.h"
	29	#include <AliESDEvent.h>
	30	#include <AliESDtrack.h>
	31	#include <AliAODEvent.h>
9c6da02e	32	#include <AliAODVZERO.h>
8eea8388	33	#include <AliVParticle.h>
	34	#include <AliTriggerAnalysis.h>
	35	#include <TH1F.h>
	36	#include <TH2F.h>
	37	#include <TList.h>
	38	#include <TString.h>
	39	#include <TCanvas.h>
	40	#include <AliPhysicsSelection.h>
	41	#include <AliMultiplicity.h>
	42
	43	ClassImp(AliNormalizationCounter)
	44
	45	//____________________________________________
	46	AliNormalizationCounter::AliNormalizationCounter():
	47	TNamed(),
	48	fCounters(),
	49	fESD(kFALSE),
	50	fMultiplicity(kFALSE),
	51	fMultiplicityEtaRange(1.0),
	52	fHistTrackFilterEvMult(0),
	53	fHistTrackAnaEvMult(0),
	54	fHistTrackFilterSpdMult(0),
	55	fHistTrackAnaSpdMult(0)
	56	{
	57	// empty constructor
	58	}
	59
	60	//__________________________________________________
	61	AliNormalizationCounter::AliNormalizationCounter(const char *name):
	62	TNamed(name,name),
	63	fCounters(name),
	64	fESD(kFALSE),
	65	fMultiplicity(kFALSE),
	66	fMultiplicityEtaRange(1.0),
	67	fHistTrackFilterEvMult(0),
	68	fHistTrackAnaEvMult(0),
	69	fHistTrackFilterSpdMult(0),
	70	fHistTrackAnaSpdMult(0)
	71	{
	72	;
	73	}
	74
	75	//______________________________________________
	76	AliNormalizationCounter::~AliNormalizationCounter()
	77	{
	78	//destructor
	79	if(fHistTrackFilterEvMult){
	80	delete fHistTrackFilterEvMult;
	81	fHistTrackFilterEvMult =0;
	82	}
	83	if(fHistTrackAnaEvMult){
	84	delete fHistTrackAnaEvMult;
	85	fHistTrackAnaEvMult=0;
	86	}
	87	if(fHistTrackFilterSpdMult){
	88	delete fHistTrackFilterSpdMult;
	89	fHistTrackFilterSpdMult=0;
	90	}
	91	if(fHistTrackAnaSpdMult){
	92	delete fHistTrackAnaSpdMult;
	93	fHistTrackAnaSpdMult=0;
	94	}
	95	}
	96
97	//______________________________________________
98	void AliNormalizationCounter::Init()
99	{
100	//variables initialization
101	fCounters.AddRubric("Event","triggered/V0AND/PileUp/PbPbC0SMH-B-NOPF-ALLNOTRD/Candles0.3/PrimaryV/countForNorm/noPrimaryV/zvtxGT10/!V0A&Candle03/!V0A&PrimaryV/Candid(Filter)/Candid(Analysis)/NCandid(Filter)/NCandid(Analysis)");
102	if(fMultiplicity) fCounters.AddRubric("Multiplicity", 5000);
103	fCounters.AddRubric("Run", 1000000);
104	fCounters.Init();
105	fHistTrackFilterEvMult=new TH2F("FiltCandidvsTracksinEv","FiltCandidvsTracksinEv",10000,-0.5,9999.5,200,-0.5,199.5);
106	fHistTrackFilterEvMult->GetYaxis()->SetTitle("NCandidates");
107	fHistTrackFilterEvMult->GetXaxis()->SetTitle("NTracksinEvent");
108	fHistTrackAnaEvMult=new TH2F("AnaCandidvsTracksinEv","AnaCandidvsTracksinEv",10000,-0.5,9999.5,100,-0.5,99.5);
109	fHistTrackAnaEvMult->GetYaxis()->SetTitle("NCandidates");
110	fHistTrackAnaEvMult->GetXaxis()->SetTitle("NTracksinEvent");
111	fHistTrackFilterSpdMult=new TH2F("FilterCandidvsSpdMult","FilterCandidvsSpdMult",5000,-0.5,4999.5,200,-0.5,199.5);
112	fHistTrackFilterSpdMult->GetYaxis()->SetTitle("NCandidates");
113	fHistTrackFilterSpdMult->GetXaxis()->SetTitle("NSPDTracklets");
114	fHistTrackAnaSpdMult=new TH2F("AnaCandidvsSpdMult","AnaCandidvsSpdMult",5000,-0.5,4999.5,100,-0.5,99.5);
115	fHistTrackAnaSpdMult->GetYaxis()->SetTitle("NCandidates");
116	fHistTrackAnaSpdMult->GetXaxis()->SetTitle("NSPDTracklets");
117	}
118
119	//______________________________________________
120	Long64_t AliNormalizationCounter::Merge(TCollection* list){
121	if (!list) return 0;
122	if (list->IsEmpty()) return 0;//(Long64_t)fCounters.Merge(list);
123
124	TIter next(list);
125	const TObject* obj = 0x0;
126	while ((obj = next())) {
127
128	// check that "obj" is an object of the class AliNormalizationCounter
129	const AliNormalizationCounter* counter = dynamic_cast<const AliNormalizationCounter*>(obj);
130	if (!counter) {
131	AliError(Form("object named %s is not AliNormalizationCounter! Skipping it.", counter->GetName()));
132	continue;
133	}
134
135	Add(counter);
136
137	}
138
139	return (Long64_t)1;//(Long64_t)fCounters->GetEntries();
140	}
141	//_______________________________________
142	void AliNormalizationCounter::Add(const AliNormalizationCounter *norm){
143	fCounters.Add(&(norm->fCounters));
144	fHistTrackFilterEvMult->Add(norm->fHistTrackFilterEvMult);
145	fHistTrackAnaEvMult->Add(norm->fHistTrackAnaEvMult);
146	fHistTrackFilterSpdMult->Add(norm->fHistTrackFilterSpdMult);
147	fHistTrackAnaSpdMult->Add(norm->fHistTrackAnaSpdMult);
148	}
149	//_______________________________________
150	/*
151	Stores the variables used for normalization as function of run number
152	returns kTRUE if the event is to be counted for normalization
153	(pass event selection cuts OR has no primary vertex)
154	*/
71d2489f	155	void AliNormalizationCounter::StoreEvent(AliVEvent event,AliRDHFCuts rdCut,Bool_t mc, Int_t multiplicity){
8eea8388	156	//
	157
	158	Bool_t isEventSelected = rdCut->IsEventSelected(event);
	159
	160	// events not passing physics selection. do nothing
	161	if(rdCut->IsEventRejectedDuePhysicsSelection()) return;
	162
	163	Bool_t v0A=kFALSE;
	164	Bool_t v0B=kFALSE;
	165	Bool_t flag03=kFALSE;
	166	Bool_t flagPV=kFALSE;
	167
	168	//Run Number
	169	Int_t runNumber = event->GetRunNumber();
	170
	171	// Evaluate the multiplicity
71d2489f	172	if(multiplicity==-9999) Multiplicity(event);
8eea8388	173
	174	//Find CINT1B
	175	AliESDEvent eventESD = (AliESDEvent)event;
	176	if(!eventESD){AliError("ESD event not available");return;}
	177	if(mc&&event->GetEventType() != 0)return;
	178	//event must be either physics or MC
	179	if(!(event->GetEventType() == 7\|\|event->GetEventType() == 0))return;
	180
	181	if(fMultiplicity)
	182	fCounters.Count(Form("Event:triggered/Run:%d/Multiplicity:%d",runNumber,multiplicity));
	183	else
	184	fCounters.Count(Form("Event:triggered/Run:%d",runNumber));
	185
	186	//Find V0AND
	187	AliTriggerAnalysis trAn; /// Trigger Analysis
9c6da02e	188	AliAODVZERO* aodV0 = (AliAODVZERO*)event->GetVZEROData();
813d4979	189	Bool_t isPP2012 = kFALSE;
	190	if(runNumber>=176326 && runNumber<=193766) isPP2012=kTRUE;
	191	if(aodV0 && !isPP2012){
9c6da02e	192	v0B = trAn.IsOfflineTriggerFired(eventESD , AliTriggerAnalysis::kV0C);
	193	v0A = trAn.IsOfflineTriggerFired(eventESD , AliTriggerAnalysis::kV0A);
	194	}
8eea8388	195	if(v0A&&v0B){
	196	if(fMultiplicity)
	197	fCounters.Count(Form("Event:V0AND/Run:%d/Multiplicity:%d",runNumber,multiplicity));
	198	else
	199	fCounters.Count(Form("Event:V0AND/Run:%d",runNumber));
	200	}
	201
	202	//FindPrimary vertex
	203	// AliVVertex vtrc = (AliVVertex)event->GetPrimaryVertex();
	204	// if(vtrc && vtrc->GetNContributors()>0){
	205	// fCounters.Count(Form("Event:PrimaryV/Run:%d",runNumber));
	206	// flagPV=kTRUE;
	207	// }
	208
	209	//trigger
	210	AliAODEvent eventAOD = (AliAODEvent)event;
	211	TString trigclass=eventAOD->GetFiredTriggerClasses();
	212	if(trigclass.Contains("C0SMH-B-NOPF-ALLNOTRD")\|\|trigclass.Contains("C0SMH-B-NOPF-ALL")){
	213	if(fMultiplicity)
	214	fCounters.Count(Form("Event:PbPbC0SMH-B-NOPF-ALLNOTRD/Run:%d/Multiplicity:%d",runNumber,multiplicity));
	215	else
	216	fCounters.Count(Form("Event:PbPbC0SMH-B-NOPF-ALLNOTRD/Run:%d",runNumber));
	217	}
	218
	219	//FindPrimary vertex
	220	if(isEventSelected){
	221	if(fMultiplicity)
	222	fCounters.Count(Form("Event:PrimaryV/Run:%d/Multiplicity:%d",runNumber,multiplicity));
	223	else
	224	fCounters.Count(Form("Event:PrimaryV/Run:%d",runNumber));
	225	flagPV=kTRUE;
	226	}else{
	227	if(rdCut->GetWhyRejection()==0){
	228	if(fMultiplicity)
	229	fCounters.Count(Form("Event:noPrimaryV/Run:%d/Multiplicity:%d",runNumber,multiplicity));
	230	else
	231	fCounters.Count(Form("Event:noPrimaryV/Run:%d",runNumber));
	232	}
	233	//find good vtx outside range
	234	if(rdCut->GetWhyRejection()==6){
	235	if(fMultiplicity) {
	236	fCounters.Count(Form("Event:zvtxGT10/Run:%d/Multiplicity:%d",runNumber,multiplicity));
	237	fCounters.Count(Form("Event:PrimaryV/Run:%d/Multiplicity:%d",runNumber,multiplicity));
	238	} else {
	239	fCounters.Count(Form("Event:zvtxGT10/Run:%d",runNumber));
	240	fCounters.Count(Form("Event:PrimaryV/Run:%d",runNumber));
	241	}
	242	flagPV=kTRUE;
	243	}
	244	if(rdCut->GetWhyRejection()==1){
	245	if(fMultiplicity)
	246	fCounters.Count(Form("Event:PileUp/Run:%d/Multiplicity:%d",runNumber,multiplicity));
	247	else
	248	fCounters.Count(Form("Event:PileUp/Run:%d",runNumber));
	249	}
	250	}
	251	//to be counted for normalization
	252	if(rdCut->CountEventForNormalization()){
	253	if(fMultiplicity)
	254	fCounters.Count(Form("Event:countForNorm/Run:%d/Multiplicity:%d",runNumber,multiplicity));
	255	else
	256	fCounters.Count(Form("Event:countForNorm/Run:%d",runNumber));
	257	}
	258
259
260	//Find Candle
261	Int_t trkEntries = (Int_t)event->GetNumberOfTracks();
262	for(Int_t i=0;i<trkEntries&&!flag03;i++){
263	AliAODTrack track=(AliAODTrack)event->GetTrack(i);
264	if((track->Pt()>0.3)&&(!flag03)){
265	if(fMultiplicity)
266	fCounters.Count(Form("Event:Candles0.3/Run:%d/Multiplicity:%d",runNumber,multiplicity));
267	else
268	fCounters.Count(Form("Event:Candles0.3/Run:%d",runNumber));
269	flag03=kTRUE;
270	break;
271	}
272	}
273
274	if(!(v0A&&v0B)&&(flag03)){
275	if(fMultiplicity)
276	fCounters.Count(Form("Event:!V0A&Candle03/Run:%d/Multiplicity:%d",runNumber,multiplicity));
277	else
278	fCounters.Count(Form("Event:!V0A&Candle03/Run:%d",runNumber));
279	}
280	if(!(v0A&&v0B)&&flagPV){
281	if(fMultiplicity)
282	fCounters.Count(Form("Event:!V0A&PrimaryV/Run:%d/Multiplicity:%d",runNumber,multiplicity));
283	else
284	fCounters.Count(Form("Event:!V0A&PrimaryV/Run:%d",runNumber));
285	}
286
287	return;
288	}
289	//_____________________________________________________________________
290	void AliNormalizationCounter::StoreCandidates(AliVEvent *event,Int_t nCand,Bool_t flagFilter){
291
292	Int_t ntracks=event->GetNumberOfTracks();
293	if(flagFilter)fHistTrackFilterEvMult->Fill(ntracks,nCand);
294	else fHistTrackAnaEvMult->Fill(ntracks,nCand);
295	Int_t nSPD=0;
296	if(fESD){
297	AliESDEvent ESDevent=(AliESDEvent)event;
298	const AliMultiplicity *alimult = ESDevent->GetMultiplicity();
299	nSPD = alimult->GetNumberOfTracklets();
300
301	}else{
302	AliAODEvent aodEvent =(AliAODEvent)event;
303	AliAODTracklets *trklets=aodEvent->GetTracklets();
304	nSPD = trklets->GetNumberOfTracklets();
305	}
306	if(flagFilter)fHistTrackFilterSpdMult->Fill(nSPD,nCand);
307	else fHistTrackAnaSpdMult->Fill(nSPD,nCand);
308
309	Int_t runNumber = event->GetRunNumber();
310	Int_t multiplicity = Multiplicity(event);
311	if(nCand==0)return;
312	if(flagFilter){
313	if(fMultiplicity)
314	fCounters.Count(Form("Event:Candid(Filter)/Run:%d/Multiplicity:%d",runNumber,multiplicity));
315	else
316	fCounters.Count(Form("Event:Candid(Filter)/Run:%d",runNumber));
317	for(Int_t i=0;i<nCand;i++){
318	if(fMultiplicity)
319	fCounters.Count(Form("Event:NCandid(Filter)/Run:%d/Multiplicity:%d",runNumber,multiplicity));
320	else
321	fCounters.Count(Form("Event:NCandid(Filter)/Run:%d",runNumber));
322	}
323	}else{
324	if(fMultiplicity)
325	fCounters.Count(Form("Event:Candid(Analysis)/Run:%d/Multiplicity:%d",runNumber,multiplicity));
326	else
327	fCounters.Count(Form("Event:Candid(Analysis)/Run:%d",runNumber));
328	for(Int_t i=0;i<nCand;i++){
329	if(fMultiplicity)
330	fCounters.Count(Form("Event:NCandid(Analysis)/Run:%d/Multiplicity:%d",runNumber,multiplicity));
331	else
332	fCounters.Count(Form("Event:NCandid(Analysis)/Run:%d",runNumber));
333	}
334	}
335	return;
336	}
337	//_______________________________________________________________________
338	TH1D* AliNormalizationCounter::DrawAgainstRuns(TString candle,Bool_t drawHist){
339	//
340	fCounters.SortRubric("Run");
341	TString selection;
342	selection.Form("event:%s",candle.Data());
343	TH1D* histoneD = fCounters.Get("run",selection.Data());
344
345	histoneD->Sumw2();
346	if(drawHist)histoneD->DrawClone();
347	return histoneD;
348	}
349	//___________________________________________________________________________
350	TH1D* AliNormalizationCounter::DrawRatio(TString candle1,TString candle2){
351	//
352	fCounters.SortRubric("Run");
353	TString name;
354
355	name.Form("%s/%s",candle1.Data(),candle2.Data());
356	TH1D* num=DrawAgainstRuns(candle1.Data(),kFALSE);
357	TH1D* den=DrawAgainstRuns(candle2.Data(),kFALSE);
358
359	den->SetTitle(candle2.Data());
360	den->SetName(candle2.Data());
361	num->Divide(num,den,1,1,"B");
362	num->SetTitle(name.Data());
363	num->SetName(name.Data());
364	num->DrawClone();
365	return num;
366	}
367	//___________________________________________________________________________
368	void AliNormalizationCounter::PrintRubrics(){
369	fCounters.PrintKeyWords();
370	}
371	//___________________________________________________________________________
372	Double_t AliNormalizationCounter::GetSum(TString candle){
373	TString selection="event:";
374	selection.Append(candle);
375	return fCounters.GetSum(selection.Data());
376	}
377	//___________________________________________________________________________
378	TH2F* AliNormalizationCounter::GetHist(Bool_t filtercuts,Bool_t spdtracklets,Bool_t drawHist){
379	if(filtercuts){
380	if(spdtracklets){
381	if(drawHist)fHistTrackFilterSpdMult->DrawCopy("LEGO2Z 0");
382	return fHistTrackFilterSpdMult;
383	}else{
384	if(drawHist)fHistTrackFilterEvMult->DrawCopy("LEGO2Z 0");
385	return fHistTrackFilterEvMult;
386	}
387	}else{
388	if(spdtracklets){
389	if(drawHist)fHistTrackAnaSpdMult->DrawCopy("LEGO2Z 0");
390	return fHistTrackAnaSpdMult;
391	}else{
392	if(drawHist)fHistTrackAnaEvMult->DrawCopy("LEGO2Z 0");
393	return fHistTrackAnaEvMult;
394	}
395	}
396	}
397	//___________________________________________________________________________
398	Double_t AliNormalizationCounter::GetNEventsForNorm(){
399	Double_t noVtxzGT10=GetSum("noPrimaryV")*GetSum("zvtxGT10")/GetSum("PrimaryV");
400	return GetSum("countForNorm")-noVtxzGT10;
401	}
402	//___________________________________________________________________________
403	Double_t AliNormalizationCounter::GetNEventsForNorm(Int_t runnumber){
404	TString listofruns = fCounters.GetKeyWords("RUN");
405	if(!listofruns.Contains(Form("%d",runnumber))){
406	printf("WARNING: %d is not a valid run number\n",runnumber);
407	fCounters.Print("Run","",kTRUE);
408	return 0.;
409	}
410	TString suffix;suffix.Form("/RUN:%d",runnumber);
411	TString zvtx;zvtx.Form("zvtxGT10%s",suffix.Data());
412	TString noPV;noPV.Form("noPrimaryV%s",suffix.Data());
413	TString pV;pV.Form("PrimaryV%s",suffix.Data());
414	TString tbc;tbc.Form("countForNorm%s",suffix.Data());
415	Double_t noVtxzGT10=GetSum(noPV.Data())*GetSum(zvtx.Data())/GetSum(pV.Data());
416	return GetSum(tbc.Data())-noVtxzGT10;
417	}
418
419	//___________________________________________________________________________
420	Double_t AliNormalizationCounter::GetNEventsForNorm(Int_t minmultiplicity, Int_t maxmultiplicity){
421
422	if(!fMultiplicity) {
423	AliInfo("Sorry, you didn't activate the multiplicity in the counter!");
424	return 0.;
425	}
426
427	TString listofruns = fCounters.GetKeyWords("Multiplicity");
428
429	Int_t nmultbins = maxmultiplicity - minmultiplicity;
430	Double_t sumnoPV=0., sumZvtx=0., sumPv=0., sumEvtNorm=0.;
431	for (Int_t ibin=0; ibin<=nmultbins; ibin++) {
432	// cout << " Looking at bin "<< ibin+minmultiplicity<<endl;
433	if(!listofruns.Contains(Form("%d",ibin+minmultiplicity))){
434	AliInfo(Form("WARNING: %d is not a valid multiplicity number. \n",ibin+minmultiplicity));
435	continue;
436	}
437	TString suffix;suffix.Form("/Multiplicity:%d",ibin+minmultiplicity);
438	TString zvtx;zvtx.Form("zvtxGT10%s",suffix.Data());
439	TString noPV;noPV.Form("noPrimaryV%s",suffix.Data());
440	TString pV;pV.Form("PrimaryV%s",suffix.Data());
441	TString tbc;tbc.Form("countForNorm%s",suffix.Data());
442	sumnoPV += GetSum(noPV.Data());
443	sumZvtx += GetSum(zvtx.Data());
444	sumPv += GetSum(pV.Data());
445	sumEvtNorm += GetSum(tbc.Data());
446	}
447	Double_t noVtxzGT10 = sumPv>0. ? sumnoPV * sumZvtx / sumPv : 0.;
448	return sumEvtNorm - noVtxzGT10;
449	}
450
451	//___________________________________________________________________________
452	TH1D* AliNormalizationCounter::DrawNEventsForNorm(Bool_t drawRatio){
453	//usare algebra histos
454	fCounters.SortRubric("Run");
455	TString selection;
456
457	selection.Form("event:noPrimaryV");
458	TH1D* hnoPrimV = fCounters.Get("run",selection.Data());
459	hnoPrimV->Sumw2();
460
461	selection.Form("event:zvtxGT10");
462	TH1D* hzvtx= fCounters.Get("run",selection.Data());
463	hzvtx->Sumw2();
464
465	selection.Form("event:PrimaryV");
466	TH1D* hPrimV = fCounters.Get("run",selection.Data());
467	hPrimV->Sumw2();
468
469	hzvtx->Multiply(hnoPrimV);
470	hzvtx->Divide(hPrimV);
471
472	selection.Form("event:countForNorm");
473	TH1D* hCountForNorm = fCounters.Get("run",selection.Data());
474	hCountForNorm->Sumw2();
475
476	hCountForNorm->Add(hzvtx,-1.);
477
478	if(drawRatio){
479	selection.Form("event:triggered");
480	TH1D* htriggered = fCounters.Get("run",selection.Data());
481	htriggered->Sumw2();
482	hCountForNorm->Divide(htriggered);
483	}
484
485	hCountForNorm->DrawClone();
486	return hCountForNorm;
487	}
488
489	//___________________________________________________________________________
490	Int_t AliNormalizationCounter::Multiplicity(AliVEvent* event){
491
492	Int_t multiplicity = 0;
493	AliAODEvent eventAOD = (AliAODEvent)event;
494	AliAODTracklets * aodTracklets = (AliAODTracklets*)eventAOD->GetTracklets();
495	Int_t ntracklets = (Int_t)aodTracklets->GetNumberOfTracklets();
496	for(Int_t i=0;i<ntracklets; i++){
497	Double_t theta = aodTracklets->GetTheta(i);
498	Double_t eta = -TMath::Log( TMath::Tan(theta/2.) ); // check the formula
499	if(TMath::Abs(eta)<fMultiplicityEtaRange){ // set the proper cut on eta
500	multiplicity++;
501	}
502	}
503
504	return multiplicity;
505	}