[u/mrichter/AliRoot.git] / EMCAL / AliEMCALCell.cxx

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

/* $Id:*/

//_________________________________________________________________________
// Top EMCAL folder which will keep all information about EMCAL itself,
// super Modules (SM), modules, towers, set of hists and so on.
//
//*-- Author: Aleksei Pavlinov (WSU, Detroit, USA) 

#include "AliEMCALCell.h"
#include "AliEMCALHistoUtilities.h"
#include "AliEMCALGeometry.h"
#include "AliEMCALFolder.h"
#include "AliEMCALSuperModule.h"
#include "AliEMCALCalibData.h"

#include "AliEMCALCalibCoefs.h"

#include <TROOT.h>
#include <TStyle.h>
#include <TList.h>
#include <TH1.h>
#include <TF1.h>
#include <TNtuple.h>
#include <TObjectSet.h>

typedef  AliEMCALHistoUtilities u;

ClassImp(AliEMCALCell)

Double_t ADCCHANNELEC = 0.0153;  // Update 24 Apr 2007: 250./16/1024 - width of one ADC channel in GeV
Double_t MPI0         = 0.13498; // mass of pi0
Double_t MPI02        = MPI0*MPI0; // mass**2

AliEMCALCell::AliEMCALCell() : 
TObjectSet(), 
fAbsId(0),fSupMod(0),fModule(0),fPhi(0),fEta(0),fPhiCell(0),fEtaCell(0),fCcIn(0),fCcOut(0),
fFun(0)
{
}

AliEMCALCell::AliEMCALCell(const Int_t absId, const char* title) : 
TObjectSet(Form("Cell%4.4i",absId)), 
fAbsId(absId),fSupMod(0),fModule(0),fPhi(0),fEta(0),fPhiCell(0),fEtaCell(0),fCcIn(0),fCcOut(0),
fFun(0)
{
  SetTitle(title);
  
  AliEMCALGeometry *g = AliEMCALGeometry::GetInstance();
  g->GetCellIndex(fAbsId, fSupMod, fModule, fPhi, fEta);
  g->GetCellPhiEtaIndexInSModule(fSupMod, fModule, fPhi, fEta, fPhiCell, fEtaCell);

} 

AliEMCALCell::~AliEMCALCell()
{
  // dtor
}

void AliEMCALCell::SetCCfromDB(AliEMCALCalibData *ccDb)
{
  if(ccDb == 0) return;
  // fADCchannelEC = fCalibData->GetADCchannel(iSupMod,ieta,iphi);
  // fetaCel-column; fPhiCell- row
  fCcIn = ccDb->GetADCchannel(fSupMod, fEtaCell, fPhiCell); // in GeV

  TH1* h = (TH1*)GetHists()->At(0);
  u::AddToNameAndTitle(h, 0, Form(", cc %5.2f MeV", fCcIn*1.e+3));
}

void AliEMCALCell::SetCCfromCCTable(AliEMCALCalibCoefs *t)
{
  if(t == 0) return;
  this->AddObject((TObject*)BookHists(), kTRUE);

  calibCoef *r = t->GetTable(fAbsId);
  if(r && r->absId == fAbsId) {
    fCcIn = r->cc;
  } else { // something wrong
    if(r) printf(" fAbsId %i : r->absId %i \n", fAbsId, r->absId);
    assert(0);
  }

  TH1* h = (TH1*)GetHists()->At(0);
  u::AddToNameAndTitle(h, 0, Form(", cc %5.2f MeV", fCcIn*1.e+3));
}

void AliEMCALCell::FillEffMass(const Double_t mgg)
{
  u::FillH1(GetHists(), 0, mgg);
}

void AliEMCALCell::FillCellNtuple(TNtuple *nt)
{
  if(nt==0) return;
  nt->Fill(fAbsId,fSupMod,fModule,fPhi,fEta,fPhiCell,fEtaCell,fCcIn,fCcOut);
}

void AliEMCALCell::FitHist(TH1* h, const char* name, const char* opt)
{
  TString optFit(""), OPT(opt);
  OPT.ToUpper();
  if(h==0) return; 
  printf("<I> AliEMCALCell::FitHist : h %p |%s| is started : opt %s\n", h, h->GetName(), opt);
  TString tit(h->GetTitle());

  TF1 *GausPol2 = 0, *g=0, *bg=0;
  if(h->GetListOfFunctions()->GetSize() == 0) {
    g = u::Gausi(name, 0.0, 0.4, h); // gaus estimation

    g->SetParLimits(0, h->Integral()/20., h->Integral());
    g->SetParLimits(1, 0.08, 0.20);
    g->SetParLimits(2, 0.001, 0.02);

    g->SetParameter(0, 1200.);
    g->SetParameter(1, h->GetBinLowEdge(h->GetMaximumBin()));
    g->SetParameter(2, 0.010);

    g->SetLineColor(kRed);
    g->SetLineWidth(2);

    optFit = "0NQ";
    h->Fit(g, optFit.Data(),"", 0.001, 0.3);

    optFit = "0NQIME";
    h->Fit(g, optFit.Data(),"", 0.001, 0.3);

    bg = new TF1(Form("bg%s",name), "pol2", 0.0, 0.3);  
    optFit = "0NQ";
    h->Fit(bg, optFit.Data(),"", 0.0, 0.3);

    GausPol2 = u::GausiPol2(name, 0.00, 0.3, g, bg);
    optFit = "0Q";
    h->Fit(GausPol2, optFit.Data(),"", 0.03, 0.28);
  // Clean up
    delete g;
    delete bg;
    optFit = "0IME+"; // no drwaing at all
    if(tit.Contains("SM") || OPT.Contains("DRAW")) optFit = "IME+";
  } else {
    GausPol2 = (TF1*)h->GetListOfFunctions()->At(0);
    optFit = "IME+";
  }
  //  optFit = "IME+";
  h->Fit(GausPol2, optFit.Data(),"", 0.01, 0.28);

  if(optFit.Contains("0") == 0) {
    gStyle->SetOptFit(111);
    u::DrawHist(h,2);
  }
  printf("<I> AliEMCALCell::FitHist : h %p |%s| is ended \n\n", h, h->GetName());
}

void AliEMCALCell::FitEffMassHist(const char* opt)
{
  AliEMCALFolder* EMCAL = (AliEMCALFolder*)(GetParent()->GetParent()->GetParent()); 
  Int_t it = EMCAL->GetIterationNumber();

  TH1* h = (TH1*)GetHists()->At(0);

  FitHist(h, GetName(), opt);

  fFun = (TF1*)h->GetListOfFunctions()->At(0);
  if(fFun) {
    Double_t mpi = fFun->GetParameter(1), mpi2 = mpi*mpi;
    Double_t ccTmp = fCcIn * MPI02 / mpi2;
    if(it<=6) {
      fCcOut = ccTmp;
    } else {
      fCcOut = (ccTmp + fCcIn)/2.;
    }
    printf(" fFun %s | %s : iet %i\n", fFun->GetName(), fFun->GetTitle(), it);
  }
  printf(" %s | fCcIn %6.5f -> % 6.5f <- fCcOut \n", GetTitle(), fCcIn , fCcOut);
}

void AliEMCALCell::Print()
{
  TH1* h = (TH1*)GetHists()->At(0);
  TF1 *f = (TF1*)h->GetListOfFunctions()->At(0);
  printf(" %s %s \n", GetName(), GetTitle());
  if(fFun) printf(" fFun : %s | %s \n", fFun->GetName(),  fFun->GetTitle());
  else fFun = f;
  if(f) f->Dump();
}

TList* AliEMCALCell::BookHists()
{
  gROOT->cd();
  TH1::AddDirectory(1);

  AliEMCALFolder* EMCAL = (AliEMCALFolder*)(GetParent()->GetParent()->GetParent()); 
  Int_t it = EMCAL->GetIterationNumber();

  new TH1F("01_EffMass", "effective mass of #gamma,#gamma(m_{#pi^{0}}=134.98 MeV) ", 60,0.0,0.3);

  TList *l = AliEMCALHistoUtilities::MoveHistsToList(Form("HistsOfCell%4.4i",fAbsId), kFALSE);
  AliEMCALHistoUtilities::AddToNameAndTitleToList(l, Form("4.4i_It%i",fAbsId,it),
						  Form(" Cell %4.4i, iter. %i",fAbsId, it));

  TH1::AddDirectory(0);
  return l;
}
Commit	Line	Data
16d3c94d	1	/**************************************************************************
	2	* Copyright(c) 1998-2007, 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	/* $Id:*/
	17
	18	//_________________________________________________________________________
	19	// Top EMCAL folder which will keep all information about EMCAL itself,
	20	// super Modules (SM), modules, towers, set of hists and so on.
	21	//
	22	//*-- Author: Aleksei Pavlinov (WSU, Detroit, USA)
	23
	24	#include "AliEMCALCell.h"
	25	#include "AliEMCALHistoUtilities.h"
	26	#include "AliEMCALGeometry.h"
	27	#include "AliEMCALFolder.h"
	28	#include "AliEMCALSuperModule.h"
	29	#include "AliEMCALCalibData.h"
	30
	31	#include "AliEMCALCalibCoefs.h"
	32
	33	#include <TROOT.h>
	34	#include <TStyle.h>
	35	#include <TList.h>
	36	#include <TH1.h>
	37	#include <TF1.h>
	38	#include <TNtuple.h>
	39	#include <TObjectSet.h>
	40
	41	typedef AliEMCALHistoUtilities u;
	42
	43	ClassImp(AliEMCALCell)
	44
	45	Double_t ADCCHANNELEC = 0.0153; // Update 24 Apr 2007: 250./16/1024 - width of one ADC channel in GeV
	46	Double_t MPI0 = 0.13498; // mass of pi0
	47	Double_t MPI02 = MPI0MPI0; // mass*2
	48
	49	AliEMCALCell::AliEMCALCell() :
	50	TObjectSet(),
	51	fAbsId(0),fSupMod(0),fModule(0),fPhi(0),fEta(0),fPhiCell(0),fEtaCell(0),fCcIn(0),fCcOut(0),
	52	fFun(0)
	53	{
	54	}
	55
	56	AliEMCALCell::AliEMCALCell(const Int_t absId, const char* title) :
	57	TObjectSet(Form("Cell%4.4i",absId)),
	58	fAbsId(absId),fSupMod(0),fModule(0),fPhi(0),fEta(0),fPhiCell(0),fEtaCell(0),fCcIn(0),fCcOut(0),
	59	fFun(0)
	60	{
	61	SetTitle(title);
	62
	63	AliEMCALGeometry *g = AliEMCALGeometry::GetInstance();
	64	g->GetCellIndex(fAbsId, fSupMod, fModule, fPhi, fEta);
65	g->GetCellPhiEtaIndexInSModule(fSupMod, fModule, fPhi, fEta, fPhiCell, fEtaCell);
66
67	}
68
69	AliEMCALCell::~AliEMCALCell()
70	{
71	// dtor
72	}
73
74	void AliEMCALCell::SetCCfromDB(AliEMCALCalibData *ccDb)
75	{
76	if(ccDb == 0) return;
77	// fADCchannelEC = fCalibData->GetADCchannel(iSupMod,ieta,iphi);
78	// fetaCel-column; fPhiCell- row
79	fCcIn = ccDb->GetADCchannel(fSupMod, fEtaCell, fPhiCell); // in GeV
80
81	TH1* h = (TH1*)GetHists()->At(0);
82	u::AddToNameAndTitle(h, 0, Form(", cc %5.2f MeV", fCcIn*1.e+3));
83	}
84
85	void AliEMCALCell::SetCCfromCCTable(AliEMCALCalibCoefs *t)
86	{
87	if(t == 0) return;
88	this->AddObject((TObject*)BookHists(), kTRUE);
89
90	calibCoef *r = t->GetTable(fAbsId);
91	if(r && r->absId == fAbsId) {
92	fCcIn = r->cc;
93	} else { // something wrong
94	if(r) printf(" fAbsId %i : r->absId %i \n", fAbsId, r->absId);
95	assert(0);
96	}
97
98	TH1* h = (TH1*)GetHists()->At(0);
99	u::AddToNameAndTitle(h, 0, Form(", cc %5.2f MeV", fCcIn*1.e+3));
100	}
101
102	void AliEMCALCell::FillEffMass(const Double_t mgg)
103	{
104	u::FillH1(GetHists(), 0, mgg);
105	}
106
107	void AliEMCALCell::FillCellNtuple(TNtuple *nt)
108	{
109	if(nt==0) return;
110	nt->Fill(fAbsId,fSupMod,fModule,fPhi,fEta,fPhiCell,fEtaCell,fCcIn,fCcOut);
111	}
112
113	void AliEMCALCell::FitHist(TH1* h, const char* name, const char* opt)
114	{
115	TString optFit(""), OPT(opt);
116	OPT.ToUpper();
117	if(h==0) return;
118	printf("<I> AliEMCALCell::FitHist : h %p \|%s\| is started : opt %s\n", h, h->GetName(), opt);
119	TString tit(h->GetTitle());
120
121	TF1 GausPol2 = 0, g=0, *bg=0;
122	if(h->GetListOfFunctions()->GetSize() == 0) {
123	g = u::Gausi(name, 0.0, 0.4, h); // gaus estimation
124
125	g->SetParLimits(0, h->Integral()/20., h->Integral());
126	g->SetParLimits(1, 0.08, 0.20);
127	g->SetParLimits(2, 0.001, 0.02);
128
129	g->SetParameter(0, 1200.);
130	g->SetParameter(1, h->GetBinLowEdge(h->GetMaximumBin()));
131	g->SetParameter(2, 0.010);
132
133	g->SetLineColor(kRed);
134	g->SetLineWidth(2);
135
136	optFit = "0NQ";
137	h->Fit(g, optFit.Data(),"", 0.001, 0.3);
138
139	optFit = "0NQIME";
140	h->Fit(g, optFit.Data(),"", 0.001, 0.3);
141
142	bg = new TF1(Form("bg%s",name), "pol2", 0.0, 0.3);
143	optFit = "0NQ";
144	h->Fit(bg, optFit.Data(),"", 0.0, 0.3);
145
146	GausPol2 = u::GausiPol2(name, 0.00, 0.3, g, bg);
147	optFit = "0Q";
148	h->Fit(GausPol2, optFit.Data(),"", 0.03, 0.28);
149	// Clean up
150	delete g;
151	delete bg;
152	optFit = "0IME+"; // no drwaing at all
153	if(tit.Contains("SM") \|\| OPT.Contains("DRAW")) optFit = "IME+";
154	} else {
155	GausPol2 = (TF1*)h->GetListOfFunctions()->At(0);
156	optFit = "IME+";
157	}
158	// optFit = "IME+";
159	h->Fit(GausPol2, optFit.Data(),"", 0.01, 0.28);
160
161	if(optFit.Contains("0") == 0) {
162	gStyle->SetOptFit(111);
163	u::DrawHist(h,2);
164	}
165	printf("<I> AliEMCALCell::FitHist : h %p \|%s\| is ended \n\n", h, h->GetName());
166	}
167
168	void AliEMCALCell::FitEffMassHist(const char* opt)
169	{
170	AliEMCALFolder* EMCAL = (AliEMCALFolder*)(GetParent()->GetParent()->GetParent());
171	Int_t it = EMCAL->GetIterationNumber();
172
173	TH1* h = (TH1*)GetHists()->At(0);
174
175	FitHist(h, GetName(), opt);
176
177	fFun = (TF1*)h->GetListOfFunctions()->At(0);
178	if(fFun) {
179	Double_t mpi = fFun->GetParameter(1), mpi2 = mpi*mpi;
180	Double_t ccTmp = fCcIn * MPI02 / mpi2;
181	if(it<=6) {
182	fCcOut = ccTmp;
183	} else {
184	fCcOut = (ccTmp + fCcIn)/2.;
185	}
186	printf(" fFun %s \| %s : iet %i\n", fFun->GetName(), fFun->GetTitle(), it);
187	}
188	printf(" %s \| fCcIn %6.5f -> % 6.5f <- fCcOut \n", GetTitle(), fCcIn , fCcOut);
189	}
190
191	void AliEMCALCell::Print()
192	{
193	TH1* h = (TH1*)GetHists()->At(0);
194	TF1 f = (TF1)h->GetListOfFunctions()->At(0);
195	printf(" %s %s \n", GetName(), GetTitle());
196	if(fFun) printf(" fFun : %s \| %s \n", fFun->GetName(), fFun->GetTitle());
197	else fFun = f;
198	if(f) f->Dump();
199	}
200
201	TList* AliEMCALCell::BookHists()
202	{
203	gROOT->cd();
204	TH1::AddDirectory(1);
205
206	AliEMCALFolder* EMCAL = (AliEMCALFolder*)(GetParent()->GetParent()->GetParent());
207	Int_t it = EMCAL->GetIterationNumber();
208
209	new TH1F("01_EffMass", "effective mass of #gamma,#gamma(m_{#pi^{0}}=134.98 MeV) ", 60,0.0,0.3);
210
211	TList *l = AliEMCALHistoUtilities::MoveHistsToList(Form("HistsOfCell%4.4i",fAbsId), kFALSE);
212	AliEMCALHistoUtilities::AddToNameAndTitleToList(l, Form("4.4i_It%i",fAbsId,it),
213	Form(" Cell %4.4i, iter. %i",fAbsId, it));
214
215	TH1::AddDirectory(0);
216	return l;
217	}