[u/mrichter/AliRoot.git] / HLT / PHOS / AliHLTPHOSPhysicsHistogramProducer.cxx


/**************************************************************************
 * This file is property of and copyright by the ALICE HLT Project        * 
 * All rights reserved.                                                   *
 *                                                                        *
 * Primary Authors: Albin Gaignette
 *                                                                        *
 * 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.                  *
 **************************************************************************/

/** 
 * @file   AliHLTPHOSPhysicsHistogramProducer
 * @author Albin Gaignette
 * @date 
 * @brief  Histogram producer for PHOS HLT 
 */

// see header file for class documentation
// or
// refer to README to build package
// or
// visit http://web.ift.uib.no/~kjeks/doc/alice-hlt

#include "AliHLTPHOSPhysicsHistogramProducer.h"
#include "AliESDCaloCluster.h"
#include "TMath.h"
#include "TClonesArray.h"
#include <iostream>
#include "TH1F.h"
#include "TH2F.h"

AliHLTPHOSPhysicsHistogramProducer::AliHLTPHOSPhysicsHistogramProducer() :
  //    AliHLTPHOSBase(),
    fHistNcls(0),
    fHistEnergy(0),
    fHistTotEnergy(0),
    fHistTwoClusterInvMass(0),
    fHistNcells(0),
    fHistNcellsPercentage(0),
    fHistCellsEnergy(0),
    fHistNclusterTotE(0),
    fHistNcellsSumCells(0),
    fHistArrayPtr(0)
{
  // See header file for documentation
  fHistArrayPtr = new TObjArray;

  fHistNcls = new TH1F("fHistNcls", "Number of clusters per event", 20, 0, 20);
  fHistNcls->GetXaxis()->SetTitle("Number of clusters per event");
  fHistNcls->GetYaxis()->SetTitle("# of counts");
  fHistArrayPtr ->AddLast(fHistNcls);
  
  fHistEnergy = new TH1F("fHistEnergy", "Energy in each cluster", 200, 0, 100);
  fHistEnergy->GetXaxis()->SetTitle("E_{T} GeV");
  fHistEnergy->GetYaxis()->SetTitle("# of counts");
  fHistEnergy->SetMarkerStyle(21);
  fHistArrayPtr->AddLast(fHistEnergy);

  fHistTotEnergy = new TH1F("fHistTotEnergy", "Total energy in each event", 200, 0, 200);
  fHistTotEnergy->GetXaxis()->SetTitle("E_{T} GeV");
  fHistTotEnergy->GetYaxis()->SetTitle("# of counts");
  fHistTotEnergy->SetMarkerStyle(21);
  fHistArrayPtr ->AddLast(fHistTotEnergy);

  fHistTwoClusterInvMass = new TH1F("fHistTwoClusterInvMass", "Invariant mass of two clusters", 200, 0, 1);
  fHistTwoClusterInvMass->GetXaxis()->SetTitle("m_{#gamma#gamma} GeV");
  fHistTwoClusterInvMass->GetYaxis()->SetTitle("Number of counts");
  fHistTwoClusterInvMass->SetMarkerStyle(21);
  fHistArrayPtr->AddLast(fHistTwoClusterInvMass);

  fHistNcells = new TH1F("fHistNcells", "Number of cells in each cluster", 200, 0, 100);
  fHistNcells->GetXaxis()->SetTitle("# of cells in the cluster");
  fHistNcells->GetYaxis()->SetTitle("# of counts");
  fHistNcells->SetMarkerStyle(21); 
  fHistArrayPtr ->AddLast(fHistNcells);

  fHistNcellsPercentage = new TH1F("fHistNcellsPercentage", "Percentage of cells in the cluster", 200, 0, 10);
  fHistNcellsPercentage->GetXaxis()->SetTitle("100* #frac{Total # of cells in the cluster}{Total # of cells in PHOS}");
  fHistNcellsPercentage->GetYaxis()->SetTitle("# of counts");
  fHistNcellsPercentage->SetMarkerStyle(21);
  fHistArrayPtr ->AddLast(fHistNcellsPercentage);

  fHistCellsEnergy = new TH2F("fHistCellsEnergy","# of cells in the cluster vs cluster Energy",100,0,60,100,1,60);
  fHistCellsEnergy->GetXaxis()->SetTitle("Cluster Energy [GeV]");
  fHistCellsEnergy->GetYaxis()->SetTitle("# of cells in the cluster");
  fHistArrayPtr ->AddLast(fHistCellsEnergy);

  fHistNclusterTotE = new TH2F("fHistNclusterTotE","# of clusters vs total Energy",100,0,200,100,1,15);
  fHistNclusterTotE->GetXaxis()->SetTitle("total Energy in the event [GeV]");
  fHistNclusterTotE->GetYaxis()->SetTitle("# of clusters in each event");
  fHistArrayPtr ->AddLast(fHistNclusterTotE);

  fHistNcellsSumCells = new TH2F("fHistNcellsSumCells","# of cells in PHOS vs sum cells in the clusters",100,0,300,100,17910,17930);
  fHistNcellsSumCells->GetXaxis()->SetTitle("# of cells in the clusters");
  fHistNcellsSumCells->GetYaxis()->SetTitle("Total # of cells in PHOS");
  fHistArrayPtr->AddLast(fHistNcellsSumCells);
}

AliHLTPHOSPhysicsHistogramProducer::~AliHLTPHOSPhysicsHistogramProducer()
{
  // See header file for documentation
  if(fHistNcls)
    {
      delete fHistNcls;
      fHistNcls = 0;
    }
  if(fHistEnergy)
    {
      delete fHistEnergy;
      fHistEnergy = 0;
    }
  if(fHistTotEnergy)
    {
      delete fHistTotEnergy;
      fHistTotEnergy = 0;
    }
  if(fHistTwoClusterInvMass)
    {
      delete fHistTwoClusterInvMass;
      fHistTwoClusterInvMass = 0;
    }
  if(fHistNcells)
    {
      delete fHistNcells;
      fHistNcells = 0;
    }
  if(fHistNcellsPercentage)
    {
      delete fHistNcellsPercentage;
      fHistNcellsPercentage = 0;
    }
  if(fHistCellsEnergy)
    {
      delete fHistCellsEnergy;
      fHistCellsEnergy = 0;
    }
  if(fHistNclusterTotE)
    {
      delete fHistNclusterTotE;
      fHistNclusterTotE = 0;
    }
  if(fHistNcellsSumCells)
    {
      delete fHistNcellsSumCells;
      fHistNcellsSumCells = 0;
    }
  if(fHistArrayPtr)
    {
      delete fHistArrayPtr;
      fHistArrayPtr = 0;
    }
}

TObjArray* AliHLTPHOSPhysicsHistogramProducer::GetHistograms()
{  
  // See header file for documentation
  return fHistArrayPtr;
}

Int_t AliHLTPHOSPhysicsHistogramProducer::AnalyseClusters(TClonesArray* clusters)
{   
  // See header file for documentation
  Int_t nPHOSModules = 3;
  Int_t totClustersAll=0;
  Int_t totClustersAllwithprob=0;
  Int_t totClustersFile=0; 
  Int_t totClustersFilewithprob=0; 
  Double_t theta=0.0; 
  Double_t m=0.0;
  Double_t TotalCells =0;
  Double_t TotalCellsPHOS= NXCOLUMNSMOD*NZROWSMOD * nPHOSModules;
  // 3584 crystals * 3 modules
  Double_t Nclusters =0;
  Double_t TotEnergy = 0.0;
     
  Int_t ncls = clusters->GetEntriesFast();

  fHistNcls->Fill(ncls);
  totClustersAll+=ncls;
  totClustersFile+=ncls;
  Double_t totE = 0;
   
  Nclusters=0;
  TotEnergy=0.0;
  TotalCells=0;
          
  Float_t** tabposition = new Float_t*[ncls];
  Double_t* tabenergy=new Double_t[ncls];
     
  for(Int_t icls = 0; icls < ncls; icls++)
    {
      AliESDCaloCluster* cluster = (AliESDCaloCluster*)clusters->At(icls);
	  	
      totClustersFilewithprob+=icls;
      totClustersAllwithprob+=icls;
      totE += cluster->E();
      fHistEnergy->Fill(cluster->E());

      tabposition[icls] = new Float_t[3];
      cluster->GetPosition(tabposition[icls]); 
      tabenergy[icls] = cluster->E(); 
  	  
      TotalCells += cluster -> GetNCells();
      fHistNcells->Fill(cluster -> GetNCells());
      fHistCellsEnergy->Fill(cluster -> GetNCells(),cluster->E());
      TotEnergy += cluster->E();
      Nclusters++;
    }
 
  fHistNcellsPercentage->Fill(100*(TotalCells/TotalCellsPHOS));
  fHistNclusterTotE->Fill(TotEnergy,Nclusters);
  fHistNcellsSumCells->Fill(TotalCells,TotalCellsPHOS);	
       
  if(totE > 0)
    fHistTotEnergy->Fill(totE);
  
  for(Int_t ipho = 0; ipho<(ncls-1) ; ipho++)
    { 
      for(Int_t jpho = ipho+1 ; jpho<ncls ; jpho++)
	{ 
	  // Calcul of the theta angle between two photons
	  theta =(2* asin(0.5*TMath::Sqrt((tabposition[ipho][0]-tabposition[jpho][0])*(tabposition[ipho][0]-tabposition[jpho][0]) +(tabposition[ipho][1]-tabposition[jpho][1])*(tabposition[ipho][1]-tabposition[jpho][1]))/460));
 	      
	  // Calcul of the mass m of the pion 
	  m =(TMath::Sqrt(2 * tabenergy[ipho]* tabenergy[jpho]*(1-TMath::Cos(theta))));
	  fHistTwoClusterInvMass->Fill(m);
	}
    }
    
  for(Int_t j=0 ; j<ncls;j++)
    {
      delete[] tabposition[j];
    }
	  
  delete[] tabposition;
  delete[] tabenergy;

  return 0;
}
Commit	Line	Data
87434909	1
	2	/**************************************************************************
	3	* This file is property of and copyright by the ALICE HLT Project *
	4	* All rights reserved. *
	5	* *
	6	* Primary Authors: Albin Gaignette
	7	* *
	8	* Permission to use, copy, modify and distribute this software and its *
	9	* documentation strictly for non-commercial purposes is hereby granted *
	10	* without fee, provided that the above copyright notice appears in all *
	11	* copies and that both the copyright notice and this permission notice *
	12	* appear in the supporting documentation. The authors make no claims *
	13	* about the suitability of this software for any purpose. It is *
	14	* provided "as is" without express or implied warranty. *
	15	**************************************************************************/
	16
	17	/**
	18	* @file AliHLTPHOSPhysicsHistogramProducer
	19	* @author Albin Gaignette
	20	* @date
	21	* @brief Histogram producer for PHOS HLT
	22	*/
	23
	24	// see header file for class documentation
	25	// or
	26	// refer to README to build package
	27	// or
	28	// visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
	29
	30	#include "AliHLTPHOSPhysicsHistogramProducer.h"
	31	#include "AliESDCaloCluster.h"
	32	#include "TMath.h"
	33	#include "TClonesArray.h"
	34	#include <iostream>
	35	#include "TH1F.h"
	36	#include "TH2F.h"
	37
	38	AliHLTPHOSPhysicsHistogramProducer::AliHLTPHOSPhysicsHistogramProducer() :
9f050726	39	// AliHLTPHOSBase(),
87434909	40	fHistNcls(0),
	41	fHistEnergy(0),
	42	fHistTotEnergy(0),
	43	fHistTwoClusterInvMass(0),
	44	fHistNcells(0),
	45	fHistNcellsPercentage(0),
	46	fHistCellsEnergy(0),
	47	fHistNclusterTotE(0),
	48	fHistNcellsSumCells(0),
	49	fHistArrayPtr(0)
	50	{
	51	// See header file for documentation
	52	fHistArrayPtr = new TObjArray;
	53
	54	fHistNcls = new TH1F("fHistNcls", "Number of clusters per event", 20, 0, 20);
	55	fHistNcls->GetXaxis()->SetTitle("Number of clusters per event");
	56	fHistNcls->GetYaxis()->SetTitle("# of counts");
	57	fHistArrayPtr ->AddLast(fHistNcls);
	58
	59	fHistEnergy = new TH1F("fHistEnergy", "Energy in each cluster", 200, 0, 100);
	60	fHistEnergy->GetXaxis()->SetTitle("E_{T} GeV");
	61	fHistEnergy->GetYaxis()->SetTitle("# of counts");
	62	fHistEnergy->SetMarkerStyle(21);
	63	fHistArrayPtr->AddLast(fHistEnergy);
	64
	65	fHistTotEnergy = new TH1F("fHistTotEnergy", "Total energy in each event", 200, 0, 200);
	66	fHistTotEnergy->GetXaxis()->SetTitle("E_{T} GeV");
	67	fHistTotEnergy->GetYaxis()->SetTitle("# of counts");
	68	fHistTotEnergy->SetMarkerStyle(21);
	69	fHistArrayPtr ->AddLast(fHistTotEnergy);
	70
	71	fHistTwoClusterInvMass = new TH1F("fHistTwoClusterInvMass", "Invariant mass of two clusters", 200, 0, 1);
	72	fHistTwoClusterInvMass->GetXaxis()->SetTitle("m_{#gamma#gamma} GeV");
	73	fHistTwoClusterInvMass->GetYaxis()->SetTitle("Number of counts");
	74	fHistTwoClusterInvMass->SetMarkerStyle(21);
	75	fHistArrayPtr->AddLast(fHistTwoClusterInvMass);
	76
	77	fHistNcells = new TH1F("fHistNcells", "Number of cells in each cluster", 200, 0, 100);
	78	fHistNcells->GetXaxis()->SetTitle("# of cells in the cluster");
	79	fHistNcells->GetYaxis()->SetTitle("# of counts");
	80	fHistNcells->SetMarkerStyle(21);
	81	fHistArrayPtr ->AddLast(fHistNcells);
	82
	83	fHistNcellsPercentage = new TH1F("fHistNcellsPercentage", "Percentage of cells in the cluster", 200, 0, 10);
	84	fHistNcellsPercentage->GetXaxis()->SetTitle("100* #frac{Total # of cells in the cluster}{Total # of cells in PHOS}");
	85	fHistNcellsPercentage->GetYaxis()->SetTitle("# of counts");
	86	fHistNcellsPercentage->SetMarkerStyle(21);
	87	fHistArrayPtr ->AddLast(fHistNcellsPercentage);
	88
	89	fHistCellsEnergy = new TH2F("fHistCellsEnergy","# of cells in the cluster vs cluster Energy",100,0,60,100,1,60);
	90	fHistCellsEnergy->GetXaxis()->SetTitle("Cluster Energy [GeV]");
	91	fHistCellsEnergy->GetYaxis()->SetTitle("# of cells in the cluster");
	92	fHistArrayPtr ->AddLast(fHistCellsEnergy);
	93
	94	fHistNclusterTotE = new TH2F("fHistNclusterTotE","# of clusters vs total Energy",100,0,200,100,1,15);
	95	fHistNclusterTotE->GetXaxis()->SetTitle("total Energy in the event [GeV]");
	96	fHistNclusterTotE->GetYaxis()->SetTitle("# of clusters in each event");
	97	fHistArrayPtr ->AddLast(fHistNclusterTotE);
	98
	99	fHistNcellsSumCells = new TH2F("fHistNcellsSumCells","# of cells in PHOS vs sum cells in the clusters",100,0,300,100,17910,17930);
	100	fHistNcellsSumCells->GetXaxis()->SetTitle("# of cells in the clusters");
	101	fHistNcellsSumCells->GetYaxis()->SetTitle("Total # of cells in PHOS");
	102	fHistArrayPtr->AddLast(fHistNcellsSumCells);
	103	}
104
105	AliHLTPHOSPhysicsHistogramProducer::~AliHLTPHOSPhysicsHistogramProducer()
106	{
107	// See header file for documentation
108	if(fHistNcls)
109	{
110	delete fHistNcls;
111	fHistNcls = 0;
112	}
113	if(fHistEnergy)
114	{
115	delete fHistEnergy;
116	fHistEnergy = 0;
117	}
118	if(fHistTotEnergy)
119	{
120	delete fHistTotEnergy;
121	fHistTotEnergy = 0;
122	}
123	if(fHistTwoClusterInvMass)
124	{
125	delete fHistTwoClusterInvMass;
126	fHistTwoClusterInvMass = 0;
127	}
128	if(fHistNcells)
129	{
130	delete fHistNcells;
131	fHistNcells = 0;
132	}
133	if(fHistNcellsPercentage)
134	{
135	delete fHistNcellsPercentage;
136	fHistNcellsPercentage = 0;
137	}
138	if(fHistCellsEnergy)
139	{
140	delete fHistCellsEnergy;
141	fHistCellsEnergy = 0;
142	}
143	if(fHistNclusterTotE)
144	{
145	delete fHistNclusterTotE;
146	fHistNclusterTotE = 0;
147	}
148	if(fHistNcellsSumCells)
149	{
150	delete fHistNcellsSumCells;
151	fHistNcellsSumCells = 0;
152	}
153	if(fHistArrayPtr)
154	{
155	delete fHistArrayPtr;
156	fHistArrayPtr = 0;
157	}
158	}
159
160	TObjArray* AliHLTPHOSPhysicsHistogramProducer::GetHistograms()
161	{
162	// See header file for documentation
163	return fHistArrayPtr;
164	}
165
166	Int_t AliHLTPHOSPhysicsHistogramProducer::AnalyseClusters(TClonesArray* clusters)
167	{
168	// See header file for documentation
169	Int_t nPHOSModules = 3;
170	Int_t totClustersAll=0;
171	Int_t totClustersAllwithprob=0;
172	Int_t totClustersFile=0;
173	Int_t totClustersFilewithprob=0;
174	Double_t theta=0.0;
175	Double_t m=0.0;
176	Double_t TotalCells =0;
177	Double_t TotalCellsPHOS= NXCOLUMNSMODNZROWSMOD nPHOSModules;
178	// 3584 crystals * 3 modules
179	Double_t Nclusters =0;
180	Double_t TotEnergy = 0.0;
181
182	Int_t ncls = clusters->GetEntriesFast();
183
184	fHistNcls->Fill(ncls);
185	totClustersAll+=ncls;
186	totClustersFile+=ncls;
187	Double_t totE = 0;
188
189	Nclusters=0;
190	TotEnergy=0.0;
191	TotalCells=0;
192
193	Float_t** tabposition = new Float_t*[ncls];
194	Double_t* tabenergy=new Double_t[ncls];
195
196	for(Int_t icls = 0; icls < ncls; icls++)
197	{
198	AliESDCaloCluster* cluster = (AliESDCaloCluster*)clusters->At(icls);
199
200	totClustersFilewithprob+=icls;
201	totClustersAllwithprob+=icls;
202	totE += cluster->E();
203	fHistEnergy->Fill(cluster->E());
204
205	tabposition[icls] = new Float_t[3];
206	cluster->GetPosition(tabposition[icls]);
207	tabenergy[icls] = cluster->E();
208
209	TotalCells += cluster -> GetNCells();
210	fHistNcells->Fill(cluster -> GetNCells());
211	fHistCellsEnergy->Fill(cluster -> GetNCells(),cluster->E());
212	TotEnergy += cluster->E();
213	Nclusters++;
214	}
215
216	fHistNcellsPercentage->Fill(100*(TotalCells/TotalCellsPHOS));
217	fHistNclusterTotE->Fill(TotEnergy,Nclusters);
218	fHistNcellsSumCells->Fill(TotalCells,TotalCellsPHOS);
219
220	if(totE > 0)
221	fHistTotEnergy->Fill(totE);
222
223	for(Int_t ipho = 0; ipho<(ncls-1) ; ipho++)
224	{
225	for(Int_t jpho = ipho+1 ; jpho<ncls ; jpho++)
226	{
227	// Calcul of the theta angle between two photons
228	theta =(2* asin(0.5TMath::Sqrt((tabposition[ipho][0]-tabposition[jpho][0])(tabposition[ipho][0]-tabposition[jpho][0]) +(tabposition[ipho][1]-tabposition[jpho][1])*(tabposition[ipho][1]-tabposition[jpho][1]))/460));
229
230	// Calcul of the mass m of the pion
231	m =(TMath::Sqrt(2 * tabenergy[ipho]* tabenergy[jpho]*(1-TMath::Cos(theta))));
232	fHistTwoClusterInvMass->Fill(m);
233	}
234	}
235
236	for(Int_t j=0 ; j<ncls;j++)
237	{
238	delete[] tabposition[j];
239	}
240
241	delete[] tabposition;
242	delete[] tabenergy;
243
244	return 0;
245	}
246
247
248
249