File indexing completed on 2025-03-13 02:31:39
0001 #include "DQM/SiStripCommissioningAnalysis/interface/PedsFullNoiseAlgorithm.h"
0002 #include "CondFormats/SiStripObjects/interface/PedsFullNoiseAnalysis.h"
0003 #include "DataFormats/SiStripCommon/interface/SiStripHistoTitle.h"
0004 #include "DataFormats/SiStripCommon/interface/SiStripEnumsAndStrings.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006
0007 #include "TProfile.h"
0008 #include "TH1.h"
0009 #include "TH2.h"
0010 #include "TF1.h"
0011 #include "TFitResult.h"
0012 #include "TMath.h"
0013 #include "Math/DistFunc.h"
0014 #include "Math/ProbFuncMathCore.h"
0015 #include "Fit/BinData.h"
0016 #include "HFitInterface.h"
0017 #include "Math/GoFTest.h"
0018 #include <iostream>
0019 #include <iomanip>
0020 #include <cmath>
0021
0022 using namespace sistrip;
0023 using namespace std;
0024
0025
0026 PedsFullNoiseAlgorithm::PedsFullNoiseAlgorithm(const edm::ParameterSet& pset, PedsFullNoiseAnalysis* const anal)
0027 : CommissioningAlgorithm(anal),
0028 hPeds_(nullptr, ""),
0029 hNoise_(nullptr, ""),
0030 hNoise2D_(nullptr, ""),
0031 maxDriftResidualCut_(pset.getParameter<double>("MaxDriftResidualCut")),
0032 minStripNoiseCut_(pset.getParameter<double>("MinStripNoiseCut")),
0033 maxStripNoiseCut_(pset.getParameter<double>("MaxStripNoiseCut")),
0034 maxStripNoiseSignificanceCut_(pset.getParameter<double>("MaxStripNoiseSignificanceCut")),
0035 adProbabCut_(pset.getParameter<double>("AdProbabCut")),
0036 ksProbabCut_(pset.getParameter<double>("KsProbabCut")),
0037 generateRandomHisto_(pset.getParameter<bool>("GenerateRandomHisto")),
0038 jbProbabCut_(pset.getParameter<double>("JbProbabCut")),
0039 chi2ProbabCut_(pset.getParameter<double>("Chi2ProbabCut")),
0040 kurtosisCut_(pset.getParameter<double>("KurtosisCut")),
0041 integralTailCut_(pset.getParameter<double>("IntegralTailCut")),
0042 integralNsigma_(pset.getParameter<int>("IntegralNsigma")),
0043 ashmanDistance_(pset.getParameter<double>("AshmanDistance")),
0044 amplitudeRatio_(pset.getParameter<double>("AmplitudeRatio")) {
0045 LogDebug(mlCommissioning_) << "[PedsFullNoiseAlgorithm::" << __func__ << "]"
0046 << " Set maximum drift of the mean value to: " << maxDriftResidualCut_
0047 << " Set minimum noise value to: " << minStripNoiseCut_
0048 << " Set maximum noise value to: " << maxStripNoiseCut_
0049 << " Set maximum noise significance value to: " << maxStripNoiseSignificanceCut_
0050 << " Set minimum Anderson-Darling p-value to: " << adProbabCut_
0051 << " Set minimum Kolmogorov-Smirnov p-value to: " << ksProbabCut_
0052 << " Set minimum Jacque-Bera p-value to: " << jbProbabCut_
0053 << " Set minimum Chi2 p-value to: " << chi2ProbabCut_
0054 << " Set N-sigma for the integral to : " << integralNsigma_
0055 << " Set maximum integral tail at N-sigma to : " << integralTailCut_
0056 << " Set maximum Kurtosis to : " << kurtosisCut_;
0057 }
0058
0059
0060
0061
0062 void PedsFullNoiseAlgorithm::extract(const std::vector<TH1*>& histos) {
0063 if (!anal()) {
0064 edm::LogWarning(mlCommissioning_) << "[PedsFullNoiseAlgorithm::" << __func__ << "]"
0065 << " NULL pointer to Analysis object!";
0066 return;
0067 }
0068
0069
0070 if (histos.size() != 3) {
0071 anal()->addErrorCode(sistrip::numberOfHistos_);
0072 }
0073
0074
0075 if (!histos.empty()) {
0076 anal()->fedKey(extractFedKey(histos.front()));
0077 }
0078
0079
0080 std::vector<TH1*>::const_iterator ihis = histos.begin();
0081 for (; ihis != histos.end(); ihis++) {
0082
0083 if (!(*ihis)) {
0084 continue;
0085 }
0086
0087 SiStripHistoTitle title((*ihis)->GetName());
0088 if (title.runType() != sistrip::PEDS_FULL_NOISE) {
0089 anal()->addErrorCode(sistrip::unexpectedTask_);
0090 continue;
0091 }
0092
0093
0094 if (title.extraInfo().find(sistrip::extrainfo::roughPedestals_) != std::string::npos) {
0095
0096 } else if (title.extraInfo().find(sistrip::extrainfo::pedestals_) != std::string::npos) {
0097 hPeds_.first = *ihis;
0098 hPeds_.second = (*ihis)->GetName();
0099 } else if (title.extraInfo().find(sistrip::extrainfo::commonMode_) != std::string::npos) {
0100
0101 } else if (title.extraInfo().find(sistrip::extrainfo::noiseProfile_) != std::string::npos) {
0102
0103 hNoise_.first = *ihis;
0104 hNoise_.second = (*ihis)->GetName();
0105 } else if (title.extraInfo().find(sistrip::extrainfo::noise2D_) != std::string::npos) {
0106 hNoise2D_.first = *ihis;
0107 hNoise2D_.second = (*ihis)->GetName();
0108 } else {
0109 anal()->addErrorCode(sistrip::unexpectedExtraInfo_);
0110 }
0111 }
0112 }
0113
0114
0115 void PedsFullNoiseAlgorithm::reset(PedsFullNoiseAnalysis* ana) {
0116 for (size_t iapv = 0; iapv < ana->peds_.size(); iapv++) {
0117 ana->pedsMean_[iapv] = 0.;
0118 ana->rawMean_[iapv] = 0.;
0119 ana->noiseMean_[iapv] = 0.;
0120 ana->pedsSpread_[iapv] = 0.;
0121 ana->noiseSpread_[iapv] = 0.;
0122 ana->rawSpread_[iapv] = 0.;
0123 ana->pedsMax_[iapv] = 0.;
0124 ana->pedsMin_[iapv] = 0.;
0125 ana->rawMax_[iapv] = 0.;
0126 ana->rawMin_[iapv] = 0.;
0127 ana->noiseMax_[iapv] = 0.;
0128 ana->noiseMin_[iapv] = 0.;
0129
0130 for (size_t istrip = 0; istrip < ana->peds_[iapv].size(); istrip++) {
0131 ana->peds_[iapv][istrip] = 0.;
0132 ana->noise_[iapv][istrip] = 0.;
0133 ana->raw_[iapv][istrip] = 0.;
0134 ana->adProbab_[iapv][istrip] = 0.;
0135 ana->ksProbab_[iapv][istrip] = 0.;
0136 ana->jbProbab_[iapv][istrip] = 0.;
0137 ana->chi2Probab_[iapv][istrip] = 0.;
0138 ana->residualRMS_[iapv][istrip] = 0.;
0139 ana->residualSigmaGaus_[iapv][istrip] = 0.;
0140 ana->noiseSignificance_[iapv][istrip] = 0.;
0141 ana->residualMean_[iapv][istrip] = 0.;
0142 ana->residualSkewness_[iapv][istrip] = 0.;
0143 ana->residualKurtosis_[iapv][istrip] = 0.;
0144 ana->residualIntegralNsigma_[iapv][istrip] = 0.;
0145 ana->residualIntegral_[iapv][istrip] = 0.;
0146 ana->deadStripBit_[iapv][istrip] = 0;
0147 ana->badStripBit_[iapv][istrip] = 0;
0148 }
0149 }
0150 }
0151
0152
0153
0154 void PedsFullNoiseAlgorithm::analyse() {
0155
0156 if (!anal()) {
0157 edm::LogWarning(mlCommissioning_) << "[PedsFullNoiseAlgorithm::" << __func__ << "]"
0158 << " NULL pointer to base Analysis object!";
0159 return;
0160 }
0161
0162 CommissioningAnalysis* tmp = const_cast<CommissioningAnalysis*>(anal());
0163 PedsFullNoiseAnalysis* ana = dynamic_cast<PedsFullNoiseAnalysis*>(tmp);
0164
0165
0166 if (!ana) {
0167 edm::LogWarning(mlCommissioning_) << "[PedsFullNoiseAlgorithm::" << __func__ << "]"
0168 << " NULL pointer to derived Analysis object!";
0169 return;
0170 }
0171
0172
0173 if (!hPeds_.first) {
0174 ana->addErrorCode(sistrip::nullPtr_);
0175 return;
0176 }
0177
0178 if (!hNoise_.first) {
0179 ana->addErrorCode(sistrip::nullPtr_);
0180 return;
0181 }
0182
0183 if (!hNoise2D_.first) {
0184 ana->addErrorCode(sistrip::nullPtr_);
0185 return;
0186 }
0187
0188
0189 TProfile* histoPeds = dynamic_cast<TProfile*>(hPeds_.first);
0190 TProfile* histoNoiseMean = dynamic_cast<TProfile*>(hNoise_.first);
0191 TH2S* histoNoise = dynamic_cast<TH2S*>(hNoise2D_.first);
0192
0193
0194 if (not histoPeds) {
0195 ana->addErrorCode(sistrip::nullPtr_);
0196 return;
0197 }
0198
0199 if (not histoNoiseMean) {
0200 ana->addErrorCode(sistrip::nullPtr_);
0201 return;
0202 }
0203
0204 if (not histoNoise) {
0205 ana->addErrorCode(sistrip::nullPtr_);
0206 return;
0207 }
0208
0209
0210 if (histoPeds->GetNbinsX() != 256) {
0211 ana->addErrorCode(sistrip::numberOfBins_);
0212 return;
0213 }
0214
0215
0216 if (histoNoiseMean->GetNbinsX() != 256) {
0217 ana->addErrorCode(sistrip::numberOfBins_);
0218 return;
0219 }
0220
0221
0222 if (histoNoise->GetNbinsY() != 256) {
0223 ana->addErrorCode(sistrip::numberOfBins_);
0224 return;
0225 }
0226
0227
0228 reset(ana);
0229
0230
0231 uint32_t apvID = -1;
0232
0233
0234 vector<float> ped_max;
0235 vector<float> ped_min;
0236 vector<float> raw_max;
0237 vector<float> raw_min;
0238 vector<float> noise_max;
0239 vector<float> noise_min;
0240
0241
0242 for (int iStrip = 0; iStrip < histoPeds->GetNbinsX(); iStrip++) {
0243 if (iStrip < histoPeds->GetNbinsX() / 2)
0244 apvID = 0;
0245 else
0246 apvID = 1;
0247
0248 int stripBin = 0;
0249 if (iStrip >= 128)
0250 stripBin = iStrip - 128;
0251 else
0252 stripBin = iStrip;
0253
0254 ana->peds_[apvID][stripBin] = histoPeds->GetBinContent(iStrip + 1);
0255 ana->noise_[apvID][stripBin] = histoNoiseMean->GetBinContent(iStrip + 1);
0256 ana->raw_[apvID][stripBin] = histoPeds->GetBinError(iStrip + 1);
0257
0258 ana->pedsMean_[apvID] += ana->peds_[apvID][stripBin];
0259 ana->rawMean_[apvID] += ana->raw_[apvID][stripBin];
0260 ana->noiseMean_[apvID] += ana->noise_[apvID][stripBin];
0261
0262
0263 if (ped_max.size() < apvID + 1)
0264 ped_max.push_back(ana->peds_[apvID][stripBin]);
0265 else {
0266 if (ana->peds_[apvID][stripBin] > ped_max.at(apvID))
0267 ped_max.at(apvID) = ana->peds_[apvID][stripBin];
0268 }
0269
0270
0271 if (ped_min.size() < apvID + 1)
0272 ped_min.push_back(ana->peds_[apvID][stripBin]);
0273 else {
0274 if (ana->peds_[apvID][stripBin] < ped_min.at(apvID))
0275 ped_min.at(apvID) = ana->peds_[apvID][stripBin];
0276 }
0277
0278
0279 if (noise_max.size() < apvID + 1)
0280 noise_max.push_back(ana->noise_[apvID][stripBin]);
0281 else {
0282 if (ana->noise_[apvID][stripBin] > noise_max.at(apvID))
0283 noise_max.at(apvID) = ana->noise_[apvID][stripBin];
0284 }
0285
0286
0287 if (noise_min.size() < apvID + 1)
0288 noise_min.push_back(ana->noise_[apvID][stripBin]);
0289 else {
0290 if (ana->noise_[apvID][stripBin] < noise_min.at(apvID))
0291 noise_min.at(apvID) = ana->noise_[apvID][stripBin];
0292 }
0293
0294
0295 if (raw_max.size() < apvID + 1)
0296 raw_max.push_back(ana->raw_[apvID][stripBin]);
0297 else {
0298 if (ana->raw_[apvID][stripBin] > raw_max.at(apvID))
0299 raw_max.at(apvID) = ana->raw_[apvID][stripBin];
0300 }
0301
0302
0303 if (raw_min.size() < apvID + 1)
0304 raw_min.push_back(ana->raw_[apvID][stripBin]);
0305 else {
0306 if (ana->raw_[apvID][stripBin] < raw_min.at(apvID))
0307 raw_min.at(apvID) = ana->raw_[apvID][stripBin];
0308 }
0309 }
0310
0311
0312 for (unsigned int iApv = 0; iApv < ana->pedsMean_.size(); iApv++) {
0313 ana->pedsMean_.at(iApv) /= (ana->peds_[iApv].size());
0314 ana->rawMean_.at(iApv) /= (ana->raw_[iApv].size());
0315 ana->noiseMean_.at(iApv) /= (ana->noise_[iApv].size());
0316 }
0317
0318
0319 for (unsigned int iApv = 0; iApv < ped_max.size(); iApv++) {
0320 if (ped_max.at(iApv) > sistrip::maximum_)
0321 ana->pedsMax_.at(iApv) = sistrip::maximum_;
0322 else if (ped_max.at(iApv) < -1. * sistrip::maximum_)
0323 ana->pedsMax_.at(iApv) = -1. * sistrip::maximum_;
0324 else
0325 ana->pedsMax_.at(iApv) = ped_max.at(iApv);
0326
0327 if (ped_min.at(iApv) > sistrip::maximum_)
0328 ana->pedsMin_.at(iApv) = sistrip::maximum_;
0329 else if (ped_min.at(iApv) < -1. * sistrip::maximum_)
0330 ana->pedsMin_.at(iApv) = -1. * sistrip::maximum_;
0331 else
0332 ana->pedsMin_.at(iApv) = ped_min.at(iApv);
0333
0334 if (noise_max.at(iApv) > sistrip::maximum_)
0335 ana->noiseMax_.at(iApv) = sistrip::maximum_;
0336 else if (noise_max.at(iApv) < -1. * sistrip::maximum_)
0337 ana->noiseMax_.at(iApv) = -1. * sistrip::maximum_;
0338 else
0339 ana->noiseMax_.at(iApv) = noise_max.at(iApv);
0340
0341 if (noise_min.at(iApv) > sistrip::maximum_)
0342 ana->noiseMin_.at(iApv) = sistrip::maximum_;
0343 else if (noise_min.at(iApv) < -1. * sistrip::maximum_)
0344 ana->noiseMin_.at(iApv) = -1. * sistrip::maximum_;
0345 else
0346 ana->noiseMin_.at(iApv) = noise_min.at(iApv);
0347
0348 if (raw_max.at(iApv) > sistrip::maximum_)
0349 ana->rawMax_.at(iApv) = sistrip::maximum_;
0350 else if (raw_max.at(iApv) < -1. * sistrip::maximum_)
0351 ana->rawMax_.at(iApv) = -1. * sistrip::maximum_;
0352 else
0353 ana->rawMax_.at(iApv) = raw_max.at(iApv);
0354
0355 if (raw_min.at(iApv) > sistrip::maximum_)
0356 ana->rawMin_.at(iApv) = sistrip::maximum_;
0357 else if (raw_min.at(iApv) < -1. * sistrip::maximum_)
0358 ana->rawMin_.at(iApv) = -1. * sistrip::maximum_;
0359 else
0360 ana->rawMin_.at(iApv) = raw_min.at(iApv);
0361 }
0362
0363
0364 for (int iStrip = 0; iStrip < histoNoiseMean->GetNbinsX(); iStrip++) {
0365 if (iStrip < histoNoiseMean->GetNbinsX() / 2)
0366 apvID = 0;
0367 else
0368 apvID = 1;
0369 ana->pedsSpread_[apvID] += pow(histoPeds->GetBinContent(iStrip + 1) - ana->pedsMean_.at(apvID), 2);
0370 ana->noiseSpread_[apvID] += pow(histoNoiseMean->GetBinContent(iStrip + 1) - ana->noiseMean_.at(apvID), 2);
0371 ana->rawSpread_[apvID] += pow(histoPeds->GetBinError(iStrip + 1) - ana->rawMean_.at(apvID), 2);
0372 }
0373
0374 for (unsigned int iApv = 0; iApv < ana->pedsSpread_.size(); iApv++) {
0375 ana->pedsSpread_[iApv] = sqrt(ana->pedsSpread_[iApv]) / sqrt(ana->peds_[iApv].size() - 1);
0376 ana->noiseSpread_[iApv] = sqrt(ana->noiseSpread_[iApv]) / sqrt(ana->noise_[iApv].size() - 1);
0377 ana->rawSpread_[iApv] = sqrt(ana->rawSpread_[iApv]) / sqrt(ana->raw_[iApv].size() - 1);
0378 }
0379
0380
0381 TH1S* histoResidualStrip = new TH1S("histoResidualStrip",
0382 "",
0383 histoNoise->GetNbinsX(),
0384 histoNoise->GetXaxis()->GetXmin(),
0385 histoNoise->GetXaxis()->GetXmax());
0386 histoResidualStrip->Sumw2();
0387 histoResidualStrip->SetDirectory(nullptr);
0388 TF1* fitFunc = new TF1("fitFunc", "gaus(0)", histoNoise->GetXaxis()->GetXmin(), histoNoise->GetXaxis()->GetXmax());
0389 TF1* fit2Gaus = nullptr;
0390 TH1F* randomHisto = nullptr;
0391 TFitResultPtr result;
0392
0393 for (int iStrip = 0; iStrip < histoNoise->GetNbinsY(); iStrip++) {
0394
0395 if (iStrip < histoNoise->GetNbinsY() / 2)
0396 apvID = 0;
0397 else
0398 apvID = 1;
0399 histoResidualStrip->Reset();
0400
0401 int stripBin = 0;
0402 if (iStrip >= 128)
0403 stripBin = iStrip - 128;
0404 else
0405 stripBin = iStrip;
0406
0407 for (int iBinX = 0; iBinX < histoNoise->GetNbinsX(); iBinX++) {
0408 histoResidualStrip->SetBinContent(iBinX + 1, histoNoise->GetBinContent(iBinX + 1, iStrip + 1));
0409 histoResidualStrip->SetBinError(iBinX + 1, histoNoise->GetBinError(iBinX + 1, iStrip + 1));
0410 }
0411
0412 if (histoResidualStrip->Integral() == 0) {
0413
0414
0415 ana->adProbab_[apvID][stripBin] = 0;
0416 ana->ksProbab_[apvID][stripBin] = 0;
0417 ana->jbProbab_[apvID][stripBin] = 0;
0418 ana->chi2Probab_[apvID][stripBin] = 0;
0419 ana->noiseSignificance_[apvID][stripBin] = 0;
0420 ana->residualMean_[apvID][stripBin] = 0;
0421 ana->residualRMS_[apvID][stripBin] = 0;
0422 ana->residualSigmaGaus_[apvID][stripBin] = 0;
0423 ana->residualSkewness_[apvID][stripBin] = 0;
0424 ana->residualKurtosis_[apvID][stripBin] = 0;
0425 ana->residualIntegralNsigma_[apvID][stripBin] = 0;
0426 ana->residualIntegral_[apvID][stripBin] = 0;
0427 ana->deadStrip_[apvID].push_back(stripBin);
0428 ana->deadStripBit_[apvID][stripBin] = 1;
0429 ana->badStripBit_[apvID][stripBin] = 0;
0430
0431 SiStripFecKey fec_key(ana->fecKey());
0432 LogTrace(mlDqmClient_) << "DeadStrip: fecCrate "
0433 << " " << fec_key.fecCrate() << " fecSlot " << fec_key.fecSlot() << " fecRing "
0434 << fec_key.fecRing() << " ccuAddr " << fec_key.ccuAddr() << " ccChan "
0435 << fec_key.ccuChan() << " lldChan " << fec_key.lldChan() << " apvID " << apvID
0436 << " stripID " << iStrip;
0437
0438 continue;
0439 }
0440
0441
0442 ana->residualMean_[apvID][stripBin] = histoResidualStrip->GetMean();
0443 ana->residualRMS_[apvID][stripBin] = histoResidualStrip->GetRMS();
0444 ana->residualSkewness_[apvID][stripBin] = histoResidualStrip->GetSkewness();
0445 ana->residualKurtosis_[apvID][stripBin] = histoResidualStrip->GetKurtosis();
0446 ana->noiseSignificance_[apvID][stripBin] =
0447 (ana->noise_[apvID][stripBin] - ana->noiseMean_[apvID]) / ana->noiseSpread_[apvID];
0448 ana->residualIntegral_[apvID][stripBin] = histoResidualStrip->Integral();
0449 ana->residualIntegralNsigma_[apvID][stripBin] =
0450 (histoResidualStrip->Integral(histoResidualStrip->FindBin(ana->residualMean_[apvID][stripBin] +
0451 ana->residualRMS_[apvID][stripBin] * integralNsigma_),
0452 histoResidualStrip->GetNbinsX() + 1) +
0453 histoResidualStrip->Integral(
0454 0,
0455 histoResidualStrip->FindBin(ana->residualMean_[apvID][stripBin] -
0456 ana->residualRMS_[apvID][stripBin] * integralNsigma_))) /
0457 ana->residualIntegral_[apvID][stripBin];
0458
0459
0460 fitFunc->SetRange(histoNoise->GetXaxis()->GetXmin(), histoNoise->GetXaxis()->GetXmax());
0461 fitFunc->SetParameters(ana->residualIntegral_[apvID][stripBin],
0462 ana->residualMean_[apvID][stripBin],
0463 ana->residualRMS_[apvID][stripBin]);
0464 result = histoResidualStrip->Fit(fitFunc, "QSRN");
0465
0466
0467 if (result.Get()) {
0468 ana->residualSigmaGaus_[apvID][stripBin] = fitFunc->GetParameter(2);
0469 ana->chi2Probab_[apvID][stripBin] = result->Prob();
0470
0471
0472 float jbVal =
0473 (ana->residualIntegral_[apvID][stripBin] / 6) *
0474 (pow(ana->residualSkewness_[apvID][stripBin], 2) + pow(ana->residualKurtosis_[apvID][stripBin], 2) / 4);
0475 ana->jbProbab_[apvID][stripBin] = ROOT::Math::chisquared_cdf_c(jbVal, 2);
0476
0477
0478 if (randomHisto == nullptr)
0479 randomHisto = (TH1F*)histoResidualStrip->Clone("randomHisto");
0480 randomHisto->Reset();
0481 randomHisto->SetDirectory(nullptr);
0482
0483 if (generateRandomHisto_) {
0484 randomHisto->FillRandom("fitFunc", histoResidualStrip->Integral());
0485 if (randomHisto->Integral() != 0) {
0486 ana->ksProbab_[apvID][stripBin] = histoResidualStrip->KolmogorovTest(randomHisto, "N");
0487 ana->adProbab_[apvID][stripBin] = histoResidualStrip->AndersonDarlingTest(randomHisto);
0488 } else {
0489 ana->ksProbab_[apvID][stripBin] = 0;
0490 ana->adProbab_[apvID][stripBin] = 0;
0491 }
0492
0493 } else {
0494 randomHisto->Add(fitFunc);
0495 if (randomHisto->Integral() != 0) {
0496 ana->ksProbab_[apvID][stripBin] = histoResidualStrip->KolmogorovTest(randomHisto, "N");
0497 ROOT::Fit::BinData data1;
0498 ROOT::Fit::BinData data2;
0499 ROOT::Fit::FillData(data1, histoResidualStrip, nullptr);
0500 data2.Initialize(randomHisto->GetNbinsX() + 1, 1);
0501 for (int ibin = 0; ibin < randomHisto->GetNbinsX(); ibin++) {
0502 if (histoResidualStrip->GetBinContent(ibin + 1) != 0 or randomHisto->GetBinContent(ibin + 1) >= 1)
0503 data2.Add(randomHisto->GetBinCenter(ibin + 1),
0504 randomHisto->GetBinContent(ibin + 1),
0505 randomHisto->GetBinError(ibin + 1));
0506 }
0507
0508 double probab, value;
0509 ROOT::Math::GoFTest::AndersonDarling2SamplesTest(data1, data2, probab, value);
0510 ana->adProbab_[apvID][stripBin] = probab;
0511 } else {
0512 ana->ksProbab_[apvID][stripBin] = 0;
0513 ana->adProbab_[apvID][stripBin] = 0;
0514 }
0515 }
0516 }
0517
0518
0519 bool badStripFlag = false;
0520 ana->deadStripBit_[apvID][stripBin] = 0;
0521
0522 if (fabs(ana->residualMean_[apvID][stripBin]) > maxDriftResidualCut_ and not badStripFlag) {
0523 ana->shiftedStrip_[apvID].push_back(stripBin);
0524 badStripFlag = true;
0525 }
0526
0527 if (ana->residualRMS_[apvID][stripBin] < minStripNoiseCut_ and not badStripFlag) {
0528 ana->lowNoiseStrip_[apvID].push_back(stripBin);
0529 badStripFlag = true;
0530 }
0531
0532 if (ana->residualRMS_[apvID][stripBin] > maxStripNoiseCut_ and not badStripFlag) {
0533 ana->largeNoiseStrip_[apvID].push_back(stripBin);
0534 badStripFlag = true;
0535 }
0536
0537 if (fabs(ana->noiseSignificance_[apvID][stripBin]) > maxStripNoiseSignificanceCut_ and
0538 not badStripFlag) {
0539 ana->largeNoiseSignificance_[apvID].push_back(stripBin);
0540 badStripFlag = true;
0541 }
0542
0543 if (result.Get() and result->Status() != 0)
0544 ana->badFitStatus_[apvID].push_back(stripBin);
0545
0546 if (ana->adProbab_[apvID][stripBin] < adProbabCut_ and not badStripFlag)
0547 ana->badADProbab_[apvID].push_back(stripBin);
0548
0549 if (ana->ksProbab_[apvID][stripBin] < ksProbabCut_ and not badStripFlag)
0550 ana->badKSProbab_[apvID].push_back(stripBin);
0551
0552 if (ana->jbProbab_[apvID][stripBin] < jbProbabCut_ and not badStripFlag)
0553 ana->badJBProbab_[apvID].push_back(stripBin);
0554
0555 if (ana->chi2Probab_[apvID][stripBin] < chi2ProbabCut_ and
0556 not badStripFlag)
0557 ana->badChi2Probab_[apvID].push_back(stripBin);
0558
0559 if (ana->adProbab_[apvID][stripBin] < adProbabCut_ and ana->ksProbab_[apvID][stripBin] < ksProbabCut_ and
0560 ana->jbProbab_[apvID][stripBin] < jbProbabCut_ and ana->chi2Probab_[apvID][stripBin] < chi2ProbabCut_)
0561 badStripFlag = true;
0562
0563 if (ana->residualKurtosis_[apvID][stripBin] > kurtosisCut_ and
0564 ana->residualIntegralNsigma_[apvID][stripBin] > integralTailCut_ and not badStripFlag) {
0565 ana->badTailStrip_[apvID].push_back(stripBin);
0566 badStripFlag = true;
0567 }
0568
0569 if (badStripFlag) {
0570
0571 fit2Gaus = new TF1("dgaus",
0572 "[0]*exp(-((x-[1])*(x-[1]))/(2*[2]*[2]))+[3]*exp(-((x-[4])*(x-[4]))/(2*[5]*[5]))",
0573 histoNoise->GetXaxis()->GetXmin(),
0574 histoNoise->GetXaxis()->GetXmax());
0575 fit2Gaus->SetParameter(0, fitFunc->GetParameter(0) / 2);
0576 fit2Gaus->SetParameter(3, fitFunc->GetParameter(0) / 2);
0577 fit2Gaus->SetParameter(1, 1.);
0578 fit2Gaus->SetParameter(4, -1.);
0579 fit2Gaus->SetParameter(2, fitFunc->GetParameter(2));
0580 fit2Gaus->SetParameter(5, fitFunc->GetParameter(2));
0581 fit2Gaus->SetParLimits(1, 0., histoNoise->GetXaxis()->GetXmax());
0582 fit2Gaus->SetParLimits(4, histoNoise->GetXaxis()->GetXmin(), 0);
0583 result = histoResidualStrip->Fit(fit2Gaus, "QSR");
0584
0585
0586 float ashman = TMath::Power(2, 0.5) * abs(fit2Gaus->GetParameter(1) - fit2Gaus->GetParameter(4)) /
0587 (sqrt(pow(fit2Gaus->GetParameter(2), 2) + pow(fit2Gaus->GetParameter(5), 2)));
0588
0589 float amplitudeRatio = std::min(fit2Gaus->GetParameter(0), fit2Gaus->GetParameter(3)) /
0590 std::max(fit2Gaus->GetParameter(0), fit2Gaus->GetParameter(3));
0591
0592 if (ashman > ashmanDistance_ and amplitudeRatio > amplitudeRatio_)
0593 ana->badDoublePeakStrip_[apvID].push_back(stripBin);
0594 }
0595
0596 if (badStripFlag) {
0597
0598 ana->badStrip_[apvID].push_back(stripBin);
0599 ana->badStripBit_[apvID][stripBin] = 1;
0600
0601 SiStripFecKey fec_key(ana->fecKey());
0602 LogTrace(mlDqmClient_) << "BadStrip: fecCrate "
0603 << " " << fec_key.fecCrate() << " fecSlot " << fec_key.fecSlot() << " fecRing "
0604 << fec_key.fecRing() << " ccuAddr " << fec_key.ccuAddr() << " ccChan "
0605 << fec_key.ccuChan() << " lldChan " << fec_key.lldChan() << " apvID " << apvID
0606 << " stripID " << stripBin;
0607
0608 } else
0609 ana->badStripBit_[apvID][stripBin] = 0;
0610 }
0611
0612 ped_max.clear();
0613 ped_min.clear();
0614 raw_max.clear();
0615 raw_min.clear();
0616 noise_max.clear();
0617 noise_min.clear();
0618 if (histoResidualStrip)
0619 delete histoResidualStrip;
0620 if (fitFunc)
0621 delete fitFunc;
0622 if (randomHisto)
0623 delete randomHisto;
0624 if (fit2Gaus)
0625 delete fit2Gaus;
0626 }