File indexing completed on 2025-03-23 15:57:17
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include "SiPixelGainCalibrationAnalysis.h"
0021 #include <sstream>
0022 #include <vector>
0023 #include <cmath>
0024 #include "TGraphErrors.h"
0025 #include "TMath.h"
0026 #include "FWCore/Utilities/interface/isFinite.h"
0027
0028 using std::cout;
0029 using std::endl;
0030
0031
0032
0033 SiPixelGainCalibrationAnalysis::SiPixelGainCalibrationAnalysis(const edm::ParameterSet &iConfig)
0034 : SiPixelOfflineCalibAnalysisBase(iConfig),
0035 conf_(iConfig),
0036 bookkeeper_(),
0037 bookkeeper_pixels_(),
0038 nfitparameters_(iConfig.getUntrackedParameter<int>("numberOfFitParameters", 2)),
0039 fitfunction_(iConfig.getUntrackedParameter<std::string>("fitFunctionRootFormula", "pol1")),
0040 listofdetids_(conf_.getUntrackedParameter<std::vector<uint32_t> >("listOfDetIDs")),
0041 ignoreMode_(conf_.getUntrackedParameter<bool>("ignoreMode", false)),
0042 reject_plateaupoints_(iConfig.getUntrackedParameter<bool>("suppressPlateauInFit", true)),
0043 reject_single_entries_(iConfig.getUntrackedParameter<bool>("suppressPointsWithOneEntryOrLess", true)),
0044 plateau_max_slope_(iConfig.getUntrackedParameter<double>("plateauSlopeMax", 1.0)),
0045 reject_first_point_(iConfig.getUntrackedParameter<bool>("rejectVCalZero", true)),
0046 reject_badpoints_frac_(iConfig.getUntrackedParameter<double>("suppressZeroAndPlateausInFitFrac", 0)),
0047 bookBIGCalibPayload_(iConfig.getUntrackedParameter<bool>("saveFullPayloads", false)),
0048 savePixelHists_(iConfig.getUntrackedParameter<bool>("savePixelLevelHists", false)),
0049 chi2Threshold_(iConfig.getUntrackedParameter<double>("minChi2NDFforHistSave", 10)),
0050 chi2ProbThreshold_(iConfig.getUntrackedParameter<double>("minChi2ProbforHistSave", 0.05)),
0051 maxGainInHist_(iConfig.getUntrackedParameter<double>("maxGainInHist", 10)),
0052 maxChi2InHist_(iConfig.getUntrackedParameter<double>("maxChi2InHist", 25)),
0053 saveALLHistograms_(iConfig.getUntrackedParameter<bool>("saveAllHistograms", false)),
0054
0055 filldb_(iConfig.getUntrackedParameter<bool>("writeDatabase", false)),
0056 writeSummary_(iConfig.getUntrackedParameter<bool>("writeSummary", true)),
0057 recordName_(conf_.getParameter<std::string>("record")),
0058
0059 appendMode_(conf_.getUntrackedParameter<bool>("appendMode", true)),
0060
0061
0062
0063
0064 gainlow_(10.),
0065 gainhi_(0.),
0066 pedlow_(255.),
0067 pedhi_(0.),
0068 useVcalHigh_(conf_.getParameter<bool>("useVCALHIGH")),
0069 scalarVcalHigh_VcalLow_(conf_.getParameter<double>("vcalHighToLowConversionFac")) {
0070 if (reject_single_entries_)
0071 min_nentries_ = 1;
0072 else
0073 min_nentries_ = 0;
0074 ::putenv((char *)"CORAL_AUTH_USER=me");
0075 ::putenv((char *)"CORAL_AUTH_PASSWORD=test");
0076 edm::LogInfo("SiPixelGainCalibrationAnalysis") << "now using fit function " << fitfunction_ << ", which has "
0077 << nfitparameters_ << " free parameters. " << std::endl;
0078 func_ = new TF1("func", fitfunction_.c_str(), 0, 256 * scalarVcalHigh_VcalLow_);
0079 graph_ = new TGraphErrors();
0080 currentDetID_ = 0;
0081 summary_.open("SummaryPerDetID.txt");
0082 statusNumbers_ = new int[10];
0083 for (int ii = 0; ii < 10; ii++)
0084 statusNumbers_[ii] = 0;
0085 }
0086
0087 SiPixelGainCalibrationAnalysis::~SiPixelGainCalibrationAnalysis() {}
0088
0089
0090
0091
0092 std::vector<float> SiPixelGainCalibrationAnalysis::CalculateAveragePerColumn(uint32_t detid, std::string label) {
0093 std::vector<float> result;
0094 int ncols = bookkeeper_[detid][label]->getNbinsX();
0095 int nrows = bookkeeper_[detid][label]->getNbinsY();
0096 for (int icol = 1; icol <= ncols; ++icol) {
0097 float val = 0;
0098 float ntimes = 0;
0099 for (int irow = 1; irow <= nrows; ++irow) {
0100 val += bookkeeper_[detid][label]->getBinContent(icol, irow);
0101 ntimes++;
0102 }
0103 val /= ntimes;
0104 result.push_back(val);
0105 }
0106 return result;
0107 }
0108
0109 bool SiPixelGainCalibrationAnalysis::checkCorrectCalibrationType() {
0110 if (calibrationMode_ == "GainCalibration")
0111 return true;
0112 else if (ignoreMode_ == true)
0113 return true;
0114 else if (calibrationMode_ == "unknown") {
0115 edm::LogInfo("SiPixelGainCalibrationAnalysis")
0116 << "calibration mode is: " << calibrationMode_ << ", continuing anyway...";
0117 return true;
0118 } else {
0119
0120 }
0121 return false;
0122 }
0123
0124 void SiPixelGainCalibrationAnalysis::calibrationSetup(const edm::EventSetup &) {}
0125
0126 void SiPixelGainCalibrationAnalysis::printSummary() {
0127 uint32_t detid = 0;
0128 for (std::map<uint32_t, std::map<std::string, MonitorElement *> >::const_iterator idet = bookkeeper_.begin();
0129 idet != bookkeeper_.end();
0130 ++idet) {
0131 if (detid == idet->first)
0132 continue;
0133 detid = idet->first;
0134 std::vector<float> gainvec = CalculateAveragePerColumn(detid, "gain_2d");
0135 std::vector<float> pedvec = CalculateAveragePerColumn(detid, "ped_2d");
0136 std::vector<float> chi2vec = CalculateAveragePerColumn(detid, "chi2_2d");
0137 std::ostringstream summarytext;
0138
0139 summarytext << "Summary for det ID " << detid << "(" << translateDetIdToString(detid) << ")\n";
0140 summarytext << "\t Following: values per column: column #, gain, pedestal, chi2\n";
0141 for (uint32_t i = 0; i < gainvec.size(); i++)
0142 summarytext << "\t " << i << " \t" << gainvec[i] << " \t" << pedvec[i] << " \t" << chi2vec[i] << "\n";
0143 summarytext << "\t list of pixels with high chi2 (chi2> " << chi2Threshold_ << "): \n";
0144
0145 for (std::map<std::string, MonitorElement *>::const_iterator ipix = bookkeeper_pixels_[detid].begin();
0146 ipix != bookkeeper_pixels_[detid].end();
0147 ++ipix)
0148 summarytext << "\t " << ipix->first << "\n";
0149 edm::LogInfo("SiPixelGainCalibrationAnalysis") << summarytext.str() << std::endl;
0150 }
0151 if (summary_.is_open()) {
0152 summary_.close();
0153 summary_.open("Summary.txt");
0154 summary_ << "Total Number of Pixel computed :" << statusNumbers_[9] << endl;
0155 summary_ << "Number of pixel tagged with status :" << endl;
0156 for (int ii = 0; ii < 9; ii++)
0157 summary_ << ii << " -> " << statusNumbers_[ii] << " ~ "
0158 << double(statusNumbers_[ii]) / double(statusNumbers_[9]) * 100. << " %" << endl;
0159
0160 summary_.close();
0161 }
0162 }
0163
0164
0165
0166 void SiPixelGainCalibrationAnalysis::calibrationEnd() {
0167 if (writeSummary_)
0168 printSummary();
0169
0170
0171 if (filldb_)
0172 fillDatabase();
0173 }
0174
0175 void SiPixelGainCalibrationAnalysis::fillDatabase() {
0176
0177
0178 edm::LogError("SiPixelGainCalibration::fillDatabase()")
0179 << "PLEASE do not fill the database directly from the gain calibration analyzer. This function is currently "
0180 "disabled and no DB payloads will be produced!"
0181 << std::endl;
0182 }
0183
0184 bool SiPixelGainCalibrationAnalysis::doFits(uint32_t detid, std::vector<SiPixelCalibDigi>::const_iterator ipix) {
0185 float lowmeanval = 255;
0186 float highmeanval = 0;
0187 bool makehistopersistent = saveALLHistograms_;
0188 std::vector<uint32_t>::const_iterator detidfinder = find(listofdetids_.begin(), listofdetids_.end(), detid);
0189 if (detidfinder != listofdetids_.end())
0190 makehistopersistent = true;
0191
0192 double xvals[257];
0193 double yvals[256];
0194 double yerrvals[256];
0195 double xvalsall[257];
0196 double yvalsall[256];
0197 double yerrvalsall[256];
0198 int npoints = 0;
0199 int nallpoints = 0;
0200 bool use_point = true;
0201 int status = 0;
0202 statusNumbers_[9]++;
0203
0204 bookkeeper_[detid]["status_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, 0);
0205 if (writeSummary_ && detid != currentDetID_) {
0206 currentDetID_ = detid;
0207 summary_ << endl << "DetId_" << currentDetID_ << endl;
0208 }
0209
0210 for (uint32_t ii = 0; ii < ipix->getnpoints() && ii < 200; ii++) {
0211
0212 nallpoints++;
0213 if (useVcalHigh_) {
0214 xvalsall[ii] = vCalValues_[ii] * scalarVcalHigh_VcalLow_;
0215 } else
0216 xvalsall[ii] = vCalValues_[ii];
0217 yerrvalsall[ii] = yvalsall[ii] = 0;
0218
0219 if (ipix->getnentries(ii) > min_nentries_) {
0220 yvalsall[ii] = ipix->getsum(ii) / (float)ipix->getnentries(ii);
0221 yerrvalsall[ii] = ipix->getsumsquares(ii) / (float)(ipix->getnentries(ii));
0222 yerrvalsall[ii] -= pow(yvalsall[ii], 2);
0223 yerrvalsall[ii] = sqrt(yerrvalsall[ii]) / sqrt(ipix->getnentries(ii));
0224
0225 if (yvalsall[ii] < lowmeanval)
0226 lowmeanval = yvalsall[ii];
0227 if (yvalsall[ii] > highmeanval)
0228 highmeanval = yvalsall[ii];
0229 }
0230 }
0231
0232
0233 double plateauval = 0;
0234 bool noPlateau = false;
0235 if (nallpoints >= 4) {
0236 for (int ii = nallpoints - 1; ii > nallpoints - 5; --ii)
0237 plateauval += yvalsall[ii];
0238 plateauval /= 4;
0239 for (int ii = nallpoints - 1; ii > nallpoints - 5; --ii) {
0240 if (fabs(yvalsall[ii] - plateauval) > 5) {
0241 plateauval = 255;
0242 noPlateau = true;
0243 continue;
0244 }
0245 }
0246
0247 int NbofPointsInPlateau = 0;
0248 for (int ii = 0; ii < nallpoints; ++ii)
0249 if (fabs(yvalsall[ii] - plateauval) < 10 || yvalsall[ii] == 0)
0250 NbofPointsInPlateau++;
0251
0252 if (NbofPointsInPlateau >= (nallpoints - 2)) {
0253 status = -2;
0254 bookkeeper_[detid]["status_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, status);
0255 if (writeSummary_) {
0256 summary_ << "row_" << ipix->row() << " col_" << ipix->col() << " status_" << status << endl;
0257 statusNumbers_[abs(status)]++;
0258 }
0259 return false;
0260 }
0261 } else
0262 plateauval = 255;
0263
0264 double maxgoodvalinfit = plateauval * (1. - reject_badpoints_frac_);
0265 npoints = 0;
0266 for (int ii = 0; ii < nallpoints; ++ii) {
0267
0268 use_point = true;
0269 if (reject_first_point_ && xvalsall[ii] < 0.1)
0270 use_point = false;
0271 if (ipix->getnentries(ii) <= min_nentries_ && reject_single_entries_)
0272 use_point = false;
0273 if (ipix->getnentries(ii) == 0 && reject_badpoints_)
0274 use_point = false;
0275 if (yvalsall[ii] > maxgoodvalinfit && !noPlateau)
0276 use_point = false;
0277 if (ii > 1 && fabs(yvalsall[ii] - yvalsall[ii - 1]) < 5. && yvalsall[ii] > 0.8 * maxgoodvalinfit &&
0278 reject_plateaupoints_) {
0279 break;
0280 }
0281
0282 if (use_point) {
0283 xvals[npoints] = xvalsall[ii];
0284 yvals[npoints] = yvalsall[ii];
0285 yerrvals[npoints] = yerrvalsall[ii];
0286 npoints++;
0287 }
0288 }
0289
0290 float chi2, slope, intercept, prob, slopeerror, intercepterror;
0291
0292
0293
0294 if (npoints < 4) {
0295 npoints = 0;
0296 for (int ii = 0; ii < nallpoints && npoints < 4 && yvalsall[ii] < plateauval * 0.97; ++ii) {
0297 if (yvalsall[ii] > 0) {
0298 if (ii > 0 && yvalsall[ii] - yvalsall[ii - 1] < 0.1)
0299 continue;
0300 xvals[npoints] = xvalsall[ii];
0301 yvals[npoints] = yvalsall[ii];
0302 yerrvals[npoints] = yerrvalsall[ii];
0303 npoints++;
0304 }
0305 }
0306 }
0307 if (npoints < 2) {
0308 status = -7;
0309 bookkeeper_[detid]["status_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, status);
0310 if (writeSummary_) {
0311 summary_ << "row_" << ipix->row() << " col_" << ipix->col() << " status_" << status << endl;
0312 statusNumbers_[abs(status)]++;
0313 }
0314 std::ostringstream pixelinfo;
0315 pixelinfo << "GainCurve_row_" << ipix->row() << "_col_" << ipix->col();
0316 std::string tempname = translateDetIdToString(detid);
0317 tempname += "_";
0318 tempname += pixelinfo.str();
0319 setDQMDirectory(detid);
0320 bookkeeper_pixels_[detid][pixelinfo.str()] = bookDQMHistogram1D(
0321 detid, pixelinfo.str(), tempname, 105 * nallpoints, xvalsall[0], xvalsall[nallpoints - 1] * 1.05);
0322 for (int ii = 0; ii < nallpoints; ++ii)
0323 bookkeeper_pixels_[detid][pixelinfo.str()]->Fill(xvalsall[ii], yvalsall[ii]);
0324 return false;
0325 }
0326
0327
0328 graph_->Set(npoints);
0329
0330 func_->SetParameter(0, 50.);
0331 func_->SetParameter(1, 0.25);
0332 for (int ipointtemp = 0; ipointtemp < npoints; ++ipointtemp) {
0333 graph_->SetPoint(ipointtemp, xvals[ipointtemp], yvals[ipointtemp]);
0334 graph_->SetPointError(ipointtemp, 0, yerrvals[ipointtemp]);
0335 }
0336 Int_t tempresult = graph_->Fit(func_, "FQ0N");
0337 slope = func_->GetParameter(1);
0338 slopeerror = func_->GetParError(1);
0339 intercept = func_->GetParameter(0);
0340 intercepterror = func_->GetParError(0);
0341 chi2 = func_->GetChisquare() / ((float)npoints - func_->GetNpar());
0342 prob = TMath::Prob(func_->GetChisquare(), npoints - func_->GetNpar());
0343 size_t ntimes = 0;
0344 while ((edm::isNotFinite(slope) || edm::isNotFinite(intercept)) && ntimes < 10) {
0345 ntimes++;
0346 makehistopersistent = true;
0347
0348 edm::LogWarning("SiPixelGainCalibrationAnalysis") << "impossible to fit values, try " << ntimes << ": ";
0349 for (int ii = 0; ii < npoints; ++ii) {
0350 edm::LogWarning("SiPixelGainCalibrationAnalysis")
0351 << "vcal " << xvals[ii] << " response: " << yvals[ii] << "+/-" << yerrvals[ii] << std::endl;
0352 }
0353 tempresult = graph_->Fit(func_, "FQ0NW");
0354 slope = func_->GetParameter(1);
0355 slopeerror = func_->GetParError(1);
0356 intercept = func_->GetParameter(0);
0357 intercepterror = func_->GetParError(0);
0358 chi2 = func_->GetChisquare() / ((float)npoints - func_->GetNpar());
0359 prob = TMath::Prob(func_->GetChisquare(), npoints - func_->GetNpar());
0360 }
0361
0362 if (tempresult == 0)
0363 status = 1;
0364 else
0365 status = 0;
0366 if (slope != 0)
0367 slope = 1. / slope;
0368 if (edm::isNotFinite(slope) || edm::isNotFinite(intercept)) {
0369 status = -6;
0370 bookkeeper_[detid]["status_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, status);
0371 if (writeSummary_) {
0372 summary_ << "row_" << ipix->row() << " col_" << ipix->col() << " status_" << status << endl;
0373 statusNumbers_[abs(status)]++;
0374 }
0375
0376 }
0377 if (chi2 > chi2Threshold_ && chi2Threshold_ >= 0)
0378 status = 5;
0379 if (prob < chi2ProbThreshold_)
0380 status = 5;
0381 if (noPlateau)
0382 status = 3;
0383 if (nallpoints < 4)
0384 status = -7;
0385 if (TMath::Abs(slope > maxGainInHist_) || slope < 0)
0386 status = -8;
0387 if (status != 1)
0388 makehistopersistent = true;
0389 statusNumbers_[abs(status)]++;
0390
0391 if (slope < gainlow_)
0392 gainlow_ = slope;
0393 if (slope > gainhi_)
0394 gainhi_ = slope;
0395 if (intercept > pedhi_)
0396 pedhi_ = intercept;
0397 if (intercept < pedlow_)
0398 pedlow_ = intercept;
0399 bookkeeper_[detid]["gain_1d"]->Fill(slope);
0400 if (slope > maxGainInHist_) {
0401 makehistopersistent = true;
0402 edm::LogWarning("SiPixelGainCalibration")
0403 << "For DETID " << detid << "pixel row,col " << ipix->row() << "," << ipix->col() << " Gain was measured to be "
0404 << slope << " which is outside the range of the summary plot (" << maxGainInHist_ << ") !!!! " << std::endl;
0405 }
0406 bookkeeper_[detid]["dynamicrange_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, highmeanval - lowmeanval);
0407 bookkeeper_[detid]["plateau_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, highmeanval);
0408 bookkeeper_[detid]["gain_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, slope);
0409 bookkeeper_[detid]["errorgain_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, slopeerror);
0410 bookkeeper_[detid]["ped_1d"]->Fill(intercept);
0411 bookkeeper_[detid]["ped_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, intercept);
0412 bookkeeper_[detid]["errorped_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, intercepterror);
0413 bookkeeper_[detid]["chi2_1d"]->Fill(chi2);
0414 bookkeeper_[detid]["chi2_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, chi2);
0415 bookkeeper_[detid]["prob_1d"]->Fill(prob);
0416 bookkeeper_[detid]["prob_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, prob);
0417 bookkeeper_[detid]["lowpoint_1d"]->Fill(xvals[0]);
0418 bookkeeper_[detid]["lowpoint_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, xvals[0]);
0419 bookkeeper_[detid]["highpoint_1d"]->Fill(xvals[npoints - 1]);
0420 bookkeeper_[detid]["highpoint_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, xvals[npoints - 1]);
0421 bookkeeper_[detid]["nfitpoints_1d"]->Fill(npoints);
0422 bookkeeper_[detid]["endpoint_1d"]->Fill((255 - intercept) * slope);
0423 bookkeeper_[detid]["status_2d"]->setBinContent(ipix->col() + 1, ipix->row() + 1, status);
0424
0425 if (!savePixelHists_)
0426 return true;
0427 if (detidfinder == listofdetids_.end() && !listofdetids_.empty())
0428 return true;
0429 if (makehistopersistent) {
0430 std::ostringstream pixelinfo;
0431 pixelinfo << "GainCurve_row_" << ipix->row() << "_col_" << ipix->col();
0432 std::string tempname = translateDetIdToString(detid);
0433 tempname += "_";
0434 tempname += pixelinfo.str();
0435
0436
0437
0438
0439 setDQMDirectory(detid);
0440 bookkeeper_pixels_[detid][pixelinfo.str()] = bookDQMHistogram1D(
0441 detid, pixelinfo.str(), tempname, 105 * nallpoints, xvalsall[0], xvalsall[nallpoints - 1] * 1.05);
0442
0443 edm::LogInfo("SiPixelGainCalibrationAnalysis")
0444 << "now saving histogram for pixel " << tempname << ", gain = " << slope << ", pedestal = " << intercept
0445 << ", chi2/NDF=" << chi2 << "(prob:" << prob << "), fit status " << status;
0446 for (int ii = 0; ii < nallpoints; ++ii) {
0447
0448 bookkeeper_pixels_[detid][pixelinfo.str()]->Fill(xvalsall[ii], yvalsall[ii]);
0449 }
0450
0451
0452
0453 if (writeSummary_) {
0454 summary_ << "row_" << ipix->row() << " col_" << ipix->col();
0455 summary_ << " status_" << status;
0456 summary_ << endl;
0457
0458
0459 }
0460 }
0461 return true;
0462 }
0463
0464 void SiPixelGainCalibrationAnalysis::newDetID(uint32_t detid) {
0465 setDQMDirectory(detid);
0466 std::string tempname = translateDetIdToString(detid);
0467 bookkeeper_[detid]["gain_1d"] = bookDQMHistogram1D(detid, "Gain1d", "gain for " + tempname, 100, 0., maxGainInHist_);
0468 bookkeeper_[detid]["gain_2d"] = bookDQMHistoPlaquetteSummary2D(detid, "Gain2d", "gain for " + tempname);
0469 bookkeeper_[detid]["errorgain_2d"] =
0470 bookDQMHistoPlaquetteSummary2D(detid, "ErrorGain2d", "error on gain for " + tempname);
0471 bookkeeper_[detid]["ped_1d"] = bookDQMHistogram1D(detid, "Pedestal1d", "pedestal for " + tempname, 256, 0., 256.0);
0472 bookkeeper_[detid]["ped_2d"] = bookDQMHistoPlaquetteSummary2D(detid, "Pedestal2d", "pedestal for " + tempname);
0473 bookkeeper_[detid]["errorped_2d"] =
0474 bookDQMHistoPlaquetteSummary2D(detid, "ErrorPedestal2d", "error on pedestal for " + tempname);
0475 bookkeeper_[detid]["chi2_1d"] =
0476 bookDQMHistogram1D(detid, "GainChi2NDF1d", "#chi^{2}/NDOF for " + tempname, 100, 0., maxChi2InHist_);
0477 bookkeeper_[detid]["chi2_2d"] =
0478 bookDQMHistoPlaquetteSummary2D(detid, "GainChi2NDF2d", "#chi^{2}/NDOF for " + tempname);
0479 bookkeeper_[detid]["prob_1d"] =
0480 bookDQMHistogram1D(detid, "GainChi2Prob1d", "P(#chi^{2},NDOF) for " + tempname, 100, 0., 1.0);
0481 bookkeeper_[detid]["prob_2d"] =
0482 bookDQMHistoPlaquetteSummary2D(detid, "GainChi2Prob2d", "P(#chi^{2},NDOF) for " + tempname);
0483 bookkeeper_[detid]["status_2d"] =
0484 bookDQMHistoPlaquetteSummary2D(detid, "GainFitResult2d", "Fit result for " + tempname);
0485 bookkeeper_[detid]["endpoint_1d"] = bookDQMHistogram1D(
0486 detid, "GainEndPoint1d", "point where fit meets ADC=255 for " + tempname, 256, 0., 256. * scalarVcalHigh_VcalLow_);
0487 bookkeeper_[detid]["lowpoint_1d"] = bookDQMHistogram1D(
0488 detid, "GainLowPoint1d", "lowest fit point for " + tempname, 256, 0., 256. * scalarVcalHigh_VcalLow_);
0489 bookkeeper_[detid]["highpoint_1d"] = bookDQMHistogram1D(
0490 detid, "GainHighPoint1d", "highest fit point for " + tempname, 256, 0., 256. * scalarVcalHigh_VcalLow_);
0491 bookkeeper_[detid]["nfitpoints_1d"] =
0492 bookDQMHistogram1D(detid, "GainNPoints1d", "number of fit point for " + tempname, 20, 0., 20);
0493 bookkeeper_[detid]["dynamicrange_2d"] = bookDQMHistoPlaquetteSummary2D(
0494 detid, "GainDynamicRange2d", "Difference lowest and highest points on gain curve for " + tempname);
0495 bookkeeper_[detid]["lowpoint_2d"] =
0496 bookDQMHistoPlaquetteSummary2D(detid, "GainLowPoint2d", "lowest fit point for " + tempname);
0497 bookkeeper_[detid]["highpoint_2d"] =
0498 bookDQMHistoPlaquetteSummary2D(detid, "GainHighPoint2d", "highest fit point for " + tempname);
0499 bookkeeper_[detid]["plateau_2d"] =
0500 bookDQMHistoPlaquetteSummary2D(detid, "GainSaturate2d", "Highest points on gain curve for " + tempname);
0501 }
0502
0503 DEFINE_FWK_MODULE(SiPixelGainCalibrationAnalysis);