CMS 3D CMS Logo
CalibrationAlgorithm.cc
Go to the documentation of this file.
1 #include "DQM/SiStripCommissioningAnalysis/interface/CalibrationAlgorithm.h"
2 #include "CondFormats/SiStripObjects/interface/CalibrationAnalysis.h"
3 #include "DataFormats/SiStripCommon/interface/SiStripHistoTitle.h"
4 #include "DataFormats/SiStripCommon/interface/SiStripEnumsAndStrings.h"
5 #include "FWCore/MessageLogger/interface/MessageLogger.h"
6 #include "DQM/SiStripCommissioningAnalysis/interface/SiStripPulseShape.h"
7 #include "TProfile.h"
8 #include "TF1.h"
9 #include "TH1.h"
10 #include "TVirtualFitter.h"
11 #include "TFitResultPtr.h"
12 #include "TFitResult.h"
13 #include <iostream>
14 #include <sstream>
15 #include <iomanip>
16 #include <cmath>
17 #include "Math/MinimizerOptions.h"
18 
19 using namespace sistrip;
20 
21 // ----------------------------------------------------------------------------
22 //
23 CalibrationAlgorithm::CalibrationAlgorithm( const edm::ParameterSet & pset, CalibrationAnalysis* const anal )
24  : CommissioningAlgorithm(anal),
25  cal_(nullptr)
26 {}
27 
28 // ----------------------------------------------------------------------------
29 //
30 void CalibrationAlgorithm::extract( const std::vector<TH1*>& histos) {
31 
32  // extract analysis object which should be already created
33  if ( !anal() ) {
34  edm::LogWarning(mlCommissioning_)
35  << "[CalibrationAlgorithm::" << __func__ << "]"
36  << " NULL pointer to base Analysis object!";
37  return;
38  }
39 
40  CommissioningAnalysis* tmp = const_cast<CommissioningAnalysis*>( anal() );
41  cal_ = dynamic_cast<CalibrationAnalysis*>( tmp );
42 
43  if ( !cal_ ) {
44  edm::LogWarning(mlCommissioning_)
45  << "[CalibrationAlgorithm::" << __func__ << "]"
46  << " NULL pointer to derived Analysis object!";
47  return;
48  }
49 
50  // Extract FED key from histo title
51  if ( !histos.empty() ) { cal_->fedKey( extractFedKey( histos.front() ) ); }
52 
53  // Extract histograms
54  std::vector<TH1*>::const_iterator ihis = histos.begin();
55  unsigned int cnt = 0;
56  for ( ; ihis != histos.end(); ihis++,cnt++ ) {
57 
58  // Check for NULL pointer
59  if ( !(*ihis) ) { continue; }
60 
61  // Check name
62  SiStripHistoTitle title( (*ihis)->GetName() );
63  if ( title.runType() != sistrip::CALIBRATION &&
64  title.runType() != sistrip::CALIBRATION_DECO
65  ) {
66  cal_->addErrorCode(sistrip::unexpectedTask_);
67  continue;
68  }
69 
71  std::vector<std::string> tokens;
72  std::string token;
73  std::istringstream tokenStream(title.extraInfo());
74  while (std::getline(tokenStream, token,'_')){
75  tokens.push_back(token);
76  }
77 
79  Histo histo_temp;
80  histo_temp.first = *ihis;
81  histo_temp.second = (*ihis)->GetTitle();
82  histo_temp.first->Sumw2();
83  histo_.push_back(histo_temp);
84  apvId_.push_back(title.channel()%2);
85  stripId_.push_back(std::stoi(tokens.at(1))*16+std::stoi(tokens.at(3)));
86  calChan_.push_back(std::stoi(tokens.at(1)));
87  }
88 }
89 
90 // ----------------------------------------------------------------------------
91 //
92 void CalibrationAlgorithm::analyse() {
93 
94  ROOT::Math::MinimizerOptions::SetDefaultMinimizer("Minuit2","Migrad");
95  ROOT::Math::MinimizerOptions::SetDefaultStrategy(0);
96 
97  if ( !cal_ ) {
98  edm::LogWarning(mlCommissioning_)
99  << "[CalibrationAlgorithm::" << __func__ << "]"
100  << " NULL pointer to derived Analysis object!";
101  return;
102  }
103 
104  float Amean[2] = {-1.,-1.};
105  float Amin[2] = {-1.,-1.};
106  float Amax[2] = {-1.,-1.};
107  float Aspread[2] = {-1.,-1.};
108  float Tmean[2] = {-1.,-1.};
109  float Tmin[2] = {-1.,-1.};
110  float Tmax[2] = {-1.,-1.};
111  float Tspread[2] = {-1.,-1.};
112  float Rmean[2] = {-1.,-1.};
113  float Rmin[2] = {-1.,-1.};
114  float Rmax[2] = {-1.,-1.};
115  float Rspread[2] = {-1.,-1.};
116  float Cmean[2] = {-1.,-1.};
117  float Cmin[2] = {-1.,-1.};
118  float Cmax[2] = {-1.,-1.};
119  float Cspread[2] = {-1.,-1.};
120  float Smean[2] = {-1.,-1.};
121  float Smin[2] = {-1.,-1.};
122  float Smax[2] = {-1.,-1.};
123  float Sspread[2] = {-1.,-1.};
124  float Kmean[2] = {-1.,-1.};
125  float Kmin[2] = {-1.,-1.};
126  float Kmax[2] = {-1.,-1.};
127  float Kspread[2] = {-1.,-1.};
128  // turnOn
129  float Omean[2] = {-1.,-1.};
130  float Omin[2] = {-1.,-1.};
131  float Omax[2] = {-1.,-1.};
132  float Ospread[2] = {-1.,-1.};
133  // maximum
134  float Mmean[2] = {-1.,-1.};
135  float Mmin[2] = {-1.,-1.};
136  float Mmax[2] = {-1.,-1.};
137  float Mspread[2] = {-1.,-1.};
138  // undershoot
139  float Umean[2] = {-1.,-1.};
140  float Umin[2] = {-1.,-1.};
141  float Umax[2] = {-1.,-1.};
142  float Uspread[2] = {-1.,-1.};
143  // baseline
144  float Bmean[2] = {-1.,-1.};
145  float Bmin[2] = {-1.,-1.};
146  float Bmax[2] = {-1.,-1.};
147  float Bspread[2] = {-1.,-1.};
148 
150  TFitResultPtr fit_result;
151  TF1* fit_function = nullptr;
152  if(cal_->deconv_){
153  fit_function = new TF1("fit_function_deco",fdeconv,0,400,7);
154  fit_function->SetParameters(4,25,25,50,250,25,0.75);
155  } else{
156  fit_function = new TF1("fit_function_peak",fpeak,0,400,6);
157  fit_function->SetParameters(4,50,50,70,250,20);
158  }
159 
161  std::vector<unsigned int> nStrips (2,0.);
162 
163  for(size_t ihist = 0; ihist < histo_.size(); ihist++){
164 
165  if ( !histo_[ihist].first ) {
166  edm::LogWarning(mlCommissioning_)
167  << " NULL pointer to histogram for: "<<histo_[ihist].second<<" !";
168  return;
169  }
170 
171  cal_->amplitude_[apvId_[ihist]][stripId_[ihist]] = 0;
172  cal_->baseline_[apvId_[ihist]][stripId_[ihist]] = 0;
173  cal_->riseTime_[apvId_[ihist]][stripId_[ihist]] = 0;
174  cal_->turnOn_[apvId_[ihist]][stripId_[ihist]] = 0;
175  cal_->peakTime_[apvId_[ihist]][stripId_[ihist]] = 0;
176  cal_->undershoot_[apvId_[ihist]][stripId_[ihist]] = 0;
177  cal_->tail_[apvId_[ihist]][stripId_[ihist]] = 0;
178  cal_->decayTime_[apvId_[ihist]][stripId_[ihist]] = 0;
179  cal_->smearing_[apvId_[ihist]][stripId_[ihist]] = 0;
180  cal_->chi2_[apvId_[ihist]][stripId_[ihist]] = 0;
181  cal_->isvalid_[apvId_[ihist]][stripId_[ihist]] = true;
182 
183  if(histo_[ihist].first->Integral() == 0){
184  cal_->isvalid_[apvId_[ihist]][stripId_[ihist]] = false;
185  continue;
186  }
187 
188  // rescale the plot and set reasonable errors
189  correctDistribution(histo_[ihist].first);
190 
191  // from NOTE2009_021 : The charge injection provided by the calibration circuit is known with a precision of 5%
192  float error = histo_[ihist].first->GetMaximum()*0.05;
193  for(int i = 1; i<=histo_[ihist].first->GetNbinsX(); ++i)
194  histo_[ihist].first->SetBinError(i,error);
195 
196  // set intial par
197  if(cal_->deconv_)
198  fit_function->SetParameters(10,15,30,10,350,50,0.75);
199  else
200  fit_function->SetParameters(6,40,40,70,350,20);
201 
202  fit_result = histo_[ihist].first->Fit(fit_function,"QS");
203 
204  // fit-result should exist and have a resonably good status
205  if(not fit_result.Get()) continue;
206 
207  float maximum_ampl = fit_function->GetMaximum();
208  float peak_time = fit_function->GetMaximumX();
209  float baseline = baseLine(fit_function);
210  float turn_on_time = turnOn(fit_function,baseline);
211  float rise_time = peak_time - turn_on_time;
212 
213  // start filling info
214  cal_->amplitude_[apvId_[ihist]][stripId_[ihist]] = maximum_ampl - baseline;
215  cal_->baseline_[apvId_[ihist]][stripId_[ihist]] = baseline;
216  cal_->riseTime_[apvId_[ihist]][stripId_[ihist]] = rise_time;
217  cal_->turnOn_[apvId_[ihist]][stripId_[ihist]] = turn_on_time;
218  cal_->peakTime_[apvId_[ihist]][stripId_[ihist]] = peak_time;
219 
220  if (cal_->deconv_ and fit_function->GetMinimumX() > peak_time) // make sure the minimum is after the peak-time
221  cal_->undershoot_[apvId_[ihist]][stripId_[ihist]] = 100*(fit_function->GetMinimum()-baseline)/(maximum_ampl - baseline);
222  else
223  cal_->undershoot_[apvId_[ihist]][stripId_[ihist]] = 0;
224 
225  // Bin related to peak + 125 ns
226  int lastBin = histo_[ihist].first->FindBin(peak_time + 125);
227  if(lastBin > histo_[ihist].first->GetNbinsX()-4)
228  lastBin = histo_[ihist].first->GetNbinsX()-4;
229 
230  // tail is the amplitude at 5 bx from the maximum
231  cal_->tail_[apvId_[ihist]][stripId_[ihist]] = 100*(histo_[ihist].first->GetBinContent(lastBin)-baseline) / (maximum_ampl - baseline);
232 
233  // reaches 1/e of the peak amplitude
234  cal_->decayTime_[apvId_[ihist]][stripId_[ihist]] = decayTime(fit_function)-peak_time;
235  cal_->smearing_[apvId_[ihist]][stripId_[ihist]] = 0;
236  cal_->chi2_[apvId_[ihist]][stripId_[ihist]] = fit_function->GetChisquare()/(histo_[ihist].first->GetNbinsX()-fit_function->GetNpar());
237 
238  // calibration channel
239  cal_->calChan_ = calChan_[ihist];
240 
241  // apply quality requirements
242  bool isvalid = true;
243  if(not cal_->deconv_){ // peak-mode
244 
245  if(cal_->amplitude_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minAmplitudeThreshold_) isvalid = false;
246  if(cal_->baseline_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minBaselineThreshold_) isvalid = false;
247  else if(cal_->baseline_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxBaselineThreshold_) isvalid = false;
248  if(cal_->decayTime_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minDecayTimeThreshold_) isvalid = false;
249  else if(cal_->decayTime_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxDecayTimeThreshold_) isvalid = false;
250  if(cal_->peakTime_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minPeakTimeThreshold_) isvalid = false;
251  else if(cal_->peakTime_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxPeakTimeThreshold_) isvalid = false;
252  if(cal_->riseTime_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minRiseTimeThreshold_) isvalid = false;
253  else if(cal_->riseTime_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxRiseTimeThreshold_) isvalid = false;
254  if(cal_->turnOn_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minTurnOnThreshold_) isvalid = false;
255  else if(cal_->turnOn_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxTurnOnThreshold_) isvalid = false;
256  if(cal_->chi2_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxChi2Threshold_) isvalid = false;
257 
258  }
259  else{
260 
261  if(fit_function->GetMinimumX() < peak_time) isvalid = false;
262  if(cal_->amplitude_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minAmplitudeThreshold_) isvalid = false;
263  if(cal_->baseline_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minBaselineThreshold_) isvalid = false;
264  if(cal_->baseline_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minBaselineThreshold_) isvalid = false;
265  else if(cal_->baseline_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxBaselineThreshold_) isvalid = false;
266  if(cal_->chi2_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxChi2Threshold_) isvalid = false;
267  if(cal_->turnOn_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minTurnOnThresholdDeco_) isvalid = false;
268  else if(cal_->turnOn_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxTurnOnThresholdDeco_) isvalid = false;
269  if(cal_->decayTime_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minDecayTimeThresholdDeco_) isvalid = false;
270  else if(cal_->decayTime_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxDecayTimeThresholdDeco_) isvalid = false;
271  if(cal_->peakTime_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minPeakTimeThresholdDeco_) isvalid = false;
272  else if(cal_->peakTime_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxPeakTimeThresholdDeco_) isvalid = false;
273  if(cal_->riseTime_[apvId_[ihist]][stripId_[ihist]] < CalibrationAnalysis::minRiseTimeThresholdDeco_) isvalid = false;
274  else if(cal_->riseTime_[apvId_[ihist]][stripId_[ihist]] > CalibrationAnalysis::maxRiseTimeThresholdDeco_) isvalid = false;
275 
276  }
277 
278  if(not isvalid){ // not valid set default to zero for all quantities
279 
280  cal_->amplitude_[apvId_[ihist]][stripId_[ihist]] = 0;
281  cal_->baseline_[apvId_[ihist]][stripId_[ihist]] = 0;
282  cal_->riseTime_[apvId_[ihist]][stripId_[ihist]] = 0;
283  cal_->turnOn_[apvId_[ihist]][stripId_[ihist]] = 0;
284  cal_->peakTime_[apvId_[ihist]][stripId_[ihist]] = 0;
285  cal_->undershoot_[apvId_[ihist]][stripId_[ihist]] = 0;
286  cal_->tail_[apvId_[ihist]][stripId_[ihist]] = 0;
287  cal_->decayTime_[apvId_[ihist]][stripId_[ihist]] = 0;
288  cal_->smearing_[apvId_[ihist]][stripId_[ihist]] = 0;
289  cal_->chi2_[apvId_[ihist]][stripId_[ihist]] = 0;
290  cal_->isvalid_[apvId_[ihist]][stripId_[ihist]] = false;
291  continue;
292 
293  }
294 
295  // in case is valid
296  nStrips[apvId_[ihist]]++;
297 
298  //compute mean, max, min, spread only for valid strips
299  Amean[apvId_[ihist]] += cal_->amplitude_[apvId_[ihist]][stripId_[ihist]];
300  Amin[apvId_[ihist]] = Amin[apvId_[ihist]]<cal_->amplitude_[apvId_[ihist]][stripId_[ihist]] ? Amin[apvId_[ihist]] : cal_->amplitude_[apvId_[ihist]][stripId_[ihist]];
301  Amax[apvId_[ihist]] = Amax[apvId_[ihist]]>cal_->amplitude_[apvId_[ihist]][stripId_[ihist]] ? Amax[apvId_[ihist]] : cal_->amplitude_[apvId_[ihist]][stripId_[ihist]];
302  Aspread[apvId_[ihist]] += cal_->amplitude_[apvId_[ihist]][stripId_[ihist]]*cal_->amplitude_[apvId_[ihist]][stripId_[ihist]];
303 
304  Tmean[apvId_[ihist]] += cal_->tail_[apvId_[ihist]][stripId_[ihist]];
305  Tmin[apvId_[ihist]] = Tmin[apvId_[ihist]]<cal_->tail_[apvId_[ihist]][stripId_[ihist]] ? Tmin[apvId_[ihist]] : cal_->tail_[apvId_[ihist]][stripId_[ihist]];
306  Tmax[apvId_[ihist]] = Tmax[apvId_[ihist]]>cal_->tail_[apvId_[ihist]][stripId_[ihist]] ? Tmax[apvId_[ihist]] : cal_->tail_[apvId_[ihist]][stripId_[ihist]];
307  Tspread[apvId_[ihist]] += cal_->tail_[apvId_[ihist]][stripId_[ihist]]*cal_->tail_[apvId_[ihist]][stripId_[ihist]];
308 
309  Rmean[apvId_[ihist]] += cal_->riseTime_[apvId_[ihist]][stripId_[ihist]];
310  Rmin[apvId_[ihist]] = Rmin[apvId_[ihist]]<cal_->riseTime_[apvId_[ihist]][stripId_[ihist]] ? Rmin[apvId_[ihist]] : cal_->riseTime_[apvId_[ihist]][stripId_[ihist]];
311  Rmax[apvId_[ihist]] = Rmax[apvId_[ihist]]>cal_->riseTime_[apvId_[ihist]][stripId_[ihist]] ? Rmax[apvId_[ihist]] : cal_->riseTime_[apvId_[ihist]][stripId_[ihist]];
312  Rspread[apvId_[ihist]] += cal_->riseTime_[apvId_[ihist]][stripId_[ihist]]*cal_->riseTime_[apvId_[ihist]][stripId_[ihist]];
313 
314  Cmean[apvId_[ihist]] += cal_->decayTime_[apvId_[ihist]][stripId_[ihist]];
315  Cmin[apvId_[ihist]] = Cmin[apvId_[ihist]]<cal_->decayTime_[apvId_[ihist]][stripId_[ihist]] ? Cmin[apvId_[ihist]] : cal_->decayTime_[apvId_[ihist]][stripId_[ihist]];
316  Cmax[apvId_[ihist]] = Cmax[apvId_[ihist]]>cal_->decayTime_[apvId_[ihist]][stripId_[ihist]] ? Cmax[apvId_[ihist]] : cal_->decayTime_[apvId_[ihist]][stripId_[ihist]];
317  Cspread[apvId_[ihist]] += cal_->decayTime_[apvId_[ihist]][stripId_[ihist]]*cal_->decayTime_[apvId_[ihist]][stripId_[ihist]];
318 
319  Smean[apvId_[ihist]] += cal_->smearing_[apvId_[ihist]][stripId_[ihist]];
320  Smin[apvId_[ihist]] = Smin[apvId_[ihist]]<cal_->smearing_[apvId_[ihist]][stripId_[ihist]] ? Smin[apvId_[ihist]] : cal_->smearing_[apvId_[ihist]][stripId_[ihist]];
321  Smax[apvId_[ihist]] = Smax[apvId_[ihist]]>cal_->smearing_[apvId_[ihist]][stripId_[ihist]] ? Smax[apvId_[ihist]] : cal_->smearing_[apvId_[ihist]][stripId_[ihist]];
322  Sspread[apvId_[ihist]] += cal_->smearing_[apvId_[ihist]][stripId_[ihist]]*cal_->smearing_[apvId_[ihist]][stripId_[ihist]];
323 
324  Kmean[apvId_[ihist]] += cal_->chi2_[apvId_[ihist]][stripId_[ihist]];
325  Kmin[apvId_[ihist]] = Kmin[apvId_[ihist]]<cal_->chi2_[apvId_[ihist]][stripId_[ihist]] ? Kmin[apvId_[ihist]] : cal_->chi2_[apvId_[ihist]][stripId_[ihist]];
326  Kmax[apvId_[ihist]] = Kmax[apvId_[ihist]]>cal_->chi2_[apvId_[ihist]][stripId_[ihist]] ? Kmax[apvId_[ihist]] : cal_->chi2_[apvId_[ihist]][stripId_[ihist]];
327  Kspread[apvId_[ihist]] += cal_->chi2_[apvId_[ihist]][stripId_[ihist]]*cal_->chi2_[apvId_[ihist]][stripId_[ihist]];
328 
329  Omean[apvId_[ihist]] += cal_->turnOn_[apvId_[ihist]][stripId_[ihist]];
330  Omin[apvId_[ihist]] = Omin[apvId_[ihist]]<cal_->turnOn_[apvId_[ihist]][stripId_[ihist]] ? Omin[apvId_[ihist]] : cal_->turnOn_[apvId_[ihist]][stripId_[ihist]];
331  Omax[apvId_[ihist]] = Omax[apvId_[ihist]]>cal_->turnOn_[apvId_[ihist]][stripId_[ihist]] ? Omax[apvId_[ihist]] : cal_->turnOn_[apvId_[ihist]][stripId_[ihist]];
332  Ospread[apvId_[ihist]] += cal_->turnOn_[apvId_[ihist]][stripId_[ihist]]*cal_->turnOn_[apvId_[ihist]][stripId_[ihist]];
333 
334  Mmean[apvId_[ihist]] += cal_->peakTime_[apvId_[ihist]][stripId_[ihist]];
335  Mmin[apvId_[ihist]] = Mmin[apvId_[ihist]]<cal_->peakTime_[apvId_[ihist]][stripId_[ihist]] ? Mmin[apvId_[ihist]] : cal_->peakTime_[apvId_[ihist]][stripId_[ihist]];
336  Mmax[apvId_[ihist]] = Mmax[apvId_[ihist]]>cal_->peakTime_[apvId_[ihist]][stripId_[ihist]] ? Mmax[apvId_[ihist]] : cal_->peakTime_[apvId_[ihist]][stripId_[ihist]];
337  Mspread[apvId_[ihist]] += cal_->peakTime_[apvId_[ihist]][stripId_[ihist]]*cal_->peakTime_[apvId_[ihist]][stripId_[ihist]];
338 
339  Umean[apvId_[ihist]] += cal_->undershoot_[apvId_[ihist]][stripId_[ihist]];
340  Umin[apvId_[ihist]] = Umin[apvId_[ihist]]<cal_->undershoot_[apvId_[ihist]][stripId_[ihist]] ? Umin[apvId_[ihist]] : cal_->undershoot_[apvId_[ihist]][stripId_[ihist]];
341  Umax[apvId_[ihist]] = Umax[apvId_[ihist]]>cal_->undershoot_[apvId_[ihist]][stripId_[ihist]] ? Umax[apvId_[ihist]] : cal_->undershoot_[apvId_[ihist]][stripId_[ihist]];
342  Uspread[apvId_[ihist]] += cal_->undershoot_[apvId_[ihist]][stripId_[ihist]]*cal_->undershoot_[apvId_[ihist]][stripId_[ihist]];
343 
344  Bmean[apvId_[ihist]] += cal_->baseline_[apvId_[ihist]][stripId_[ihist]];
345  Bmin[apvId_[ihist]] = Bmin[apvId_[ihist]]<cal_->baseline_[apvId_[ihist]][stripId_[ihist]] ? Bmin[apvId_[ihist]] : cal_->baseline_[apvId_[ihist]][stripId_[ihist]];
346  Bmax[apvId_[ihist]] = Bmax[apvId_[ihist]]>cal_->baseline_[apvId_[ihist]][stripId_[ihist]] ? Bmax[apvId_[ihist]] : cal_->baseline_[apvId_[ihist]][stripId_[ihist]];
347  Bspread[apvId_[ihist]] += cal_->baseline_[apvId_[ihist]][stripId_[ihist]]*cal_->baseline_[apvId_[ihist]][stripId_[ihist]];
348  }
349 
350 
351  // make mean values
352  for(int i=0; i < 2; i++){
353  if(nStrips[i] != 0){
354  Amean[i] = Amean[i]/nStrips[i];
355  Tmean[i] = Tmean[i]/nStrips[i];
356  Rmean[i] = Rmean[i]/nStrips[i];
357  Cmean[i] = Cmean[i]/nStrips[i];
358  Omean[i] = Omean[i]/nStrips[i];
359  Mmean[i] = Mmean[i]/nStrips[i];
360  Umean[i] = Umean[i]/nStrips[i];
361  Bmean[i] = Bmean[i]/nStrips[i];
362  Smean[i] = Smean[i]/nStrips[i];
363  Kmean[i] = Kmean[i]/nStrips[i];
364 
365  Aspread[i] = Aspread[i]/nStrips[i];
366  Tspread[i] = Tspread[i]/nStrips[i];
367  Rspread[i] = Rspread[i]/nStrips[i];
368  Cspread[i] = Cspread[i]/nStrips[i];
369  Ospread[i] = Ospread[i]/nStrips[i];
370  Mspread[i] = Mspread[i]/nStrips[i];
371  Uspread[i] = Uspread[i]/nStrips[i];
372  Bspread[i] = Bspread[i]/nStrips[i];
373  Sspread[i] = Sspread[i]/nStrips[i];
374  Kspread[i] = Kspread[i]/nStrips[i];
375  }
376  }
377 
378  // fill the mean, max, min, spread, ... histograms.
379  for(int i=0; i < 2; ++i) {
380 
381  cal_->mean_amplitude_[i] = Amean[i];
382  cal_->mean_tail_[i] = Tmean[i];
383  cal_->mean_riseTime_[i] = Rmean[i];
384  cal_->mean_decayTime_[i] = Cmean[i];
385  cal_->mean_turnOn_[i] = Omean[i];
386  cal_->mean_peakTime_[i] = Mmean[i];
387  cal_->mean_undershoot_[i] = Umean[i];
388  cal_->mean_baseline_[i] = Bmean[i];
389  cal_->mean_smearing_[i] = Smean[i];
390  cal_->mean_chi2_[i] = Kmean[i];
391 
392  cal_->min_amplitude_[i] = Amin[i];
393  cal_->min_tail_[i] = Tmin[i];
394  cal_->min_riseTime_[i] = Rmin[i];
395  cal_->min_decayTime_[i] = Cmin[i];
396  cal_->min_turnOn_[i] = Omin[i];
397  cal_->min_peakTime_[i] = Mmin[i];
398  cal_->min_undershoot_[i] = Umin[i];
399  cal_->min_baseline_[i] = Bmin[i];
400  cal_->min_smearing_[i] = Smin[i];
401  cal_->min_chi2_[i] = Kmin[i];
402 
403  cal_->max_amplitude_[i] = Amax[i];
404  cal_->max_tail_[i] = Tmax[i];
405  cal_->max_riseTime_[i] = Rmax[i];
406  cal_->max_decayTime_[i] = Cmax[i];
407  cal_->max_turnOn_[i] = Omax[i];
408  cal_->max_peakTime_[i] = Mmax[i];
409  cal_->max_undershoot_[i] = Umax[i];
410  cal_->max_baseline_[i] = Bmax[i];
411  cal_->max_smearing_[i] = Smax[i];
412  cal_->max_chi2_[i] = Kmax[i];
413 
414  cal_->spread_amplitude_[i] = sqrt(fabs(Aspread[i]-Amean[i]*Amean[i]));
415  cal_->spread_tail_[i] = sqrt(fabs(Tspread[i]-Tmean[i]*Tmean[i]));
416  cal_->spread_riseTime_[i] = sqrt(fabs(Rspread[i]-Rmean[i]*Rmean[i]));
417  cal_->spread_decayTime_[i] = sqrt(fabs(Cspread[i]-Cmean[i]*Cmean[i]));
418  cal_->spread_turnOn_[i] = sqrt(fabs(Ospread[i]-Omean[i]*Omean[i]));
419  cal_->spread_peakTime_[i] = sqrt(fabs(Mspread[i]-Mmean[i]*Mmean[i]));
420  cal_->spread_undershoot_[i] = sqrt(fabs(Uspread[i]-Umean[i]*Umean[i]));
421  cal_->spread_baseline_[i] = sqrt(fabs(Bspread[i]-Bmean[i]*Bmean[i]));
422  cal_->spread_smearing_[i] = sqrt(fabs(Sspread[i]-Smean[i]*Smean[i]));
423  cal_->spread_chi2_[i] = sqrt(fabs(Kspread[i]-Kmean[i]*Kmean[i]));
424  }
425 
426  if(fit_function) delete fit_function;
427 }
428 
429 // ------
430 void CalibrationAlgorithm::correctDistribution( TH1* histo ) const {
431  // 20 events per point in the TM loop --> divide by 20 to have the amplitude of a single event readout
432  for(int iBin = 0; iBin < histo->GetNbinsX(); iBin++){
433  histo->SetBinContent(iBin+1,-histo->GetBinContent(iBin+1)/20.);
434  }
435 }
436 
437 // ----------------------------------------------------------------------------
438 float CalibrationAlgorithm::baseLine(TF1* f){
439  float xmax = 10;
440  float baseline = 0;
441  int npoints = 0;
442  float x = f->GetXmin();
443  for( ; x < xmax; x += 0.1) {
444  baseline += f->Eval(x);
445  npoints ++;
446  }
447  return baseline/npoints;
448 }
449 
450 
451 // ----------------------------------------------------------------------------
452 float CalibrationAlgorithm::turnOn(TF1* f, const float & baseline){ // should happen within 100 ns in both deco and peak modes
453  float max_amplitude = f->GetMaximum();
454  float time = 10.;
455  for( ; time < 100 && (f->Eval(time) - baseline) < 0.05 * (max_amplitude - baseline); time += 0.1) {} // flucutation higher than 5% of the pulse height
456  return time;
457 }
458 
459 // ----------------------------------------------------------------------------
460 float CalibrationAlgorithm::decayTime(TF1* f){ // if we approximate the decay to an exp(-t/tau), in one constant unit, the amplited is reduced by e^{-1}
461  float xval = f->GetMaximumX();
462  float max_amplitude = f->GetMaximum();
463  float x = xval;
464  for(; x < 1000; x = x+0.1){
465  if(f->Eval(x) < max_amplitude*exp(-1))
466  break;
467  }
468  return x;
469 }
sistrip::unexpectedTask_
static const char unexpectedTask_[]
Definition: ConstantsForCommissioningAnalysis.h:21
CalibrationAnalysis::minRiseTimeThresholdDeco_
static const float minRiseTimeThresholdDeco_
Definition: CalibrationAnalysis.h:122
CalibrationAnalysis::maxTurnOnThreshold_
static const float maxTurnOnThreshold_
Definition: CalibrationAnalysis.h:116
CalibrationAnalysis::maxRiseTimeThresholdDeco_
static const float maxRiseTimeThresholdDeco_
Definition: CalibrationAnalysis.h:123
CalibrationAnalysis::minPeakTimeThreshold_
static const float minPeakTimeThreshold_
Definition: CalibrationAnalysis.h:111
CalibrationAnalysis::maxDecayTimeThresholdDeco_
static const float maxDecayTimeThresholdDeco_
Definition: CalibrationAnalysis.h:119
mps_fire.i
i
Definition: mps_fire.py:338
CommissioningAnalysis::fedKey
const uint32_t & fedKey() const
Definition: CommissioningAnalysis.h:141
SiStripHistoTitle
Utility class that holds histogram title.
Definition: SiStripHistoTitle.h:20
CalibrationAnalysis::maxRiseTimeThreshold_
static const float maxRiseTimeThreshold_
Definition: CalibrationAnalysis.h:114
CommissioningAlgorithm::Histo
std::pair< TH1 *, std::string > Histo
Definition: CommissioningAlgorithm.h:27
AlCaHLTBitMon_QueryRunRegistry.string
string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
nullptr
#define nullptr
Definition: GCC11Compatibility.h:37
CalibrationAlgorithm::extract
void extract(const std::vector< TH1 * > &) override
Definition: CalibrationAlgorithm.cc:30
muonDTDigis_cfi.pset
pset
Definition: muonDTDigis_cfi.py:27
CalibrationAnalysis::maxBaselineThreshold_
static const float maxBaselineThreshold_
Definition: CalibrationAnalysis.h:106
CalibrationAlgorithm::calChan_
std::vector< int > calChan_
Definition: CalibrationAlgorithm.h:47
sistrip
sistrip classes
Definition: EnsembleCalibrationLA.cc:11
CalibrationAlgorithm::cal_
CalibrationAnalysis * cal_
Definition: CalibrationAlgorithm.h:51
CalibrationAnalysis::minDecayTimeThresholdDeco_
static const float minDecayTimeThresholdDeco_
Definition: CalibrationAnalysis.h:118
CalibrationAnalysis::minAmplitudeThreshold_
static const float minAmplitudeThreshold_
Definition: CalibrationAnalysis.h:104
CalibrationAlgorithm::histo
const Histo & histo(int &i)
Definition: CalibrationAlgorithm.h:26
CalibrationAnalysis::maxDecayTimeThreshold_
static const float maxDecayTimeThreshold_
Definition: CalibrationAnalysis.h:110
CalibrationAnalysis
Analysis for calibration runs.
Definition: CalibrationAnalysis.h:16
CalibrationAnalysis::minTurnOnThresholdDeco_
static const float minTurnOnThresholdDeco_
Definition: CalibrationAnalysis.h:124
CalibrationAnalysis::minDecayTimeThreshold_
static const float minDecayTimeThreshold_
Definition: CalibrationAnalysis.h:109
sistrip::mlCommissioning_
static const char mlCommissioning_[]
Definition: ConstantsForLogger.h:15
CalibrationAlgorithm::apvId_
std::vector< int > apvId_
Definition: CalibrationAlgorithm.h:48
OfflinePixel3DPrimaryVertices_cfi.Tmin
Tmin
Definition: OfflinePixel3DPrimaryVertices_cfi.py:25
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:18
CommissioningAlgorithm::extractFedKey
uint32_t extractFedKey(const TH1 *const )
Definition: CommissioningAlgorithm.cc:31
CalibrationAnalysis::minPeakTimeThresholdDeco_
static const float minPeakTimeThresholdDeco_
Definition: CalibrationAnalysis.h:120
CommissioningAnalysis::addErrorCode
virtual void addErrorCode(const std::string &error)
Definition: CommissioningAnalysis.h:155
f
double f[11][100]
Definition: MuScleFitUtils.cc:78
npoints
static const int npoints
Definition: NuclearInteractionFTFSimulator.h:38
CalibrationAnalysis::minBaselineThreshold_
static const float minBaselineThreshold_
Definition: CalibrationAnalysis.h:105
CalibrationAlgorithm::histo_
std::vector< Histo > histo_
Definition: CalibrationAlgorithm.h:45
CalibrationAnalysis::maxTurnOnThresholdDeco_
static const float maxTurnOnThresholdDeco_
Definition: CalibrationAnalysis.h:125
CalibrationAlgorithm::stripId_
std::vector< int > stripId_
Definition: CalibrationAlgorithm.h:46
CalibrationAnalysis::minRiseTimeThreshold_
static const float minRiseTimeThreshold_
Definition: CalibrationAnalysis.h:113
TrackerOfflineValidation_Dqm_cff.xmax
xmax
Definition: TrackerOfflineValidation_Dqm_cff.py:11
fdeconv
double fdeconv(double *x, double *par)
Definition: SiStripPulseShape.cc:107
CalibrationAlgorithm::correctDistribution
void correctDistribution(TH1 *) const
Definition: CalibrationAlgorithm.cc:430
CalibrationAnalysis::minTurnOnThreshold_
static const float minTurnOnThreshold_
Definition: CalibrationAnalysis.h:115
CalibrationAnalysis::maxChi2Threshold_
static const float maxChi2Threshold_
Definition: CalibrationAnalysis.h:107
tmp
std::vector< std::vector< double > > tmp
Definition: MVATrainer.cc:100
fftjetcommon_cfi.title
title
Definition: fftjetcommon_cfi.py:32
CalibrationAnalysis::maxPeakTimeThresholdDeco_
static const float maxPeakTimeThresholdDeco_
Definition: CalibrationAnalysis.h:121
fpeak
double fpeak(double *x, double *par)
Definition: SiStripPulseShape.cc:69
CalibrationAlgorithm::turnOn
float turnOn(TF1 *, const float &)
Definition: CalibrationAlgorithm.cc:452
CalibrationAnalysis::maxPeakTimeThreshold_
static const float maxPeakTimeThreshold_
Definition: CalibrationAnalysis.h:112
Gflash::Rmax
const double Rmax[kNumberCalorimeter]
Definition: GflashNameSpace.h:31
CommissioningAnalysis
Abstract base for derived classes that provide analysis of commissioning histograms.
Definition: CommissioningAnalysis.h:18
JetChargeProducer_cfi.exp
exp
Definition: JetChargeProducer_cfi.py:6
Gflash::Rmin
const double Rmin[kNumberCalorimeter]
Definition: GflashNameSpace.h:30
CommissioningAlgorithm::anal
CommissioningAnalysis *const anal() const
Definition: CommissioningAlgorithm.h:55
Generated for CMSSW Reference Manual by doxygen 1.8.11