ATLAS Offline Software
MuonDQAFitFunc.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
3 */
4 
6 // Package : MuonDQAFitFunc
7 // Authors: N. Benekos(Illinois) - A. Cortes (Illinois)
8 // Aug. 2008
10 
12 
13 #include "TF1.h"
14 #include "TH1.h"
15 
16 namespace {
17  constexpr double Zmass = 91.1876;
18 }
19 
20 namespace Muon {
21 
22  MuonDQAFitFunc::MuonDQAFitFunc(const std::string& ty,const std::string& na,const IInterface* pa)
23  : AthAlgTool(ty,na,pa),
24  m_minMuonEffWindow(0.8),
25  m_maxMuonEffWindow(1.05)
26  {
27  declareProperty("ZmumuPDGmass", m_ZmumuPDGmass=Zmass);
28  declareProperty("minMuonEffWindow", m_minMuonEffWindow);
29  declareProperty("maxMuonEffWindow", m_maxMuonEffWindow);
30  declareInterface<MuonDQAFitFunc>(this);
31  }
32 
33  //================================================================
34  void MuonDQAFitFunc::ZmumuFitHistograms(TH1F* hmass, TH1F* hwidth, TH1F* h1[], int nbins)
35  {
36  double mass=0.;
37  double mass_error=0.;
38  double width=0.;
39  double width_error=0.;
40  double maxmass=0;
41  double maxwidth=0;
42  TF1 *ZmumFit = new TF1("ZmumFit","gaus",89.,93.);
43  for (int i=0; i<nbins; i++) {
44  if (h1[i]->GetEntries()>50) {
45  ZmumFit->SetParameter(1,3.097);
46  ZmumFit->SetParameter(2,0.05);
47  h1[i]->Fit("ZmumFit","RQNM");
48  mass=ZmumFit->GetParameter(1)-m_ZmumuPDGmass;
49  mass_error=ZmumFit->GetParError(1);
50  width=ZmumFit->GetParameter(2);
51  width_error=ZmumFit->GetParError(2);
52  hmass->SetBinContent(i+1,mass);
53  hmass->SetBinError(i+1,mass_error);
54  hwidth->SetBinContent(i+1,width);
55  hwidth->SetBinError(i+1,width_error);
56  if( (std::abs(mass)+mass_error)>maxmass ) maxmass=std::abs(mass)+mass_error;
57  if( (std::abs(width)+width_error)>maxwidth ) maxwidth=std::abs(width)+width_error;
58  }
59  }
60  if (maxmass>0.1) maxmass=0.1;
61  if (maxwidth>0.2) maxwidth=0.2;
62  hmass->SetAxisRange(-maxmass,maxmass,"Y");
63  hwidth->SetAxisRange(0.,maxwidth,"Y");
64 
65  delete ZmumFit;
66  }
67  //================================================================
68  void MuonDQAFitFunc::MuonEffHisto1D(TH1F* h_Num, TH1F* h_Denom, TProfile* h_Eff) const
69  {
70  int nBins;
71  nBins = h_Num->GetNbinsX();
72  for (int bin=0; bin!=nBins; ++bin) {
73  int nPass = int(h_Num->GetBinContent(bin+1));
74  int nFail = int(h_Denom->GetBinContent(bin+1)) - nPass;
75  double x = h_Denom->GetBinCenter(bin+1);
76  for (int pass=0; pass!=nPass; ++pass) {
77  h_Eff->Fill(x,1.);
78  }
79  for (int fail=0; fail!=nFail; ++fail) {
80  h_Eff->Fill(x,0.);
81  }
82  }
84 
85  }
86 
87  //================================================================
88  void MuonDQAFitFunc::MuonEffHisto2D(TH2F* h_Num, TH2F* h_Denom, TProfile2D* h_Eff) const
89  {
90  int nBins;
91  int nBinsY;
92  nBins = h_Num->GetNbinsX();
93  nBinsY = h_Num->GetNbinsY();
94  for (int bin=0; bin!=nBins; ++bin) {
95  for (int binY=0; binY!=nBinsY; ++binY) {
96  int nPass = int(h_Num->GetBinContent(bin+1,binY+1));
97  int nFail = int(h_Denom->GetBinContent(bin+1,binY+1)) - nPass;
98  double x = h_Denom->GetXaxis()->GetBinCenter(bin+1);
99  double y = h_Denom->GetYaxis()->GetBinCenter(binY+1);
100  for (int pass=0; pass!=nPass; ++pass) {
101  h_Eff->Fill(x,y,1.);
102  }
103  for (int fail=0; fail!=nFail; ++fail) {
104  h_Eff->Fill(x,y,0.);
105  }
106  }
107  }
109 
110  }
111 
112  //================================================================
113  // Establishes a minimim window for the TProfile
114  //================================================================
115  void MuonDQAFitFunc::MinWindow1Set_from_TProf(TProfile* hProf, float windowMin, float windowMax) const
116  {
117  float min=hProf->GetMinimum();
118  float max=hProf->GetMaximum();
119  float margin=0.05*(max-min);
120  if(min > windowMin) min=windowMin-margin;
121  if(max < windowMax) max=windowMax+margin;
122  hProf->SetMinimum(min);
123  hProf->SetMaximum(max);
124 
125  }
126 
127  //================================================================
128  void MuonDQAFitFunc::MinWindow2Set_from_TProf(TProfile2D* hProf, float windowMin, float windowMax) const
129  {
130  float min=hProf->GetMinimum();
131  float max=hProf->GetMaximum();
132  float margin=0.05*(max-min);
133  if(min > windowMin) min=windowMin-margin;
134  if(max < windowMax) max=windowMax+margin;
135  hProf->SetMinimum(min);
136  hProf->SetMaximum(max);
137  }
138 
139  //================================================================
140  void MuonDQAFitFunc::FillMeanRmsProj(TH2F* h2d, TH1F* h1d,int MeanRms) const
141  {
142 
143  // fills a 1-D histogram with either the mean or RMS of the residual distribution for each layer in the barrel.
144  // Does this by projecting a 2-D histo of residual vs layer.
145 
146  int nBins_2d = h2d->GetNbinsX();
147  int nBins_1d = h1d->GetNbinsX();
148 
149  if(nBins_2d!=nBins_1d) ATH_MSG_DEBUG("Mean/RMS Histograms not set up correctly - nBins mismatch");
150 
151  //calling this means that the histogram bin content is flagged
152  //as being an average and so adding histos from different jobs
153  //will produce weighted mean
154 
155  h1d->SetBit(TH1::kIsAverage);
156 
157  for(int j = 1; j!=nBins_2d+1; j++){
158 
159  TH1F* hproj = (TH1F*)h2d->ProjectionY("Proj",j,j,"e");
160 
161  //do not fill if there are no entries in the bin
162  if(hproj->GetEntries()<=0) {
163  delete hproj;
164  continue;
165  }
166 
167  if(MeanRms==0){
168  h1d->SetBinContent(j,hproj->GetMean());
169  h1d->SetBinError(j,hproj->GetMeanError());
170 
171  }
172  else if(MeanRms==1){
173  h1d->SetBinContent(j,hproj->GetRMS());
174  h1d->SetBinError(j,hproj->GetRMSError());
175  }
176  else ATH_MSG_DEBUG("Incorrect switch in MeanRMSProjectionsBarrel()");
177 
178  delete hproj;
179  }
180 
181  return;
182  }
183  //================================================================
184 
186  {
187  //Uses the spread information in an already filled TProfile
188  //to fill a second TH1F with the (correctly weigthed by N events) RMS of each bin
189 
190  int nBins = hProf->GetNbinsX();
191 
192  for(int j = 1; j!=nBins+1; j++){
193 
194  float bincentre = (float)hProf->GetBinCenter(j);
195 
196  //By default this is RMS/sqrt(N)
197  float binerr = (float)hProf->GetBinError(j);
198 
199  //getting binentries
200  float binentries = (float)hProf->GetBinEntries(j);
201 
202  //calculating the RMS from the above quantities
203  float rms = (float)binerr*(std::sqrt(binentries));
204 
205  hRms->Fill(bincentre,rms);
206  } // Loop over bins
207 
208  return;
209  }
210 
211  //================================================================
212 
213  void MuonDQAFitFunc::FillRMSFromCharProfile(TProfile* hProf_char, TH1F* hRms_char) const
214  {
215 
216  //Uses the spread information in an already filled TProfile (with character bins)
217  //to fill a second TH1F with the (correctly weigthed by N events) RMS of each bin
218 
219  int nBins = hProf_char->GetNbinsX();
220 
221  for(int j = 1; j!=nBins+1; j++){
222 
223  //float bincentre = (float)hProf->GetBinCenter(j);
224  const char* BinLabel_char = hProf_char->GetXaxis()->GetBinLabel(j);
225 
226  //By default this is RMS/sqrt(N)
227  float binerr = (float)hProf_char->GetBinError(j);
228 
229  //getting binentries
230  float binentries = (float)hProf_char->GetBinEntries(j);
231 
232  //calculating the RMS from the above quantities
233  float rms = (float)binerr*(std::sqrt(binentries));
234 
235  hRms_char->SetCanExtend(TH1::kAllAxes);
236  hRms_char->Fill(BinLabel_char,rms);
237 
238  } // Loop over bins
239 
240  hRms_char->LabelsDeflate("X");
241  if (hRms_char->GetEntries() > 0) hRms_char->LabelsOption("v");
242  return;
243  }
244 
245  //================================================================
246  void MuonDQAFitFunc::FillGausMeanOrWidth(TH2F* h2d, TH1F* h1d, float MinFit, float MaxFit, int iopt) const
247  {
248 
249 
250  //Makes a Gaussian fit to each bin of a TH2F and fills a TH1F with
251  //the either the mean or width of this Gaussian
252 
253  //calling this means that the histogram bin content is flagged
254  //as being an average and so adding histos from different jobs
255  //will produce weighted mean
256  h1d->SetBit(TH1::kIsAverage);
257 
258  int nBins_2d = h2d->GetNbinsX();
259  int nBins_1d = h1d->GetNbinsX();
260 
261  if(nBins_2d!=nBins_1d) ATH_MSG_DEBUG("Mean/Width Histograms not set up correctly - nBins mismatch");
262 
263  for(int i = 1; i!=nBins_2d+1; i++){
264 
265  TH1F* hProj = (TH1F*)h2d->ProjectionY("Proj",i,i,"e");
266 
267  //do not fill if there are no entries in the bin
268  if(hProj->GetEntries()<=0) {
269  delete hProj;
270  continue;
271  }
272 
273  TF1 *fit = new TF1("fit","gaus",MinFit,MaxFit);
274 
275  hProj->Fit("fit","RQNM");
276  float Mean = fit->GetParameter(1);
277  float MeanSigma = fit->GetParError(1);
278  float Width = fit->GetParameter(2);
279  float WidthSigma = fit->GetParError(2);
280 
281  if(iopt==0){
282  h1d->SetBinContent(i,Mean);
283  h1d->SetBinError(i,MeanSigma);
284  }
285  else if(iopt==1){
286  h1d->SetBinContent(i,Width);
287  h1d->SetBinError(i,WidthSigma);
288  }
289  else ATH_MSG_DEBUG("Incorrect switch in FillGausMeanOrWidth");
290 
291  delete hProj;
292  delete fit;
293  }
294 
295  return;
296 
297  }
298 
299 
300  //================================================================
301  void MuonDQAFitFunc::FillGaussMeanAndWidth(TH2F* h2d, TH1F* h1d, float MinFit, float MaxFit) const
302  {
303 
304 
305  //Makes a Gaussian fit to each bin of a TH2F and fills a TH1F,
306  // the mean of this Gaussian, and the errors the width.
307 
308  //calling this means that the histogram bin content is flagged
309  //as being an average and so adding histos from different jobs
310  //will produce weighted mean
311  h1d->SetBit(TH1::kIsAverage);
312 
313  int nBins_2d = h2d->GetNbinsX();
314  int nBins_1d = h1d->GetNbinsX();
315 
316  if(nBins_2d!=nBins_1d) ATH_MSG_DEBUG("Mean/Width Histograms not set up correctly - nBins mismatch");
317 
318  for(int i = 1; i!=nBins_2d+1; i++){
319 
320  TH1F* hProj = (TH1F*)h2d->ProjectionY("Proj",i,i,"e");
321 
322  int numNonEmptyBins = 0;
323  for (int ii = 0; ii < (int)hProj->GetNbinsX()+1; ++ii)
324  {
325  if (hProj->GetBinContent(ii)>0.) ++numNonEmptyBins;
326  }
327 
328  //do not fill if there are no entries in the bin
329  //if(hProj->GetEntries()<=0) {
330  if((hProj->GetEntries()<=0)||(numNonEmptyBins<2)) {
331  delete hProj;
332  continue;
333  }
334 
335  TF1 *fit = new TF1("fit","gaus",MinFit,MaxFit);
336 
337  hProj->Fit("fit","RQNM");
338  float Mean = fit->GetParameter(1);
339  // float MeanSigma = fit->GetParError(1);
340  float Width = fit->GetParameter(2);
341  // float WidthSigma = fit->GetParError(2);
342 
343  h1d->SetBinContent(i,Mean);
344  h1d->SetBinError(i, Width);
345 
346  delete hProj;
347  delete fit;
348  }
349 
350  h1d->SetEntries(h2d->GetEntries());
351  return;
352 
353  }
354 
355 
356 
357 
358  //================================================================
359  void MuonDQAFitFunc::FillGaussMeanAndWidth(TH1F* h1, TH1F* h2, int nbin, float MinFit, float MaxFit) const
360  {
361  // Makes a Gaussian fit to a TH1F and fills the given bin of another TH1F
362  // with the mean of this Gaussian, and the errors with the width.
363 
364  //calling this means that the histogram bin content is flagged
365  //as being an average and so adding histos from different jobs
366  //will produce weighted mean
367  h2->SetBit(TH1::kIsAverage);
368 
369  //do not fill if there are no entries in h1
370  if(h1->GetEntries()<=0) return;
371 
372  TF1 *fit = new TF1("fit","gaus",MinFit,MaxFit);
373  h1->Fit("fit","RQNM");
374  float Mean = fit->GetParameter(1);
375  // float MeanSigma = fit->GetParError(1);
376  float Width = fit->GetParameter(2);
377  // float WidthSigma = fit->GetParError(2);
378 
379  h2->SetBinContent(nbin,Mean);
380  h2->SetBinError(nbin,Width);
381  delete fit;
382  h2->SetEntries(h1->GetEntries());
383  return;
384  }
385 
386 } // namespace
Muon::MuonDQAFitFunc::MuonEffHisto2D
void MuonEffHisto2D(TH2F *h_Num, TH2F *h_Denom, TProfile2D *h_Eff) const
Definition: MuonDQAFitFunc.cxx:88
Muon::MuonDQAFitFunc::FillRMSFromCharProfile
void FillRMSFromCharProfile(TProfile *, TH1F *) const
Fills a histogram with the RMS values of a TProfile, which labeled bins->keeping the labeling for the...
Definition: MuonDQAFitFunc.cxx:213
CaloCellPos2Ntuple.int
int
Definition: CaloCellPos2Ntuple.py:24
Base_Fragment.mass
mass
Definition: Sherpa_i/share/common/Base_Fragment.py:59
AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T > &t)
Definition: AthCommonDataStore.h:145
max
constexpr double max()
Definition: ap_fixedTest.cxx:33
Muon::MuonDQAFitFunc::FillGausMeanOrWidth
void FillGausMeanOrWidth(TH2F *, TH1F *, float, float, int) const
Definition: MuonDQAFitFunc.cxx:246
min
constexpr double min()
Definition: ap_fixedTest.cxx:26
Muon::MuonDQAFitFunc::ZmumuFitHistograms
void ZmumuFitHistograms(TH1F *hmass, TH1F *hwidth, TH1F *h1[], int nbins)
Definition: MuonDQAFitFunc.cxx:34
bin
Definition: BinsDiffFromStripMedian.h:43
GetEntries
TGraphErrors * GetEntries(TH2F *histo)
Definition: TRTCalib_makeplots.cxx:4019
Muon
NRpcCablingAlg reads raw condition data and writes derived condition data to the condition store.
Definition: TrackSystemController.h:45
read_hist_ntuple.h1
h1
Definition: read_hist_ntuple.py:21
python.TrigEgammaMonitorHelper.TH2F
def TH2F(name, title, nxbins, bins_par2, bins_par3, bins_par4, bins_par5=None, bins_par6=None, path='', **kwargs)
Definition: TrigEgammaMonitorHelper.py:45
x
#define x
Trk::binY
@ binY
Definition: BinningType.h:48
lumiFormat.i
int i
Definition: lumiFormat.py:85
Muon::MuonDQAFitFunc::m_ZmumuPDGmass
double m_ZmumuPDGmass
Definition: MuonDQAFitFunc.h:71
MuonDQAFitFunc.h
python.TrigEgammaMonitorHelper.TProfile
def TProfile(*args, **kwargs)
Definition: TrigEgammaMonitorHelper.py:81
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
Muon::MuonDQAFitFunc::FillMeanRmsProj
void FillMeanRmsProj(TH2F *, TH1F *, int) const
Definition: MuonDQAFitFunc.cxx:140
Muon::MuonDQAFitFunc::MinWindow2Set_from_TProf
void MinWindow2Set_from_TProf(TProfile2D *hProf, float windowMin, float windowMax) const
Definition: MuonDQAFitFunc.cxx:128
Muon::MuonDQAFitFunc::m_maxMuonEffWindow
float m_maxMuonEffWindow
Definition: MuonDQAFitFunc.h:68
MuonValidation_CreateResolutionProfiles.fit
def fit(h, emin, emax)
Definition: MuonValidation_CreateResolutionProfiles.py:69
Muon::MuonDQAFitFunc::MuonDQAFitFunc
MuonDQAFitFunc(const std::string &type, const std::string &name, const IInterface *parent)
constructor
Definition: MuonDQAFitFunc.cxx:22
dumpTgcDigiJitter.nBins
list nBins
Definition: dumpTgcDigiJitter.py:29
plotBeamSpotVxVal.bin
int bin
Definition: plotBeamSpotVxVal.py:83
SCT_CalibAlgs::nbins
@ nbins
Definition: SCT_CalibNumbers.h:10
Muon::MuonDQAFitFunc::FillGaussMeanAndWidth
void FillGaussMeanAndWidth(TH2F *, TH1F *, float, float) const
Definition: MuonDQAFitFunc.cxx:301
Muon::MuonDQAFitFunc::FillRMSFromProfile
void FillRMSFromProfile(TProfile *, TH1F *) const
Fills a histogram with the RMS values of a TProfile.
Definition: MuonDQAFitFunc.cxx:185
y
#define y
Base_Fragment.width
width
Definition: Sherpa_i/share/common/Base_Fragment.py:59
Muon::MuonDQAFitFunc::MuonEffHisto1D
void MuonEffHisto1D(TH1F *h_Num, TH1F *h_Denom, TProfile *h_Eff) const
Definition: MuonDQAFitFunc.cxx:68
beamspotnt.rms
rms
Definition: bin/beamspotnt.py:1266
Muon::MuonDQAFitFunc::m_minMuonEffWindow
float m_minMuonEffWindow
Definition: MuonDQAFitFunc.h:67
python.TrigEgammaMonitorHelper.TH1F
def TH1F(name, title, nxbins, bins_par2, bins_par3=None, path='', **kwargs)
Definition: TrigEgammaMonitorHelper.py:24
Muon::MuonDQAFitFunc::MinWindow1Set_from_TProf
void MinWindow1Set_from_TProf(TProfile *hProf, float windowMin, float windowMax) const
================================================================
Definition: MuonDQAFitFunc.cxx:115
AthAlgTool
Definition: AthAlgTool.h:26
readCCLHist.float
float
Definition: readCCLHist.py:83
LArMonTransforms.Mean
def Mean(inputs)
Definition: LArMonTransforms.py:438
beamspotman.fail
def fail(message)
Definition: beamspotman.py:201