File indexing completed on 2024-12-20 03:13:39
0001 #include "DQM/SiStripCommissioningAnalysis/interface/CalibrationScanAlgorithm.h"
0002 #include "CondFormats/SiStripObjects/interface/CalibrationScanAnalysis.h"
0003 #include "DataFormats/SiStripCommon/interface/SiStripHistoTitle.h"
0004 #include "DataFormats/SiStripCommon/interface/SiStripEnumsAndStrings.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include "DQM/SiStripCommissioningAnalysis/interface/SiStripPulseShape.h"
0007 #include "TProfile.h"
0008 #include "TF1.h"
0009 #include "TH1.h"
0010 #include "TVirtualFitter.h"
0011 #include "TFitResult.h"
0012 #include "TMath.h"
0013 #include "TGraph2D.h"
0014 #include "TH2F.h"
0015 #include "TCanvas.h"
0016 #include "TROOT.h"
0017 #include <iostream>
0018 #include <sstream>
0019 #include <iomanip>
0020 #include <cmath>
0021 #include "Math/MinimizerOptions.h"
0022
0023 using namespace sistrip;
0024
0025
0026
0027 CalibrationScanAlgorithm::CalibrationScanAlgorithm(const edm::ParameterSet& pset, CalibrationScanAnalysis* const anal)
0028 : CommissioningAlgorithm(anal), cal_(nullptr) {}
0029
0030
0031
0032 void CalibrationScanAlgorithm::extract(const std::vector<TH1*>& histos) {
0033 if (!anal()) {
0034 edm::LogWarning(mlCommissioning_) << "[CalibrationScanAlgorithm::" << __func__ << "]"
0035 << " NULL pointer to base Analysis object!";
0036 return;
0037 }
0038
0039 CommissioningAnalysis* tmp = const_cast<CommissioningAnalysis*>(anal());
0040 cal_ = dynamic_cast<CalibrationScanAnalysis*>(tmp);
0041 if (!cal_) {
0042 edm::LogWarning(mlCommissioning_) << "[CalibrationScanAlgorithm::" << __func__ << "]"
0043 << " NULL pointer to derived Analysis object!";
0044 return;
0045 }
0046
0047
0048 if (!histos.empty()) {
0049 cal_->fedKey(extractFedKey(histos.front()));
0050 }
0051
0052
0053 std::vector<TH1*>::const_iterator ihis = histos.begin();
0054 unsigned int cnt = 0;
0055 for (; ihis != histos.end(); ihis++, cnt++) {
0056
0057 if (!(*ihis)) {
0058 continue;
0059 }
0060
0061
0062 SiStripHistoTitle title((*ihis)->GetName());
0063 if (title.runType() != sistrip::CALIBRATION_SCAN && title.runType() != sistrip::CALIBRATION_SCAN_DECO) {
0064 cal_->addErrorCode(sistrip::unexpectedTask_);
0065 continue;
0066 }
0067
0068
0069 Histo histo_temp;
0070 histo_temp.first = *ihis;
0071 histo_temp.first->Sumw2();
0072 histo_temp.second = (*ihis)->GetTitle();
0073 histo_[title.extraInfo()].resize(2);
0074 if (title.channel() % 2 == 0)
0075 histo_[title.extraInfo()][0] = histo_temp;
0076 else
0077 histo_[title.extraInfo()][1] = histo_temp;
0078 }
0079 }
0080
0081
0082 void CalibrationScanAlgorithm::analyse() {
0083 ROOT::Math::MinimizerOptions::SetDefaultMinimizer("Minuit2", "Migrad");
0084 ROOT::Math::MinimizerOptions::SetDefaultStrategy(0);
0085
0086 if (!cal_) {
0087 edm::LogWarning(mlCommissioning_) << "[CalibrationScanAlgorithm::" << __func__ << "]"
0088 << " NULL pointer to derived Analysis object!";
0089 return;
0090 }
0091
0092
0093 TFitResultPtr fit_result;
0094 TF1* fit_function_turnOn = nullptr;
0095 TF1* fit_function_decay = nullptr;
0096 TF1* fit_function_deco = nullptr;
0097 if (cal_->deconv_) {
0098 fit_function_deco = new TF1("fit_function_deco", fdeconv, 0, 400, 7);
0099 fit_function_deco->SetParameters(4, 25, 25, 50, 250, 25, 0.75);
0100 } else {
0101 fit_function_turnOn = new TF1("fit_function_turnOn", fturnOn, 0, 400, 4);
0102 fit_function_decay = new TF1("fit_function_decay", fdecay, 0, 400, 3);
0103 fit_function_turnOn->SetParameters(50, 50, 40, 20);
0104 fit_function_decay->SetParameters(-150, -0.01, -0.1);
0105 }
0106
0107
0108 for (auto map_element : histo_) {
0109
0110 cal_->addOneCalibrationPoint(map_element.first);
0111
0112
0113 std::vector<std::string> tokens;
0114 std::string token;
0115 std::istringstream tokenStream(map_element.first);
0116 while (std::getline(tokenStream, token, '_')) {
0117 tokens.push_back(token);
0118 }
0119
0120 scanned_isha_.push_back(std::stoi(tokens.at(1)));
0121 scanned_vfs_.push_back(std::stoi(tokens.at(3)));
0122
0123
0124 for (size_t iapv = 0; iapv < 2; iapv++) {
0125 if (!map_element.second[iapv].first) {
0126 edm::LogWarning(mlCommissioning_)
0127 << " NULL pointer to histogram for: " << map_element.second[iapv].second << " !";
0128 return;
0129 }
0130
0131 cal_->amplitude_[map_element.first][iapv] = 0;
0132 cal_->baseline_[map_element.first][iapv] = 0;
0133 cal_->riseTime_[map_element.first][iapv] = 0;
0134 cal_->turnOn_[map_element.first][iapv] = 0;
0135 cal_->peakTime_[map_element.first][iapv] = 0;
0136 cal_->undershoot_[map_element.first][iapv] = 0;
0137 cal_->tail_[map_element.first][iapv] = 0;
0138 cal_->decayTime_[map_element.first][iapv] = 0;
0139 cal_->smearing_[map_element.first][iapv] = 0;
0140 cal_->chi2_[map_element.first][iapv] = 0;
0141 cal_->isvalid_[map_element.first][iapv] = true;
0142
0143 if (map_element.second[iapv].first->Integral() == 0) {
0144 cal_->isvalid_[map_element.first][iapv] = false;
0145 continue;
0146 }
0147
0148
0149 correctDistribution(map_element.second[iapv].first, false);
0150
0151
0152 float error = (map_element.second[iapv].first->GetMaximum() * 0.05);
0153 for (int i = 1; i <= map_element.second[iapv].first->GetNbinsX(); ++i)
0154 map_element.second[iapv].first->SetBinError(i, error);
0155
0156
0157 if (cal_->deconv_) {
0158 fit_function_deco->SetParameters(4, 25, 25, 50, 250, 25, 0.75);
0159 fit_result = map_element.second[iapv].first->Fit(fit_function_deco, "QRS");
0160
0161 if (not fit_result.Get()) {
0162 cal_->isvalid_[map_element.first][iapv] = false;
0163 continue;
0164 }
0165
0166
0167 float maximum_ampl = fit_function_deco->GetMaximum();
0168 float peak_time = fit_function_deco->GetMaximumX();
0169 float baseline = baseLine(fit_function_deco);
0170 float turn_on_time = turnOn(fit_function_deco, baseline);
0171 float rise_time = peak_time - turn_on_time;
0172
0173
0174 cal_->amplitude_[map_element.first][iapv] = maximum_ampl - baseline;
0175 cal_->baseline_[map_element.first][iapv] = baseline;
0176 cal_->riseTime_[map_element.first][iapv] = rise_time;
0177 cal_->turnOn_[map_element.first][iapv] = turn_on_time;
0178 cal_->peakTime_[map_element.first][iapv] = peak_time;
0179 if (fit_function_deco->GetMinimumX() > rise_time)
0180 cal_->undershoot_[map_element.first][iapv] =
0181 100 * (fit_function_deco->GetMinimum() - baseline) / (maximum_ampl - baseline);
0182 else
0183 cal_->undershoot_[map_element.first][iapv] = 0;
0184
0185
0186 int lastBin = map_element.second[iapv].first->FindBin(peak_time + 125);
0187 if (lastBin > map_element.second[iapv].first->GetNbinsX() - 4)
0188 lastBin = map_element.second[iapv].first->GetNbinsX() - 4;
0189
0190
0191 cal_->tail_[map_element.first][iapv] =
0192 100 * (map_element.second[iapv].first->GetBinContent(lastBin) - baseline) / (maximum_ampl - baseline);
0193
0194
0195 cal_->decayTime_[map_element.first][iapv] = decayTime(fit_function_deco) - peak_time;
0196 cal_->smearing_[map_element.first][iapv] = 0;
0197 cal_->chi2_[map_element.first][iapv] =
0198 fit_function_deco->GetChisquare() /
0199 (map_element.second[iapv].first->GetNbinsX() - fit_function_deco->GetNpar());
0200
0201 } else {
0202
0203 fit_function_turnOn->SetParameters(50, 50, 40, 20);
0204 fit_function_turnOn->SetRange(fit_function_turnOn->GetXmin(),
0205 map_element.second[iapv].first->GetBinCenter(
0206 map_element.second[iapv].first->GetMaximumBin()));
0207 fit_result = map_element.second[iapv].first->Fit(fit_function_turnOn, "QSR");
0208 if (not fit_result.Get()) {
0209 cal_->isvalid_[map_element.first][iapv] = false;
0210 continue;
0211 }
0212
0213
0214 float maximum_ampl = fit_function_turnOn->GetMaximum();
0215 float peak_time = fit_function_turnOn->GetMaximumX();
0216 float baseline = baseLine(fit_function_turnOn);
0217 float turn_on_time = turnOn(fit_function_turnOn, baseline);
0218 float rise_time = peak_time - turn_on_time;
0219
0220
0221 cal_->amplitude_[map_element.first][iapv] = maximum_ampl - baseline;
0222 cal_->baseline_[map_element.first][iapv] = baseline;
0223 cal_->riseTime_[map_element.first][iapv] = rise_time;
0224 cal_->turnOn_[map_element.first][iapv] = turn_on_time;
0225 cal_->peakTime_[map_element.first][iapv] = peak_time;
0226
0227 fit_function_decay->SetParameters(-150, -0.01, -0.1);
0228 fit_function_decay->SetRange(
0229 map_element.second[iapv].first->GetBinCenter(map_element.second[iapv].first->GetMaximumBin()) + 10.,
0230 fit_function_decay->GetXmax());
0231 fit_result = map_element.second[iapv].first->Fit(fit_function_decay, "QSR+");
0232
0233 if (fit_result.Get() and fit_result->Status() >= 4) {
0234 cal_->isvalid_[map_element.first][iapv] = false;
0235 continue;
0236 }
0237
0238 cal_->undershoot_[map_element.first][iapv] = 0;
0239
0240
0241 int lastBin = map_element.second[iapv].first->FindBin(peak_time + 125);
0242 if (lastBin > map_element.second[iapv].first->GetNbinsX() - 4)
0243 lastBin = map_element.second[iapv].first->GetNbinsX() - 4;
0244
0245
0246 cal_->tail_[map_element.first][iapv] =
0247 100 * (map_element.second[iapv].first->GetBinContent(lastBin) - baseline) / (maximum_ampl - baseline);
0248
0249
0250 cal_->decayTime_[map_element.first][iapv] = decayTime(fit_function_decay) - peak_time;
0251 cal_->smearing_[map_element.first][iapv] = 0;
0252 cal_->chi2_[map_element.first][iapv] =
0253 (fit_function_turnOn->GetChisquare() + fit_function_decay->GetChisquare()) /
0254 (map_element.second[iapv].first->GetNbinsX() - fit_function_turnOn->GetNpar() -
0255 fit_function_decay->GetNpar());
0256
0257
0258 bool isvalid = true;
0259 if (cal_->amplitude_[map_element.first][iapv] < CalibrationScanAnalysis::minAmplitudeThreshold_)
0260 isvalid = false;
0261 if (cal_->baseline_[map_element.first][iapv] < CalibrationScanAnalysis::minBaselineThreshold_)
0262 isvalid = false;
0263 else if (cal_->baseline_[map_element.first][iapv] > CalibrationScanAnalysis::maxBaselineThreshold_)
0264 isvalid = false;
0265 if (cal_->decayTime_[map_element.first][iapv] < CalibrationScanAnalysis::minDecayTimeThreshold_)
0266 isvalid = false;
0267 else if (cal_->decayTime_[map_element.first][iapv] > CalibrationScanAnalysis::maxDecayTimeThreshold_)
0268 isvalid = false;
0269 if (cal_->peakTime_[map_element.first][iapv] < CalibrationScanAnalysis::minPeakTimeThreshold_)
0270 isvalid = false;
0271 else if (cal_->peakTime_[map_element.first][iapv] > CalibrationScanAnalysis::maxPeakTimeThreshold_)
0272 isvalid = false;
0273 if (cal_->riseTime_[map_element.first][iapv] < CalibrationScanAnalysis::minRiseTimeThreshold_)
0274 isvalid = false;
0275 else if (cal_->riseTime_[map_element.first][iapv] > CalibrationScanAnalysis::maxRiseTimeThreshold_)
0276 isvalid = false;
0277 if (cal_->turnOn_[map_element.first][iapv] < CalibrationScanAnalysis::minTurnOnThreshold_)
0278 isvalid = false;
0279 else if (cal_->turnOn_[map_element.first][iapv] > CalibrationScanAnalysis::maxTurnOnThreshold_)
0280 isvalid = false;
0281 if (cal_->chi2_[map_element.first][iapv] > CalibrationScanAnalysis::maxChi2Threshold_)
0282 isvalid = false;
0283
0284 if (not isvalid) {
0285 cal_->amplitude_[map_element.first][iapv] = 0;
0286 cal_->baseline_[map_element.first][iapv] = 0;
0287 cal_->riseTime_[map_element.first][iapv] = 0;
0288 cal_->turnOn_[map_element.first][iapv] = 0;
0289 cal_->peakTime_[map_element.first][iapv] = 0;
0290 cal_->undershoot_[map_element.first][iapv] = 0;
0291 cal_->tail_[map_element.first][iapv] = 0;
0292 cal_->decayTime_[map_element.first][iapv] = 0;
0293 cal_->smearing_[map_element.first][iapv] = 0;
0294 cal_->chi2_[map_element.first][iapv] = 0;
0295 cal_->isvalid_[map_element.first][iapv] = false;
0296 }
0297 }
0298 }
0299 }
0300
0301 if (fit_function_deco)
0302 delete fit_function_deco;
0303 if (fit_function_decay)
0304 delete fit_function_decay;
0305 if (fit_function_turnOn)
0306 delete fit_function_turnOn;
0307 }
0308
0309
0310 void CalibrationScanAlgorithm::correctDistribution(TH1* histo, const bool& isShape) const {
0311
0312
0313 if (not isShape) {
0314 for (int iBin = 0; iBin < histo->GetNbinsX(); iBin++) {
0315 histo->SetBinContent(iBin + 1, -histo->GetBinContent(iBin + 1) / 16.);
0316 histo->SetBinContent(iBin + 1, histo->GetBinContent(iBin + 1) / 5.);
0317 }
0318 } else
0319 histo->Scale(1. / histo->Integral());
0320 }
0321
0322
0323 float CalibrationScanAlgorithm::baseLine(TF1* f) {
0324 float x = f->GetXmin();
0325 float xmax = 10;
0326 float baseline = 0;
0327 int npoints = 0;
0328 for (; x < xmax; x += 0.1) {
0329 baseline += f->Eval(x);
0330 npoints++;
0331 }
0332 return baseline / npoints;
0333 }
0334
0335
0336 float CalibrationScanAlgorithm::turnOn(
0337 TF1* f, const float& baseline) {
0338 float max_amplitude = f->GetMaximum();
0339 float time = 10.;
0340 for (; time < 100 && (f->Eval(time) - baseline) < 0.05 * (max_amplitude - baseline); time += 0.1) {
0341 }
0342 return time;
0343 }
0344
0345
0346 float CalibrationScanAlgorithm::decayTime(
0347 TF1* f) {
0348 float xval = std::max(f->GetXmin(), f->GetMaximumX());
0349 float max_amplitude = f->GetMaximum();
0350 float x = xval;
0351 for (; x < 1000;
0352 x = x +
0353 0.1) {
0354 if (f->Eval(x) < max_amplitude * exp(-1))
0355 break;
0356 }
0357 return x;
0358 }
0359
0360
0361 void CalibrationScanAlgorithm::tuneIndependently(const int& iapv,
0362 const float& targetRiseTime,
0363 const float& targetDecayTime) {
0364 std::map<int, std::vector<float> > decayTime_vs_vfs;
0365 TString name;
0366 int imap = 0;
0367
0368 for (auto map_element : histo_) {
0369
0370 if (scanned_isha_.at(imap) <= CalibrationScanAnalysis::minISHAforVFSTune_ or
0371 scanned_isha_.at(imap) >= CalibrationScanAnalysis::maxISHAforVFSTune_) {
0372 imap++;
0373 continue;
0374 }
0375
0376 if (cal_->isValid(map_element.first)[iapv])
0377 decayTime_vs_vfs[scanned_vfs_.at(imap)].push_back(cal_->decayTime(map_element.first)[iapv]);
0378
0379 if (name == "") {
0380 name = Form("%s", map_element.second[iapv].first->GetName());
0381 name.ReplaceAll("_" + map_element.first, "");
0382 }
0383 imap++;
0384 }
0385
0386
0387 for (auto iter : decayTime_vs_vfs)
0388 sort(iter.second.begin(), iter.second.end());
0389
0390 name.ReplaceAll("ExpertHisto_", "");
0391
0392
0393 cal_->decayTime_vs_vfs_.push_back(new TGraph());
0394 cal_->decayTime_vs_vfs_.back()->SetName(Form("decayTime_%s", name.Data()));
0395
0396
0397 cal_->riseTime_vs_isha_.push_back(new TGraph());
0398 cal_->riseTime_vs_isha_.back()->SetName(Form("riseTime_%s", name.Data()));
0399
0400 if (!decayTime_vs_vfs.empty()) {
0401 int ipoint = 0;
0402 for (auto map_element : decayTime_vs_vfs) {
0403 if (!map_element.second.empty()) {
0404 cal_->decayTime_vs_vfs_.at(iapv)->SetPoint(
0405 ipoint, map_element.second.at(round(map_element.second.size() / 2)), map_element.first);
0406 ipoint++;
0407 }
0408 }
0409
0410 double max_apv =
0411 TMath::MaxElement(cal_->decayTime_vs_vfs_.at(iapv)->GetN(), cal_->decayTime_vs_vfs_.at(iapv)->GetY());
0412 double min_apv =
0413 TMath::MinElement(cal_->decayTime_vs_vfs_.at(iapv)->GetN(), cal_->decayTime_vs_vfs_.at(iapv)->GetY());
0414
0415 cal_->vfs_[iapv] = cal_->decayTime_vs_vfs_.at(iapv)->Eval(targetDecayTime);
0416
0417
0418 if (cal_->vfs_[iapv] < min_apv)
0419 cal_->vfs_[iapv] = min_apv;
0420 else if (cal_->vfs_[iapv] > max_apv)
0421 cal_->vfs_[iapv] = max_apv;
0422
0423
0424 std::map<int, std::vector<float> > riseTime_vs_isha;
0425 imap = 0;
0426
0427 for (auto map_element : histo_) {
0428 if (std::abs(scanned_vfs_.at(imap) - cal_->vfs_[iapv]) < CalibrationScanAnalysis::VFSrange_ and
0429 cal_->isValid(map_element.first)[iapv])
0430 riseTime_vs_isha[scanned_isha_.at(imap)].push_back(cal_->riseTime(map_element.first)[iapv]);
0431 if (name == "") {
0432 name = Form("%s", map_element.second[iapv].first->GetName());
0433 name.ReplaceAll("_" + map_element.first, "");
0434 }
0435 imap++;
0436 }
0437
0438
0439 for (auto iter : riseTime_vs_isha)
0440 sort(iter.second.begin(), iter.second.end());
0441 name.ReplaceAll("ExpertHisto_", "");
0442
0443
0444 if (!riseTime_vs_isha.empty()) {
0445 int ipoint = 0;
0446 for (auto map_element : riseTime_vs_isha) {
0447 if (!map_element.second.empty()) {
0448 cal_->riseTime_vs_isha_.at(iapv)->SetPoint(
0449 ipoint, map_element.second.at(round(map_element.second.size() / 2)), map_element.first);
0450 ipoint++;
0451 }
0452 }
0453
0454 double max_apv =
0455 TMath::MaxElement(cal_->riseTime_vs_isha_.at(iapv)->GetN(), cal_->riseTime_vs_isha_.at(iapv)->GetY());
0456 double min_apv =
0457 TMath::MinElement(cal_->riseTime_vs_isha_.at(iapv)->GetN(), cal_->riseTime_vs_isha_.at(iapv)->GetY());
0458
0459 cal_->isha_[iapv] = cal_->riseTime_vs_isha_.at(iapv)->Eval(targetRiseTime);
0460
0461 if (cal_->isha_[iapv] < min_apv)
0462 cal_->isha_[iapv] = min_apv;
0463 else if (cal_->isha_[iapv] > max_apv)
0464 cal_->isha_[iapv] = max_apv;
0465 } else
0466 cal_->isha_[iapv] = -1;
0467 }
0468 }
0469
0470
0471 void CalibrationScanAlgorithm::tuneSimultaneously(const int& iapv,
0472 const float& targetRiseTime,
0473 const float& targetDecayTime) {
0474
0475 cal_->decayTime_vs_isha_vfs_.push_back(new TGraph2D());
0476 cal_->riseTime_vs_isha_vfs_.push_back(new TGraph2D());
0477
0478
0479 TString name_apv;
0480 int ipoint_apv = 0;
0481 int imap = 0;
0482
0483 for (auto map_element : histo_) {
0484 if (cal_->isValid(map_element.first)[iapv]) {
0485 cal_->decayTime_vs_isha_vfs_.at(iapv)->SetPoint(
0486 ipoint_apv, scanned_isha_.at(imap), scanned_vfs_.at(imap), cal_->decayTime(map_element.first)[iapv]);
0487 cal_->riseTime_vs_isha_vfs_.at(iapv)->SetPoint(
0488 ipoint_apv, scanned_isha_.at(imap), scanned_vfs_.at(imap), cal_->riseTime(map_element.first)[iapv]);
0489 ipoint_apv++;
0490 }
0491 if (name_apv == "") {
0492 name_apv = Form("%s", map_element.second[iapv].first->GetName());
0493 name_apv.ReplaceAll("_" + map_element.first, "");
0494 }
0495 imap++;
0496 }
0497
0498 name_apv.ReplaceAll("ExpertHisto_", "");
0499
0500 cal_->decayTime_vs_isha_vfs_.at(iapv)->SetName(Form("decayTime_%s", name_apv.Data()));
0501 cal_->riseTime_vs_isha_vfs_.at(iapv)->SetName(Form("riseTime_%s", name_apv.Data()));
0502
0503
0504 TH2F* hist_decay_apv = new TH2F("hist_decay_apv",
0505 "hist_decay_apv",
0506 500,
0507 *min_element(scanned_isha_.begin(), scanned_isha_.end()),
0508 *max_element(scanned_isha_.begin(), scanned_isha_.end()),
0509 500,
0510 *min_element(scanned_vfs_.begin(), scanned_vfs_.end()),
0511 *max_element(scanned_vfs_.begin(), scanned_vfs_.end()));
0512
0513 TH2F* hist_rise_apv = (TH2F*)hist_decay_apv->Clone();
0514 hist_rise_apv->SetName("hist_rise_apv");
0515 hist_rise_apv->Reset();
0516
0517 TH2F* hist_distance = (TH2F*)hist_decay_apv->Clone();
0518 hist_distance->SetName("hist_distance");
0519 hist_distance->Reset();
0520
0521 for (int iBin = 1; iBin <= hist_decay_apv->GetNbinsX(); iBin++) {
0522 for (int jBin = 1; jBin <= hist_decay_apv->GetNbinsY(); jBin++) {
0523 if (ipoint_apv != 0) {
0524 if (cal_->decayTime_vs_isha_vfs_.at(iapv)->GetN() > 10)
0525 hist_decay_apv->SetBinContent(
0526 iBin,
0527 jBin,
0528 cal_->decayTime_vs_isha_vfs_.at(iapv)->Interpolate(hist_decay_apv->GetXaxis()->GetBinCenter(iBin),
0529 hist_decay_apv->GetYaxis()->GetBinCenter(jBin)));
0530 if (cal_->riseTime_vs_isha_vfs_.at(iapv)->GetN() > 10)
0531 hist_rise_apv->SetBinContent(
0532 iBin,
0533 jBin,
0534 cal_->riseTime_vs_isha_vfs_.at(iapv)->Interpolate(hist_rise_apv->GetXaxis()->GetBinCenter(iBin),
0535 hist_rise_apv->GetYaxis()->GetBinCenter(jBin)));
0536 }
0537 }
0538 }
0539
0540
0541 hist_decay_apv->Smooth();
0542 hist_rise_apv->Smooth();
0543
0544 for (int iBin = 1; iBin <= hist_decay_apv->GetNbinsX(); iBin++) {
0545 for (int jBin = 1; jBin <= hist_decay_apv->GetNbinsY(); jBin++) {
0546 hist_distance->SetBinContent(
0547 iBin,
0548 jBin,
0549 sqrt(pow((hist_decay_apv->GetBinContent(iBin, jBin) - targetDecayTime) / targetDecayTime, 2) +
0550 pow((hist_rise_apv->GetBinContent(iBin, jBin) - targetRiseTime) / targetRiseTime, 2)));
0551 }
0552 }
0553
0554 int minx, miny, minz;
0555 hist_distance->GetMinimumBin(minx, miny, minz);
0556
0557 cal_->isha_[iapv] = round(hist_distance->GetXaxis()->GetBinCenter(minx));
0558 cal_->vfs_[iapv] = round(hist_distance->GetYaxis()->GetBinCenter(miny));
0559
0560 delete hist_decay_apv;
0561 delete hist_rise_apv;
0562 delete hist_distance;
0563 }
0564
0565 void CalibrationScanAlgorithm::fillTunedObservables(const int& apvid) {
0566
0567 int distance_apv = 10000;
0568
0569
0570 for (size_t i = 0; i < scanned_isha_.size(); i++) {
0571 if (std::abs(scanned_isha_.at(i) - cal_->bestISHA().at(apvid)) < distance_apv) {
0572 distance_apv = std::abs(scanned_isha_.at(i) - cal_->bestISHA().at(apvid));
0573 cal_->tunedISHA_.at(apvid) = scanned_isha_.at(i);
0574 }
0575 }
0576
0577 distance_apv = 10000;
0578
0579
0580 for (size_t i = 0; i < scanned_vfs_.size(); i++) {
0581 if (std::abs(scanned_vfs_.at(i) - cal_->bestVFS().at(apvid)) < distance_apv) {
0582 distance_apv = std::abs(scanned_vfs_.at(i) - cal_->bestVFS().at(apvid));
0583 cal_->tunedVFS_.at(apvid) = scanned_vfs_.at(i);
0584 }
0585 }
0586
0587
0588 std::string key_apv = std::string(Form("isha_%d_vfs_%d", cal_->tunedISHA().at(apvid), cal_->tunedVFS().at(apvid)));
0589 if (!cal_->amplitude(key_apv).empty()) {
0590 cal_->tunedAmplitude_[apvid] = cal_->amplitude(key_apv)[apvid];
0591 cal_->tunedTail_[apvid] = cal_->tail(key_apv)[apvid];
0592 cal_->tunedRiseTime_[apvid] = cal_->riseTime(key_apv)[apvid];
0593 cal_->tunedDecayTime_[apvid] = cal_->decayTime(key_apv)[apvid];
0594 cal_->tunedTurnOn_[apvid] = cal_->turnOn(key_apv)[apvid];
0595 cal_->tunedPeakTime_[apvid] = cal_->peakTime(key_apv)[apvid];
0596 cal_->tunedUndershoot_[apvid] = cal_->undershoot(key_apv)[apvid];
0597 cal_->tunedBaseline_[apvid] = cal_->baseline(key_apv)[apvid];
0598 cal_->tunedSmearing_[apvid] = cal_->smearing(key_apv)[apvid];
0599 cal_->tunedChi2_[apvid] = cal_->chi2(key_apv)[apvid];
0600 } else {
0601 cal_->tunedAmplitude_[apvid] = 0;
0602 cal_->tunedTail_[apvid] = 0;
0603 cal_->tunedRiseTime_[apvid] = 0;
0604 cal_->tunedDecayTime_[apvid] = 0;
0605 cal_->tunedTurnOn_[apvid] = 0;
0606 cal_->tunedPeakTime_[apvid] = 0;
0607 cal_->tunedUndershoot_[apvid] = 0;
0608 cal_->tunedBaseline_[apvid] = 0;
0609 cal_->tunedSmearing_[apvid] = 0;
0610 cal_->tunedChi2_[apvid] = 0;
0611 }
0612 }