New calibration classes. They depend on TTable, so libTable.so is added to the list...
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALCell.cxx
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        = MPI0*MPI0; // 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 }