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