[u/mrichter/AliRoot.git] / PWG2 / SPECTRA / esdV0.C

#define esdV0_cxx
// The class definition in esdV0.h has been generated automatically
// by the ROOT utility TTree::MakeSelector(). This class is derived
// from the ROOT class TSelector. For more information on the TSelector
// framework see $ROOTSYS/README/README.SELECTOR or the ROOT User Manual.

// The following methods are defined in this file:
//    Begin():        called everytime a loop on the tree starts,
//                    a convenient place to create your histograms.
//    SlaveBegin():   called after Begin(), when on PROOF called only on the
//                    slave servers.
//    Process():      called for each event, in this function you decide what
//                    to read and fill your histograms.
//    SlaveTerminate: called at the end of the loop on the tree, when on PROOF
//                    called only on the slave servers.
//    Terminate():    called at the end of the loop on the tree,
//                    a convenient place to draw/fit your histograms.
//
// To use this file, try the following session on your Tree T:
//
// Root > T->Process("esdV0.C")
// Root > T->Process("esdV0.C","some options")
// Root > T->Process("esdV0.C+")
//

#include "esdV0.h"
#include <TPDGCode.h>
#include <TStyle.h>
#include <TCanvas.h>
#include "AliPID.h"

esdV0::esdV0(TTree *) :
  TSelector(),
  fChain(0),
  fESD(0),
  fHistMassK0(0),
  fHistMassLambda(0),
  fHistMassAntiLambda(0),
  fHistMassLambdaVsProb(0),
  fHistMassLambdaCut(0),
  fHistMassAntiLambdaCut(0)
{
  // Constructor. Initialization of pointers
}

esdV0::~esdV0() {
  // Remove all pointers

  delete fESD;

  // histograms are in the output list and deleted when the output
  // list is deleted by the TSelector dtor
}

void esdV0::Begin(TTree *)
{
  // The Begin() function is called at the start of the query.
  // When running with PROOF Begin() is only called on the client.
  // The tree argument is deprecated (on PROOF 0 is passed).

  TString option = GetOption();
}

void esdV0::SlaveBegin(TTree * tree)
{
  // The SlaveBegin() function is called after the Begin() function.
  // When running with PROOF SlaveBegin() is called on each slave server.
  // The tree argument is deprecated (on PROOF 0 is passed).

  Init(tree);
  
  TString option = GetOption();

  // Create histograms on each slave server
  fHistMassK0 = new TH1F("hMassK0", "K^{0} candidates", 100, 0.4, 0.6);
  fHistMassK0->GetXaxis()->SetTitle("M(#pi^{+}#pi^{-}) [GeV/c^{2}]");
  fHistMassK0->SetOption("E");
  fHistMassK0->SetMarkerStyle(kFullCircle);
  fHistMassK0->SetStats(kFALSE);

  fHistMassLambda = new TH1F("hMassLambda", "#Lambda^{0} candidates", 75, 1.05, 1.2);
  fHistMassLambda->GetXaxis()->SetTitle("M(p#pi^{-}) [GeV/c^{2}]");
  fHistMassLambda->SetOption("E");
  fHistMassLambda->SetMarkerStyle(kFullCircle);
  fHistMassLambda->SetMarkerColor(kRed);
  fHistMassLambda->SetFillColor(4);
  fHistMassLambda->SetStats(kFALSE);

  fHistMassAntiLambda = new TH1F("hMassAntiLambda", "#bar{#Lambda}^{0} candidates", 75, 1.05, 1.2);
  fHistMassAntiLambda->GetXaxis()->SetTitle("M(#bar{p}#pi^{+}) [GeV/c^{2}]");
  fHistMassAntiLambda->SetOption("E");
  fHistMassAntiLambda->SetMarkerStyle(kFullCircle);
  fHistMassAntiLambda->SetMarkerColor(kRed);
  fHistMassAntiLambda->SetFillColor(4);
  fHistMassAntiLambda->SetStats(kFALSE);

  fHistMassLambdaVsProb = new TH2F("hMassLambdaVsProb", "#Lambda^{0} and #bar{#Lambda}^{0} candidates", 21, -2.5, 102.5, 75, 1.05, 1.2);
  fHistMassLambdaVsProb->GetXaxis()->SetTitle("prob(p) [%]");
  fHistMassLambdaVsProb->GetYaxis()->SetTitle("M(p#pi) [GeV/c^{2}]");
  fHistMassLambdaVsProb->SetOption("BOX");
  fHistMassLambdaVsProb->SetFillStyle(0);
  fHistMassLambdaVsProb->SetStats(kFALSE);

  fHistMassLambdaCut = new TH1F("hMassLambdaCut", "#Lambda^{0} candidates", 75, 1.05, 1.2);
  fHistMassLambdaCut->GetXaxis()->SetTitle("M(p#pi^{-}) [GeV/c^{2}]");
  fHistMassLambdaCut->SetOption("E");
  fHistMassLambdaCut->SetMarkerStyle(kFullCircle);
  fHistMassLambdaCut->SetMarkerColor(kRed);
  fHistMassLambdaCut->SetStats(kFALSE);

  fHistMassAntiLambdaCut = new TH1F("hMassAntiLambdaCut", "#bar{#Lambda}^{0} candidates", 75, 1.05, 1.2);
  fHistMassAntiLambdaCut->GetXaxis()->SetTitle("M(#bar{p}#pi^{+}) [GeV/c^{2}]");
  fHistMassAntiLambdaCut->SetOption("E");
  fHistMassAntiLambdaCut->SetMarkerStyle(kFullCircle);
  fHistMassAntiLambdaCut->SetMarkerColor(kRed);
  fHistMassAntiLambdaCut->SetStats(kFALSE);
}

void esdV0::Init(TTree *tree)
{
  // The Init() function is called when the selector needs to initialize
  // a new tree or chain. Typically here the branch addresses of the tree
  // will be set. It is normaly not necessary to make changes to the
  // generated code, but the routine can be extended by the user if needed.
  // Init() will be called many times when running with PROOF.

  // Set branch addresses
  if (tree == 0) return;
  fChain = tree;

  fChain->SetBranchAddress("ESD",&fESD);
}

Bool_t esdV0::Notify()
{
  // The Notify() function is called when a new file is opened. This
  // can be either for a new TTree in a TChain or when when a new TTree
  // is started when using PROOF. Typically here the branch pointers
  // will be retrieved. It is normaly not necessary to make changes
  // to the generated code, but the routine can be extended by the
  // user if needed.
  
  return kTRUE;
}


Bool_t esdV0::Process(Long64_t entry)
{
  // The Process() function is called for each entry in the tree (or possibly
  // keyed object in the case of PROOF) to be processed. The entry argument
  // specifies which entry in the currently loaded tree is to be processed.
  // It can be passed to either TTree::GetEntry() or TBranch::GetEntry()
  // to read either all or the required parts of the data. When processing
  // keyed objects with PROOF, the object is already loaded and is available
  // via the fObject pointer.
  //
  // This function should contain the "body" of the analysis. It can contain
  // simple or elaborate selection criteria, run algorithms on the data
  // of the event and typically fill histograms.

  // WARNING when a selector is used with a TChain, you must use
  //  the pointer to the current TTree to call GetEntry(entry).
  //  The entry is always the local entry number in the current tree.
  //  Assuming that fChain is the pointer to the TChain being processed,
  //  use fChain->GetEntry(entry).

  fChain->GetEntry(entry);

  if (!fESD) return kFALSE;

  Int_t nv0s = fESD->GetNumberOfV0s();

  Double_t mass;
  Double_t pPrior[5] = {0, 0, 0.85, 0.1, 0.05};
  Double_t p[5];
  Double_t pSum, pNorm;
  for (Int_t iV0 = 0; iV0 < fESD->GetNumberOfV0s(); iV0++) {
    AliESDv0* v0 = fESD->GetV0(iV0);
    if (!v0) continue;

    fHistMassK0->Fill(v0->GetEffMass());

    v0->ChangeMassHypothesis(kLambda0);
    mass = v0->GetEffMass();
    fESD->GetTrack(v0->GetPindex())->GetESDpid(p);
    pSum = 0;
    for (Int_t i = 0; i < AliPID::kSPECIES; i++) pSum += p[i] * pPrior[i];
    if (pSum <= 0) pSum = 1.;
    pNorm = p[AliPID::kProton] * pPrior[AliPID::kProton] / pSum;
    fHistMassLambdaVsProb->Fill(100.*pNorm, mass); 
    if (pNorm > 0.1) fHistMassLambda->Fill(mass);
    
    v0->ChangeMassHypothesis(kLambda0Bar);
    mass = v0->GetEffMass();
    fESD->GetTrack(v0->GetNindex())->GetESDpid(p);
    pSum = 0;
    for (Int_t i = 0; i < AliPID::kSPECIES; i++) pSum += p[i] * pPrior[i];
    if (pSum <= 0) pSum = 1.;
    pNorm = p[AliPID::kProton] * pPrior[AliPID::kProton] / pSum;
    fHistMassLambdaVsProb->Fill(100.*pNorm, mass); 
    if (pNorm > 0.1) fHistMassAntiLambda->Fill(mass);
   
  }//loop v0s
  
  return kTRUE;
}


void esdV0::SlaveTerminate()
{
  // The SlaveTerminate() function is called after all entries or objects
  // have been processed. When running with PROOF SlaveTerminate() is called
  // on each slave server.

  // Add the histograms to the output on each slave server
  fOutput->Add(fHistMassK0);
  fOutput->Add(fHistMassLambda);
  fOutput->Add(fHistMassAntiLambda);
  fOutput->Add(fHistMassLambdaVsProb);
  fOutput->Add(fHistMassLambdaCut);
  fOutput->Add(fHistMassAntiLambdaCut);
}

void esdV0::Terminate()
{
  // The Terminate() function is the last function to be called during
  // a query. It always runs on the client, it can be used to present
  // the results graphically or save the results to file.

  fHistMassK0 = dynamic_cast<TH1F*>(fOutput->FindObject("hMassK0"));
  fHistMassLambda = dynamic_cast<TH1F*>(fOutput->FindObject("hMassLambda"));
  fHistMassAntiLambda = dynamic_cast<TH1F*>(fOutput->FindObject("hMassAntiLambda"));
  fHistMassLambdaVsProb = dynamic_cast<TH2F*>(fOutput->FindObject("hMassLambdaVsProb"));
  fHistMassLambdaCut = dynamic_cast<TH1F*>(fOutput->FindObject("hMassLambdaCut"));
  fHistMassAntiLambdaCut = dynamic_cast<TH1F*>(fOutput->FindObject("hMassAntiLambdaCut"));

  TFile* file = TFile::Open("V0.root", "RECREATE");
  fHistMassK0->Write();
  fHistMassLambda->Write();
  fHistMassAntiLambda->Write();
  fHistMassLambdaVsProb->Write();
  fHistMassLambdaCut->Write();
  fHistMassAntiLambdaCut->Write();
  file->Close();
  delete file;

  if (!gROOT->IsBatch()) 
    {
      TCanvas *c1 = new TCanvas("c1","V0s",10,10,610,610);
      c1->SetFillColor(10);
      c1->SetHighLightColor(10);
      c1->Divide(2,1);
 
      c1->cd(1);
      fHistMassK0->DrawCopy("E");
      fHistMassK0->Fit("gaus","q","",0.49,0.51);
      TVirtualPad* pad = (TVirtualPad*)c1->cd(2);
      pad->Divide(1,2);
      pad->cd(1);
      fHistMassLambda->DrawCopy("E");
      pad->cd(2);
      fHistMassAntiLambda->DrawCopy("E");
    }
}
Commit	Line	Data
bef34687	1	#define esdV0_cxx
	2	// The class definition in esdV0.h has been generated automatically
	3	// by the ROOT utility TTree::MakeSelector(). This class is derived
	4	// from the ROOT class TSelector. For more information on the TSelector
	5	// framework see $ROOTSYS/README/README.SELECTOR or the ROOT User Manual.
	6
	7	// The following methods are defined in this file:
	8	// Begin(): called everytime a loop on the tree starts,
	9	// a convenient place to create your histograms.
	10	// SlaveBegin(): called after Begin(), when on PROOF called only on the
	11	// slave servers.
	12	// Process(): called for each event, in this function you decide what
	13	// to read and fill your histograms.
	14	// SlaveTerminate: called at the end of the loop on the tree, when on PROOF
	15	// called only on the slave servers.
	16	// Terminate(): called at the end of the loop on the tree,
	17	// a convenient place to draw/fit your histograms.
	18	//
	19	// To use this file, try the following session on your Tree T:
	20	//
	21	// Root > T->Process("esdV0.C")
	22	// Root > T->Process("esdV0.C","some options")
	23	// Root > T->Process("esdV0.C+")
	24	//
	25
	26	#include "esdV0.h"
	27	#include <TPDGCode.h>
	28	#include <TStyle.h>
	29	#include <TCanvas.h>
	30	#include "AliPID.h"
	31
	32	esdV0::esdV0(TTree *) :
	33	TSelector(),
	34	fChain(0),
	35	fESD(0),
	36	fHistMassK0(0),
	37	fHistMassLambda(0),
	38	fHistMassAntiLambda(0),
	39	fHistMassLambdaVsProb(0),
	40	fHistMassLambdaCut(0),
	41	fHistMassAntiLambdaCut(0)
	42	{
	43	// Constructor. Initialization of pointers
	44	}
	45
	46	esdV0::~esdV0() {
	47	// Remove all pointers
	48
	49	delete fESD;
	50
	51	// histograms are in the output list and deleted when the output
	52	// list is deleted by the TSelector dtor
	53	}
	54
	55	void esdV0::Begin(TTree *)
	56	{
	57	// The Begin() function is called at the start of the query.
	58	// When running with PROOF Begin() is only called on the client.
	59	// The tree argument is deprecated (on PROOF 0 is passed).
	60
	61	TString option = GetOption();
	62	}
	63
	64	void esdV0::SlaveBegin(TTree * tree)
65	{
66	// The SlaveBegin() function is called after the Begin() function.
67	// When running with PROOF SlaveBegin() is called on each slave server.
68	// The tree argument is deprecated (on PROOF 0 is passed).
69
70	Init(tree);
71
72	TString option = GetOption();
73
74	// Create histograms on each slave server
75	fHistMassK0 = new TH1F("hMassK0", "K^{0} candidates", 100, 0.4, 0.6);
76	fHistMassK0->GetXaxis()->SetTitle("M(#pi^{+}#pi^{-}) [GeV/c^{2}]");
77	fHistMassK0->SetOption("E");
78	fHistMassK0->SetMarkerStyle(kFullCircle);
79	fHistMassK0->SetStats(kFALSE);
80
81	fHistMassLambda = new TH1F("hMassLambda", "#Lambda^{0} candidates", 75, 1.05, 1.2);
82	fHistMassLambda->GetXaxis()->SetTitle("M(p#pi^{-}) [GeV/c^{2}]");
83	fHistMassLambda->SetOption("E");
84	fHistMassLambda->SetMarkerStyle(kFullCircle);
85	fHistMassLambda->SetMarkerColor(kRed);
86	fHistMassLambda->SetFillColor(4);
87	fHistMassLambda->SetStats(kFALSE);
88
89	fHistMassAntiLambda = new TH1F("hMassAntiLambda", "#bar{#Lambda}^{0} candidates", 75, 1.05, 1.2);
90	fHistMassAntiLambda->GetXaxis()->SetTitle("M(#bar{p}#pi^{+}) [GeV/c^{2}]");
91	fHistMassAntiLambda->SetOption("E");
92	fHistMassAntiLambda->SetMarkerStyle(kFullCircle);
93	fHistMassAntiLambda->SetMarkerColor(kRed);
94	fHistMassAntiLambda->SetFillColor(4);
95	fHistMassAntiLambda->SetStats(kFALSE);
96
97	fHistMassLambdaVsProb = new TH2F("hMassLambdaVsProb", "#Lambda^{0} and #bar{#Lambda}^{0} candidates", 21, -2.5, 102.5, 75, 1.05, 1.2);
98	fHistMassLambdaVsProb->GetXaxis()->SetTitle("prob(p) [%]");
99	fHistMassLambdaVsProb->GetYaxis()->SetTitle("M(p#pi) [GeV/c^{2}]");
100	fHistMassLambdaVsProb->SetOption("BOX");
101	fHistMassLambdaVsProb->SetFillStyle(0);
102	fHistMassLambdaVsProb->SetStats(kFALSE);
103
104	fHistMassLambdaCut = new TH1F("hMassLambdaCut", "#Lambda^{0} candidates", 75, 1.05, 1.2);
105	fHistMassLambdaCut->GetXaxis()->SetTitle("M(p#pi^{-}) [GeV/c^{2}]");
106	fHistMassLambdaCut->SetOption("E");
107	fHistMassLambdaCut->SetMarkerStyle(kFullCircle);
108	fHistMassLambdaCut->SetMarkerColor(kRed);
109	fHistMassLambdaCut->SetStats(kFALSE);
110
111	fHistMassAntiLambdaCut = new TH1F("hMassAntiLambdaCut", "#bar{#Lambda}^{0} candidates", 75, 1.05, 1.2);
112	fHistMassAntiLambdaCut->GetXaxis()->SetTitle("M(#bar{p}#pi^{+}) [GeV/c^{2}]");
113	fHistMassAntiLambdaCut->SetOption("E");
114	fHistMassAntiLambdaCut->SetMarkerStyle(kFullCircle);
115	fHistMassAntiLambdaCut->SetMarkerColor(kRed);
116	fHistMassAntiLambdaCut->SetStats(kFALSE);
117	}
118
119	void esdV0::Init(TTree *tree)
120	{
121	// The Init() function is called when the selector needs to initialize
122	// a new tree or chain. Typically here the branch addresses of the tree
123	// will be set. It is normaly not necessary to make changes to the
124	// generated code, but the routine can be extended by the user if needed.
125	// Init() will be called many times when running with PROOF.
126
127	// Set branch addresses
128	if (tree == 0) return;
129	fChain = tree;
130
131	fChain->SetBranchAddress("ESD",&fESD);
132	}
133
134	Bool_t esdV0::Notify()
135	{
136	// The Notify() function is called when a new file is opened. This
137	// can be either for a new TTree in a TChain or when when a new TTree
138	// is started when using PROOF. Typically here the branch pointers
139	// will be retrieved. It is normaly not necessary to make changes
140	// to the generated code, but the routine can be extended by the
141	// user if needed.
142
143	return kTRUE;
144	}
145
146
147	Bool_t esdV0::Process(Long64_t entry)
148	{
149	// The Process() function is called for each entry in the tree (or possibly
150	// keyed object in the case of PROOF) to be processed. The entry argument
151	// specifies which entry in the currently loaded tree is to be processed.
152	// It can be passed to either TTree::GetEntry() or TBranch::GetEntry()
153	// to read either all or the required parts of the data. When processing
154	// keyed objects with PROOF, the object is already loaded and is available
155	// via the fObject pointer.
156	//
157	// This function should contain the "body" of the analysis. It can contain
158	// simple or elaborate selection criteria, run algorithms on the data
159	// of the event and typically fill histograms.
160
161	// WARNING when a selector is used with a TChain, you must use
162	// the pointer to the current TTree to call GetEntry(entry).
163	// The entry is always the local entry number in the current tree.
164	// Assuming that fChain is the pointer to the TChain being processed,
165	// use fChain->GetEntry(entry).
166
167	fChain->GetEntry(entry);
168
169	if (!fESD) return kFALSE;
170
171	Int_t nv0s = fESD->GetNumberOfV0s();
172
173	Double_t mass;
174	Double_t pPrior[5] = {0, 0, 0.85, 0.1, 0.05};
175	Double_t p[5];
176	Double_t pSum, pNorm;
177	for (Int_t iV0 = 0; iV0 < fESD->GetNumberOfV0s(); iV0++) {
178	AliESDv0* v0 = fESD->GetV0(iV0);
179	if (!v0) continue;
180
181	fHistMassK0->Fill(v0->GetEffMass());
182
183	v0->ChangeMassHypothesis(kLambda0);
184	mass = v0->GetEffMass();
185	fESD->GetTrack(v0->GetPindex())->GetESDpid(p);
186	pSum = 0;
187	for (Int_t i = 0; i < AliPID::kSPECIES; i++) pSum += p[i] * pPrior[i];
188	if (pSum <= 0) pSum = 1.;
189	pNorm = p[AliPID::kProton] * pPrior[AliPID::kProton] / pSum;
190	fHistMassLambdaVsProb->Fill(100.*pNorm, mass);
191	if (pNorm > 0.1) fHistMassLambda->Fill(mass);
192
193	v0->ChangeMassHypothesis(kLambda0Bar);
194	mass = v0->GetEffMass();
195	fESD->GetTrack(v0->GetNindex())->GetESDpid(p);
196	pSum = 0;
197	for (Int_t i = 0; i < AliPID::kSPECIES; i++) pSum += p[i] * pPrior[i];
198	if (pSum <= 0) pSum = 1.;
199	pNorm = p[AliPID::kProton] * pPrior[AliPID::kProton] / pSum;
200	fHistMassLambdaVsProb->Fill(100.*pNorm, mass);
201	if (pNorm > 0.1) fHistMassAntiLambda->Fill(mass);
202
203	}//loop v0s
204
205	return kTRUE;
206	}
207
208
209
210	void esdV0::SlaveTerminate()
211	{
212	// The SlaveTerminate() function is called after all entries or objects
213	// have been processed. When running with PROOF SlaveTerminate() is called
214	// on each slave server.
215
216	// Add the histograms to the output on each slave server
217	fOutput->Add(fHistMassK0);
218	fOutput->Add(fHistMassLambda);
219	fOutput->Add(fHistMassAntiLambda);
220	fOutput->Add(fHistMassLambdaVsProb);
221	fOutput->Add(fHistMassLambdaCut);
222	fOutput->Add(fHistMassAntiLambdaCut);
223	}
224
225	void esdV0::Terminate()
226	{
227	// The Terminate() function is the last function to be called during
228	// a query. It always runs on the client, it can be used to present
229	// the results graphically or save the results to file.
230
231	fHistMassK0 = dynamic_cast<TH1F*>(fOutput->FindObject("hMassK0"));
232	fHistMassLambda = dynamic_cast<TH1F*>(fOutput->FindObject("hMassLambda"));
233	fHistMassAntiLambda = dynamic_cast<TH1F*>(fOutput->FindObject("hMassAntiLambda"));
234	fHistMassLambdaVsProb = dynamic_cast<TH2F*>(fOutput->FindObject("hMassLambdaVsProb"));
235	fHistMassLambdaCut = dynamic_cast<TH1F*>(fOutput->FindObject("hMassLambdaCut"));
236	fHistMassAntiLambdaCut = dynamic_cast<TH1F*>(fOutput->FindObject("hMassAntiLambdaCut"));
237
238	TFile* file = TFile::Open("V0.root", "RECREATE");
239	fHistMassK0->Write();
240	fHistMassLambda->Write();
241	fHistMassAntiLambda->Write();
242	fHistMassLambdaVsProb->Write();
243	fHistMassLambdaCut->Write();
244	fHistMassAntiLambdaCut->Write();
245	file->Close();
246	delete file;
247
248	if (!gROOT->IsBatch())
249	{
250	TCanvas *c1 = new TCanvas("c1","V0s",10,10,610,610);
251	c1->SetFillColor(10);
252	c1->SetHighLightColor(10);
253	c1->Divide(2,1);
254
255	c1->cd(1);
256	fHistMassK0->DrawCopy("E");
257	fHistMassK0->Fit("gaus","q","",0.49,0.51);
258	TVirtualPad* pad = (TVirtualPad*)c1->cd(2);
259	pad->Divide(1,2);
260	pad->cd(1);
261	fHistMassLambda->DrawCopy("E");
262	pad->cd(2);
263	fHistMassAntiLambda->DrawCopy("E");
264	}
265	}