[u/mrichter/AliRoot.git] / PWG2 / FEMTOSCOPY / AliFemtoUser / AliFemtoCutMonitorParticleEtCorr.cxx

////////////////////////////////////////////////////////////////////////////////
//                                                                            //
// AliFemtoCutMonitorParticleEtCorr - the cut monitor for particles           //
// which saves particles' et histogram and makes the bin-by-bin correlation   //
//                                                                            //
// Author: Adam.Kisiel@cern.ch                                                //
//                                                                            //
////////////////////////////////////////////////////////////////////////////////
#include "AliFemtoCutMonitorParticleEtCorr.h"
#include "AliFemtoModelHiddenInfo.h"
#include <TH1D.h>
#include <TH2D.h>
#include <TList.h>
#include <TMath.h>

AliFemtoCutMonitorParticleEtCorr::AliFemtoCutMonitorParticleEtCorr():
  AliFemtoCutMonitor(),
  fPhiBins(60),
  fPtPerPhi(0),     
  fPtCovPerPhi(0),
  fPtMultPerPhi(0),
  fNEventsProcessed(0)
{
  // Default constructor
  for (int ib=0; ib<200; ib++) {
    fPtSumEvent[ib] = 0; 
    fMultSumEvent[ib] = 0;
  }
}

AliFemtoCutMonitorParticleEtCorr::AliFemtoCutMonitorParticleEtCorr(const char *aName, int aPhiBins):
  AliFemtoCutMonitor(),
  fPhiBins(aPhiBins),
  fPtPerPhi(0),     
  fPtCovPerPhi(0),
  fPtMultPerPhi(0),
  fNEventsProcessed(0)
{
  // Normal constructor
  char name[200];
  snprintf(name, 200, "EtCorrAvgPt%s", aName);
  fPtPerPhi = new TH1D(name, "Average Pt Per Phi", aPhiBins, -0.5, aPhiBins-0.5);
  snprintf(name, 200, "EtCorrMult%s", aName);
  fPtMultPerPhi = new TH2D(name, "Multiplicity Per Phi", aPhiBins, -0.5, aPhiBins-0.5, aPhiBins, -0.5, aPhiBins-0.5);
  snprintf(name, 200, "EtCorrAvgPtCov%s", aName);
  fPtCovPerPhi = new TH2D(name, "Covariance of Average Pt Per Phi", aPhiBins, -0.5, aPhiBins-0.5, aPhiBins, -0.5, aPhiBins-0.5);

  fPtPerPhi->Sumw2();
  fPtCovPerPhi->Sumw2();
  fPtMultPerPhi->Sumw2();
  fPhiBins = aPhiBins;

  for (int ib=0; ib<200; ib++) {
    fPtSumEvent[ib] = 0; 
    fMultSumEvent[ib] = 0;
  }

}

AliFemtoCutMonitorParticleEtCorr::AliFemtoCutMonitorParticleEtCorr(const AliFemtoCutMonitorParticleEtCorr &aCut):
  AliFemtoCutMonitor(),
  fPhiBins(0),
  fPtPerPhi(0),     
  fPtCovPerPhi(0),
  fPtMultPerPhi(0),
  fNEventsProcessed(0)
{
  // copy constructor
  if (fPtCovPerPhi) delete fPtCovPerPhi;
  fPtCovPerPhi = new TH2D(*aCut.fPtCovPerPhi);
  if (fPtPerPhi) delete fPtPerPhi;
  fPtPerPhi = new TH1D(*aCut.fPtPerPhi);
  if (fPtMultPerPhi) delete fPtMultPerPhi;
  fPtMultPerPhi = new TH2D(*aCut.fPtMultPerPhi);
  fPhiBins = aCut.fPhiBins;
  fNEventsProcessed = aCut.fNEventsProcessed;

  for (int ib=0; ib<200; ib++) {
    fPtSumEvent[ib] = 0; 
    fMultSumEvent[ib] = 0;
  }
}

AliFemtoCutMonitorParticleEtCorr::~AliFemtoCutMonitorParticleEtCorr()
{
  // Destructor
  delete fPtPerPhi;
  delete fPtMultPerPhi;
  delete fPtCovPerPhi;
}

AliFemtoCutMonitorParticleEtCorr& AliFemtoCutMonitorParticleEtCorr::operator=(const AliFemtoCutMonitorParticleEtCorr& aCut)
{
  // assignment operator
  if (this == &aCut) 
    return *this;

  if (fPtCovPerPhi) delete fPtCovPerPhi;
  fPtCovPerPhi = new TH2D(*aCut.fPtCovPerPhi);
  if (fPtPerPhi) delete fPtPerPhi;
  fPtPerPhi = new TH1D(*aCut.fPtPerPhi);
  if (fPtMultPerPhi) delete fPtMultPerPhi;
  fPtMultPerPhi = new TH2D(*aCut.fPtMultPerPhi);
  fPhiBins = aCut.fPhiBins;
  fNEventsProcessed = aCut.fNEventsProcessed;
  
  return *this;
}

AliFemtoString AliFemtoCutMonitorParticleEtCorr::Report(){ 
  // Prepare report from the execution
  string stemp = "*** AliFemtoCutMonitorParticleEtCorr report"; 
  AliFemtoString returnThis = stemp;
  return returnThis; 
}

void AliFemtoCutMonitorParticleEtCorr::Fill(const AliFemtoTrack* aTrack)
{
  // Fill in the monitor histograms with the values from the current track
  //  float tEnergy = ::sqrt(aTrack->P().mag2()+fMass*fMass);
  //  float tRapidity = 0.5*::log((tEnergy+aTrack->P().z())/(tEnergy-aTrack->P().z()));
  float tPt = ::sqrt((aTrack->P().x())*(aTrack->P().x())+(aTrack->P().y())*(aTrack->P().y()));
  //  float tEta = -TMath::Log(TMath::Tan(aTrack->P().theta()/2.0));
  float tPhi = aTrack->P().Phi();
  Double_t tPiTwo = TMath::Pi()*2;

  while (tPhi > tPiTwo) tPhi -= tPiTwo;
  while (tPhi < 0) tPhi += tPiTwo;

  int nbin = (int) floor(tPhi * fPhiBins / tPiTwo);
  fPtSumEvent[nbin] += tPt;
  fMultSumEvent[nbin] += 1;
}

void AliFemtoCutMonitorParticleEtCorr::Write()
{
  // Write out the relevant histograms
  fPtPerPhi->Write();
  fPtCovPerPhi->Write();
  fPtMultPerPhi->Write();
}

TList *AliFemtoCutMonitorParticleEtCorr::GetOutputList()
{
  TList *tOutputList = new TList();
  tOutputList->Add(fPtPerPhi);
  tOutputList->Add(fPtCovPerPhi);
  tOutputList->Add(fPtMultPerPhi);

  return tOutputList;
}

void AliFemtoCutMonitorParticleEtCorr::EventBegin(const AliFemtoEvent* aEvent)
{
  if (aEvent)
    for (int iter=0; iter<fPhiBins; iter++) {
      fPtSumEvent[iter] = 0;
      fMultSumEvent[iter] = 0;
    }
}

void AliFemtoCutMonitorParticleEtCorr::EventEnd(const AliFemtoEvent* aEvent)
{
  if (aEvent) {
    for (int ispt=0; ispt<fPhiBins; ispt++) {
      fPtPerPhi->Fill(ispt, fPtSumEvent[ispt]);
      for (int ispt2=0; ispt2<fPhiBins; ispt2++) {
	fPtCovPerPhi->Fill(ispt, ispt2, fPtSumEvent[ispt]*fPtSumEvent[ispt2]);
	fPtMultPerPhi->Fill(ispt, ispt2, fMultSumEvent[ispt]*fMultSumEvent[ispt2]);
      }
    }
    fNEventsProcessed++;
  }
}
Commit	Line	Data
95150c4b	1	////////////////////////////////////////////////////////////////////////////////
	2	// //
	3	// AliFemtoCutMonitorParticleEtCorr - the cut monitor for particles //
	4	// which saves particles' et histogram and makes the bin-by-bin correlation //
	5	// //
	6	// Author: Adam.Kisiel@cern.ch //
	7	// //
	8	////////////////////////////////////////////////////////////////////////////////
	9	#include "AliFemtoCutMonitorParticleEtCorr.h"
	10	#include "AliFemtoModelHiddenInfo.h"
	11	#include <TH1D.h>
	12	#include <TH2D.h>
	13	#include <TList.h>
	14	#include <TMath.h>
	15
	16	AliFemtoCutMonitorParticleEtCorr::AliFemtoCutMonitorParticleEtCorr():
	17	AliFemtoCutMonitor(),
	18	fPhiBins(60),
	19	fPtPerPhi(0),
	20	fPtCovPerPhi(0),
	21	fPtMultPerPhi(0),
	22	fNEventsProcessed(0)
	23	{
	24	// Default constructor
0002d348	25	for (int ib=0; ib<200; ib++) {
	26	fPtSumEvent[ib] = 0;
	27	fMultSumEvent[ib] = 0;
	28	}
95150c4b	29	}
	30
	31	AliFemtoCutMonitorParticleEtCorr::AliFemtoCutMonitorParticleEtCorr(const char *aName, int aPhiBins):
	32	AliFemtoCutMonitor(),
	33	fPhiBins(aPhiBins),
	34	fPtPerPhi(0),
	35	fPtCovPerPhi(0),
	36	fPtMultPerPhi(0),
	37	fNEventsProcessed(0)
	38	{
	39	// Normal constructor
	40	char name[200];
	41	snprintf(name, 200, "EtCorrAvgPt%s", aName);
	42	fPtPerPhi = new TH1D(name, "Average Pt Per Phi", aPhiBins, -0.5, aPhiBins-0.5);
	43	snprintf(name, 200, "EtCorrMult%s", aName);
	44	fPtMultPerPhi = new TH2D(name, "Multiplicity Per Phi", aPhiBins, -0.5, aPhiBins-0.5, aPhiBins, -0.5, aPhiBins-0.5);
	45	snprintf(name, 200, "EtCorrAvgPtCov%s", aName);
	46	fPtCovPerPhi = new TH2D(name, "Covariance of Average Pt Per Phi", aPhiBins, -0.5, aPhiBins-0.5, aPhiBins, -0.5, aPhiBins-0.5);
	47
	48	fPtPerPhi->Sumw2();
	49	fPtCovPerPhi->Sumw2();
	50	fPtMultPerPhi->Sumw2();
	51	fPhiBins = aPhiBins;
0002d348	52
	53	for (int ib=0; ib<200; ib++) {
	54	fPtSumEvent[ib] = 0;
	55	fMultSumEvent[ib] = 0;
	56	}
	57
95150c4b	58	}
	59
	60	AliFemtoCutMonitorParticleEtCorr::AliFemtoCutMonitorParticleEtCorr(const AliFemtoCutMonitorParticleEtCorr &aCut):
	61	AliFemtoCutMonitor(),
	62	fPhiBins(0),
	63	fPtPerPhi(0),
	64	fPtCovPerPhi(0),
	65	fPtMultPerPhi(0),
	66	fNEventsProcessed(0)
	67	{
	68	// copy constructor
	69	if (fPtCovPerPhi) delete fPtCovPerPhi;
	70	fPtCovPerPhi = new TH2D(*aCut.fPtCovPerPhi);
	71	if (fPtPerPhi) delete fPtPerPhi;
	72	fPtPerPhi = new TH1D(*aCut.fPtPerPhi);
	73	if (fPtMultPerPhi) delete fPtMultPerPhi;
	74	fPtMultPerPhi = new TH2D(*aCut.fPtMultPerPhi);
	75	fPhiBins = aCut.fPhiBins;
	76	fNEventsProcessed = aCut.fNEventsProcessed;
0002d348	77
	78	for (int ib=0; ib<200; ib++) {
	79	fPtSumEvent[ib] = 0;
	80	fMultSumEvent[ib] = 0;
	81	}
95150c4b	82	}
	83
	84	AliFemtoCutMonitorParticleEtCorr::~AliFemtoCutMonitorParticleEtCorr()
	85	{
	86	// Destructor
	87	delete fPtPerPhi;
	88	delete fPtMultPerPhi;
	89	delete fPtCovPerPhi;
	90	}
	91
	92	AliFemtoCutMonitorParticleEtCorr& AliFemtoCutMonitorParticleEtCorr::operator=(const AliFemtoCutMonitorParticleEtCorr& aCut)
	93	{
	94	// assignment operator
	95	if (this == &aCut)
	96	return *this;
	97
	98	if (fPtCovPerPhi) delete fPtCovPerPhi;
	99	fPtCovPerPhi = new TH2D(*aCut.fPtCovPerPhi);
	100	if (fPtPerPhi) delete fPtPerPhi;
	101	fPtPerPhi = new TH1D(*aCut.fPtPerPhi);
	102	if (fPtMultPerPhi) delete fPtMultPerPhi;
	103	fPtMultPerPhi = new TH2D(*aCut.fPtMultPerPhi);
	104	fPhiBins = aCut.fPhiBins;
	105	fNEventsProcessed = aCut.fNEventsProcessed;
	106
	107	return *this;
	108	}
	109
	110	AliFemtoString AliFemtoCutMonitorParticleEtCorr::Report(){
	111	// Prepare report from the execution
	112	string stemp = "*** AliFemtoCutMonitorParticleEtCorr report";
	113	AliFemtoString returnThis = stemp;
	114	return returnThis;
	115	}
	116
	117	void AliFemtoCutMonitorParticleEtCorr::Fill(const AliFemtoTrack* aTrack)
	118	{
	119	// Fill in the monitor histograms with the values from the current track
	120	// float tEnergy = ::sqrt(aTrack->P().mag2()+fMass*fMass);
	121	// float tRapidity = 0.5*::log((tEnergy+aTrack->P().z())/(tEnergy-aTrack->P().z()));
	122	float tPt = ::sqrt((aTrack->P().x())(aTrack->P().x())+(aTrack->P().y())(aTrack->P().y()));
	123	// float tEta = -TMath::Log(TMath::Tan(aTrack->P().theta()/2.0));
f6944668	124	float tPhi = aTrack->P().Phi();
95150c4b	125	Double_t tPiTwo = TMath::Pi()*2;
	126
	127	while (tPhi > tPiTwo) tPhi -= tPiTwo;
	128	while (tPhi < 0) tPhi += tPiTwo;
	129
	130	int nbin = (int) floor(tPhi * fPhiBins / tPiTwo);
	131	fPtSumEvent[nbin] += tPt;
	132	fMultSumEvent[nbin] += 1;
	133	}
	134
	135	void AliFemtoCutMonitorParticleEtCorr::Write()
	136	{
	137	// Write out the relevant histograms
	138	fPtPerPhi->Write();
	139	fPtCovPerPhi->Write();
	140	fPtMultPerPhi->Write();
	141	}
	142
	143	TList *AliFemtoCutMonitorParticleEtCorr::GetOutputList()
	144	{
	145	TList *tOutputList = new TList();
	146	tOutputList->Add(fPtPerPhi);
	147	tOutputList->Add(fPtCovPerPhi);
	148	tOutputList->Add(fPtMultPerPhi);
	149
	150	return tOutputList;
	151	}
	152
	153	void AliFemtoCutMonitorParticleEtCorr::EventBegin(const AliFemtoEvent* aEvent)
	154	{
3ccc556f	155	if (aEvent)
	156	for (int iter=0; iter<fPhiBins; iter++) {
	157	fPtSumEvent[iter] = 0;
	158	fMultSumEvent[iter] = 0;
	159	}
95150c4b	160	}
	161
	162	void AliFemtoCutMonitorParticleEtCorr::EventEnd(const AliFemtoEvent* aEvent)
	163	{
3ccc556f	164	if (aEvent) {
	165	for (int ispt=0; ispt<fPhiBins; ispt++) {
	166	fPtPerPhi->Fill(ispt, fPtSumEvent[ispt]);
	167	for (int ispt2=0; ispt2<fPhiBins; ispt2++) {
	168	fPtCovPerPhi->Fill(ispt, ispt2, fPtSumEvent[ispt]*fPtSumEvent[ispt2]);
	169	fPtMultPerPhi->Fill(ispt, ispt2, fMultSumEvent[ispt]*fMultSumEvent[ispt2]);
	170	}
95150c4b	171	}
3ccc556f	172	fNEventsProcessed++;
95150c4b	173	}
95150c4b	174	}
95150c4b	175