ATLAS Offline Software
dependence.cxx
Go to the documentation of this file.
1 
10 #include <string>
11 #include <iostream>
12 #include <cstdlib>
13 #include <cmath>
14 
16 #include "ReadCards.h"
17 
18 #include "TH1F.h"
19 #include "TH1D.h"
20 #include "TF1.h"
21 #include "TFile.h"
22 #include "TPad.h"
23 #include "TCanvas.h"
24 #include "TStyle.h"
25 
26 #include "AtlasStyle.h"
27 #include "DrawLabel.h"
28 
29 
30 
31 int usage( const std::string& err_msg="", int status=0 ) {
32 
33  std::string name = "dependence";
34 
35  if ( err_msg!="" ) std::cerr << err_msg << std::endl;
36 
37  std::ostream& s = std::cout;
38 
39  s << "Usage: " << name << "\t [OPTIONS] \n\n";
40  s << "\t" << " plots comparison histograms";
41  s << " - compiled on " << __DATE__ << " at " << __TIME__ << "\n\n";
42  s << "Options: \n";
43  s << " -c, --config value \t configure which histograms to plot from config file,\n\n";
44  s << " -f, --fit \t fit a pol2 to the efficiencies\n\n";
45  s << " -lx, --logx \t drawwith log x axis\n";
46  s << " -as, --atlasstyle \t usethe atlas style\n\n";
47 
48  s << " -h, --help \t this help\n";
49 
50  s << std::endl;
51 
52  return status;
53 }
54 
55 
56 
57 bool doFit = false;
58 bool doLogx = false;
59 
60 double ptmax = 0;
61 
62 void add_to_bin( double d, TH1* h, TH1* h0 ) {
63 
64  int ibin = h->FindBin(d);
65  if ( ibin<1 || ibin>h->GetNbinsX() ) return;
66 
67  if ( ptmax ) {
68  double v = 0;
69  for ( int i=1 ; i<=h0->GetNbinsX() ; i++ ) if ( h0->GetBinCenter(i)<ptmax ) v += h0->GetBinContent(i);
70  h->SetBinContent( ibin, v );
71  }
72  else {
73  double v = 0;
74  for ( int i=1 ; i<=h0->GetNbinsX() ; i++ ) v += h0->GetBinContent(i);
75  h->SetBinContent( ibin, v );
76  }
77 }
78 
79 
80 std::string runfit( TH1F* h ) {
81 
82  double chi2dof = 0;
83  double chi2 = 0;
84  int ndof = 0;
85 
86  std::string fitf = "";
87 
88  TF1* hfit = new TF1( "hyper", "[0]-[3]*( (x-[1]) + sqrt( (x-[1])*(x-[1]) + 4*[2]*[2] ) )/(2*[2])", 0, 5 );
89 
90  hfit->SetParameter(0, 0.95);
91  hfit->SetParameter(1, 1);
92  hfit->SetParameter(2, 1);
93  hfit->SetParameter(3, 1);
94 
95  std::string fitname[4] = { "pol0", "pol1", "hyper", "pol2" };
96 
97  for ( int ifit=0 ; ifit<4 ; ifit++ ) {
98 
99  h->Fit( fitname[ifit].c_str() );
100  TF1* f1 = (TF1*)h->GetListOfFunctions()->FindObject( fitname[ifit].c_str() );
101 
102  f1->SetLineWidth(1);
103  f1->SetLineColor(kRed);
104 
105  chi2 = f1->GetChisquare();
106  ndof = f1->GetNDF();
107 
108  chi2dof = chi2/ndof;
109 
110  fitf = fitname[ifit];
111 
112  char flabel[256];
113 
114  std::sprintf( flabel, "%s: #chi^{2}/dof = %6.2lf/%d", fitf.c_str(), chi2, ndof );
115 
116  std::string chi2label = flabel;
117 
118  if ( chi2dof < 1.5 ) return chi2label;
119 
120  }
121 
122  return "";
123 
124 }
125 
126 
127 
128 void efficiency( std::vector<double>& bins, std::vector<double>& values, const std::vector<std::string>& files,
129  const std::string& histname, const std::string& tplotname, const std::string& label="" ) {
130 
131 
132  std::string plotname = tplotname;
133  std::string labels = ";value;mean";
134 
135  size_t pos = plotname.find(';');
136 
137  if ( pos!=std::string::npos ) {
138  labels = plotname.substr( pos, plotname.size() );
139  plotname.resize(pos);// 'pos' guaranteed to be <= string length
140  std::cout << "plotname: " << plotname << "\tlabels: " << labels << std::endl;
141  }
142 
143 
144  if ( values.size() != files.size() ) {
145  std::cerr << "number of values (" << values.size() << ") and files (" << files.size() << ") do not match" << std::endl;
146  std::exit(-1);
147  }
148 
149  TH1F* hd = 0;
150  TH1F* hn = 0;
151 
152  if ( bins.size()==3 ) {
153  hd = new TH1F( "denominator", labels.c_str(), int(bins[0]), bins[1], bins[2] );
154  hn = new TH1F( "numerator", labels.c_str(), int(bins[0]), bins[1], bins[2] );
155  }
156  else if ( bins.size()>3 ) {
157  hd = new TH1F( "denominator", labels.c_str(), bins.size()-1, &bins[0] );
158  hn = new TH1F( "numerator", labels.c_str(), bins.size()-1, &bins[0] );
159  }
160  else {
161  std::exit( usage( "unusable binning", -1) );
162  }
163 
164  std::cout << "looping over files ..." << std::endl;
165 
166  for ( size_t i=0 ; i<files.size() ; i++ ) {
167 
168  const std::string& s = files[i];
169 
170  std::cout << "file: "<< s << "\tvalue: "<< values[i] << std::endl;
171 
172  TFile* f = new TFile( s.c_str() );
173 
174  if ( f==0 ) {
175  std::cerr << "could not open file" << s << std::endl;
176  std::exit( -1 );
177  }
178 
181 
182  TH1F* he_n = (TH1F*)f->Get( (histname+"_n").c_str() );
183  TH1F* he_d = (TH1F*)f->Get( (histname+"_d").c_str() );
184 
185  if ( he_n==0 ) {
186  std::cerr << "histogram " << histname+"_n" << "could not be retrieved" << std::endl;
187  continue;
188  }
189 
190  if ( he_d==0 ) {
191  std::cerr << "histogram " << histname+"_d" << "could not be retrieved" << std::endl;
192  continue;
193  }
194 
195  add_to_bin( values[i], hn, he_n );
196  add_to_bin( values[i], hd, he_d );
197 
198  delete he_n;
199  delete he_d;
200 
201  delete f;
202  }
203 
204  double scale_eff = 1;
205  Efficiency1D e( hn, hd, "", scale_eff );
206 
207  TGraphAsymmErrors* tgtest = e.Bayes(scale_eff);
208 
209  TCanvas* c2 = new TCanvas( plotname.c_str(), "eff", 10, 10, 700, 500 );
210  c2->cd();
211 
212  TH1F* h = e.Hist();
213 
214  h->SetMarkerStyle(20);
215  h->SetLineColor(h->GetMarkerColor());
216 
217  h->GetYaxis()->SetRangeUser(0.85,1.02);
218  h->GetXaxis()->SetRangeUser(0.4,1.7);
219 
220  tgtest->SetMarkerStyle(20);
221  tgtest->SetMarkerColor(h->GetMarkerColor());
222  tgtest->SetLineColor(h->GetMarkerColor());
223 
224 
225  std::string chi2label = "";
226 
227  if ( doFit ) chi2label = runfit( h );
228 
229  h->GetXaxis()->SetMoreLogLabels(true);
230 
231  h->Draw("e1");
232 
233  tgtest->Draw("samep");
234 
235  if ( doLogx ) gPad->SetLogx(true);
236 
237  DrawLabel( 0.2, 0.25, label );
238  if ( chi2label!="" ) DrawLabel( 0.2, 0.21, chi2label );
239 
240 
241  gPad->Print( plotname.c_str() );
242 
243  delete hd;
244  delete hn;
245  delete c2;
246 
247 
248 }
249 
250 
251 
252 
253 
254 void mean( std::vector<double>& bins, std::vector<double>& values, const std::vector<std::string>& files,
255  const std::string& histname, const std::string& tplotname, const std::string& label="" ) {
256 
257  if ( values.size() != files.size() ) {
258  std::cerr << "number of values (" << values.size() << ") and files (" << files.size() << ") do not match" << std::endl;
259  std::exit(-1);
260  }
261 
262  TH1F* h = 0;
263 
264 
265  std::string plotname = tplotname;
266  std::string labels = ";value;mean";
267 
268  size_t pos = plotname.find(';');
269 
270  if ( pos!=std::string::npos ) {
271  labels = plotname.substr( pos, plotname.size() );
272  plotname.resize(pos); // 'pos' guaranteed to be <= string length
273  std::cout << "plotname: " << plotname << "\tlabels: " << labels << std::endl;
274  }
275 
276 
277  if ( bins.size()==3 ) {
278  h = new TH1F( "mean", labels.c_str(), int(bins[0]), bins[1], bins[2] );
279  }
280  else if ( bins.size()>3 ) {
281  h = new TH1F( "mean", labels.c_str(), bins.size()-1, &bins[0] );
282  }
283  else {
284  std::exit( usage( "unusable binning", -1) );
285  }
286 
287  std::cout << "looping over files ..." << std::endl;
288 
289  for ( size_t i=0 ; i<files.size() ; i++ ) {
290 
291  const std::string& s = files[i];
292 
293  TFile* f = new TFile( s.c_str() );
294 
295  if ( f==0 ) {
296  std::cerr << "could not open file" << s << std::endl;
297  std::exit( -1 );
298  }
299 
301 
302  TH1F* he = (TH1F*)f->Get( histname.c_str() );
303 
304  if ( he==0 ) {
305  std::cerr << "histogram " << histname << "could not be retrieved" << std::endl;
306  continue;
307  }
308 
309  std::cout << "\t he: " << he << std::endl;
310 
311  int ibin = h->FindBin(values[i]);
312  if ( ibin<1 || ibin>h->GetNbinsX() ) continue;
313 
314 
315  std::cout << "bin: " << ibin << std::endl;
316  std::cout << "mean: " << he->GetMean() << " +- " << he->GetMeanError() << std::endl;
317 
318 
319  h->SetBinContent( ibin, he->GetMean() );
320  h->SetBinError( ibin, he->GetMeanError() );
321 
325  delete he;
326  delete f;
327  }
328 
329  // TCanvas* c = new TCanvas( plotname.c_str(), "eff", 10, 10, 1000, 500 );
330  TCanvas* c = new TCanvas( plotname.c_str(), "eff", 10, 10, 700, 500 );
331  c->cd();
332 
333  h->SetMarkerStyle(20);
334  h->SetLineColor(h->GetMarkerColor());
335 
336  std::string chi2label = "";
337 
338  if ( doFit ) chi2label = runfit( h );
339 
340  h->GetXaxis()->SetMoreLogLabels(true);
341 
342  h->Draw("e1");
343 
344  if ( doLogx ) gPad->SetLogx(true);
345 
346  DrawLabel( 0.2, 0.25, label );
347  if ( chi2label!="" ) DrawLabel( 0.2, 0.21, chi2label );
348 
349  gPad->Print( plotname.c_str() );
350 
351  delete h;
352  delete c;
353 
354 }
355 
356 
357 
358 
359 
360 int main( int argc, char** argv ) {
361 
363 
364  std::string config_file = "";
365 
366  bool atlasstyle = false;
367 
369 
370  for ( int i=1 ; i<argc ; i++ ) {
371  std::string arg = argv[i];
372 
373  if ( arg=="-h" || arg=="--help" ) return usage();
374  else if ( arg=="-f" || arg=="--fit" ) doFit = true;
375  else if ( arg=="-as" || arg=="--atlasstyle" ) atlasstyle = true;
376  else if ( arg=="-lx" || arg=="--logx" ) doLogx = true;
377  else if ( arg=="-c" ) {
378  if ( ++i<argc ) config_file = argv[i];
379  else return usage( "no config file provided", -1 );
380  }
381  }
382 
383  if ( config_file == "" ) return usage( "no config file provided", -1 );
384 
385  if ( atlasstyle ) SetAtlasStyle();
386 
387  gStyle->SetPadTopMargin(0.05);
388  gStyle->SetPadRightMargin(0.05);
389 
390  gStyle->SetOptStat(0);
391  gStyle->SetErrorX(0);
392 
394 
395  ReadCards config(config_file);
396 
397  bool do_mean = false;
398  bool do_efficiency = false;
399 
400  if ( config.isTagDefined("pTmax") ) ptmax = config.GetValue( "pTmax" );
401 
402  if ( config.isTagDefined("mean") ) do_mean = ( config.GetValue( "mean")>0 ? true : false );
403  if ( config.isTagDefined("efficiency") ) do_efficiency = ( config.GetValue( "efficiency")>0 ? true : false );
404 
405  // void (*function)( std::vector<double>& bins, std::vector<double>& values, const std::vector<std::string>& files,
406  // const std::string& histname, const std::string& plotname ) = 0;
407 
408  std::vector<void (*)( std::vector<double>& bins, std::vector<double>& values, const std::vector<std::string>& files,
409  const std::string& histname, const std::string& plotname, const std::string& label ) > function;
410 
411  std::vector<std::vector<std::string> > plots;
412 
413  std::vector<std::vector<std::string> > files;
414 
415  if ( do_efficiency ) {
416  function.push_back( efficiency );
417  files.push_back( config.GetStringVector("efiles") );
418  plots.push_back( config.GetStringVector("eplots") );
419  }
420 
421  if ( do_mean ) {
422  function.push_back( mean );
423  files.push_back( config.GetStringVector("mfiles") );
424  plots.push_back( config.GetStringVector("mplots") );
425  }
426 
427  if ( function.empty() ) return 0;
428 
430 
431 
433 
434  std::vector<double> bins = config.GetVector("bins");
435 
436  std::vector<double> values = config.GetVector("values");
437 
439 
440 
441 
442 
443  // output plot name
444 
445  // std::string plotname = config.GetString("plotname");
446 
447 
449 
450  std::cout << "looping over plots ..." << std::endl;
451 
452  for ( size_t j=function.size() ; j-- ; ) {
453 
454  for ( size_t i=0 ; i<plots[j].size() ; i+=3 ) {
455 
456  std::string hist = plots[j][i];
457 
458  std::string plotname = plots[j][i+1];
459 
460  std::string label = plots[j][i+2];
461 
462  std::cout << "plot: " << plotname << " :\t" << hist << "\tlabel: " << label << std::endl;
463 
464  function[j]( bins, values, files[j], hist, plotname, label );
465 
466  }
467 
468  }
469 
470  return 0;
471 
472 }
IDTPM::ndof
float ndof(const U &p)
Definition: TrackParametersHelper.h:127
python.CaloRecoConfig.f
f
Definition: CaloRecoConfig.py:127
TileDCSDataPlotter.h0
h0
Definition: TileDCSDataPlotter.py:873
python.SystemOfUnits.s
int s
Definition: SystemOfUnits.py:131
SetAtlasStyle
void SetAtlasStyle()
Definition: InnerDetector/InDetCalibAlgs/PixelCalibAlgs/Macro/AtlasStyle.h:17
mean
void mean(std::vector< double > &bins, std::vector< double > &values, const std::vector< std::string > &files, const std::string &histname, const std::string &tplotname, const std::string &label="")
Definition: dependence.cxx:254
ptmax
double ptmax
Definition: dependence.cxx:60
doFit
bool doFit
Definition: dependence.cxx:57
python.App.bins
bins
Definition: App.py:410
hist_file_dump.d
d
Definition: hist_file_dump.py:137
PlotCalibFromCool.label
label
Definition: PlotCalibFromCool.py:78
plotmaker.hist
hist
Definition: plotmaker.py:148
LArCellConditions.argv
argv
Definition: LArCellConditions.py:112
empty
bool empty(TH1 *h)
Definition: computils.cxx:294
add_to_bin
void add_to_bin(double d, TH1 *h, TH1 *h0)
Definition: dependence.cxx:62
DrawLabel
int DrawLabel(float xstart, float ystart, string label)
Definition: GraphToolKit.h:13
config
Definition: PhysicsAnalysis/AnalysisCommon/AssociationUtils/python/config.py:1
python.Bindings.values
values
Definition: Control/AthenaPython/python/Bindings.py:797
python.setupRTTAlg.size
int size
Definition: setupRTTAlg.py:39
efficiency
void efficiency(std::vector< double > &bins, std::vector< double > &values, const std::vector< std::string > &files, const std::string &histname, const std::string &tplotname, const std::string &label="")
Definition: dependence.cxx:128
beamspotnt.labels
list labels
Definition: bin/beamspotnt.py:1447
Efficiency1D
Definition: Efficiency1D.h:19
Efficiency1D.h
PixelAthClusterMonAlgCfg.histname
histname
Definition: PixelAthClusterMonAlgCfg.py:106
ReadCards
Get tag-value pairs from a file.
Definition: ReadCards.h:50
ReadCards.h
lumiFormat.i
int i
Definition: lumiFormat.py:92
usage
int usage(const std::string &err_msg="", int status=0)
Definition: dependence.cxx:31
extractSporadic.h
list h
Definition: extractSporadic.py:97
generateReferenceFile.files
files
Definition: generateReferenceFile.py:12
chi2
double chi2(TH1 *h0, TH1 *h1)
Definition: comparitor.cxx:522
runfit
std::string runfit(TH1F *h)
Definition: dependence.cxx:80
DQHistogramMergeRegExp.argc
argc
Definition: DQHistogramMergeRegExp.py:20
calibdata.exit
exit
Definition: calibdata.py:236
create_dcsc_inputs_sqlite.arg
list arg
Definition: create_dcsc_inputs_sqlite.py:48
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:192
compileRPVLLRates.c2
c2
Definition: compileRPVLLRates.py:361
TRT_PAI_physicsConstants::he
const double he
same in ev
Definition: TRT_PAI_physicsConstants.h:22
python.LumiBlobConversion.pos
pos
Definition: LumiBlobConversion.py:18
AtlasStyle.h
ATLAS Style, based on a style file from BaBar.
python.PyAthena.v
v
Definition: PyAthena.py:157
DiTauMassTools::MaxHistStrategyV2::e
e
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:26
main
int main(int argc, char **argv)
Definition: dependence.cxx:360
h
TH1F
Definition: rootspy.cxx:320
TH1
Definition: rootspy.cxx:268
config
std::vector< std::string > config
Definition: fbtTestBasics.cxx:72
covarianceTool.plots
plots
Definition: covarianceTool.py:698
merge.status
status
Definition: merge.py:17
DrawLabel.h
python.TrigEgammaMonitorHelper.TH1F
def TH1F(name, title, nxbins, bins_par2, bins_par3=None, path='', **kwargs)
Definition: TrigEgammaMonitorHelper.py:24
python.compressB64.c
def c
Definition: compressB64.py:93
read_hist_ntuple.f1
f1
Definition: read_hist_ntuple.py:4
doLogx
bool doLogx
Definition: dependence.cxx:58