eb90e2bd8c057de414288d1785d7c45797c39b29
[u/mrichter/AliRoot.git] / HLT / global / physics / AliHLTCaloHistoInvMass.cxx
1 //-*- Mode: C++ -*-
2 /**************************************************************************
3  * This file is property of and copyright by the ALICE HLT Project        * 
4  * All rights reserved.                                                   *
5  *                                                                        *
6  * Primary Authors: Svein Lindal                                          *
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   AliHLTCaloHistoInvMass
19  * @author Svein Lindal <slindal@fys.uio.no>
20  * @date 
21  * @brief  Produces plots of invariant mass of two clusters. 
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 "AliHLTCaloHistoInvMass.h"
31 #include "AliHLTCaloClusterDataStruct.h"
32 #include "TObjArray.h"
33 #include "AliESDEvent.h"
34 #include "TRefArray.h"
35 #include "TH1F.h"
36 #include "TString.h"
37 #include "AliESDCaloCluster.h"
38 #include "TVector3.h"
39 #include "TLorentzVector.h"
40
41 AliHLTCaloHistoInvMass::AliHLTCaloHistoInvMass(TString det) :
42   fHistTwoClusterInvMass0(NULL),
43   fHistTwoClusterInvMass1(NULL),
44   fHistTwoClusterInvMass2(NULL),
45   fHistTwoClusterInvMass3(NULL),
46   fHistTwoClusterInvMass4(NULL)
47 {
48   // See header file for documentation
49   fHistTwoClusterInvMass0 = new TH1F(Form("%s fHistTwoClusterInvMass0", det.Data()), Form("Invariant mass of two clusters in %s, 1.5 GeV < p_{T} < 2.0", det.Data()), 200, 0, 1);
50   fHistTwoClusterInvMass0->GetXaxis()->SetTitle("m_{#gamma#gamma} GeV");
51   fHistTwoClusterInvMass0->GetYaxis()->SetTitle("Counts");
52   fHistTwoClusterInvMass0->SetMarkerStyle(21);
53   fHistArray->AddLast(fHistTwoClusterInvMass0);
54
55   fHistTwoClusterInvMass1 = new TH1F(Form("%s fHistTwoClusterInvMass1", det.Data()), Form("Invariant mass of two clusters in %s, 2.0 GeV < p_{T} < 2.5", det.Data()), 200, 0, 1);
56   fHistTwoClusterInvMass1->GetXaxis()->SetTitle("m_{#gamma#gamma} GeV");
57   fHistTwoClusterInvMass1->GetYaxis()->SetTitle("Counts");
58   fHistTwoClusterInvMass1->SetMarkerStyle(21);
59   fHistArray->AddLast(fHistTwoClusterInvMass1);
60
61   fHistTwoClusterInvMass2 = new TH1F(Form("%s fHistTwoClusterInvMass2", det.Data()), Form("Invariant mass of two clusters in %s, 2.5 GeV < p_{T} < 3.0", det.Data()), 200, 0, 1);
62   fHistTwoClusterInvMass2->GetXaxis()->SetTitle("m_{#gamma#gamma} GeV");
63   fHistTwoClusterInvMass2->GetYaxis()->SetTitle("Counts");
64   fHistTwoClusterInvMass2->SetMarkerStyle(21);
65   fHistArray->AddLast(fHistTwoClusterInvMass2);
66
67   fHistTwoClusterInvMass3 = new TH1F(Form("%s fHistTwoClusterInvMass3", det.Data()), Form("Invariant mass of two clusters in %s, 3.0 GeV < p_{T} < 5.0", det.Data()), 200, 0, 1);
68   fHistTwoClusterInvMass3->GetXaxis()->SetTitle("m_{#gamma#gamma} GeV");
69   fHistTwoClusterInvMass3->GetYaxis()->SetTitle("Counts");
70   fHistTwoClusterInvMass3->SetMarkerStyle(21);
71   fHistArray->AddLast(fHistTwoClusterInvMass3);
72
73   fHistTwoClusterInvMass4 = new TH1F(Form("%s fHistTwoClusterInvMass4", det.Data()), Form("Invariant mass of two clusters in %s p_{T} > 5.0 GeV", det.Data()), 200, 0, 1);
74   fHistTwoClusterInvMass4->GetXaxis()->SetTitle("m_{#gamma#gamma} GeV");
75   fHistTwoClusterInvMass4->GetYaxis()->SetTitle("Counts");
76   fHistTwoClusterInvMass4->SetMarkerStyle(21);
77   fHistArray->AddLast(fHistTwoClusterInvMass4);
78
79
80 }
81
82 AliHLTCaloHistoInvMass::~AliHLTCaloHistoInvMass()
83 {
84
85   if(fHistTwoClusterInvMass0)
86     delete fHistTwoClusterInvMass0;
87   fHistTwoClusterInvMass0 = NULL;
88
89   if(fHistTwoClusterInvMass1)
90     delete fHistTwoClusterInvMass1;
91   fHistTwoClusterInvMass1 = NULL;
92
93   if(fHistTwoClusterInvMass2)
94     delete fHistTwoClusterInvMass2;
95   fHistTwoClusterInvMass2 = NULL;
96
97   if(fHistTwoClusterInvMass3)
98     delete fHistTwoClusterInvMass3;
99   fHistTwoClusterInvMass3 = NULL;
100
101   if(fHistTwoClusterInvMass4)
102     delete fHistTwoClusterInvMass4;
103   fHistTwoClusterInvMass4 = NULL;
104
105
106 }
107
108
109 Int_t AliHLTCaloHistoInvMass::FillHistograms(Int_t nc, vector<AliHLTCaloClusterDataStruct*> &cVec) {
110   //See header file for documentation
111   
112   Float_t cPos[nc][3];
113   Float_t cEnergy[nc];
114
115   for(int ic = 0; ic < nc; ic++) {
116     AliHLTCaloClusterDataStruct * cluster = cVec.at(ic);
117     cluster->GetPosition(cPos[ic]);
118     cEnergy[ic] = cluster->E();
119   }
120
121   return FillInvariantMassHistograms(nc, cPos, cEnergy);
122 }
123
124
125
126 Int_t AliHLTCaloHistoInvMass::FillHistograms(const Int_t nc, TRefArray * clusterArray) {
127   //See header file for documentation
128   
129   Float_t cPos[nc][3];
130   Float_t cEnergy[nc];
131
132   for(int ic = 0; ic < nc; ic++) {
133     AliESDCaloCluster * cluster = static_cast<AliESDCaloCluster*>(clusterArray->At(ic));
134     cluster->GetPosition(cPos[ic]);
135     cEnergy[ic] = cluster->E();
136   }
137
138   return FillInvariantMassHistograms(nc, cPos, cEnergy);
139
140 }
141
142
143 Int_t AliHLTCaloHistoInvMass::FillInvariantMassHistograms(const Int_t nc, Float_t cPos[][3], Float_t  cEnergy[]){ 
144
145   Int_t iResult = 0;
146
147   for(Int_t ic = 0; ic<(nc-1); ic++) { 
148      
149     //BALLE hardcoded variable
150     if(cEnergy[ic] < 0.5)
151       continue;
152
153     //Get momentum vector
154     TVector3 iVec(cPos[ic]);
155     iVec = iVec.Unit();
156     iVec = cEnergy[ic] * iVec;
157     
158     
159     for(Int_t jc = ic+1; jc<nc; jc++) { 
160      
161     //BALLE hardcoded variable
162       if(cEnergy[jc] < 0.5)
163         continue;
164
165       //Get second momentum vector
166       TVector3 jVec(cPos[jc]);
167       jVec = jVec.Unit();
168       jVec = cEnergy[jc] * jVec;
169       
170       //Calculate inv mass
171       Double_t m2 = 2 *(cEnergy[ic]* cEnergy[jc] - iVec.Dot(jVec));
172       Double_t m = 0.0;
173       if(m2 > 0) { 
174         m = TMath::Sqrt(m2);
175       }
176       
177       //Fill histograms
178       
179       Float_t sum = cEnergy[ic]+cEnergy[jc];
180       if(sum > 1.5)
181       {
182          if(sum > 2.0)
183          {
184             if(sum > 2.5)
185             {
186                if(sum > 3.0)
187                {
188                   if(sum > 5.0)
189                   {
190                      fHistTwoClusterInvMass4->Fill(m);
191                   }
192                   else
193                   {
194                      fHistTwoClusterInvMass3->Fill(m);
195                   }     
196                }
197                else
198                {
199                   fHistTwoClusterInvMass2->Fill(m);
200                }
201             }
202             else
203             {
204                fHistTwoClusterInvMass1->Fill(m);
205             }
206          }
207          else
208          {
209             fHistTwoClusterInvMass0->Fill(m);
210          }
211       }
212     }
213   }
214
215   return iResult;
216
217 }
218
219 // template <class T>
220 // Int_t AliHLTCaloHistoInvMass::FillHistograms(Int_t nc, vector<T*> clusterVec ) {
221 //   Float_t cPos[nc][3];
222 //   Float_t cEnergy[nc];
223   
224 //   for(int ic = 0; ic < nc; ic++) {
225 //     T * cluster = cVec.at(ic);
226 //     cluster->GetPosition(cPos[ic]);
227 //     cEnergy[ic] = cluster->E();
228 //   }
229   
230
231 //   for(Int_t ic = 0; ic<(nc-1); ic++) { 
232
233 //     //Get momentum vector
234 //     TVector3 iVec(cPos[ic]);
235 //     iVec = iVec.Unit();
236 //     iVec = cEnergy[ic] * iVec;
237
238     
239 //     for(Int_t jc = ic+1; jc<nc; jc++) { 
240
241 //       //Get second momentum vector
242 //       TVector3 jVec(cPos[jc]);
243 //       jVec = jVec.Unit();
244 //       jVec = cEnergy[jc] * jVec;
245
246 //       //Calculate inv mass
247 //       Double_t m = TMath::Sqrt( 2 *(cEnergy[ic]* cEnergy[jc] - iVec.Dot(jVec) ) );
248
249 //       //Fill histogram
250 //       fHistTwoClusterInvMass->Fill(m);
251 //     }
252 //   }
253 // }