[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis2 / AliFMDMCCorrector.cxx

// 
// This class calculates the exclusive charged particle density
// in each for the 5 FMD rings. 
//
// Input:
//   - 5 RingHistos objects - each with a number of vertex dependent 
//     2D histograms of the inclusive charge particle density 
//
// Output:
//   - 5 RingHistos objects - each with a number of vertex dependent 
//     2D histograms of the exclusive charge particle density 
// 
// Corrections used: 
//   - AliFMDCorrSecondaryMap;
//   - AliFMDCorrVertexBias
//   - AliFMDCorrMergingEfficiency
//
#include "AliFMDMCCorrector.h"
#include <AliESDFMD.h>
#include <TAxis.h>
#include <TList.h>
#include <TMath.h>
#include "AliForwardCorrectionManager.h"
#include "AliLog.h"
#include <TH2D.h>
#include <TROOT.h>
#include <TProfile2D.h>

ClassImp(AliFMDMCCorrector)
#if 0
; // For Emacs
#endif 


//____________________________________________________________________
AliFMDMCCorrector::~AliFMDMCCorrector()
{
  // 
  // Destructor 
  //
  if (fComps) fComps->Clear();
  if (fFMD1i) delete fFMD1i;
  if (fFMD2i) delete fFMD2i;
  if (fFMD2o) delete fFMD2o;
  if (fFMD3i) delete fFMD3i;
  if (fFMD3o) delete fFMD3o;
}

//____________________________________________________________________
AliFMDMCCorrector&
AliFMDMCCorrector::operator=(const AliFMDMCCorrector& o)
{
  // 
  // Assignement operator
  // 
  // Parameters:
  //    o Object to assign from 
  // 
  // Return:
  //    Reference to this object
  //
  AliFMDCorrector::operator=(o);

  return *this;
}

//____________________________________________________________________
Bool_t
AliFMDMCCorrector::CorrectMC(AliForwardUtil::Histos& hists,
			       UShort_t                vtxbin)
{
  // 
  // Do the calculations 
  // 
  // Parameters:
  //    hists    Cache of histograms 
  //    vtxBin   Vertex bin 
  // 
  // Return:
  //    true on successs 
  //
  AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();

  UShort_t uvb = vtxbin;
  for (UShort_t d=1; d<=3; d++) { 
    UShort_t nr = (d == 1 ? 1 : 2);
    for (UShort_t q=0; q<nr; q++) { 
      Char_t      r  = (q == 0 ? 'I' : 'O');
      TH2D*       h  = hists.Get(d,r);


      if (fUseSecondaryMap) {
        TH2D*       bg = fcm.GetSecondaryMap()->GetCorrection(d,r,uvb);
        if (!bg) {
          AliWarning(Form("No secondary correction for FMDM%d%c in vertex bin %d",
              d, r, uvb));
          continue;
        }
        // Divide by primary/total ratio
        h->Divide(bg);
      }
      if (fUseVertexBias) {
        TH2D*       ef = fcm.GetVertexBias()->GetCorrection(r, uvb);
        if (!ef) {
          AliWarning(Form("No event vertex bias correction in vertex bin %d",
			uvb));
          continue;
        }
        // Divide by the event selection efficiency
        h->Divide(ef);
      }
    }
  }
  
  return kTRUE;
}

//____________________________________________________________________
void
AliFMDMCCorrector::Init(const TAxis& eAxis)
{
  // 
  // Initialize this object 
  // 
  // Parameters:
  //    etaAxis Eta axis to use 
  //
  AliFMDCorrector::Init(eAxis);

  fFMD1i = Make(1,'I',eAxis);
  fFMD2i = Make(2,'I',eAxis);
  fFMD2o = Make(2,'O',eAxis);
  fFMD3i = Make(3,'I',eAxis);
  fFMD3o = Make(3,'O',eAxis);

  fComps->Add(fFMD1i);
  fComps->Add(fFMD2i);
  fComps->Add(fFMD2o);
  fComps->Add(fFMD3i);
  fComps->Add(fFMD3o);
}

//____________________________________________________________________
TProfile2D*
AliFMDMCCorrector::Make(UShort_t d, Char_t r, 
				const TAxis& axis) const
{
  // 
  // MAke comparison profiles
  // 
  // Parameters:
  //    d     Detector 
  //    r     Ring 
  //    axis  Eta axis 
  // 
  // Return:
  //    Newly allocated profile object
  //
  TProfile2D* ret = new TProfile2D(Form("FMD%d%c_esd_vs_mc", d, r),
				   Form("ESD/MC signal for FMD%d%c", d, r),
				   axis.GetNbins(), 
				   axis.GetXmin(),
				   axis.GetXmax(), 
				   (r == 'I' || r == 'i') ? 20 : 40,
				   0, 2*TMath::Pi());
  ret->GetXaxis()->SetTitle("#eta");
  ret->GetYaxis()->SetTitle("#varphi [degrees]");
  ret->GetZaxis()->SetTitle("#LT primary density ESD/MC#GT");
  ret->SetDirectory(0);
  return ret;
}
//____________________________________________________________________
void
AliFMDMCCorrector::Fill(UShort_t d, Char_t r, TH2* esd, TH2* mc)
{
  // 
  // Fill comparison profiles
  // 
  // Parameters:
  //    d    Detector 
  //    r    Ring 
  //    esd  ESD histogram
  //    mc   MC histogram
  //
  if (!esd || !mc) return;
  TProfile2D* p = 0;
  switch (d) { 
  case 1:  p = fFMD1i;                                   break;
  case 2:  p = (r == 'I' || r == 'i' ? fFMD2i : fFMD2o); break;
  case 3:  p = (r == 'I' || r == 'i' ? fFMD3i : fFMD3o); break;
  }
  if (!p) return;

  for (Int_t iEta = 1; iEta <= esd->GetNbinsX(); iEta++) { 
    Double_t eta = esd->GetXaxis()->GetBinCenter(iEta);
    for (Int_t iPhi = 1; iPhi <= esd->GetNbinsY(); iPhi++) { 
      Double_t phi  = esd->GetYaxis()->GetBinCenter(iPhi);
      Double_t mEsd = esd->GetBinContent(iEta,iPhi);
      Double_t mMc  = mc->GetBinContent(iEta,iPhi);
      
      p->Fill(eta, phi, (mMc > 0 ? mEsd / mMc : 0));
    }
  }
}

//____________________________________________________________________
Bool_t
AliFMDMCCorrector::CompareResults(AliForwardUtil::Histos& esd,
					  AliForwardUtil::Histos& mc)
{
  // 
  // Compare the result of analysing the ESD for 
  // the inclusive charged particle density to analysing 
  // MC truth 
  // 
  // Parameters:
  //    esd 
  //    mc 
  // 
  // Return:
  //   true 
  //
  Fill(1, 'I', esd.Get(1,'I'), mc.Get(1,'I'));
  Fill(2, 'I', esd.Get(2,'I'), mc.Get(2,'I'));
  Fill(2, 'O', esd.Get(2,'O'), mc.Get(2,'O'));
  Fill(3, 'I', esd.Get(3,'I'), mc.Get(3,'I'));
  Fill(3, 'O', esd.Get(3,'O'), mc.Get(3,'O'));

  return kTRUE;
}

//____________________________________________________________________
void
AliFMDMCCorrector::DefineOutput(TList* dir)
{
  // 
  // Output diagnostic histograms to directory 
  // 
  // Parameters:
  //    dir List to write in
  //  
  AliFMDCorrector::DefineOutput(dir);
  TList* d = static_cast<TList*>(dir->FindObject(GetName()));

  fComps = new TList;
  fComps->SetName("esd_mc_comparison");
  d->Add(fComps);
}

//____________________________________________________________________
//
// EOF
//
Commit	Line	Data
7984e5f7	1	//
	2	// This class calculates the exclusive charged particle density
	3	// in each for the 5 FMD rings.
	4	//
	5	// Input:
	6	// - 5 RingHistos objects - each with a number of vertex dependent
	7	// 2D histograms of the inclusive charge particle density
	8	//
	9	// Output:
	10	// - 5 RingHistos objects - each with a number of vertex dependent
	11	// 2D histograms of the exclusive charge particle density
	12	//
	13	// Corrections used:
	14	// - AliFMDCorrSecondaryMap;
	15	// - AliFMDCorrVertexBias
	16	// - AliFMDCorrMergingEfficiency
	17	//
72cc12cd	18	#include "AliFMDMCCorrector.h"
886dd429	19	#include <AliESDFMD.h>
	20	#include <TAxis.h>
	21	#include <TList.h>
	22	#include <TMath.h>
	23	#include "AliForwardCorrectionManager.h"
886dd429	24	#include "AliLog.h"
	25	#include <TH2D.h>
	26	#include <TROOT.h>
	27	#include <TProfile2D.h>
886dd429	28
72cc12cd	29	ClassImp(AliFMDMCCorrector)
886dd429	30	#if 0
	31	; // For Emacs
	32	#endif
	33
	34
	35	//____________________________________________________________________
72cc12cd	36	AliFMDMCCorrector::~AliFMDMCCorrector()
886dd429	37	{
7984e5f7	38	//
	39	// Destructor
	40	//
886dd429	41	if (fComps) fComps->Clear();
	42	if (fFMD1i) delete fFMD1i;
	43	if (fFMD2i) delete fFMD2i;
	44	if (fFMD2o) delete fFMD2o;
	45	if (fFMD3i) delete fFMD3i;
	46	if (fFMD3o) delete fFMD3o;
	47	}
	48
	49	//____________________________________________________________________
72cc12cd	50	AliFMDMCCorrector&
72cc12cd	51	AliFMDMCCorrector::operator=(const AliFMDMCCorrector& o)
886dd429	52	{
7984e5f7	53	//
	54	// Assignement operator
	55	//
	56	// Parameters:
	57	// o Object to assign from
	58	//
	59	// Return:
	60	// Reference to this object
	61	//
72cc12cd	62	AliFMDCorrector::operator=(o);
886dd429	63
	64	return *this;
	65	}
	66
	67	//____________________________________________________________________
	68	Bool_t
72cc12cd	69	AliFMDMCCorrector::CorrectMC(AliForwardUtil::Histos& hists,
886dd429	70	UShort_t vtxbin)
886dd429	71	{
7984e5f7	72	//
	73	// Do the calculations
	74	//
	75	// Parameters:
	76	// hists Cache of histograms
	77	// vtxBin Vertex bin
	78	//
	79	// Return:
	80	// true on successs
	81	//
886dd429	82	AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
	83
	84	UShort_t uvb = vtxbin;
	85	for (UShort_t d=1; d<=3; d++) {
	86	UShort_t nr = (d == 1 ? 1 : 2);
	87	for (UShort_t q=0; q<nr; q++) {
	88	Char_t r = (q == 0 ? 'I' : 'O');
	89	TH2D* h = hists.Get(d,r);
7ec4d843	90
	91
	92	if (fUseSecondaryMap) {
	93	TH2D* bg = fcm.GetSecondaryMap()->GetCorrection(d,r,uvb);
	94	if (!bg) {
	95	AliWarning(Form("No secondary correction for FMDM%d%c in vertex bin %d",
	96	d, r, uvb));
	97	continue;
	98	}
	99	// Divide by primary/total ratio
	100	h->Divide(bg);
886dd429	101	}
7ec4d843	102	if (fUseVertexBias) {
	103	TH2D* ef = fcm.GetVertexBias()->GetCorrection(r, uvb);
	104	if (!ef) {
	105	AliWarning(Form("No event vertex bias correction in vertex bin %d",
886dd429	106	uvb));
7ec4d843	107	continue;
	108	}
	109	// Divide by the event selection efficiency
	110	h->Divide(ef);
886dd429	111	}
886dd429	112	}
	113	}
	114
	115	return kTRUE;
	116	}
	117
	118	//____________________________________________________________________
	119	void
72cc12cd	120	AliFMDMCCorrector::Init(const TAxis& eAxis)
886dd429	121	{
7984e5f7	122	//
	123	// Initialize this object
	124	//
	125	// Parameters:
	126	// etaAxis Eta axis to use
	127	//
7ec4d843	128	AliFMDCorrector::Init(eAxis);
7ec4d843	129
886dd429	130	fFMD1i = Make(1,'I',eAxis);
	131	fFMD2i = Make(2,'I',eAxis);
	132	fFMD2o = Make(2,'O',eAxis);
	133	fFMD3i = Make(3,'I',eAxis);
	134	fFMD3o = Make(3,'O',eAxis);
	135
	136	fComps->Add(fFMD1i);
	137	fComps->Add(fFMD2i);
	138	fComps->Add(fFMD2o);
	139	fComps->Add(fFMD3i);
	140	fComps->Add(fFMD3o);
	141	}
	142
	143	//____________________________________________________________________
	144	TProfile2D*
72cc12cd	145	AliFMDMCCorrector::Make(UShort_t d, Char_t r,
886dd429	146	const TAxis& axis) const
886dd429	147	{
7984e5f7	148	//
	149	// MAke comparison profiles
	150	//
	151	// Parameters:
	152	// d Detector
	153	// r Ring
	154	// axis Eta axis
	155	//
	156	// Return:
	157	// Newly allocated profile object
	158	//
886dd429	159	TProfile2D* ret = new TProfile2D(Form("FMD%d%c_esd_vs_mc", d, r),
	160	Form("ESD/MC signal for FMD%d%c", d, r),
	161	axis.GetNbins(),
	162	axis.GetXmin(),
	163	axis.GetXmax(),
	164	(r == 'I' \|\| r == 'i') ? 20 : 40,
	165	0, 2*TMath::Pi());
	166	ret->GetXaxis()->SetTitle("#eta");
	167	ret->GetYaxis()->SetTitle("#varphi [degrees]");
	168	ret->GetZaxis()->SetTitle("#LT primary density ESD/MC#GT");
	169	ret->SetDirectory(0);
	170	return ret;
	171	}
	172	//____________________________________________________________________
	173	void
72cc12cd	174	AliFMDMCCorrector::Fill(UShort_t d, Char_t r, TH2* esd, TH2* mc)
886dd429	175	{
7984e5f7	176	//
	177	// Fill comparison profiles
	178	//
	179	// Parameters:
	180	// d Detector
	181	// r Ring
	182	// esd ESD histogram
	183	// mc MC histogram
	184	//
886dd429	185	if (!esd \|\| !mc) return;
	186	TProfile2D* p = 0;
	187	switch (d) {
	188	case 1: p = fFMD1i; break;
	189	case 2: p = (r == 'I' \|\| r == 'i' ? fFMD2i : fFMD2o); break;
	190	case 3: p = (r == 'I' \|\| r == 'i' ? fFMD3i : fFMD3o); break;
	191	}
	192	if (!p) return;
	193
	194	for (Int_t iEta = 1; iEta <= esd->GetNbinsX(); iEta++) {
	195	Double_t eta = esd->GetXaxis()->GetBinCenter(iEta);
	196	for (Int_t iPhi = 1; iPhi <= esd->GetNbinsY(); iPhi++) {
	197	Double_t phi = esd->GetYaxis()->GetBinCenter(iPhi);
	198	Double_t mEsd = esd->GetBinContent(iEta,iPhi);
	199	Double_t mMc = mc->GetBinContent(iEta,iPhi);
	200
	201	p->Fill(eta, phi, (mMc > 0 ? mEsd / mMc : 0));
	202	}
	203	}
	204	}
	205
	206	//____________________________________________________________________
	207	Bool_t
72cc12cd	208	AliFMDMCCorrector::CompareResults(AliForwardUtil::Histos& esd,
886dd429	209	AliForwardUtil::Histos& mc)
886dd429	210	{
7984e5f7	211	//
	212	// Compare the result of analysing the ESD for
	213	// the inclusive charged particle density to analysing
	214	// MC truth
	215	//
	216	// Parameters:
	217	// esd
	218	// mc
	219	//
	220	// Return:
	221	// true
	222	//
886dd429	223	Fill(1, 'I', esd.Get(1,'I'), mc.Get(1,'I'));
	224	Fill(2, 'I', esd.Get(2,'I'), mc.Get(2,'I'));
	225	Fill(2, 'O', esd.Get(2,'O'), mc.Get(2,'O'));
	226	Fill(3, 'I', esd.Get(3,'I'), mc.Get(3,'I'));
	227	Fill(3, 'O', esd.Get(3,'O'), mc.Get(3,'O'));
	228
	229	return kTRUE;
	230	}
	231
	232	//____________________________________________________________________
	233	void
72cc12cd	234	AliFMDMCCorrector::DefineOutput(TList* dir)
886dd429	235	{
7984e5f7	236	//
	237	// Output diagnostic histograms to directory
	238	//
	239	// Parameters:
	240	// dir List to write in
	241	//
72cc12cd	242	AliFMDCorrector::DefineOutput(dir);
886dd429	243	TList* d = static_cast<TList*>(dir->FindObject(GetName()));
	244
	245	fComps = new TList;
	246	fComps->SetName("esd_mc_comparison");
	247	d->Add(fComps);
	248	}
	249
	250	//____________________________________________________________________
	251	//
	252	// EOF
	253	//
	254
	255
	256