[u/mrichter/AliRoot.git] / PWG3 / muon / AliMuonsHFHeader.cxx

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

/////////////////////////////////////////////////////////////
//
// class used to extract and store info at event level
//
// Author: X-M. Zhang, zhang@clermont.in2p3.fr
//                     zhangxm@iopp.ccnu.edu.cn
/////////////////////////////////////////////////////////////

#include <TMath.h>
#include <TH1.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TList.h>
#include <TObjArray.h>
#include <TObjString.h>

#include "AliTriggerAnalysis.h"
#include "AliBackgroundSelection.h"
#include "AliAODVertex.h"
#include "AliAODEvent.h"
#include "AliESDEvent.h"
#include "AliMultiplicity.h"
#include "AliMuonInfoStoreRD.h"
#include "AliMuonInfoStoreMC.h"
#include "AliDimuInfoStoreRD.h"
#include "AliDimuInfoStoreMC.h"
#include "AliMuonsHFHeader.h"

class TNamed;
class AliESDVertex;

ClassImp(AliMuonsHFHeader)

const TString AliMuonsHFHeader::fgkStdBranchName("MuEvsH");
Int_t         AliMuonsHFHeader::fgAnaMode         = 0;
Bool_t        AliMuonsHFHeader::fgIsMC            = kFALSE;
Bool_t        AliMuonsHFHeader::fgIsEventSelected = kFALSE;
Double_t      AliMuonsHFHeader::fgCuts[3] = { -999999., 999999., 999999.};

//_____________________________________________________________________________
AliMuonsHFHeader::AliMuonsHFHeader() :
TNamed(),
fTriggerMask(0),
fFiredTrigger(0),
fNFiredTrigger(0),
fIsPhysicsTriggered(kFALSE),
fIsPhysicsAccepted(kFALSE),
fEventType(0),
fUnrecoVertex(kFALSE),
fNContributors(0),
fUnrecoVtxSPD(kFALSE),
fNContributorsSPD(0),
fNTrackletsSPD(0),
fCentrality(0.)
{
  //
  // default constructor
  //
  for (Int_t i=3; i--;) fVtx[i] = 0.;
  for (Int_t i=3; i--;) fVtxSPD[i] = 0.;
}

//_____________________________________________________________________________
AliMuonsHFHeader::AliMuonsHFHeader(const AliMuonsHFHeader &src) :
TNamed(),
fTriggerMask(src.fTriggerMask),
fFiredTrigger(src.fFiredTrigger),
fNFiredTrigger(src.fNFiredTrigger),
fIsPhysicsTriggered(src.fIsPhysicsTriggered),
fIsPhysicsAccepted(src.fIsPhysicsAccepted),
fEventType(src.fEventType),
fUnrecoVertex(src.fUnrecoVertex),
fNContributors(src.fNContributors),
fUnrecoVtxSPD(src.fUnrecoVtxSPD),
fNContributorsSPD(src.fNContributorsSPD),
fNTrackletsSPD(src.fNTrackletsSPD),
fCentrality(src.fCentrality)
{
  //
  // copy constructor
  //
  for (Int_t i=3; i--;) fVtx[i]    = src.fVtx[i];
  for (Int_t i=3; i--;) fVtxSPD[i] = src.fVtxSPD[i];
}

//_____________________________________________________________________________
AliMuonsHFHeader& AliMuonsHFHeader::operator=(const AliMuonsHFHeader &src)
{
  //
  // assignment constructor
  //
  fTriggerMask        = src.fTriggerMask;
  fFiredTrigger       = src.fFiredTrigger;
  fNFiredTrigger      = src.fNFiredTrigger; 
  fIsPhysicsTriggered = src.fIsPhysicsTriggered;
  fIsPhysicsAccepted  = src.fIsPhysicsAccepted;
  fEventType          = src.fEventType;
  fUnrecoVertex       = src.fUnrecoVertex;
  fNContributors      = src.fNContributors;
  fUnrecoVtxSPD       = src.fUnrecoVtxSPD;
  fNContributorsSPD   = src.fNContributorsSPD;
  fNTrackletsSPD      = src.fNTrackletsSPD;
  fCentrality         = src.fCentrality;

  for (Int_t i=3; i--;) fVtx[i]      = src.fVtx[i];
  for (Int_t i=3; i--;) fVtxSPD[i]   = src.fVtxSPD[i];

  return *this;
}

//_____________________________________________________________________________
AliMuonsHFHeader::~AliMuonsHFHeader()
{
  //
  // default destructor
  //
}

//_____________________________________________________________________________
void AliMuonsHFHeader::SetEvent(AliAODEvent *event)
{
  // extract event info from AOD event

  fTriggerMask = event->GetTriggerMask();

  AliAODVertex *vertex = event->GetPrimaryVertex(); 
  vertex->GetXYZ(fVtx);
  fNContributors = vertex->GetNContributors();
  fUnrecoVertex = (TMath::Abs(fVtx[0])<1e-6 && TMath::Abs(fVtx[1])<1e-6 &&
                   TMath::Abs(fVtx[2])<1e-6);
  this->SetTitle(vertex->GetTitle());

  this->SetFiredTrigger(event->GetFiredTriggerClasses());
  this->EventSelection();
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::SetEvent(AliESDEvent *event)
{
  // extract event info from ESD event
  // SPD vertex and Physics event selection are implimented

  fTriggerMask = event->GetTriggerMask();

  const AliESDVertex *vertex = event->GetPrimaryVertex(); 
  vertex->GetXYZ(fVtx);
  fNContributors = vertex->GetNContributors();
  fUnrecoVertex = (TMath::Abs(fVtx[0])<1e-6 && TMath::Abs(fVtx[1])<1e-6 &&
                   TMath::Abs(fVtx[2])<1e-6);
  this->SetTitle(vertex->GetTitle());

  const AliESDVertex *vtxSPD = event->GetPrimaryVertexSPD();
  vtxSPD->GetXYZ(fVtxSPD);
  fNContributorsSPD = vtxSPD->GetNContributors();
  fUnrecoVtxSPD = (TMath::Abs(fVtxSPD[0])<1e-6 && TMath::Abs(fVtxSPD[1])<1e-6 &&
                   TMath::Abs(fVtxSPD[2])<1e-6);
  fNTrackletsSPD = event->GetMultiplicity()->GetNumberOfTracklets();

  this->SetFiredTrigger(event->GetFiredTriggerClasses());
  this->PhysicsTriggerAna(event);
  this->EventSelection();
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::PhysicsTriggerAna(const AliESDEvent *esd)
{
  // ESD event trigger condition analysis
  // according to the method in $ALICE_ROOT/ANALYSIS/AliPhysicsSelection.cxx

  fEventType = esd->GetHeader()->GetEventType();
  fIsPhysicsTriggered = kFALSE;
  fIsPhysicsAccepted  = kFALSE;

  AliTriggerAnalysis *triggerAna = new AliTriggerAnalysis();
  triggerAna->SetAnalyzeMC(fgIsMC);
  triggerAna->SetSPDGFOThreshhold(1);

  Int_t  triggerHW  = triggerAna->SPDFiredChips(esd, 1);
  Bool_t isFiredv0A = triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0A);
  Bool_t isFiredv0C = triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0C);
  if (triggerHW==0 && !isFiredv0A && !isFiredv0C) {
    delete triggerAna;
    triggerAna = 0;
    return;
  }

  Bool_t triggerBG = (triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0ABG) ||
                      triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0CBG));
  Int_t  isFiredSPD = triggerAna->SPDFiredChips(esd, 0);
  Bool_t triggerFD = ((isFiredSPD>1) || (isFiredSPD>0 && (isFiredv0A || isFiredv0C)) || (isFiredv0A && isFiredv0C));
  if ((!triggerBG) && triggerFD) fIsPhysicsTriggered = kTRUE; 
  delete triggerAna;
  triggerAna = 0;

  if (fIsPhysicsTriggered) {
    AliBackgroundSelection *bkgId = new AliBackgroundSelection();
    fIsPhysicsAccepted = bkgId->IsSelected(const_cast<AliESDEvent*>(esd));
    delete bkgId; bkgId=0;
  }
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::SetFiredTrigger(TString trigger)
{
  // get info of fired trigger classes of event

  fFiredTrigger = trigger;
  TObjArray *tokens = trigger.Tokenize(" ");
  fNFiredTrigger = tokens->GetEntries();
  return;
}

//_____________________________________________________________________________
Bool_t AliMuonsHFHeader::IsTriggerFired(TString trigger)
{
  // check whether the trigger class "trigger" is fired in this event

  if (fNFiredTrigger<=0) return kFALSE;

  TObjArray *tokens = fFiredTrigger.Tokenize(" "); 
  for (Int_t i=fNFiredTrigger; i--;) {
    TString fired = ((TObjString*)tokens->At(i))->String(); 
    if (fired.CompareTo(trigger)==0) return kTRUE;
  }
  return kFALSE;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::EventSelection(TString triggerName)
{
  // select event according to the "triggerName" & event selection cuts

  fgIsEventSelected = kFALSE;
  if (!this->IsPhysicsAccepted()) return;
  if (!fgIsMC && !this->IsTriggerFired(triggerName)) return;
  this->EventSelection();
  return; 
}

//_____________________________________________________________________________
void AliMuonsHFHeader::EventSelection()
{
  // select event according to the event selection cuts

  fgIsEventSelected = kFALSE;
  if (this->NVtxContributorsSPD()<fgCuts[0]) return;
  if (TMath::Abs(this->VzSPD())>fgCuts[1]) return;
  if (this->VtSPD()>fgCuts[2]) return;
  fgIsEventSelected = kTRUE;
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::CreateHistograms(TList *listEvent, TList *listMuon, TList *listDimu)
{
  // create output histos of muon analysis according to the analysis mode & MC flag

  this->CreateHistosEventH(listEvent);
  if (fgIsMC) {
    if (fgAnaMode!=2) this->CreateHistosMuonMC(listMuon);
    if (fgAnaMode!=1) this->CreateHistosDimuMC(listDimu);
  } else {
    if (fgAnaMode!=2) this->CreateHistosMuonRD(listMuon);
    if (fgAnaMode!=1) this->CreateHistosDimuRD(listDimu);
  }
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::CreateHistosEventH(TList *list)
{
  // create histograms at event level

  if (!list) list = new TList();
  list->SetOwner();
  Bool_t oldStatus = TH1::AddDirectoryStatus();
  TH1::AddDirectory(kFALSE);

  const Int_t nHistos = 5;
  char    *name[nHistos] = { "hVztSPD", "hVzxSPD", "hVzySPD", "hVzNcontr", "hVtNcontr"  };
  Int_t  nbinsX[nHistos] = {      800 ,      800 ,      800 ,       800  ,        400   };
  Double_t xlow[nHistos] = {      -40.,      -40.,      -40.,       -40. ,          0.  };
  Double_t  xup[nHistos] = {       40.,       40.,       40.,        40. ,          4.  };
  Int_t  nbinsY[nHistos] = {      400 ,      600 ,      600 ,       202  ,        202   };
  Double_t ylow[nHistos] = {        0.,       -3.,       -3.,        -2.5,         -2.5 };
  Double_t  yup[nHistos] = {        4.,        3.,        3.,       199.5,        199.5 };

  TH2F *histo = 0;
  for (Int_t i=0; i<nHistos; i++) {
    histo = new TH2F(name[i], name[i], nbinsX[i], xlow[i], xup[i], nbinsY[i], ylow[i], yup[i]);
    histo->Sumw2(); list->AddAt(histo, i); histo = 0;
  }

  TH1::AddDirectory(oldStatus);
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::CreateHistosMuonRD(TList *list)
{
  // create histograms for single muon

  if (!list) list = new TList();
  list->SetOwner();
  Bool_t oldStatus = TH1::AddDirectoryStatus();
  TH1::AddDirectory(kFALSE);
  
  const Int_t nHistos = 8; 
  char    *name[nHistos] = {  "hP", "hPt", "hEta", "hDCA", "hTrg", "hCharge", "hUnfVtx" , "hEtaTimesDCA"};
  Int_t   nbins[nHistos] = { 1500 ,  300 ,   100 ,   500 ,    4  ,       3  ,       80  ,         10000 };
  Double_t xlow[nHistos] = {    0.,    0.,   -10.,     0.,   -0.5,      -1.5,      -40. ,             0.};
  Double_t  xup[nHistos] = {  150.,   30.,     0.,   500.,    3.5,       1.5,       40. ,          1000.};

  TH1F *histo = 0;
  for (Int_t i=0; i<nHistos; i++) {
    histo = new TH1F(name[i], name[i], nbins[i], xlow[i], xup[i]);
    histo->Sumw2(); list->AddAt(histo, i); histo = 0;
  }

  TH1::AddDirectory(oldStatus);
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::CreateHistosDimuRD(TList *list)
{
  // create histograms for dimuon

  if (!list) list = new TList();
  list->SetOwner();
  Bool_t oldStatus = TH1::AddDirectoryStatus();
  TH1::AddDirectory(kFALSE);

  TH1F *histo = 0;
  const Int_t nHistos = 3;
  char    *name[nHistos] = {  "hP", "hPt", "hInvM" };
  Int_t   nbins[nHistos] = { 1500 ,  300 ,    300  };
  Double_t xlow[nHistos] = {    0.,    0.,      0. };
  Double_t  xup[nHistos] = {  150.,   30.,     30. };
  char *dimuName[3] = {"DimuNN", "DimuNP", "DimuPP"};
  for (Int_t i=0; i<3; i++) {
    for (Int_t j=0; j<nHistos; j++) {
      histo = new TH1F(Form("%s_%s",name[j],dimuName[i]), name[j], nbins[j], xlow[j], xup[j]);
      histo->Sumw2(); list->AddAt(histo,i*nHistos+j); histo = 0;
    }
  }

  TH1::AddDirectory(oldStatus);
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::CreateHistosMuonMC(TList *list)
{
  // create histograms for single muon with MC info

  if (!list) list = new TList[AliMuonInfoStoreMC::NSources()];
  for (Int_t i=AliMuonInfoStoreMC::NSources(); i--;) this->CreateHistosMuonRD(&list[i]);
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::CreateHistosDimuMC(TList *list)
{
  // create histograms for dimuon with MC info

  if (!list) list = new TList[AliDimuInfoStoreMC::NSources()];
  for (Int_t i=AliDimuInfoStoreMC::NSources(); i--;) this->CreateHistosDimuRD(&list[i]);
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::FillHistosEventH(TList *list)
{
  // fill histograms at event level according to event selection cuts

  if (!list) return;
  if (!AliMuonsHFHeader::IsEventSelected()) return;

  const Int_t nHistos = 5;
  Double_t vz = this->VzSPD();
  Double_t vt = this->VtSPD();
  Int_t    nc = this->NVtxContributorsSPD();
  Double_t distX[nHistos] = { vz,            vz,            vz, vz, vt };
  Double_t distY[nHistos] = { vt, this->VxSPD(), this->VySPD(), nc, nc };
  for (Int_t i=nHistos; i--;) ((TH2F*)list->At(i))->Fill(distX[i], distY[i]);
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::FillHistosMuonRD(TList *list, AliMuonInfoStoreRD* const muonStoreRD)
{
  // fill histograms for single muon according to event & muon track selection cuts

  if (!list) return;
  if (!AliMuonsHFHeader::IsEventSelected()) return;
  if (!muonStoreRD->MuonSelection()) return;

  const Int_t nHistos = 8;
  Double_t dist[nHistos] = {muonStoreRD->Momentum().Mag(),
                            muonStoreRD->Momentum().Pt(),
                            muonStoreRD->Momentum().Eta(),
                            muonStoreRD->DCA(),
                            muonStoreRD->MatchTrigger(),
                            muonStoreRD->Charge(),
                            this->VzSPD(),
                            muonStoreRD->DCA()*TMath::Exp(-2.*muonStoreRD->Momentum().Eta())/1000.};
  for (Int_t i=nHistos; i--;) ((TH1F*)list->At(i))->Fill(dist[i]);
  return; 
}

//_____________________________________________________________________________
void AliMuonsHFHeader::FillHistosDimuRD(TList *list, AliDimuInfoStoreRD* const dimuStoreRD)
{
  // fill histograms for dimuon according to evnet & dimuon candidates selection cuts

  if (!list) return;
  if (!AliMuonsHFHeader::IsEventSelected()) return;
  if (!dimuStoreRD->DimuSelection()) return;

  Int_t theDimu = 0;
  if (dimuStoreRD->Charge()==0)     theDimu = 1;
  else if (dimuStoreRD->Charge()>0) theDimu = 2;

  const Int_t nHistos = 3;
  Double_t dist[nHistos] = {dimuStoreRD->Momentum().Mag(),
                            dimuStoreRD->Momentum().Pt(),
                            dimuStoreRD->InvM()};
  for (Int_t i=nHistos; i--;) ((TH1F*)list->At(theDimu*nHistos+i))->Fill(dist[i]);
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::FillHistosMuonMC(TList *list, AliMuonInfoStoreMC* const muonStoreMC)
{
  // fill histograms for single muon with MC info

  Int_t srcMuon = muonStoreMC->MuonSource(); 
  this->FillHistosMuonRD(&list[6], (AliMuonInfoStoreRD*)muonStoreMC);
  if (srcMuon<0) this->FillHistosMuonRD(&list[5], (AliMuonInfoStoreRD*)muonStoreMC);
  else this->FillHistosMuonRD(&list[srcMuon], (AliMuonInfoStoreRD*)muonStoreMC);
  return;
}

//_____________________________________________________________________________
void AliMuonsHFHeader::FillHistosDimuMC(TList *list, AliDimuInfoStoreMC* const dimuStoreMC)
{
  // fill histograms for dimuon with MC info

  Int_t srcDimu = dimuStoreMC->DimuSource();
  this->FillHistosDimuRD(&list[6], (AliDimuInfoStoreRD*)dimuStoreMC);
  this->FillHistosDimuRD(&list[srcDimu], (AliDimuInfoStoreRD*)dimuStoreMC);
  return;
}
Commit	Line	Data
fd1d0cb9	1	/**************************************************************************
	2	* Copyright(c) 1998-2006, 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	/////////////////////////////////////////////////////////////
	17	//
	18	// class used to extract and store info at event level
	19	//
	20	// Author: X-M. Zhang, zhang@clermont.in2p3.fr
	21	// zhangxm@iopp.ccnu.edu.cn
	22	/////////////////////////////////////////////////////////////
	23
4292b3b6	24	#include <TMath.h>
fd1d0cb9	25	#include <TH1.h>
	26	#include <TH1F.h>
	27	#include <TH2F.h>
	28	#include <TList.h>
	29	#include <TObjArray.h>
	30	#include <TObjString.h>
4292b3b6	31
fd1d0cb9	32	#include "AliTriggerAnalysis.h"
	33	#include "AliBackgroundSelection.h"
	34	#include "AliAODVertex.h"
	35	#include "AliAODEvent.h"
	36	#include "AliESDEvent.h"
	37	#include "AliMultiplicity.h"
	38	#include "AliMuonInfoStoreRD.h"
	39	#include "AliMuonInfoStoreMC.h"
	40	#include "AliDimuInfoStoreRD.h"
	41	#include "AliDimuInfoStoreMC.h"
4292b3b6	42	#include "AliMuonsHFHeader.h"
4292b3b6	43
fd1d0cb9	44	class TNamed;
	45	class AliESDVertex;
	46
4292b3b6	47	ClassImp(AliMuonsHFHeader)
4292b3b6	48
fd1d0cb9	49	const TString AliMuonsHFHeader::fgkStdBranchName("MuEvsH");
	50	Int_t AliMuonsHFHeader::fgAnaMode = 0;
	51	Bool_t AliMuonsHFHeader::fgIsMC = kFALSE;
	52	Bool_t AliMuonsHFHeader::fgIsEventSelected = kFALSE;
	53	Double_t AliMuonsHFHeader::fgCuts[3] = { -999999., 999999., 999999.};
	54
4292b3b6	55	//_____________________________________________________________________________
	56	AliMuonsHFHeader::AliMuonsHFHeader() :
	57	TNamed(),
	58	fTriggerMask(0),
fd1d0cb9	59	fFiredTrigger(0),
	60	fNFiredTrigger(0),
	61	fIsPhysicsTriggered(kFALSE),
	62	fIsPhysicsAccepted(kFALSE),
	63	fEventType(0),
	64	fUnrecoVertex(kFALSE),
4292b3b6	65	fNContributors(0),
fd1d0cb9	66	fUnrecoVtxSPD(kFALSE),
	67	fNContributorsSPD(0),
	68	fNTrackletsSPD(0),
	69	fCentrality(0.)
4292b3b6	70	{
	71	//
	72	// default constructor
	73	//
fd1d0cb9	74	for (Int_t i=3; i--;) fVtx[i] = 0.;
	75	for (Int_t i=3; i--;) fVtxSPD[i] = 0.;
	76	}
	77
	78	//_____________________________________________________________________________
	79	AliMuonsHFHeader::AliMuonsHFHeader(const AliMuonsHFHeader &src) :
	80	TNamed(),
	81	fTriggerMask(src.fTriggerMask),
	82	fFiredTrigger(src.fFiredTrigger),
	83	fNFiredTrigger(src.fNFiredTrigger),
	84	fIsPhysicsTriggered(src.fIsPhysicsTriggered),
	85	fIsPhysicsAccepted(src.fIsPhysicsAccepted),
	86	fEventType(src.fEventType),
	87	fUnrecoVertex(src.fUnrecoVertex),
	88	fNContributors(src.fNContributors),
	89	fUnrecoVtxSPD(src.fUnrecoVtxSPD),
	90	fNContributorsSPD(src.fNContributorsSPD),
	91	fNTrackletsSPD(src.fNTrackletsSPD),
	92	fCentrality(src.fCentrality)
	93	{
	94	//
	95	// copy constructor
	96	//
	97	for (Int_t i=3; i--;) fVtx[i] = src.fVtx[i];
	98	for (Int_t i=3; i--;) fVtxSPD[i] = src.fVtxSPD[i];
	99	}
	100
	101	//_____________________________________________________________________________
	102	AliMuonsHFHeader& AliMuonsHFHeader::operator=(const AliMuonsHFHeader &src)
	103	{
	104	//
	105	// assignment constructor
	106	//
	107	fTriggerMask = src.fTriggerMask;
	108	fFiredTrigger = src.fFiredTrigger;
	109	fNFiredTrigger = src.fNFiredTrigger;
	110	fIsPhysicsTriggered = src.fIsPhysicsTriggered;
	111	fIsPhysicsAccepted = src.fIsPhysicsAccepted;
	112	fEventType = src.fEventType;
	113	fUnrecoVertex = src.fUnrecoVertex;
	114	fNContributors = src.fNContributors;
	115	fUnrecoVtxSPD = src.fUnrecoVtxSPD;
	116	fNContributorsSPD = src.fNContributorsSPD;
	117	fNTrackletsSPD = src.fNTrackletsSPD;
	118	fCentrality = src.fCentrality;
	119
	120	for (Int_t i=3; i--;) fVtx[i] = src.fVtx[i];
	121	for (Int_t i=3; i--;) fVtxSPD[i] = src.fVtxSPD[i];
	122
	123	return *this;
4292b3b6	124	}
	125
	126	//_____________________________________________________________________________
	127	AliMuonsHFHeader::~AliMuonsHFHeader()
	128	{
	129	//
	130	// default destructor
	131	//
	132	}
	133
	134	//_____________________________________________________________________________
fd1d0cb9	135	void AliMuonsHFHeader::SetEvent(AliAODEvent *event)
4292b3b6	136	{
fd1d0cb9	137	// extract event info from AOD event
	138
	139	fTriggerMask = event->GetTriggerMask();
	140
	141	AliAODVertex *vertex = event->GetPrimaryVertex();
	142	vertex->GetXYZ(fVtx);
	143	fNContributors = vertex->GetNContributors();
	144	fUnrecoVertex = (TMath::Abs(fVtx[0])<1e-6 && TMath::Abs(fVtx[1])<1e-6 &&
	145	TMath::Abs(fVtx[2])<1e-6);
	146	this->SetTitle(vertex->GetTitle());
	147
	148	this->SetFiredTrigger(event->GetFiredTriggerClasses());
	149	this->EventSelection();
4292b3b6	150	return;
	151	}
	152
	153	//_____________________________________________________________________________
fd1d0cb9	154	void AliMuonsHFHeader::SetEvent(AliESDEvent *event)
4292b3b6	155	{
fd1d0cb9	156	// extract event info from ESD event
	157	// SPD vertex and Physics event selection are implimented
	158
	159	fTriggerMask = event->GetTriggerMask();
	160
	161	const AliESDVertex *vertex = event->GetPrimaryVertex();
	162	vertex->GetXYZ(fVtx);
4292b3b6	163	fNContributors = vertex->GetNContributors();
fd1d0cb9	164	fUnrecoVertex = (TMath::Abs(fVtx[0])<1e-6 && TMath::Abs(fVtx[1])<1e-6 &&
	165	TMath::Abs(fVtx[2])<1e-6);
	166	this->SetTitle(vertex->GetTitle());
	167
	168	const AliESDVertex *vtxSPD = event->GetPrimaryVertexSPD();
	169	vtxSPD->GetXYZ(fVtxSPD);
	170	fNContributorsSPD = vtxSPD->GetNContributors();
	171	fUnrecoVtxSPD = (TMath::Abs(fVtxSPD[0])<1e-6 && TMath::Abs(fVtxSPD[1])<1e-6 &&
	172	TMath::Abs(fVtxSPD[2])<1e-6);
	173	fNTrackletsSPD = event->GetMultiplicity()->GetNumberOfTracklets();
	174
	175	this->SetFiredTrigger(event->GetFiredTriggerClasses());
	176	this->PhysicsTriggerAna(event);
	177	this->EventSelection();
	178	return;
	179	}
	180
	181	//_____________________________________________________________________________
	182	void AliMuonsHFHeader::PhysicsTriggerAna(const AliESDEvent *esd)
	183	{
	184	// ESD event trigger condition analysis
	185	// according to the method in $ALICE_ROOT/ANALYSIS/AliPhysicsSelection.cxx
	186
	187	fEventType = esd->GetHeader()->GetEventType();
	188	fIsPhysicsTriggered = kFALSE;
	189	fIsPhysicsAccepted = kFALSE;
	190
	191	AliTriggerAnalysis *triggerAna = new AliTriggerAnalysis();
	192	triggerAna->SetAnalyzeMC(fgIsMC);
	193	triggerAna->SetSPDGFOThreshhold(1);
	194
	195	Int_t triggerHW = triggerAna->SPDFiredChips(esd, 1);
	196	Bool_t isFiredv0A = triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0A);
	197	Bool_t isFiredv0C = triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0C);
	198	if (triggerHW==0 && !isFiredv0A && !isFiredv0C) {
	199	delete triggerAna;
	200	triggerAna = 0;
	201	return;
	202	}
	203
	204	Bool_t triggerBG = (triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0ABG) \|\|
	205	triggerAna->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0CBG));
	206	Int_t isFiredSPD = triggerAna->SPDFiredChips(esd, 0);
	207	Bool_t triggerFD = ((isFiredSPD>1) \|\| (isFiredSPD>0 && (isFiredv0A \|\| isFiredv0C)) \|\| (isFiredv0A && isFiredv0C));
	208	if ((!triggerBG) && triggerFD) fIsPhysicsTriggered = kTRUE;
	209	delete triggerAna;
	210	triggerAna = 0;
	211
	212	if (fIsPhysicsTriggered) {
	213	AliBackgroundSelection *bkgId = new AliBackgroundSelection();
	214	fIsPhysicsAccepted = bkgId->IsSelected(const_cast<AliESDEvent*>(esd));
	215	delete bkgId; bkgId=0;
	216	}
	217	return;
	218	}
	219
	220	//_____________________________________________________________________________
	221	void AliMuonsHFHeader::SetFiredTrigger(TString trigger)
	222	{
	223	// get info of fired trigger classes of event
	224
	225	fFiredTrigger = trigger;
	226	TObjArray *tokens = trigger.Tokenize(" ");
	227	fNFiredTrigger = tokens->GetEntries();
228	return;
229	}
230
231	//_____________________________________________________________________________
232	Bool_t AliMuonsHFHeader::IsTriggerFired(TString trigger)
233	{
234	// check whether the trigger class "trigger" is fired in this event
235
236	if (fNFiredTrigger<=0) return kFALSE;
237
238	TObjArray *tokens = fFiredTrigger.Tokenize(" ");
239	for (Int_t i=fNFiredTrigger; i--;) {
240	TString fired = ((TObjString*)tokens->At(i))->String();
241	if (fired.CompareTo(trigger)==0) return kTRUE;
242	}
243	return kFALSE;
244	}
245
246	//_____________________________________________________________________________
247	void AliMuonsHFHeader::EventSelection(TString triggerName)
248	{
249	// select event according to the "triggerName" & event selection cuts
250
251	fgIsEventSelected = kFALSE;
252	if (!this->IsPhysicsAccepted()) return;
253	if (!fgIsMC && !this->IsTriggerFired(triggerName)) return;
254	this->EventSelection();
255	return;
256	}
257
258	//_____________________________________________________________________________
259	void AliMuonsHFHeader::EventSelection()
260	{
261	// select event according to the event selection cuts
262
263	fgIsEventSelected = kFALSE;
264	if (this->NVtxContributorsSPD()<fgCuts[0]) return;
265	if (TMath::Abs(this->VzSPD())>fgCuts[1]) return;
266	if (this->VtSPD()>fgCuts[2]) return;
267	fgIsEventSelected = kTRUE;
268	return;
269	}
270
271	//_____________________________________________________________________________
272	void AliMuonsHFHeader::CreateHistograms(TList listEvent, TList listMuon, TList *listDimu)
273	{
274	// create output histos of muon analysis according to the analysis mode & MC flag
275
276	this->CreateHistosEventH(listEvent);
277	if (fgIsMC) {
278	if (fgAnaMode!=2) this->CreateHistosMuonMC(listMuon);
279	if (fgAnaMode!=1) this->CreateHistosDimuMC(listDimu);
280	} else {
281	if (fgAnaMode!=2) this->CreateHistosMuonRD(listMuon);
282	if (fgAnaMode!=1) this->CreateHistosDimuRD(listDimu);
283	}
284	return;
285	}
286
287	//_____________________________________________________________________________
288	void AliMuonsHFHeader::CreateHistosEventH(TList *list)
289	{
290	// create histograms at event level
291
292	if (!list) list = new TList();
293	list->SetOwner();
294	Bool_t oldStatus = TH1::AddDirectoryStatus();
295	TH1::AddDirectory(kFALSE);
296
297	const Int_t nHistos = 5;
298	char *name[nHistos] = { "hVztSPD", "hVzxSPD", "hVzySPD", "hVzNcontr", "hVtNcontr" };
299	Int_t nbinsX[nHistos] = { 800 , 800 , 800 , 800 , 400 };
300	Double_t xlow[nHistos] = { -40., -40., -40., -40. , 0. };
301	Double_t xup[nHistos] = { 40., 40., 40., 40. , 4. };
302	Int_t nbinsY[nHistos] = { 400 , 600 , 600 , 202 , 202 };
303	Double_t ylow[nHistos] = { 0., -3., -3., -2.5, -2.5 };
304	Double_t yup[nHistos] = { 4., 3., 3., 199.5, 199.5 };
305
306	TH2F *histo = 0;
307	for (Int_t i=0; i<nHistos; i++) {
308	histo = new TH2F(name[i], name[i], nbinsX[i], xlow[i], xup[i], nbinsY[i], ylow[i], yup[i]);
309	histo->Sumw2(); list->AddAt(histo, i); histo = 0;
310	}
311
312	TH1::AddDirectory(oldStatus);
313	return;
314	}
315
316	//_____________________________________________________________________________
317	void AliMuonsHFHeader::CreateHistosMuonRD(TList *list)
318	{
319	// create histograms for single muon
320
321	if (!list) list = new TList();
322	list->SetOwner();
323	Bool_t oldStatus = TH1::AddDirectoryStatus();
324	TH1::AddDirectory(kFALSE);
325
326	const Int_t nHistos = 8;
327	char *name[nHistos] = { "hP", "hPt", "hEta", "hDCA", "hTrg", "hCharge", "hUnfVtx" , "hEtaTimesDCA"};
328	Int_t nbins[nHistos] = { 1500 , 300 , 100 , 500 , 4 , 3 , 80 , 10000 };
329	Double_t xlow[nHistos] = { 0., 0., -10., 0., -0.5, -1.5, -40. , 0.};
330	Double_t xup[nHistos] = { 150., 30., 0., 500., 3.5, 1.5, 40. , 1000.};
331
332	TH1F *histo = 0;
333	for (Int_t i=0; i<nHistos; i++) {
334	histo = new TH1F(name[i], name[i], nbins[i], xlow[i], xup[i]);
335	histo->Sumw2(); list->AddAt(histo, i); histo = 0;
336	}
337
338	TH1::AddDirectory(oldStatus);
339	return;
340	}
341
342	//_____________________________________________________________________________
343	void AliMuonsHFHeader::CreateHistosDimuRD(TList *list)
344	{
345	// create histograms for dimuon
346
347	if (!list) list = new TList();
348	list->SetOwner();
349	Bool_t oldStatus = TH1::AddDirectoryStatus();
350	TH1::AddDirectory(kFALSE);
351
352	TH1F *histo = 0;
353	const Int_t nHistos = 3;
354	char *name[nHistos] = { "hP", "hPt", "hInvM" };
355	Int_t nbins[nHistos] = { 1500 , 300 , 300 };
356	Double_t xlow[nHistos] = { 0., 0., 0. };
357	Double_t xup[nHistos] = { 150., 30., 30. };
358	char *dimuName[3] = {"DimuNN", "DimuNP", "DimuPP"};
359	for (Int_t i=0; i<3; i++) {
360	for (Int_t j=0; j<nHistos; j++) {
361	histo = new TH1F(Form("%s_%s",name[j],dimuName[i]), name[j], nbins[j], xlow[j], xup[j]);
362	histo->Sumw2(); list->AddAt(histo,i*nHistos+j); histo = 0;
363	}
364	}
365
366	TH1::AddDirectory(oldStatus);
367	return;
368	}
369
370	//_____________________________________________________________________________
371	void AliMuonsHFHeader::CreateHistosMuonMC(TList *list)
372	{
373	// create histograms for single muon with MC info
374
375	if (!list) list = new TList[AliMuonInfoStoreMC::NSources()];
376	for (Int_t i=AliMuonInfoStoreMC::NSources(); i--;) this->CreateHistosMuonRD(&list[i]);
377	return;
378	}
379
380	//_____________________________________________________________________________
381	void AliMuonsHFHeader::CreateHistosDimuMC(TList *list)
382	{
383	// create histograms for dimuon with MC info
384
385	if (!list) list = new TList[AliDimuInfoStoreMC::NSources()];
386	for (Int_t i=AliDimuInfoStoreMC::NSources(); i--;) this->CreateHistosDimuRD(&list[i]);
387	return;
388	}
389
390	//_____________________________________________________________________________
391	void AliMuonsHFHeader::FillHistosEventH(TList *list)
392	{
393	// fill histograms at event level according to event selection cuts
394
395	if (!list) return;
396	if (!AliMuonsHFHeader::IsEventSelected()) return;
397
398	const Int_t nHistos = 5;
399	Double_t vz = this->VzSPD();
400	Double_t vt = this->VtSPD();
401	Int_t nc = this->NVtxContributorsSPD();
402	Double_t distX[nHistos] = { vz, vz, vz, vz, vt };
403	Double_t distY[nHistos] = { vt, this->VxSPD(), this->VySPD(), nc, nc };
404	for (Int_t i=nHistos; i--;) ((TH2F*)list->At(i))->Fill(distX[i], distY[i]);
405	return;
406	}
407
408	//_____________________________________________________________________________
409	void AliMuonsHFHeader::FillHistosMuonRD(TList list, AliMuonInfoStoreRD const muonStoreRD)
410	{
411	// fill histograms for single muon according to event & muon track selection cuts
412
413	if (!list) return;
414	if (!AliMuonsHFHeader::IsEventSelected()) return;
415	if (!muonStoreRD->MuonSelection()) return;
416
417	const Int_t nHistos = 8;
418	Double_t dist[nHistos] = {muonStoreRD->Momentum().Mag(),
419	muonStoreRD->Momentum().Pt(),
420	muonStoreRD->Momentum().Eta(),
421	muonStoreRD->DCA(),
422	muonStoreRD->MatchTrigger(),
423	muonStoreRD->Charge(),
424	this->VzSPD(),
425	muonStoreRD->DCA()TMath::Exp(-2.muonStoreRD->Momentum().Eta())/1000.};
426	for (Int_t i=nHistos; i--;) ((TH1F*)list->At(i))->Fill(dist[i]);
427	return;
428	}
429
430	//_____________________________________________________________________________
431	void AliMuonsHFHeader::FillHistosDimuRD(TList list, AliDimuInfoStoreRD const dimuStoreRD)
432	{
433	// fill histograms for dimuon according to evnet & dimuon candidates selection cuts
434
435	if (!list) return;
436	if (!AliMuonsHFHeader::IsEventSelected()) return;
437	if (!dimuStoreRD->DimuSelection()) return;
438
439	Int_t theDimu = 0;
440	if (dimuStoreRD->Charge()==0) theDimu = 1;
441	else if (dimuStoreRD->Charge()>0) theDimu = 2;
442
443	const Int_t nHistos = 3;
444	Double_t dist[nHistos] = {dimuStoreRD->Momentum().Mag(),
445	dimuStoreRD->Momentum().Pt(),
446	dimuStoreRD->InvM()};
447	for (Int_t i=nHistos; i--;) ((TH1F)list->At(theDimunHistos+i))->Fill(dist[i]);
448	return;
449	}
450
451	//_____________________________________________________________________________
452	void AliMuonsHFHeader::FillHistosMuonMC(TList list, AliMuonInfoStoreMC const muonStoreMC)
453	{
454	// fill histograms for single muon with MC info
455
456	Int_t srcMuon = muonStoreMC->MuonSource();
457	this->FillHistosMuonRD(&list[6], (AliMuonInfoStoreRD*)muonStoreMC);
458	if (srcMuon<0) this->FillHistosMuonRD(&list[5], (AliMuonInfoStoreRD*)muonStoreMC);
459	else this->FillHistosMuonRD(&list[srcMuon], (AliMuonInfoStoreRD*)muonStoreMC);
460	return;
461	}
462
463	//_____________________________________________________________________________
464	void AliMuonsHFHeader::FillHistosDimuMC(TList list, AliDimuInfoStoreMC const dimuStoreMC)
465	{
466	// fill histograms for dimuon with MC info
467
468	Int_t srcDimu = dimuStoreMC->DimuSource();
469	this->FillHistosDimuRD(&list[6], (AliDimuInfoStoreRD*)dimuStoreMC);
470	this->FillHistosDimuRD(&list[srcDimu], (AliDimuInfoStoreRD*)dimuStoreMC);
4292b3b6	471	return;
4292b3b6	472	}