]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/CaloTrackCorrelations/AliAnaInsideClusterInvariantMass.cxx
Correct the calculation of cluster local maxima and cell neighbours plus coding violation
[u/mrichter/AliRoot.git] / PWGGA / CaloTrackCorrelations / AliAnaInsideClusterInvariantMass.cxx
CommitLineData
992b14a7 1/**************************************************************************
2 * Copyright(c) 1998-1999, 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//_________________________________________________________________________
17//
18// Split clusters with some criteria and calculate invariant mass
19// to identify them as pi0 or conversion
20//
21//
22//-- Author: Gustavo Conesa (LPSC-Grenoble)
23//_________________________________________________________________________
24
25//////////////////////////////////////////////////////////////////////////////
26
27
28// --- ROOT system ---
29#include <TList.h>
30#include <TClonesArray.h>
31#include <TObjString.h>
32#include <TH3F.h>
29ca9cad 33#include <TCanvas.h>
34#include <TStyle.h>
35#include <TPad.h>
992b14a7 36
37// --- Analysis system ---
38#include "AliAnaInsideClusterInvariantMass.h"
39#include "AliCaloTrackReader.h"
40#include "AliMCAnalysisUtils.h"
41#include "AliStack.h"
42#include "AliFiducialCut.h"
43#include "TParticle.h"
44#include "AliVCluster.h"
45#include "AliAODEvent.h"
46#include "AliAODMCParticle.h"
47#include "AliEMCALGeoParams.h"
48
c5693f62 49// --- Detectors ---
50//#include "AliPHOSGeoUtils.h"
51#include "AliEMCALGeometry.h"
52
992b14a7 53ClassImp(AliAnaInsideClusterInvariantMass)
54
55//__________________________________________________________________
56AliAnaInsideClusterInvariantMass::AliAnaInsideClusterInvariantMass() :
745913ae 57 AliAnaCaloTrackCorrBaseClass(),
29ca9cad 58 fCalorimeter(""),
71e3889f 59 fM02MaxCut(0), fM02MinCut(0),
60 fMinNCells(0),
243c2909 61 fMassEtaMin(0), fMassEtaMax(0),
62 fMassPi0Min(0), fMassPi0Max(0),
71e3889f 63 fMassConMin(0), fMassConMax(0),
64 fPlotCluster(0)
992b14a7 65{
66 //default ctor
67
68 // Init array of histograms
5c46c992 69 for(Int_t i = 0; i < 7; i++)
70 {
71 for(Int_t j = 0; j < 2; j++)
72 {
73
74 fhMassNLocMax1[i][j] = 0;
75 fhMassNLocMax2[i][j] = 0;
76 fhMassNLocMaxN[i][j] = 0;
77 fhNLocMax[i][j] = 0;
71e3889f 78 fhNLocMaxNLabel[i][j] = 0;
5c46c992 79 fhNLocMaxEMax[i][j] = 0;
80 fhNLocMaxEFrac[i][j] = 0;
81 fhNLocMaxM02Cut[i][j] = 0;
82 fhM02NLocMax1[i][j] = 0;
83 fhM02NLocMax2[i][j] = 0;
84 fhM02NLocMaxN[i][j] = 0;
85 fhNCellNLocMax1[i][j] = 0;
86 fhNCellNLocMax2[i][j] = 0;
87 fhNCellNLocMaxN[i][j] = 0;
88 fhM02Pi0LocMax1[i][j] = 0;
89 fhM02EtaLocMax1[i][j] = 0;
90 fhM02ConLocMax1[i][j] = 0;
91 fhM02Pi0LocMax2[i][j] = 0;
92 fhM02EtaLocMax2[i][j] = 0;
93 fhM02ConLocMax2[i][j] = 0;
94 fhM02Pi0LocMaxN[i][j] = 0;
95 fhM02EtaLocMaxN[i][j] = 0;
96 fhM02ConLocMaxN[i][j] = 0;
97
0137016b 98 fhMassM02NLocMax1[i][j]= 0;
99 fhMassM02NLocMax2[i][j]= 0;
100 fhMassM02NLocMaxN[i][j]= 0;
101
5c46c992 102 }
103
104 fhTrackMatchedDEtaLocMax1[i] = 0;
105 fhTrackMatchedDPhiLocMax1[i] = 0;
106 fhTrackMatchedDEtaLocMax2[i] = 0;
107 fhTrackMatchedDPhiLocMax2[i] = 0;
108 fhTrackMatchedDEtaLocMaxN[i] = 0;
109 fhTrackMatchedDPhiLocMaxN[i] = 0;
243c2909 110
992b14a7 111 }
112
5c46c992 113 for(Int_t i = 0; i < 2; i++)
114 {
115 fhAnglePairLocMax1 [i] = 0;
116 fhAnglePairLocMax2 [i] = 0;
117 fhAnglePairLocMaxN [i] = 0;
118 fhAnglePairMassLocMax1[i] = 0;
119 fhAnglePairMassLocMax2[i] = 0;
120 fhAnglePairMassLocMaxN[i] = 0;
121 }
122
7b686344 123 for(Int_t i = 0; i < 4; i++)
124 {
125 fhMassM02NLocMax1Ebin[i] = 0 ;
126 fhMassM02NLocMax2Ebin[i] = 0 ;
127 fhMassM02NLocMaxNEbin[i] = 0 ;
128 }
129
992b14a7 130 InitParameters();
131
132}
133
134//_______________________________________________________________
135TObjString * AliAnaInsideClusterInvariantMass::GetAnalysisCuts()
136{
137 //Save parameters used for analysis
138 TString parList ; //this will be list of parameters used for this analysis.
139 const Int_t buffersize = 255;
140 char onePar[buffersize] ;
141
142 snprintf(onePar,buffersize,"--- AliAnaInsideClusterInvariantMass ---\n") ;
143 parList+=onePar ;
144
243c2909 145 snprintf(onePar,buffersize,"Calorimeter: %s\n", fCalorimeter.Data()) ;
992b14a7 146 parList+=onePar ;
dbba06ca 147 snprintf(onePar,buffersize,"fLocMaxCutE =%2.2f \n", GetCaloUtils()->GetLocalMaximaCutE()) ;
29ca9cad 148 parList+=onePar ;
dbba06ca 149 snprintf(onePar,buffersize,"fLocMaxCutEDiff =%2.2f \n",GetCaloUtils()->GetLocalMaximaCutEDiff()) ;
29ca9cad 150 parList+=onePar ;
71e3889f 151 snprintf(onePar,buffersize,"%2.2f< M02 < %2.2f \n", fM02MinCut, fM02MaxCut) ;
992b14a7 152 parList+=onePar ;
e23a0471 153 snprintf(onePar,buffersize,"fMinNCells =%d \n", fMinNCells) ;
992b14a7 154 parList+=onePar ;
243c2909 155 snprintf(onePar,buffersize,"pi0 : %2.1f < m <%2.1f\n", fMassPi0Min,fMassPi0Max);
156 parList+=onePar ;
157 snprintf(onePar,buffersize,"eta : %2.1f < m <%2.1f\n", fMassEtaMin,fMassEtaMax);
158 parList+=onePar ;
159 snprintf(onePar,buffersize,"conv: %2.1f < m <%2.1f\n", fMassConMin,fMassConMax);
160 parList+=onePar ;
161
992b14a7 162 return new TObjString(parList) ;
163
164}
165
992b14a7 166
167//_____________________________________________________________________________________
168TLorentzVector AliAnaInsideClusterInvariantMass::GetCellMomentum(const Int_t absId,
243c2909 169 Float_t en,
992b14a7 170 AliVCaloCells * cells)
171{
172
173 // Cell momentum calculation for invariant mass
174
175 Double_t cellpos[] = {0, 0, 0};
176 GetEMCALGeometry()->GetGlobal(absId, cellpos);
177
178 if(GetVertex(0)){//calculate direction from vertex
179 cellpos[0]-=GetVertex(0)[0];
180 cellpos[1]-=GetVertex(0)[1];
181 cellpos[2]-=GetVertex(0)[2];
182 }
183
243c2909 184 Double_t r = TMath::Sqrt(cellpos[0]*cellpos[0]+cellpos[1]*cellpos[1]+cellpos[2]*cellpos[2] ) ;
185
186 //If not calculated before, get the energy
29ca9cad 187 if(en <=0 )
243c2909 188 {
189 en = cells->GetCellAmplitude(absId);
dbba06ca 190 GetCaloUtils()->RecalibrateCellAmplitude(en,fCalorimeter,absId);
243c2909 191 }
992b14a7 192
992b14a7 193 TLorentzVector cellMom ;
194 cellMom.SetPxPyPzE( en*cellpos[0]/r, en*cellpos[1]/r, en*cellpos[2]/r, en) ;
195
196 return cellMom;
197
198}
199
200//________________________________________________________________
201TList * AliAnaInsideClusterInvariantMass::GetCreateOutputObjects()
202{
203 // Create histograms to be saved in output file and
204 // store them in outputContainer
205 TList * outputContainer = new TList() ;
206 outputContainer->SetName("InsideClusterHistos") ;
207
745913ae 208 Int_t nptbins = GetHistogramRanges()->GetHistoPtBins(); Float_t ptmax = GetHistogramRanges()->GetHistoPtMax(); Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
209 Int_t ssbins = GetHistogramRanges()->GetHistoShowerShapeBins(); Float_t ssmax = GetHistogramRanges()->GetHistoShowerShapeMax(); Float_t ssmin = GetHistogramRanges()->GetHistoShowerShapeMin();
210 Int_t mbins = GetHistogramRanges()->GetHistoMassBins(); Float_t mmax = GetHistogramRanges()->GetHistoMassMax(); Float_t mmin = GetHistogramRanges()->GetHistoMassMin();
211 Int_t ncbins = GetHistogramRanges()->GetHistoNClusterCellBins(); Int_t ncmax = GetHistogramRanges()->GetHistoNClusterCellMax(); Int_t ncmin = GetHistogramRanges()->GetHistoNClusterCellMin();
992b14a7 212
5c46c992 213 Int_t nresetabins = GetHistogramRanges()->GetHistoTrackResidualEtaBins();
214 Float_t resetamax = GetHistogramRanges()->GetHistoTrackResidualEtaMax();
215 Float_t resetamin = GetHistogramRanges()->GetHistoTrackResidualEtaMin();
216 Int_t nresphibins = GetHistogramRanges()->GetHistoTrackResidualPhiBins();
217 Float_t resphimax = GetHistogramRanges()->GetHistoTrackResidualPhiMax();
218 Float_t resphimin = GetHistogramRanges()->GetHistoTrackResidualPhiMin();
219
992b14a7 220 TString ptype[] ={"","#gamma","#gamma->e^{#pm}","#pi^{0}","#eta","e^{#pm}", "hadron"};
221 TString pname[] ={"","Photon","Conversion", "Pi0", "Eta", "Electron","Hadron"};
222
223 Int_t n = 1;
224
225 if(IsDataMC()) n = 7;
226
243c2909 227 Int_t nMaxBins = 10;
228
5c46c992 229 TString sMatched[] = {"","Matched"};
230
231 for(Int_t i = 0; i < n; i++)
232 {
233
234 for(Int_t j = 0; j < 2; j++)
235 {
236
237 fhMassNLocMax1[i][j] = new TH2F(Form("hMassNLocMax1%s%s",pname[i].Data(),sMatched[j].Data()),
0137016b 238 Form("Invariant mass of 2 highest energy cells vs E, %s %s",ptype[i].Data(),sMatched[j].Data()),
5c46c992 239 nptbins,ptmin,ptmax,mbins,mmin,mmax);
240 fhMassNLocMax1[i][j]->SetYTitle("M (GeV/c^{2})");
241 fhMassNLocMax1[i][j]->SetXTitle("E (GeV)");
242 outputContainer->Add(fhMassNLocMax1[i][j]) ;
243
244 fhMassNLocMax2[i][j] = new TH2F(Form("hMassNLocMax2%s%s",pname[i].Data(),sMatched[j].Data()),
0137016b 245 Form("Invariant mass of 2 local maxima cells vs E,%s %s",ptype[i].Data(),sMatched[j].Data()),
5c46c992 246 nptbins,ptmin,ptmax,mbins,mmin,mmax);
247 fhMassNLocMax2[i][j]->SetYTitle("M (GeV/c^{2})");
248 fhMassNLocMax2[i][j]->SetXTitle("E (GeV)");
249 outputContainer->Add(fhMassNLocMax2[i][j]) ;
250
251 fhMassNLocMaxN[i][j] = new TH2F(Form("hMassNLocMaxN%s%s",pname[i].Data(),sMatched[j].Data()),
0137016b 252 Form("Invariant mass of N>2 local maxima cells vs E, %s %s",ptype[i].Data(),sMatched[j].Data()),
5c46c992 253 nptbins,ptmin,ptmax,mbins,mmin,mmax);
254 fhMassNLocMaxN[i][j]->SetYTitle("M (GeV/c^{2})");
255 fhMassNLocMaxN[i][j]->SetXTitle("E (GeV)");
256 outputContainer->Add(fhMassNLocMaxN[i][j]) ;
257
0137016b 258 fhMassM02NLocMax1[i][j] = new TH2F(Form("hMassM02NLocMax1%s%s",pname[i].Data(),sMatched[j].Data()),
259 Form("Invariant mass of 2 highest energy cells #lambda_{0}^{2}, E > 7 GeV,%s %s",ptype[i].Data(),sMatched[j].Data()),
260 ssbins,ssmin,ssmax,mbins,mmin,mmax);
261 fhMassM02NLocMax1[i][j]->SetYTitle("M (GeV/c^{2})");
262 fhMassM02NLocMax1[i][j]->SetXTitle("#lambda_{0}^{2}");
263 outputContainer->Add(fhMassM02NLocMax1[i][j]) ;
264
265 fhMassM02NLocMax2[i][j] = new TH2F(Form("hMassM02NLocMax2%s%s",pname[i].Data(),sMatched[j].Data()),
266 Form("Invariant mass of 2 local maxima cells #lambda_{0}^{2}, E > 7 GeV, %s %s",ptype[i].Data(),sMatched[j].Data()),
267 ssbins,ssmin,ssmax,mbins,mmin,mmax);
268 fhMassM02NLocMax2[i][j]->SetYTitle("M (GeV/c^{2})");
269 fhMassM02NLocMax2[i][j]->SetXTitle("#lambda_{0}^{2}");
270 outputContainer->Add(fhMassM02NLocMax2[i][j]) ;
271
272 fhMassM02NLocMaxN[i][j] = new TH2F(Form("hMassM02NLocMaxN%s%s",pname[i].Data(),sMatched[j].Data()),
273 Form("Invariant mass of N>2 local maxima cells vs #lambda_{0}^{2}, %s %s",ptype[i].Data(),sMatched[j].Data()),
274 ssbins,ssmin,ssmax,mbins,mmin,mmax);
275 fhMassM02NLocMaxN[i][j]->SetYTitle("M (GeV/c^{2})");
276 fhMassM02NLocMaxN[i][j]->SetXTitle("#lambda_{0}^{2}");
277 outputContainer->Add(fhMassM02NLocMaxN[i][j]) ;
278
5c46c992 279
280 fhNLocMax[i][j] = new TH2F(Form("hNLocMax%s%s",pname[i].Data(),sMatched[j].Data()),
281 Form("Number of local maxima in cluster %s %s",ptype[i].Data(),sMatched[j].Data()),
282 nptbins,ptmin,ptmax,nMaxBins,0,nMaxBins);
283 fhNLocMax[i][j] ->SetYTitle("N maxima");
284 fhNLocMax[i][j] ->SetXTitle("E (GeV)");
285 outputContainer->Add(fhNLocMax[i][j]) ;
286
71e3889f 287 if(IsDataMC())
288 {
289 fhNLocMaxNLabel[i][j] = new TH2F(Form("hNLocMaxNLabel%s%s",pname[i].Data(),sMatched[j].Data()),
290 Form("Number of local maxima in cluster vs number of MC labels %s %s",ptype[i].Data(),sMatched[j].Data()),
291 nMaxBins,0,nMaxBins,nMaxBins,0,nMaxBins);
292 fhNLocMaxNLabel[i][j] ->SetYTitle("N maxima");
293 fhNLocMaxNLabel[i][j] ->SetXTitle("N MC labels");
294 outputContainer->Add(fhNLocMaxNLabel[i][j]) ;
295 }
296
5c46c992 297 fhNLocMaxEMax[i][j] = new TH2F(Form("hNLocMaxEMax%s%s",pname[i].Data(),sMatched[j].Data()),
298 Form("Number of local maxima in cluster vs energy of maxima %s %s",ptype[i].Data(),sMatched[j].Data()),
299 nptbins*10,ptmin,ptmax,nMaxBins,0,nMaxBins);
300 fhNLocMaxEMax[i][j] ->SetYTitle("N maxima");
301 fhNLocMaxEMax[i][j] ->SetXTitle("E of maxima (GeV)");
302 outputContainer->Add(fhNLocMaxEMax[i][j]) ;
303
304 fhNLocMaxEFrac[i][j] = new TH2F(Form("hNLocMaxEFrac%s%s",pname[i].Data(),sMatched[j].Data()),
305 Form("Number of local maxima in cluster vs fraction of cluster energy of maxima %s %s",ptype[i].Data(),sMatched[j].Data()),
306 100,0,1,nMaxBins,0,nMaxBins);
307 fhNLocMaxEFrac[i][j] ->SetYTitle("N maxima");
308 fhNLocMaxEFrac[i][j] ->SetXTitle("E maxima / E cluster");
309 outputContainer->Add(fhNLocMaxEFrac[i][j]) ;
310
311 fhNLocMaxM02Cut[i][j] = new TH2F(Form("hNLocMaxM02Cut%s%s",pname[i].Data(),sMatched[j].Data()),
71e3889f 312 Form("Number of local maxima in cluster %s for %2.2f < M02 < %2.2f",ptype[i].Data(),fM02MinCut,fM02MaxCut),
5c46c992 313 nptbins,ptmin,ptmax,nMaxBins,0,nMaxBins);
314 fhNLocMaxM02Cut[i][j]->SetYTitle("N maxima");
315 fhNLocMaxM02Cut[i][j]->SetXTitle("E (GeV)");
316 outputContainer->Add(fhNLocMaxM02Cut[i][j]) ;
317
318
319 fhM02NLocMax1[i][j] = new TH2F(Form("hM02NLocMax1%s%s",pname[i].Data(),sMatched[j].Data()),
320 Form("#lambda_{0}^{2} vs E for N max = 1 %s %s",ptype[i].Data(),sMatched[j].Data()),
321 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
322 fhM02NLocMax1[i][j] ->SetYTitle("#lambda_{0}^{2}");
323 fhM02NLocMax1[i][j] ->SetXTitle("E (GeV)");
324 outputContainer->Add(fhM02NLocMax1[i][j]) ;
325
326 fhM02NLocMax2[i][j] = new TH2F(Form("hM02NLocMax2%s%s",pname[i].Data(),sMatched[j].Data()),
327 Form("#lambda_{0}^{2} vs E for N max = 2 %s %s",ptype[i].Data(),sMatched[j].Data()),
328 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
329 fhM02NLocMax2[i][j] ->SetYTitle("#lambda_{0}^{2}");
330 fhM02NLocMax2[i][j] ->SetXTitle("E (GeV)");
331 outputContainer->Add(fhM02NLocMax2[i][j]) ;
332
333
334 fhM02NLocMaxN[i][j] = new TH2F(Form("hM02NLocMaxN%s%s",pname[i].Data(),sMatched[j].Data()),
335 Form("#lambda_{0}^{2} vs E for N max > 2 %s %s",ptype[i].Data(),sMatched[j].Data()),
336 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
337 fhM02NLocMaxN[i][j] ->SetYTitle("#lambda_{0}^{2}");
338 fhM02NLocMaxN[i][j] ->SetXTitle("E (GeV)");
339 outputContainer->Add(fhM02NLocMaxN[i][j]) ;
340
341
342 fhNCellNLocMax1[i][j] = new TH2F(Form("hNCellNLocMax1%s%s",pname[i].Data(),sMatched[j].Data()),
343 Form("#lambda_{0}^{2} vs E for N max = 1 %s %s",ptype[i].Data(),sMatched[j].Data()),
344 nptbins,ptmin,ptmax,ncbins,ncmin,ncmax);
345 fhNCellNLocMax1[i][j] ->SetYTitle("N cells");
346 fhNCellNLocMax1[i][j] ->SetXTitle("E (GeV)");
347 outputContainer->Add(fhNCellNLocMax1[i][j]) ;
348
349 fhNCellNLocMax2[i][j] = new TH2F(Form("hNCellNLocMax2%s%s",pname[i].Data(),sMatched[j].Data()),
350 Form("#lambda_{0}^{2} vs E for N max = 2 %s %s",ptype[i].Data(),sMatched[j].Data()),
351 nptbins,ptmin,ptmax,ncbins,ncmin,ncmax);
352 fhNCellNLocMax2[i][j] ->SetYTitle("N cells");
353 fhNCellNLocMax2[i][j] ->SetXTitle("E (GeV)");
354 outputContainer->Add(fhNCellNLocMax2[i][j]) ;
355
356
357 fhNCellNLocMaxN[i][j] = new TH2F(Form("hNCellNLocMaxN%s%s",pname[i].Data(),sMatched[j].Data()),
358 Form("#lambda_{0}^{2} vs E for N max > 2 %s %s",ptype[i].Data(),sMatched[j].Data()),
359 nptbins,ptmin,ptmax,ncbins,ncmin,ncmax);
360 fhNCellNLocMaxN[i][j] ->SetYTitle("N cells");
361 fhNCellNLocMaxN[i][j] ->SetXTitle("E (GeV)");
362 outputContainer->Add(fhNCellNLocMaxN[i][j]) ;
363
364
365 fhM02Pi0LocMax1[i][j] = new TH2F(Form("hM02Pi0LocMax1%s%s",pname[i].Data(),sMatched[j].Data()),
366 Form("#lambda_{0}^{2} vs E for mass range [%2.2f-%2.2f] MeV/c^{2} %s, for N Local max = 1",fMassPi0Min,fMassPi0Max,ptype[i].Data()),
367 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
368 fhM02Pi0LocMax1[i][j] ->SetYTitle("#lambda_{0}^{2}");
369 fhM02Pi0LocMax1[i][j] ->SetXTitle("E (GeV)");
370 outputContainer->Add(fhM02Pi0LocMax1[i][j]) ;
371
372 fhM02EtaLocMax1[i][j] = new TH2F(Form("hM02EtaLocMax1%s%s",pname[i].Data(),sMatched[j].Data()),
373 Form("#lambda_{0}^{2} vs E for mass range [%2.2f-%2.2f] MeV/c^{2}, %s, for N Local max = 1",fMassEtaMin,fMassEtaMax,ptype[i].Data()),
374 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
375 fhM02EtaLocMax1[i][j] ->SetYTitle("#lambda_{0}^{2}");
376 fhM02EtaLocMax1[i][j] ->SetXTitle("E (GeV)");
377 outputContainer->Add(fhM02EtaLocMax1[i][j]) ;
378
379 fhM02ConLocMax1[i][j] = new TH2F(Form("hM02ConLocMax1%s%s",pname[i].Data(),sMatched[j].Data()),
380 Form("#lambda_{0}^{2} vs E for mass range [%2.2f-%2.2f] MeV/c^{2}, %s, for N Local max = 1",fMassConMin,fMassConMax,ptype[i].Data()),
381 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
382 fhM02ConLocMax1[i][j] ->SetYTitle("#lambda_{0}^{2}");
383 fhM02ConLocMax1[i][j] ->SetXTitle("E (GeV)");
384 outputContainer->Add(fhM02ConLocMax1[i][j]) ;
385
386 fhM02Pi0LocMax2[i][j] = new TH2F(Form("hM02Pi0LocMax2%s%s",pname[i].Data(),sMatched[j].Data()),
387 Form("#lambda_{0}^{2} vs E for mass range [%2.2f-%2.2f] MeV/c^{2} %s, for N Local max = 2",fMassPi0Min,fMassPi0Max,ptype[i].Data()),
388 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
389 fhM02Pi0LocMax2[i][j] ->SetYTitle("#lambda_{0}^{2}");
390 fhM02Pi0LocMax2[i][j] ->SetXTitle("E (GeV)");
391 outputContainer->Add(fhM02Pi0LocMax2[i][j]) ;
392
393 fhM02EtaLocMax2[i][j] = new TH2F(Form("hM02EtaLocMax2%s%s",pname[i].Data(),sMatched[j].Data()),
394 Form("#lambda_{0}^{2} vs E for mass range [%2.2f-%2.2f] MeV/c^{2}, %s, for N Local max = 2",fMassEtaMin,fMassEtaMax,ptype[i].Data()),
395 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
396 fhM02EtaLocMax2[i][j] ->SetYTitle("#lambda_{0}^{2}");
397 fhM02EtaLocMax2[i][j] ->SetXTitle("E (GeV)");
398 outputContainer->Add(fhM02EtaLocMax2[i][j]) ;
399
400 fhM02ConLocMax2[i][j] = new TH2F(Form("hM02ConLocMax2%s%s",pname[i].Data(),sMatched[j].Data()),
401 Form("#lambda_{0}^{2} vs E for mass range [%2.2f-%2.2f] MeV/c^{2}, %s, for N Local max = 2",fMassConMin,fMassConMax,ptype[i].Data()),
402 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
403 fhM02ConLocMax2[i][j] ->SetYTitle("#lambda_{0}^{2}");
404 fhM02ConLocMax2[i][j] ->SetXTitle("E (GeV)");
405 outputContainer->Add(fhM02ConLocMax2[i][j]) ;
406
407 fhM02Pi0LocMaxN[i][j] = new TH2F(Form("hM02Pi0LocMaxN%s%s",pname[i].Data(),sMatched[j].Data()),
408 Form("#lambda_{0}^{2} vs E for mass range [%2.2f-%2.2f] MeV/c^{2} %s, for N Local max > 2",fMassPi0Min,fMassPi0Max,ptype[i].Data()),
409 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
410 fhM02Pi0LocMaxN[i][j] ->SetYTitle("#lambda_{0}^{2}");
411 fhM02Pi0LocMaxN[i][j] ->SetXTitle("E (GeV)");
412 outputContainer->Add(fhM02Pi0LocMaxN[i][j]) ;
413
414 fhM02EtaLocMaxN[i][j] = new TH2F(Form("hM02EtaLocMaxN%s%s",pname[i].Data(),sMatched[j].Data()),
415 Form("#lambda_{0}^{2} vs E for mass range [%2.2f-%2.2f] MeV/c^{2}, %s, for N Local max > 2", fMassEtaMin,fMassEtaMax,ptype[i].Data()),
416 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
417 fhM02EtaLocMaxN[i][j] ->SetYTitle("#lambda_{0}^{2}");
418 fhM02EtaLocMaxN[i][j] ->SetXTitle("E (GeV)");
419 outputContainer->Add(fhM02EtaLocMaxN[i][j]) ;
420
421 fhM02ConLocMaxN[i][j] = new TH2F(Form("hM02ConLocMaxN%s%s",pname[i].Data(),sMatched[j].Data()),
422 Form("#lambda_{0}^{2} vs E for mass range [%2.2f-%2.2f], %s, for N Local max > 2",fMassConMin,fMassConMax,ptype[i].Data()),
423 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
424 fhM02ConLocMaxN[i][j] ->SetYTitle("#lambda_{0}^{2}");
425 fhM02ConLocMaxN[i][j] ->SetXTitle("E (GeV)");
426 outputContainer->Add(fhM02ConLocMaxN[i][j]) ;
427
428 } // matched, not matched
429
7b686344 430
5c46c992 431 } // MC particle list
432
7b686344 433 for(Int_t i = 0; i < 4; i++)
434 {
435 fhMassM02NLocMax1Ebin[i] = new TH2F(Form("hMassM02NLocMax1Ebin%d",i),
436 Form("Invariant mass of 2 highest energy cells #lambda_{0}^{2}, E bin %d",i),
437 ssbins,ssmin,ssmax,mbins,mmin,mmax);
438 fhMassM02NLocMax1Ebin[i]->SetYTitle("M (GeV/c^{2})");
439 fhMassM02NLocMax1Ebin[i]->SetXTitle("#lambda_{0}^{2}");
440 outputContainer->Add(fhMassM02NLocMax1Ebin[i]) ;
441
442 fhMassM02NLocMax2Ebin[i] = new TH2F(Form("hMassM02NLocMax2Ebin%d",i),
443 Form("Invariant mass of 2 local maxima cells #lambda_{0}^{2}, E bin %d",i),
444 ssbins,ssmin,ssmax,mbins,mmin,mmax);
445 fhMassM02NLocMax2Ebin[i]->SetYTitle("M (GeV/c^{2})");
446 fhMassM02NLocMax2Ebin[i]->SetXTitle("#lambda_{0}^{2}");
447 outputContainer->Add(fhMassM02NLocMax2Ebin[i]) ;
448
449 fhMassM02NLocMaxNEbin[i] = new TH2F(Form("hMassM02NLocMaxNEbin%d",i),
450 Form("Invariant mass of N>2 local maxima cells vs #lambda_{0}^{2}, E bin %d",i),
451 ssbins,ssmin,ssmax,mbins,mmin,mmax);
452 fhMassM02NLocMaxNEbin[i]->SetYTitle("M (GeV/c^{2})");
453 fhMassM02NLocMaxNEbin[i]->SetXTitle("#lambda_{0}^{2}");
454 outputContainer->Add(fhMassM02NLocMaxNEbin[i]) ;
455 }
456
5c46c992 457 for(Int_t i = 0; i < n; i++)
458 {
459 fhTrackMatchedDEtaLocMax1[i] = new TH2F
460 (Form("hTrackMatchedDEtaLocMax1%s",pname[i].Data()),
461 Form("d#eta of cluster-track vs cluster energy, 1 Local Maxima, %s",ptype[i].Data()),
462 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
463 fhTrackMatchedDEtaLocMax1[i]->SetYTitle("d#eta");
464 fhTrackMatchedDEtaLocMax1[i]->SetXTitle("E_{cluster} (GeV)");
465
466 fhTrackMatchedDPhiLocMax1[i] = new TH2F
467 (Form("hTrackMatchedDPhiLocMax1%s",pname[i].Data()),
468 Form("d#phi of cluster-track vs cluster energy, 1 Local Maxima, %s",ptype[i].Data()),
469 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
470 fhTrackMatchedDPhiLocMax1[i]->SetYTitle("d#phi (rad)");
471 fhTrackMatchedDPhiLocMax1[i]->SetXTitle("E_{cluster} (GeV)");
472
473 outputContainer->Add(fhTrackMatchedDEtaLocMax1[i]) ;
474 outputContainer->Add(fhTrackMatchedDPhiLocMax1[i]) ;
475
476 fhTrackMatchedDEtaLocMax2[i] = new TH2F
477 (Form("hTrackMatchedDEtaLocMax2%s",pname[i].Data()),
478 Form("d#eta of cluster-track vs cluster energy, 2 Local Maxima, %s",ptype[i].Data()),
479 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
480 fhTrackMatchedDEtaLocMax2[i]->SetYTitle("d#eta");
481 fhTrackMatchedDEtaLocMax2[i]->SetXTitle("E_{cluster} (GeV)");
482
483 fhTrackMatchedDPhiLocMax2[i] = new TH2F
484 (Form("hTrackMatchedDPhiLocMax2%s",pname[i].Data()),
485 Form("d#phi of cluster-track vs cluster energy, 2 Local Maxima, %s",ptype[i].Data()),
486 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
487 fhTrackMatchedDPhiLocMax2[i]->SetYTitle("d#phi (rad)");
488 fhTrackMatchedDPhiLocMax2[i]->SetXTitle("E_{cluster} (GeV)");
489
490 outputContainer->Add(fhTrackMatchedDEtaLocMax2[i]) ;
491 outputContainer->Add(fhTrackMatchedDPhiLocMax2[i]) ;
492
493 fhTrackMatchedDEtaLocMaxN[i] = new TH2F
494 (Form("hTrackMatchedDEtaLocMaxN%s",pname[i].Data()),
495 Form("d#eta of cluster-track vs cluster energy, N>2 Local Maxima, %s",ptype[i].Data()),
496 nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
497 fhTrackMatchedDEtaLocMaxN[i]->SetYTitle("d#eta");
498 fhTrackMatchedDEtaLocMaxN[i]->SetXTitle("E_{cluster} (GeV)");
499
500 fhTrackMatchedDPhiLocMaxN[i] = new TH2F
501 (Form("hTrackMatchedDPhiLocMaxN%s",pname[i].Data()),
502 Form("d#phi of cluster-track vs cluster energy, N>2 Local Maxima, %s",ptype[i].Data()),
503 nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
504 fhTrackMatchedDPhiLocMaxN[i]->SetYTitle("d#phi (rad)");
505 fhTrackMatchedDPhiLocMaxN[i]->SetXTitle("E_{cluster} (GeV)");
506
507 outputContainer->Add(fhTrackMatchedDEtaLocMaxN[i]) ;
508 outputContainer->Add(fhTrackMatchedDPhiLocMaxN[i]) ;
509
510 }
511
512 for(Int_t j = 0; j < 2; j++)
513 {
514
515 fhAnglePairLocMax1[j] = new TH2F(Form("hAnglePairLocMax1%s",sMatched[j].Data()),
516 Form("Opening angle of 2 highest energy cells vs pair Energy, %s",sMatched[j].Data()),
517 nptbins,ptmin,ptmax,200,0,0.2);
518 fhAnglePairLocMax1[j]->SetYTitle("#alpha (rad)");
519 fhAnglePairLocMax1[j]->SetXTitle("E (GeV)");
520 outputContainer->Add(fhAnglePairLocMax1[j]) ;
521
522 fhAnglePairLocMax2[j] = new TH2F(Form("hAnglePairLocMax2%s",sMatched[j].Data()),
523 Form("Opening angle of 2 local maxima cells vs Energy, %s",sMatched[j].Data()),
524 nptbins,ptmin,ptmax,200,0,0.2);
525 fhAnglePairLocMax2[j]->SetYTitle("#alpha (rad)");
526 fhAnglePairLocMax2[j]->SetXTitle("E (GeV)");
527 outputContainer->Add(fhAnglePairLocMax2[j]) ;
528
529 fhAnglePairLocMaxN[j] = new TH2F(Form("hAnglePairLocMaxN%s",sMatched[j].Data()),
530 Form("Opening angle of N>2 local maxima cells vs Energy, %s",sMatched[j].Data()),
531 nptbins,ptmin,ptmax,200,0,0.2);
532 fhAnglePairLocMaxN[j]->SetYTitle("#alpha (rad)");
533 fhAnglePairLocMaxN[j]->SetXTitle("E (GeV)");
534 outputContainer->Add(fhAnglePairLocMaxN[j]) ;
535
536 fhAnglePairMassLocMax1[j] = new TH2F(Form("hAnglePairMassLocMax1%s",sMatched[j].Data()),
0137016b 537 Form("Opening angle of 2 highest energy cells vs Mass for E > 7 GeV, %s",sMatched[j].Data()),
5c46c992 538 mbins,mmin,mmax,200,0,0.2);
539 fhAnglePairMassLocMax1[j]->SetXTitle("M (GeV/c^{2})");
540 fhAnglePairMassLocMax1[j]->SetYTitle("#alpha (rad)");
541 outputContainer->Add(fhAnglePairMassLocMax1[j]) ;
542
543 fhAnglePairMassLocMax2[j] = new TH2F(Form("hAnglePairMassLocMax2%s",sMatched[j].Data()),
0137016b 544 Form("Opening angle of 2 local maxima cells vs Mass for E > 7 GeV, %s",sMatched[j].Data()),
5c46c992 545 mbins,mmin,mmax,200,0,0.2);
546 fhAnglePairMassLocMax2[j]->SetXTitle("M (GeV/c^{2})");
547 fhAnglePairMassLocMax2[j]->SetYTitle("#alpha (rad)");
548 outputContainer->Add(fhAnglePairMassLocMax2[j]) ;
549
550 fhAnglePairMassLocMaxN[j] = new TH2F(Form("hAnglePairMassLocMaxN%s",sMatched[j].Data()),
0137016b 551 Form("Opening angle of N>2 local maxima cells vs Mass for E > 7 GeV, %s",sMatched[j].Data()),
5c46c992 552 mbins,mmin,mmax,200,0,0.2);
553 fhAnglePairMassLocMaxN[j]->SetXTitle("M (GeV/c^{2})");
554 fhAnglePairMassLocMaxN[j]->SetYTitle("#alpha (rad)");
555 outputContainer->Add(fhAnglePairMassLocMaxN[j]) ;
992b14a7 556
557 }
558
559 return outputContainer ;
560
561}
562
992b14a7 563//___________________________________________
564void AliAnaInsideClusterInvariantMass::Init()
565{
566 //Init
567 //Do some checks
568 if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD()){
569 printf("AliAnaInsideClusterInvariantMass::Init() - !!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!\n");
570 abort();
571 }
572 else if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD()){
573 printf("AliAnaInsideClusterInvariantMass::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");
574 abort();
575 }
576
577 if( GetReader()->GetDataType() == AliCaloTrackReader::kMC ){
578 printf("AliAnaInsideClusterInvariantMass::Init() - !!STOP: You want to use pure MC data!!\n");
579 abort();
580
581 }
582
583}
584
585//_____________________________________________________
586void AliAnaInsideClusterInvariantMass::InitParameters()
587{
588 //Initialize the parameters of the analysis.
589 AddToHistogramsName("AnaPi0InsideClusterInvariantMass_");
590
591 fCalorimeter = "EMCAL" ;
29ca9cad 592
71e3889f 593 fM02MinCut = 0.26 ;
594 fM02MaxCut = 10 ;
595
992b14a7 596 fMinNCells = 4 ;
243c2909 597
29ca9cad 598 fMassEtaMin = 0.4;
599 fMassEtaMax = 0.6;
243c2909 600
29ca9cad 601 fMassPi0Min = 0.08;
602 fMassPi0Max = 0.20;
243c2909 603
29ca9cad 604 fMassConMin = 0.0;
605 fMassConMax = 0.05;
243c2909 606
992b14a7 607}
608
609
610//__________________________________________________________________
611void AliAnaInsideClusterInvariantMass::MakeAnalysisFillHistograms()
612{
613 //Search for pi0 in fCalorimeter with shower shape analysis
614
615 TObjArray * pl = 0x0;
616 AliVCaloCells* cells = 0x0;
617
618 //Select the Calorimeter of the photon
619 if(fCalorimeter == "PHOS"){
620 pl = GetPHOSClusters();
621 cells = GetPHOSCells();
622 }
623 else if (fCalorimeter == "EMCAL"){
624 pl = GetEMCALClusters();
625 cells = GetEMCALCells();
626 }
627
628 if(!pl || !cells) {
629 Info("MakeAnalysisFillHistograms","TObjArray with %s clusters is NULL!\n",fCalorimeter.Data());
630 return;
631 }
632
633 if(fCalorimeter == "PHOS") return; // Not implemented for PHOS yet
634
635 for(Int_t icluster = 0; icluster < pl->GetEntriesFast(); icluster++){
636 AliVCluster * cluster = (AliVCluster*) (pl->At(icluster));
637
638 // Study clusters with large shape parameter
639 Float_t en = cluster->E();
640 Float_t l0 = cluster->GetM02();
641 Int_t nc = cluster->GetNCells();
243c2909 642
992b14a7 643 //If too small or big E or low number of cells, skip it
e23a0471 644 if( en < GetMinEnergy() || en > GetMaxEnergy() || nc < fMinNCells) continue ;
992b14a7 645
646 Int_t absId1 = -1; Int_t absId2 = -1;
647 Int_t *absIdList = new Int_t [nc];
648 Float_t *maxEList = new Float_t[nc];
dbba06ca 649 Int_t nMax = GetCaloUtils()->GetNumberOfLocalMaxima(cluster, cells, absIdList, maxEList) ;
5c46c992 650 Bool_t matched = IsTrackMatched(cluster,GetReader()->GetInputEvent());
992b14a7 651
5c46c992 652 if (nMax <= 0)
653 {
29ca9cad 654 printf("AliAnaInsideClusterInvariantMass::MakeAnalysisFillHistograms() - No local maximum found!\n");
5c46c992 655
cab316a8 656 /*
29ca9cad 657 for(Int_t iDigit = 0; iDigit < cluster->GetNCells(); iDigit++ ) {
29ca9cad 658 Float_t ec = cells->GetCellAmplitude(cluster->GetCellsAbsId()[iDigit]);
dbba06ca 659 GetCaloUtils()->RecalibrateCellAmplitude(ec,fCalorimeter,cluster->GetCellsAbsId()[iDigit]);
29ca9cad 660 printf("iDigit %d, absId %d, Ecell %f\n",iDigit,cluster->GetCellsAbsId()[iDigit], ec);
29ca9cad 661 }
cab316a8 662 */
5c46c992 663
c5693f62 664 delete [] absIdList ;
665 delete [] maxEList ;
992b14a7 666 return;
667 }
668
5c46c992 669 fhNLocMax[0][matched]->Fill(en,nMax);
29ca9cad 670 for(Int_t imax = 0; imax < nMax; imax++)
671 {
5c46c992 672 fhNLocMaxEMax [0][matched]->Fill(maxEList[imax] ,nMax);
673 fhNLocMaxEFrac[0][matched]->Fill(maxEList[imax]/en,nMax);
29ca9cad 674 }
675
992b14a7 676
5c46c992 677 if ( nMax == 1 ) { fhM02NLocMax1[0][matched]->Fill(en,l0) ; fhNCellNLocMax1[0][matched]->Fill(en,nc) ; }
678 else if( nMax == 2 ) { fhM02NLocMax2[0][matched]->Fill(en,l0) ; fhNCellNLocMax2[0][matched]->Fill(en,nc) ; }
679 else if( nMax >= 3 ) { fhM02NLocMaxN[0][matched]->Fill(en,l0) ; fhNCellNLocMaxN[0][matched]->Fill(en,nc) ; }
992b14a7 680 else printf("N max smaller than 1 -> %d \n",nMax);
5c46c992 681
682 Float_t dZ = cluster->GetTrackDz();
683 Float_t dR = cluster->GetTrackDx();
684
685 if(cluster->IsEMCAL() && GetCaloUtils()->IsRecalculationOfClusterTrackMatchingOn())
686 {
687 dR = 2000., dZ = 2000.;
688 GetCaloUtils()->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dZ,dR);
689 }
690 //printf("Pi0EbE: dPhi %f, dEta %f\n",dR,dZ);
691
692 if(TMath::Abs(dR) < 999)
693 {
694 if ( nMax == 1 ) { fhTrackMatchedDEtaLocMax1[0]->Fill(en,dZ); fhTrackMatchedDPhiLocMax1[0]->Fill(en,dR); }
695 else if( nMax == 2 ) { fhTrackMatchedDEtaLocMax2[0]->Fill(en,dZ); fhTrackMatchedDPhiLocMax2[0]->Fill(en,dR); }
696 else if( nMax >= 3 ) { fhTrackMatchedDEtaLocMaxN[0]->Fill(en,dZ); fhTrackMatchedDPhiLocMaxN[0]->Fill(en,dR); }
697 }
992b14a7 698
699 // Play with the MC stack if available
700 // Check origin of the candidates
701 Int_t mcindex = -1;
702 if(IsDataMC()){
703
704 Int_t tag = GetMCAnalysisUtils()->CheckOrigin(cluster->GetLabels(),cluster->GetNLabels(), GetReader(), 0);
705
c5693f62 706 if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0) ) mcindex = kmcPi0;
707 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) ) mcindex = kmcEta;
992b14a7 708 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) &&
c5693f62 709 !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) mcindex = kmcPhoton;
992b14a7 710 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton) &&
c5693f62 711 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) mcindex = kmcConversion;
712 else if ( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCElectron)) mcindex = kmcElectron;
713 else mcindex = kmcHadron;
992b14a7 714
29ca9cad 715 //GetMCAnalysisUtils()->PrintMCTag(tag);
31ae6d59 716 //printf("\t MC index Assigned %d \n",mcindex);
71e3889f 717 fhNLocMaxNLabel[0] [matched]->Fill(en,nMax);
718 fhNLocMaxNLabel[mcindex][matched]->Fill(cluster->GetNLabels(),nMax);
719
5c46c992 720 fhNLocMax[mcindex][matched]->Fill(en,nMax);
29ca9cad 721 for(Int_t imax = 0; imax < nMax; imax++)
722 {
5c46c992 723 fhNLocMaxEMax [mcindex][matched]->Fill(maxEList[imax] ,nMax);
724 fhNLocMaxEFrac[mcindex][matched]->Fill(maxEList[imax]/en,nMax);
725 }
726 if (nMax == 1 ) { fhM02NLocMax1[mcindex][matched]->Fill(en,l0) ; fhNCellNLocMax1[mcindex][matched]->Fill(en,nc) ; }
727 else if(nMax == 2 ) { fhM02NLocMax2[mcindex][matched]->Fill(en,l0) ; fhNCellNLocMax2[mcindex][matched]->Fill(en,nc) ; }
728 else if(nMax >= 3 ) { fhM02NLocMaxN[mcindex][matched]->Fill(en,l0) ; fhNCellNLocMaxN[mcindex][matched]->Fill(en,nc) ; }
729
730 if(TMath::Abs(dR) < 999)
731 {
732 if ( nMax == 1 ) { fhTrackMatchedDEtaLocMax1[mcindex]->Fill(en,dZ); fhTrackMatchedDPhiLocMax1[mcindex]->Fill(en,dR); }
733 else if( nMax == 2 ) { fhTrackMatchedDEtaLocMax2[mcindex]->Fill(en,dZ); fhTrackMatchedDPhiLocMax2[mcindex]->Fill(en,dR); }
734 else if( nMax >= 3 ) { fhTrackMatchedDEtaLocMaxN[mcindex]->Fill(en,dZ); fhTrackMatchedDPhiLocMaxN[mcindex]->Fill(en,dR); }
29ca9cad 735 }
992b14a7 736
737 }
738
29ca9cad 739 //---------------------------------------------------------------------
740 // From here only if M02 is large, fill histograms or split the cluster
741 //---------------------------------------------------------------------
742
71e3889f 743 if( l0 < fM02MinCut || l0 > fM02MaxCut )
29ca9cad 744 {
37f53c3e 745 delete [] absIdList ;
746 delete [] maxEList ;
747 continue;
748 }
992b14a7 749
5c46c992 750 fhNLocMaxM02Cut[0][matched]->Fill(en,nMax);
751 if(IsDataMC()) fhNLocMaxM02Cut[mcindex][matched]->Fill(en,nMax);
243c2909 752
29ca9cad 753 //---------------------------------------------------------------------
992b14a7 754 // Get the 2 max indeces and do inv mass
29ca9cad 755 //---------------------------------------------------------------------
756
243c2909 757 if ( nMax == 2 ) {
758 absId1 = absIdList[0];
759 absId2 = absIdList[1];
760 }
761 else if( nMax == 1 )
762 {
763
764 absId1 = absIdList[0];
992b14a7 765
992b14a7 766 //Find second highest energy cell
243c2909 767
04c5e581 768 Float_t enmax = 0 ;
992b14a7 769 for(Int_t iDigit = 0 ; iDigit < cluster->GetNCells() ; iDigit++){
770 Int_t absId = cluster->GetCellsAbsId()[iDigit];
771 if( absId == absId1 ) continue ;
772 Float_t endig = cells->GetCellAmplitude(absId);
dbba06ca 773 GetCaloUtils()->RecalibrateCellAmplitude(endig,fCalorimeter,absId);
992b14a7 774 if(endig > enmax) {
775 enmax = endig ;
776 absId2 = absId ;
777 }
992b14a7 778 }// cell loop
243c2909 779 }// 1 maxima
780 else { // n max > 2
781 // loop on maxima, find 2 highest
992b14a7 782
243c2909 783 // First max
784 Float_t enmax = 0 ;
992b14a7 785 for(Int_t iDigit = 0 ; iDigit < nMax ; iDigit++){
786 Float_t endig = maxEList[iDigit];
787 if(endig > enmax) {
243c2909 788 enmax = endig ;
789 absId1 = absIdList[iDigit];
790 }
791 }// first maxima loop
792
793 // Second max
794 Float_t enmax2 = 0;
795 for(Int_t iDigit = 0 ; iDigit < nMax ; iDigit++){
796 if(absIdList[iDigit]==absId1) continue;
797 Float_t endig = maxEList[iDigit];
798 if(endig > enmax2) {
799 enmax2 = endig ;
992b14a7 800 absId2 = absIdList[iDigit];
801 }
243c2909 802 }// second maxima loop
803
804 } // n local maxima > 2
29ca9cad 805
806 //---------------------------------------------------------------------
807 // Split the cluster energy in 2, around the highest 2 local maxima
808 //---------------------------------------------------------------------
809
5c46c992 810 //Float_t en1 = 0, en2 = 0;
811 //SplitEnergy(absId1,absId2,cluster, cells, en1, en2, nMax /*absIdList, maxEList,*/);
812
813 AliAODCaloCluster *cluster1 = new AliAODCaloCluster(0, 0,NULL,0.,NULL,NULL,1,0);
814 AliAODCaloCluster *cluster2 = new AliAODCaloCluster(1, 0,NULL,0.,NULL,NULL,1,0);
815
816 SplitEnergy(absId1,absId2,cluster, cells, cluster1, cluster2, nMax /*absIdList, maxEList,*/);
817
29ca9cad 818 //---------------------------------------------------------------------
819 // Get mass of pair of clusters
820 //---------------------------------------------------------------------
821
822 // First set position of cluster as local maxima position,
823 // assign splitted energy to calculate momentum
992b14a7 824
5c46c992 825 //TLorentzVector cellMom1 = GetCellMomentum(absId1, en1, cells);
826 //TLorentzVector cellMom2 = GetCellMomentum(absId2, en2, cells);
29ca9cad 827
5c46c992 828 TLorentzVector cellMom1;
829 TLorentzVector cellMom2;
830
831 cluster1->GetMomentum(cellMom1,GetVertex(0));
832 cluster2->GetMomentum(cellMom2,GetVertex(0));
992b14a7 833
5c46c992 834 Float_t mass = (cellMom1+cellMom2).M();
835 Double_t angle = cellMom2.Angle(cellMom1.Vect());
836
243c2909 837 if (nMax==1)
838 {
5c46c992 839 fhAnglePairLocMax1[matched]->Fill(en,angle);
0137016b 840 fhMassNLocMax1[0][matched] ->Fill(en,mass );
841
842 if( en > 7 )
843 {
5c46c992 844 fhAnglePairMassLocMax1[matched]->Fill(mass,angle);
0137016b 845 fhMassM02NLocMax1[0] [matched]->Fill(l0, mass );
846 }
5c46c992 847
7b686344 848 if(en > 6 && en <= 10) fhMassM02NLocMax1Ebin[0]->Fill(l0, mass );
849 if(en > 10 && en <= 15) fhMassM02NLocMax1Ebin[1]->Fill(l0, mass );
850 if(en > 15 && en <= 20) fhMassM02NLocMax1Ebin[2]->Fill(l0, mass );
851 if(en > 20) fhMassM02NLocMax1Ebin[3]->Fill(l0, mass );
852
5c46c992 853 if (mass < fMassConMax && mass > fMassConMin) fhM02ConLocMax1[0][matched]->Fill(en,l0);
854 else if(mass < fMassPi0Max && mass > fMassPi0Min) fhM02Pi0LocMax1[0][matched]->Fill(en,l0);
855 else if(mass < fMassEtaMax && mass > fMassEtaMin) fhM02EtaLocMax1[0][matched]->Fill(en,l0);
243c2909 856 }
5c46c992 857 else if(nMax==2)
858 {
859 fhAnglePairLocMax2[matched]->Fill(en,angle);
0137016b 860 fhMassNLocMax2[0] [matched]->Fill(en,mass );
861
862 if( en > 7 )
863 {
5c46c992 864 fhAnglePairMassLocMax2[matched]->Fill(mass,angle);
0137016b 865 fhMassM02NLocMax2[0][matched] ->Fill(l0, mass );
866
867 }
5c46c992 868
7b686344 869 if(en > 6 && en <= 10) fhMassM02NLocMax2Ebin[0]->Fill(l0, mass );
870 if(en > 10 && en <= 15) fhMassM02NLocMax2Ebin[1]->Fill(l0, mass );
871 if(en > 15 && en <= 20) fhMassM02NLocMax2Ebin[2]->Fill(l0, mass );
872 if(en > 20) fhMassM02NLocMax2Ebin[3]->Fill(l0, mass );
873
5c46c992 874 if (mass < fMassConMax && mass > fMassConMin) fhM02ConLocMax2[0][matched]->Fill(en,l0);
875 else if(mass < fMassPi0Max && mass > fMassPi0Min) fhM02Pi0LocMax2[0][matched]->Fill(en,l0);
876 else if(mass < fMassEtaMax && mass > fMassEtaMin) fhM02EtaLocMax2[0][matched]->Fill(en,l0);
243c2909 877 }
5c46c992 878 else if(nMax >2)
879 {
880 fhAnglePairLocMaxN[matched]->Fill(en,angle);
0137016b 881 fhMassNLocMaxN[0] [matched]->Fill(en,mass );
882
883 if( en > 7 )
884 {
5c46c992 885 fhAnglePairMassLocMaxN[matched]->Fill(mass,angle);
0137016b 886 fhMassM02NLocMaxN[0] [matched]->Fill(l0 ,mass );
887 }
5c46c992 888
7b686344 889 if(en > 6 && en <= 10) fhMassM02NLocMaxNEbin[0]->Fill(l0, mass );
890 if(en > 10 && en <= 15) fhMassM02NLocMaxNEbin[1]->Fill(l0, mass );
891 if(en > 15 && en <= 20) fhMassM02NLocMaxNEbin[2]->Fill(l0, mass );
892 if(en > 20) fhMassM02NLocMaxNEbin[3]->Fill(l0, mass );
893
5c46c992 894 if (mass < fMassConMax && mass > fMassConMin) fhM02ConLocMaxN[0][matched]->Fill(en,l0);
895 else if(mass < fMassPi0Max && mass > fMassPi0Min) fhM02Pi0LocMaxN[0][matched]->Fill(en,l0);
896 else if(mass < fMassEtaMax && mass > fMassEtaMin) fhM02EtaLocMaxN[0][matched]->Fill(en,l0);
243c2909 897 }
898
992b14a7 899 if(IsDataMC()){
900
243c2909 901 if (nMax==1)
902 {
5c46c992 903 fhMassNLocMax1[mcindex][matched]->Fill(en,mass);
0137016b 904 if( en > 7 ) fhMassM02NLocMax1[mcindex][matched]->Fill(l0, mass );
5c46c992 905 if (mass < fMassConMax && mass > fMassConMin) fhM02ConLocMax1[mcindex][matched]->Fill(en,l0);
906 else if(mass < fMassPi0Max && mass > fMassPi0Min) fhM02Pi0LocMax1[mcindex][matched]->Fill(en,l0);
907 else if(mass < fMassEtaMax && mass > fMassEtaMin) fhM02EtaLocMax1[mcindex][matched]->Fill(en,l0);
243c2909 908 }
0137016b 909 else if(nMax==2)
910 {
5c46c992 911 fhMassNLocMax2[mcindex][matched]->Fill(en,mass);
0137016b 912 if( en > 7 ) fhMassM02NLocMax2[mcindex][matched]->Fill(l0, mass );
5c46c992 913 if (mass < fMassConMax && mass > fMassConMin) fhM02ConLocMax2[mcindex][matched]->Fill(en,l0);
914 else if(mass < fMassPi0Max && mass > fMassPi0Min) fhM02Pi0LocMax2[mcindex][matched]->Fill(en,l0);
915 else if(mass < fMassEtaMax && mass > fMassEtaMin) fhM02EtaLocMax2[mcindex][matched]->Fill(en,l0);
243c2909 916 }
0137016b 917 else if(nMax >2)
918 {
5c46c992 919 fhMassNLocMaxN[mcindex][matched]->Fill(en,mass);
0137016b 920 if( en > 7 ) fhMassM02NLocMaxN[mcindex][matched]->Fill(l0, mass );
5c46c992 921 if (mass < fMassConMax && mass > fMassConMin) fhM02ConLocMaxN[mcindex][matched]->Fill(en,l0);
922 else if(mass < fMassPi0Max && mass > fMassPi0Min) fhM02Pi0LocMaxN[mcindex][matched]->Fill(en,l0);
923 else if(mass < fMassEtaMax && mass > fMassEtaMin) fhM02EtaLocMaxN[mcindex][matched]->Fill(en,l0);
243c2909 924 }
992b14a7 925
926 }//Work with MC truth first
5c46c992 927
928 delete cluster1 ;
929 delete cluster2 ;
c5693f62 930 delete [] absIdList ;
931 delete [] maxEList ;
992b14a7 932
933 }//loop
934
935 if(GetDebug() > 1) printf("AliAnaInsideClusterInvariantMass::MakeAnalysisFillHistograms() - END \n");
e23a0471 936
992b14a7 937}
938
939//______________________________________________________________________
940void AliAnaInsideClusterInvariantMass::Print(const Option_t * opt) const
941{
942 //Print some relevant parameters set for the analysis
943 if(! opt)
944 return;
945
946 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
745913ae 947 AliAnaCaloTrackCorrBaseClass::Print("");
243c2909 948 printf("Calorimeter = %s\n", fCalorimeter.Data()) ;
dbba06ca 949 printf("Loc. Max. E > %2.2f\n", GetCaloUtils()->GetLocalMaximaCutE());
950 printf("Loc. Max. E Diff > %2.2f\n", GetCaloUtils()->GetLocalMaximaCutEDiff());
71e3889f 951 printf("%2.2f < lambda_0^2 <%2.2f \n",fM02MinCut,fM02MaxCut);
243c2909 952 printf("pi0 : %2.2f<m<%2.2f \n", fMassPi0Min,fMassPi0Max);
953 printf("eta : %2.2f<m<%2.2f \n", fMassEtaMin,fMassEtaMax);
954 printf("conv: %2.2f<m<%2.2f \n", fMassConMin,fMassConMax);
955
992b14a7 956 printf(" \n") ;
957
958}
959
dbba06ca 960
992b14a7 961
c5693f62 962//________________________________________________________________________________________
963void AliAnaInsideClusterInvariantMass::SplitEnergy(const Int_t absId1, const Int_t absId2,
243c2909 964 AliVCluster* cluster,
965 AliVCaloCells* cells,
5c46c992 966 //Float_t & e1, Float_t & e2,
967 AliAODCaloCluster* cluster1,
968 AliAODCaloCluster* cluster2,
29ca9cad 969 const Int_t nMax)
992b14a7 970{
5c46c992 971
29ca9cad 972 // Split energy of cluster between the 2 local maxima, sum energy on 3x3, and if the 2
973 // maxima are too close and have common cells, split the energy between the 2
974
71e3889f 975 TH2F* hClusterMap = 0 ;
976 TH2F* hClusterLocMax = 0 ;
977 TH2F* hCluster1 = 0 ;
978 TH2F* hCluster2 = 0 ;
979
980 if(fPlotCluster)
981 {
982 hClusterMap = new TH2F("hClusterMap","Cluster Map",48,0,48,24,0,24);
983 hClusterLocMax = new TH2F("hClusterLocMax","Cluster 2 highest local maxima",48,0,48,24,0,24);
984 hCluster1 = new TH2F("hCluster1","Cluster 1",48,0,48,24,0,24);
985 hCluster2 = new TH2F("hCluster2","Cluster 2",48,0,48,24,0,24);
986 hClusterMap ->SetXTitle("column");
987 hClusterMap ->SetYTitle("row");
988 hClusterLocMax ->SetXTitle("column");
989 hClusterLocMax ->SetYTitle("row");
990 hCluster1 ->SetXTitle("column");
991 hCluster1 ->SetYTitle("row");
992 hCluster2 ->SetXTitle("column");
993 hCluster2 ->SetYTitle("row");
994 }
29ca9cad 995
243c2909 996 const Int_t ncells = cluster->GetNCells();
29ca9cad 997 Int_t absIdList[ncells];
998
5c46c992 999 Float_t e1 = 0, e2 = 0 ;
71e3889f 1000 Int_t icol = -1, irow = -1, iRCU = -1, sm = -1;
29ca9cad 1001 Float_t eCluster = 0;
71e3889f 1002 Float_t minCol = 100, minRow = 100, maxCol = -1, maxRow = -1;
243c2909 1003 for(Int_t iDigit = 0; iDigit < ncells; iDigit++ ) {
1004 absIdList[iDigit] = cluster->GetCellsAbsId()[iDigit];
243c2909 1005
71e3889f 1006
29ca9cad 1007 Float_t ec = cells->GetCellAmplitude(absIdList[iDigit]);
dbba06ca 1008 GetCaloUtils()->RecalibrateCellAmplitude(ec,fCalorimeter, absIdList[iDigit]);
29ca9cad 1009 eCluster+=ec;
1010
71e3889f 1011 if(fPlotCluster)
1012 {
1013 //printf("iDigit %d, absId %d, Ecell %f\n",iDigit,absIdList[iDigit], cells->GetCellAmplitude(absIdList[iDigit]));
1014 sm = GetCaloUtils()->GetModuleNumberCellIndexes(absIdList[iDigit], fCalorimeter, icol, irow, iRCU) ;
1015 if(icol > maxCol) maxCol = icol;
1016 if(icol < minCol) minCol = icol;
1017 if(irow > maxRow) maxRow = irow;
1018 if(irow < minRow) minRow = irow;
1019 hClusterMap->Fill(icol,irow,ec);
1020 }
29ca9cad 1021
1022 }
1023
243c2909 1024 // Init counters and variables
5c46c992 1025 Int_t ncells1 = 1 ;
1026 UShort_t absIdList1[9] ;
1027 Double_t fracList1 [9] ;
1028 absIdList1[0] = absId1 ;
1029 fracList1 [0] = 1. ;
243c2909 1030
1031 Float_t ecell1 = cells->GetCellAmplitude(absId1);
dbba06ca 1032 GetCaloUtils()->RecalibrateCellAmplitude(ecell1, fCalorimeter, absId1);
243c2909 1033 e1 = ecell1;
1034
5c46c992 1035 Int_t ncells2 = 1 ;
1036 UShort_t absIdList2[9] ;
1037 Double_t fracList2 [9] ;
1038 absIdList2[0] = absId2 ;
1039 fracList2 [0] = 1. ;
1040
243c2909 1041 Float_t ecell2 = cells->GetCellAmplitude(absId2);
dbba06ca 1042 GetCaloUtils()->RecalibrateCellAmplitude(ecell2, fCalorimeter, absId2);
243c2909 1043 e2 = ecell2;
29ca9cad 1044
71e3889f 1045 if(fPlotCluster)
1046 {
1047 Int_t icol1 = -1, irow1 = -1, icol2 = -1, irow2 = -1;
1048 sm = GetCaloUtils()->GetModuleNumberCellIndexes(absId1, fCalorimeter, icol1, irow1, iRCU) ;
1049 hClusterLocMax->Fill(icol1,irow1,ecell1);
1050 sm = GetCaloUtils()->GetModuleNumberCellIndexes(absId2, fCalorimeter, icol2, irow2, iRCU) ;
1051 hClusterLocMax->Fill(icol2,irow2,ecell2);
1052 }
29ca9cad 1053
1054 // Very rough way to share the cluster energy
1055 Float_t eRemain = (eCluster-ecell1-ecell2)/2;
1056 Float_t shareFraction1 = ecell1/eCluster+eRemain/eCluster;
1057 Float_t shareFraction2 = ecell2/eCluster+eRemain/eCluster;
1058
1059 for(Int_t iDigit = 0; iDigit < ncells; iDigit++){
1060 Int_t absId = absIdList[iDigit];
1061
1062 if(absId==absId1 || absId==absId2 || absId < 0) continue;
1063
1064 Float_t ecell = cells->GetCellAmplitude(absId);
dbba06ca 1065 GetCaloUtils()->RecalibrateCellAmplitude(ecell, fCalorimeter, absId);
29ca9cad 1066
5c46c992 1067 if(GetCaloUtils()->AreNeighbours(fCalorimeter, absId1,absId ))
1068 {
1069 absIdList1[ncells1]= absId;
29ca9cad 1070
5c46c992 1071 if(GetCaloUtils()->AreNeighbours(fCalorimeter, absId2,absId ))
1072 {
1073 fracList1[ncells1] = shareFraction1;
1074 e1 += ecell*shareFraction1;
1075 }
1076 else
1077 {
1078 fracList1[ncells1] = 1.;
1079 e1 += ecell;
1080 }
1081
1082 ncells1++;
1083
29ca9cad 1084 } // neigbour to cell1
243c2909 1085
5c46c992 1086 if(GetCaloUtils()->AreNeighbours(fCalorimeter, absId2,absId ))
1087 {
1088 absIdList2[ncells2]= absId;
29ca9cad 1089
5c46c992 1090 if(GetCaloUtils()->AreNeighbours(fCalorimeter, absId1,absId ))
1091 {
1092 fracList2[ncells2] = shareFraction2;
29ca9cad 1093 e2 += ecell*shareFraction2;
5c46c992 1094 }
1095 else
1096 {
1097 fracList2[ncells2] = 1.;
29ca9cad 1098 e2 += ecell;
5c46c992 1099 }
1100
1101 ncells2++;
1102
29ca9cad 1103 } // neigbour to cell2
5c46c992 1104
29ca9cad 1105 }
1106
1107 if(GetDebug() > 1) printf("AliAnaInsideClusterInvariantMass::SplitEnergy() - n Local Max %d, Cluster energy = %f, Ecell1 = %f, Ecell2 = %f, Enew1 = %f, Enew2 = %f, Remain %f, \n ncells %d, ncells1 %d, ncells2 %d, f1 %f, f2 %f, sum f12 = %f \n",
1108 nMax, eCluster,ecell1,ecell2,e1,e2,eCluster-e1-e2,ncells,ncells1,ncells2,shareFraction1,shareFraction2,shareFraction1+shareFraction2);
5c46c992 1109
1110 cluster1->SetE(e1);
1111 cluster2->SetE(e2);
1112
1113 cluster1->SetNCells(ncells1);
1114 cluster2->SetNCells(ncells2);
1115
1116 cluster1->SetCellsAbsId(absIdList1);
1117 cluster2->SetCellsAbsId(absIdList2);
1118
1119 cluster1->SetCellsAmplitudeFraction(fracList1);
1120 cluster2->SetCellsAmplitudeFraction(fracList2);
1121
1122 GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterPosition(GetEMCALGeometry(), cells, cluster1);
1123 GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterPosition(GetEMCALGeometry(), cells, cluster2);
1124
29ca9cad 1125
71e3889f 1126 if(fPlotCluster)
243c2909 1127 {
71e3889f 1128 //printf("Cells of cluster1: ");
1129 for(Int_t iDigit = 0; iDigit < ncells1; iDigit++ )
1130 {
1131 //printf(" %d ",absIdList1[iDigit]);
1132
1133 sm = GetCaloUtils()->GetModuleNumberCellIndexes(absIdList1[iDigit], fCalorimeter, icol, irow, iRCU) ;
1134
1135 if( GetCaloUtils()->AreNeighbours(fCalorimeter, absId2,absIdList1[iDigit]) )
1136 hCluster1->Fill(icol,irow,cells->GetCellAmplitude(absIdList1[iDigit])*shareFraction1);
1137 else
1138 hCluster1->Fill(icol,irow,cells->GetCellAmplitude(absIdList1[iDigit]));
1139 }
1140
1141 //printf(" \n ");
1142 //printf("Cells of cluster2: ");
1143
1144 for(Int_t iDigit = 0; iDigit < ncells2; iDigit++ )
1145 {
1146 //printf(" %d ",absIdList2[iDigit]);
1147
1148 sm = GetCaloUtils()->GetModuleNumberCellIndexes(absIdList2[iDigit], fCalorimeter, icol, irow, iRCU) ;
1149 if( GetCaloUtils()->AreNeighbours(fCalorimeter, absId1,absIdList2[iDigit]) )
1150 hCluster2->Fill(icol,irow,cells->GetCellAmplitude(absIdList2[iDigit])*shareFraction2);
1151 else
1152 hCluster2->Fill(icol,irow,cells->GetCellAmplitude(absIdList2[iDigit]));
1153
1154 }
1155 //printf(" \n ");
243c2909 1156
71e3889f 1157 gStyle->SetPadRightMargin(0.1);
1158 gStyle->SetPadLeftMargin(0.1);
1159 gStyle->SetOptStat(0);
1160 gStyle->SetOptFit(000000);
243c2909 1161
71e3889f 1162 if(maxCol-minCol > maxRow-minRow)
1163 {
1164 maxRow+= maxCol-minCol;
1165 }
29ca9cad 1166 else
71e3889f 1167 {
1168 maxCol+= maxRow-minRow;
1169 }
29ca9cad 1170
71e3889f 1171 TCanvas * c= new TCanvas("canvas", "canvas", 4000, 4000) ;
1172 c->Divide(2,2);
1173 c->cd(1);
1174 gPad->SetGridy();
1175 gPad->SetGridx();
1176 hClusterMap ->SetAxisRange(minCol, maxCol,"X");
1177 hClusterMap ->SetAxisRange(minRow, maxRow,"Y");
1178 hClusterMap ->Draw("colz");
1179 c->cd(2);
1180 gPad->SetGridy();
1181 gPad->SetGridx();
1182 hClusterLocMax ->SetAxisRange(minCol, maxCol,"X");
1183 hClusterLocMax ->SetAxisRange(minRow, maxRow,"Y");
1184 hClusterLocMax ->Draw("colz");
1185 c->cd(3);
1186 gPad->SetGridy();
1187 gPad->SetGridx();
1188 hCluster1 ->SetAxisRange(minCol, maxCol,"X");
1189 hCluster1 ->SetAxisRange(minRow, maxRow,"Y");
1190 hCluster1 ->Draw("colz");
1191 c->cd(4);
1192 gPad->SetGridy();
1193 gPad->SetGridx();
1194 hCluster2 ->SetAxisRange(minCol, maxCol,"X");
1195 hCluster2 ->SetAxisRange(minRow, maxRow,"Y");
1196 hCluster2 ->Draw("colz");
1197
1198 if(eCluster > 6 )c->Print(Form("clusterFigures/Event%d_E%1.0f_nMax%d_NCell1_%d_NCell2_%d.eps",GetEventNumber(),cluster->E(),nMax,ncells1,ncells2));
1199
1200 delete c;
1201 delete hClusterMap;
1202 delete hClusterLocMax;
1203 delete hCluster1;
1204 delete hCluster2;
29ca9cad 1205 }
992b14a7 1206}
1207