ATLAS Offline Software
JetCaloVariables.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 // #include "GaudiKernel/MsgStream.h"
6 // #include "GaudiKernel/Algorithm.h"
7 // #include "GaudiKernel/ListItem.h"
8 // #include "GaudiKernel/IToolSvc.h"
9 // #include "GaudiKernel/PropertyMgr.h"
10 // #include "GaudiKernel/IMessageSvc.h"
11 // #include "GaudiKernel/StatusCode.h"
12 // #include "StoreGate/StoreGate.h"
13 
14 // #include "PathResolver/PathResolver.h"
15 
16 // #include "EventKernel/SignalStateHelper.h"
17 // #include "JetEvent/Jet.h"
18 
19 // //#include "JetUtils/JetCaloHelper.h"
20 
21 // #include "JetUtils/JetCaloVariables.h"
22 
23 // #include "TAxis.h"
24 // #include "TFile.h"
25 // #include "TH1.h"
26 // #include "TH2.h"
27 // #include <cmath>
28 // #include <functional>
29 
30 
31 // /////////////////////////////////////////////////////////////////////////////////////
32 // JetCaloVariables::JetCaloVariables(const std::string& m_JetContainerName){
33 // cosmicPdf_EMFraction.clear();
34 // dijetPdf_EMFraction.clear();
35 // cosmicPdf_RconeANDRatio.clear();
36 // dijetPdf_RconeANDRatio.clear();
37 // ptbin_pdf.clear();
38 // filePDF=NULL;
39 // JetContainerName=m_JetContainerName;
40 // }
41 
42 // JetCaloVariables::~JetCaloVariables(){
43 // for(unsigned int i=0; i<cosmicPdf_EMFraction.size(); i++){
44 // if(cosmicPdf_EMFraction[i]!=NULL) delete cosmicPdf_EMFraction[i];
45 // if(dijetPdf_EMFraction[i]!=NULL) delete dijetPdf_EMFraction[i];
46 // if(cosmicPdf_RconeANDRatio[i]!=NULL) delete cosmicPdf_RconeANDRatio[i];
47 // if(dijetPdf_RconeANDRatio[i]!=NULL) delete dijetPdf_RconeANDRatio[i];
48 // }
49 // if(filePDF!=NULL) delete filePDF;
50 // }
51 
52 
53 // StatusCode JetCaloVariables::getpdf(){
54 // Athena::IMessageSvcHolder ims;
55 // MsgStream mLog(ims, "JetCaloVariables");
56 // //////////////////////////////////////////////////////////////////////////////
57 // ptbin_pdf.clear();
58 // ptbin_pdf.push_back(40);
59 // ptbin_pdf.push_back(80);
60 // ptbin_pdf.push_back(160);
61 // //////////////////////////////////////////////////////////////////////////////
62 // cosmicPdf_EMFraction.clear();
63 // dijetPdf_EMFraction.clear();
64 // cosmicPdf_RconeANDRatio.clear();
65 // dijetPdf_RconeANDRatio.clear();
66 
67 // // Find the path to inputFile
68 // std::string fname = PathResolver::find_file("cos_LLR_pdf.root", "DATAPATH");
69 // if(fname==""){
70 // mLog << MSG::ERROR << "Could not find input file " << fname << endmsg;
71 // return StatusCode::FAILURE;
72 // }
73 
74 // filePDF = new TFile(fname.c_str());
75 // mLog << MSG::INFO << "The input file is " << fname.c_str() << endmsg;
76 // if(filePDF->IsOpen()){
77 // if(JetContainerName.compare("Cone4H1TopoJets")==0){
78 // cosmicPdf_EMFraction.push_back((TH1D*)filePDF->Get("cone4topo_JetEMf_bin0_s_pdf_more"));
79 // cosmicPdf_EMFraction.push_back((TH1D*)filePDF->Get("cone4topo_JetEMf_bin1_s_pdf_more"));
80 // cosmicPdf_EMFraction.push_back((TH1D*)filePDF->Get("cone4topo_JetEMf_bin2_s_pdf_more"));
81 // cosmicPdf_EMFraction.push_back((TH1D*)filePDF->Get("cone4topo_JetEMf_bin3_s_pdf_more"));
82 
83 // cosmicPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("cone4topo_Rcone_Ratio_bin0_s_pdf_more"));
84 // cosmicPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("cone4topo_Rcone_Ratio_bin1_s_pdf_more"));
85 // cosmicPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("cone4topo_Rcone_Ratio_bin2_s_pdf_more"));
86 // cosmicPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("cone4topo_Rcone_Ratio_bin3_s_pdf_more"));
87 
88 // dijetPdf_EMFraction.push_back((TH1D*)filePDF->Get("cone4topo_JetEMf_bin0_b_pdf_more"));
89 // dijetPdf_EMFraction.push_back((TH1D*)filePDF->Get("cone4topo_JetEMf_bin1_b_pdf_more"));
90 // dijetPdf_EMFraction.push_back((TH1D*)filePDF->Get("cone4topo_JetEMf_bin2_b_pdf_more"));
91 // dijetPdf_EMFraction.push_back((TH1D*)filePDF->Get("cone4topo_JetEMf_bin3_b_pdf_more"));
92 
93 // dijetPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("cone4topo_Rcone_Ratio_bin0_b_pdf_more"));
94 // dijetPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("cone4topo_Rcone_Ratio_bin1_b_pdf_more"));
95 // dijetPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("cone4topo_Rcone_Ratio_bin2_b_pdf_more"));
96 // dijetPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("cone4topo_Rcone_Ratio_bin3_b_pdf_more"));
97 // }else{
98 // // we temporary use the pdf obtained with AntiKt4H1TopoJets for all other kinds of jets.//
99 // // These pdf should be updated later. //
100 // cosmicPdf_EMFraction.push_back((TH1D*)filePDF->Get("antikt4topo_JetEMf_bin0_s_pdf_more"));
101 // cosmicPdf_EMFraction.push_back((TH1D*)filePDF->Get("antikt4topo_JetEMf_bin1_s_pdf_more"));
102 // cosmicPdf_EMFraction.push_back((TH1D*)filePDF->Get("antikt4topo_JetEMf_bin2_s_pdf_more"));
103 // cosmicPdf_EMFraction.push_back((TH1D*)filePDF->Get("antikt4topo_JetEMf_bin3_s_pdf_more"));
104 
105 // cosmicPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("antikt4topo_Rcone_Ratio_bin0_s_pdf_more"));
106 // cosmicPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("antikt4topo_Rcone_Ratio_bin1_s_pdf_more"));
107 // cosmicPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("antikt4topo_Rcone_Ratio_bin2_s_pdf_more"));
108 // cosmicPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("antikt4topo_Rcone_Ratio_bin3_s_pdf_more"));
109 
110 // dijetPdf_EMFraction.push_back((TH1D*)filePDF->Get("antikt4topo_JetEMf_bin0_b_pdf_more"));
111 // dijetPdf_EMFraction.push_back((TH1D*)filePDF->Get("antikt4topo_JetEMf_bin1_b_pdf_more"));
112 // dijetPdf_EMFraction.push_back((TH1D*)filePDF->Get("antikt4topo_JetEMf_bin2_b_pdf_more"));
113 // dijetPdf_EMFraction.push_back((TH1D*)filePDF->Get("antikt4topo_JetEMf_bin3_b_pdf_more"));
114 
115 // dijetPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("antikt4topo_Rcone_Ratio_bin0_b_pdf_more"));
116 // dijetPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("antikt4topo_Rcone_Ratio_bin1_b_pdf_more"));
117 // dijetPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("antikt4topo_Rcone_Ratio_bin2_b_pdf_more"));
118 // dijetPdf_RconeANDRatio.push_back((TH2D*)filePDF->Get("antikt4topo_Rcone_Ratio_bin3_b_pdf_more"));
119 // }
120 // }else{
121 // mLog << MSG::ERROR << "Could not open file " << fname << endmsg;
122 // return StatusCode::FAILURE;
123 // }
124 // /////////////////////////////////////////////////////////////////////////
125 // //for(int i=0; i<int(cosmicPdf_EMFraction.size()); i++)
126 // //mLog << MSG::INFO <<cosmicPdf_EMFraction[i]->GetName()<<":"<<cosmicPdf_EMFraction[i]->GetSumOfWeights()<<endmsg;
127 
128 // //for(int i=0; i<int(cosmicPdf_RconeANDRatio.size()); i++)
129 // //mLog << MSG::INFO <<cosmicPdf_RconeANDRatio[i]->GetName()<<":"<<cosmicPdf_RconeANDRatio[i]->GetSumOfWeights()<<endmsg;
130 
131 // //for(int i=0; i<int(dijetPdf_EMFraction.size()); i++)
132 // //mLog << MSG::INFO <<dijetPdf_EMFraction[i]->GetName()<<":"<<dijetPdf_EMFraction[i]->GetSumOfWeights()<<endmsg;
133 
134 // //for(int i=0; i<int(dijetPdf_RconeANDRatio.size()); i++)
135 // //mLog << MSG::INFO <<dijetPdf_RconeANDRatio[i]->GetName()<<":"<<dijetPdf_RconeANDRatio[i]->GetSumOfWeights()<<endmsg;
136 
137 // return StatusCode::SUCCESS;
138 // }
139 
140 
141 // double JetCaloVariables::compute_RCone(const Jet* jet){
142 // // double sumE_cells=0;
143 // // double rcone=0;
144 
145 
146 // // SignalStateHelper emscale_h(jet, P4SignalState::JETEMSCALE);//Set the jet to EMSCALE
147 
148 // // double em_eta=jet->eta();
149 // // double em_phi=jet->phi();
150 // // emscale_h.resetSignalState(); //Set it back to the original state
151 
152 // // NavigationToken<CaloCell,double> cellToken;
153 // // jet->fillToken(cellToken,double(1.));
154 // // if ( cellToken.size() > 0 ){
155 // // NavigationToken<CaloCell,double>::const_iterator fCell = cellToken.begin();
156 // // NavigationToken<CaloCell,double>::const_iterator lCell = cellToken.end();
157 // // for( ; fCell != lCell; fCell++){
158 // // if( (*fCell)->e()>0 ){
159 // // sumE_cells+=(*fCell)->e();
160 // // double deta = ((*fCell)->eta()-em_eta);
161 // // double dphi = fabs((*fCell)->phi()-em_phi);
162 // // if(dphi>3.1415)dphi=2*3.1415-dphi;
163 // // double dr = sqrt(deta*deta+dphi*dphi);
164 // // rcone+=dr*fabs((*fCell)->e());
165 // // }
166 // // }
167 // // if(sumE_cells!=0) rcone=rcone/sumE_cells;
168 // // }
169 // // return rcone;
170 // return 0;
171 // }
172 
173 
174 // double JetCaloVariables::compute_RatioLeadingCells(const Jet* jet){
175 // //------------------------------------------
176 // // for the definition of Ratio_LeadingCells
177 // //------------------------------------------
178 // // double ratio_leading_cells=0;
179 // // double sumE_cells=0;
180 // // std::vector<double> cell_energies_em; cell_energies_em.clear();
181 // // std::vector<double> cell_energies_had; cell_energies_had.clear();
182 // // NavigationToken<CaloCell,double> cellToken;
183 // // jet->fillToken(cellToken,double(1.));
184 // // if ( cellToken.size() > 0 ){
185 // // NavigationToken<CaloCell,double>::const_iterator fCell = cellToken.begin();
186 // // NavigationToken<CaloCell,double>::const_iterator lCell = cellToken.end();
187 // // for( ; fCell != lCell; fCell++){
188 // // if( (*fCell)->e()>0 ){
189 // // sumE_cells+=(*fCell)->e();
190 // // if( (*fCell)->caloDDE()->getSubCalo()==0 ){ //0=EM
191 // // cell_energies_em.push_back((*fCell)->e());
192 // // }else{
193 // // cell_energies_had.push_back((*fCell)->e());
194 // // }
195 // // }
196 // // }
197 // // std::sort(cell_energies_had.rbegin(),cell_energies_had.rend());
198 // // for(unsigned int i=0; i<cell_energies_had.size()&&i<2; i++){
199 // // ratio_leading_cells+=cell_energies_had[i];
200 // // }
201 // // std::sort(cell_energies_em.rbegin(),cell_energies_em.rend());
202 // // for(unsigned int i=0; i<cell_energies_em.size()&&i<32; i++){
203 // // ratio_leading_cells+=cell_energies_em[i];
204 // // }
205 // // if(sumE_cells!=0) ratio_leading_cells=ratio_leading_cells/sumE_cells;
206 // // }
207 // // return ratio_leading_cells;
208 // return 0;
209 // }
210 
211 
212 // double JetCaloVariables::compute_LLREmFraction(const Jet* jet){
213 // // Athena::IMessageSvcHolder ims;
214 // // MsgStream mLog(ims, "JetCaloVariables");
215 
216 // // //--- to get pdf ---
217 // // if(filePDF==NULL){
218 // // StatusCode sc = getpdf();
219 // // if(sc.isFailure()) mLog << MSG::WARNING << "Some problems with getpdf()! " << endmsg;
220 // // }
221 
222 // // //--- to get emfraction and em_pt ---
223 // // double emfraction=JetCaloHelper::jetEMFraction(jet);
224 
225 // // SignalStateHelper emscale_h(jet, P4SignalState::JETEMSCALE);//Set the jet to EMSCALE
226 // // double em_pt=jet->pt();
227 // // double em_eta=jet->eta();
228 // // emscale_h.resetSignalState(); //Set it back to the original state
229 
230 // // //mLog << MSG::DEBUG << "compute_LLREmFraction: em_pt = " <<em_pt<< endmsg;
231 
232 // // if(em_pt/1000.<20. || fabs(em_eta)>2.5) return -20.;
233 
234 // // int iptbin=0;
235 // // //em_pt: MeV, ptbin_pdf[iptbin] : GeV
236 // // while(ptbin_pdf[iptbin]<em_pt/1000.&&iptbin<int(ptbin_pdf.size())){
237 // // iptbin++;
238 // // }
239 // // //mLog << MSG::INFO << "compute_LLREmFraction iptbin = " <<iptbin<< endmsg;
240 // // TH1* hist_cos=cosmicPdf_EMFraction[iptbin];
241 // // TH1* hist_dij=dijetPdf_EMFraction[iptbin];
242 
243 // // double xmax=hist_cos->GetXaxis()->GetXmax();
244 // // double xmin=hist_cos->GetXaxis()->GetXmin();
245 // // double emfra=emfraction;
246 // // if(emfra>xmax){
247 // // emfra=xmax-0.001;
248 // // }else if(emfra<xmin){
249 // // emfra=xmin+0.001;
250 // // }
251 // // double sig = hist_cos->Interpolate(emfra);
252 // // double bkg = hist_dij->Interpolate(emfra);
253 
254 // // if(sig!=0&&bkg!=0){
255 // // return log(sig/bkg);
256 // // }else if(sig==0){
257 // // return -20;
258 // // }else{
259 // // return 20;
260 // // }
261 // return 0;
262 // }
263 
264 
265 // double JetCaloVariables::compute_LLRRconeANDRatio(const Jet* jet){
266 // // Athena::IMessageSvcHolder ims;
267 // // MsgStream mLog(ims, "JetCaloVariables");
268 
269 // // //--- to get pdf ---
270 // // if(filePDF==NULL){
271 // // StatusCode sc = getpdf();
272 // // if(sc.isFailure()) mLog << MSG::WARNING << "Some problems with getpdf()! " << endmsg;
273 // // }
274 
275 // // //--- to get rcone, ratio_leadingcells and em_pt ---
276 // // double rcone=compute_RCone(jet);
277 // // double ratio=compute_RatioLeadingCells(jet);
278 
279 // // SignalStateHelper emscale_h(jet, P4SignalState::JETEMSCALE);//Set the jet to EMSCALE
280 // // double em_pt=jet->pt();
281 // // double em_eta=jet->eta();
282 // // emscale_h.resetSignalState(); //Set it back to the original state
283 // // //mLog << MSG::DEBUG << "compute_LLRRconeANDRatio: em_pt = " <<em_pt<< endmsg;
284 
285 // // if(em_pt/1000.<20. || fabs(em_eta)>2.5) return -20.;
286 
287 // // int iptbin=0;
288 // // //pt: MeV, ptbin_pdf[iptbin] : GeV
289 // // while(ptbin_pdf[iptbin]<em_pt/1000.&&iptbin<int(ptbin_pdf.size())){
290 // // iptbin++;
291 // // }
292 
293 // // //mLog << MSG::INFO << "compute_LLRLeading iptbin = " <<iptbin<< endmsg;
294 // // TH2* hist_cos=cosmicPdf_RconeANDRatio[iptbin];
295 // // TH2* hist_dij=dijetPdf_RconeANDRatio[iptbin];
296 // // double xmax=hist_cos->GetXaxis()->GetXmax();
297 // // double ymax=hist_cos->GetYaxis()->GetXmax();
298 // // double xmin=hist_cos->GetXaxis()->GetXmin();
299 // // double ymin=hist_cos->GetYaxis()->GetXmin();
300 // // double rconetmp=rcone;
301 // // double ratiotmp=ratio;
302 // // if(rconetmp>xmax){
303 // // rconetmp=xmax-0.001;
304 // // }else if(rconetmp<xmin){
305 // // rconetmp=xmin+0.001;
306 // // }
307 // // if(ratiotmp>ymax){
308 // // ratiotmp=ymax-0.001;
309 // // }else if(ratiotmp<ymin){
310 // // ratiotmp=ymin+0.001;
311 // // }
312 
313 // // double sig=hist_cos->Interpolate(rconetmp,ratiotmp);
314 // // double bkg=hist_dij->Interpolate(rconetmp,ratiotmp);
315 // // //mLog << MSG::INFO <<"LLRLeading = "<<log(sig/bkg)<< endmsg;
316 
317 // // if(sig!=0&&bkg!=0){
318 // // return log(sig/bkg);
319 // // }else if(sig==0){
320 // // return -20;
321 // // }else{
322 // // return 100;
323 // // }
324 // return 0;
325 // }
326 
327 
328 // double JetCaloVariables::compute_jetTime(const Jet* jet){
329 // // double sumE_cells=0;
330 // // double jettime=0;
331 // // NavigationToken<CaloCell,double> cellToken;
332 // // jet->fillToken(cellToken,double(1.));
333 // // if ( cellToken.size() > 0 ){
334 // // NavigationToken<CaloCell,double>::const_iterator fCell = cellToken.begin();
335 // // NavigationToken<CaloCell,double>::const_iterator lCell = cellToken.end();
336 // // for( ; fCell != lCell; fCell++){
337 // // sumE_cells+=(*fCell)->e();
338 // // double tmptime=(*fCell)->time();
339 // // jettime+=tmptime*((*fCell)->e());
340 // // }
341 // // if(sumE_cells!=0) jettime=jettime/sumE_cells;
342 // // }else{
343 // // jettime=-900.;
344 // // }
345 // // return jettime;
346 // return 0.;
347 // }
Generated on Thu Jul 4 2024 21:18:36 for ATLAS Offline Software by   doxygen 1.8.18