File indexing completed on 2024-04-06 11:59:40
0001
0002
0003
0004
0005
0006 #include "CalibTracker/SiStripChannelGain/plugins/SiStripGainCosmicCalculator.h"
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0010 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0011 #include "Geometry/CommonTopologies/interface/StripTopology.h"
0012 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
0013 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetType.h"
0014 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
0015 #include "CLHEP/Random/RandFlat.h"
0016 #include "CLHEP/Random/RandGauss.h"
0017 #include "Geometry/CommonDetUnit/interface/PixelGeomDetType.h"
0018 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0019 #include "DataFormats/TrackReco/interface/Track.h"
0020 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0021 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
0022 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
0023 #include "DataFormats/TrackerCommon/interface/SiStripSubStructure.h"
0024 #include "DataFormats/DetId/interface/DetId.h"
0025 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0026 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0027
0028
0029
0030 SiStripGainCosmicCalculator::SiStripGainCosmicCalculator(const edm::ParameterSet& iConfig)
0031 : ConditionDBWriter<SiStripApvGain>(iConfig) {
0032 edm::LogInfo("SiStripGainCosmicCalculator::SiStripGainCosmicCalculator");
0033 ExpectedChargeDeposition = 200.;
0034 edm::LogInfo("SiStripApvGainCalculator::SiStripApvGainCalculator")
0035 << "ExpectedChargeDeposition=" << ExpectedChargeDeposition;
0036
0037 TrackProducer = iConfig.getParameter<std::string>("TrackProducer");
0038 TrackLabel = iConfig.getParameter<std::string>("TrackLabel");
0039
0040 detModulesToBeExcluded.clear();
0041 detModulesToBeExcluded = iConfig.getParameter<std::vector<unsigned> >("detModulesToBeExcluded");
0042 MinNrEntries = iConfig.getUntrackedParameter<unsigned>("minNrEntries", 20);
0043 MaxChi2OverNDF = iConfig.getUntrackedParameter<double>("maxChi2OverNDF", 5.);
0044
0045 outputHistogramsInRootFile = iConfig.getParameter<bool>("OutputHistogramsInRootFile");
0046 outputFileName = iConfig.getParameter<std::string>("OutputFileName");
0047
0048 edm::LogInfo("SiStripApvGainCalculator")
0049 << "Clusters from " << detModulesToBeExcluded.size() << " modules will be ignored in the calibration:";
0050 edm::LogInfo("SiStripApvGainCalculator") << "The calibration for these DetIds will be set to a default value";
0051 for (std::vector<uint32_t>::const_iterator imod = detModulesToBeExcluded.begin();
0052 imod != detModulesToBeExcluded.end();
0053 imod++) {
0054 edm::LogInfo("SiStripApvGainCalculator") << "exclude detid = " << *imod;
0055 }
0056
0057 printdebug_ = iConfig.getUntrackedParameter<bool>("printDebug", false);
0058
0059 tTopoToken_ = esConsumes<edm::Transition::BeginRun>();
0060 detCablingToken_ = esConsumes<edm::Transition::BeginRun>();
0061 tkGeomToken_ = esConsumes<edm::Transition::BeginRun>();
0062 }
0063
0064 SiStripGainCosmicCalculator::~SiStripGainCosmicCalculator() {
0065 edm::LogInfo("SiStripGainCosmicCalculator::~SiStripGainCosmicCalculator");
0066 }
0067
0068 void SiStripGainCosmicCalculator::algoEndJob() {}
0069
0070 void SiStripGainCosmicCalculator::algoBeginJob(const edm::EventSetup& iSetup) {
0071 tTopo_ = &iSetup.getData(tTopoToken_);
0072 siStripDetCabling_ = &iSetup.getData(detCablingToken_);
0073 tkGeom_ = &iSetup.getData(tkGeomToken_);
0074
0075 std::cout << "SiStripGainCosmicCalculator::algoBeginJob called" << std::endl;
0076 total_nr_of_events = 0;
0077 HlistAPVPairs = new TObjArray();
0078 HlistOtherHistos = new TObjArray();
0079
0080 HlistOtherHistos->Add(new TH1F(Form("APVPairCorrections"), Form("APVPairCorrections"), 50, -1., 4.));
0081 HlistOtherHistos->Add(new TH1F(Form("APVPairCorrectionsTIB1mono"), Form("APVPairCorrectionsTIB1mono"), 50, -1., 4.));
0082 HlistOtherHistos->Add(
0083 new TH1F(Form("APVPairCorrectionsTIB1stereo"), Form("APVPairCorrectionsTIB1stereo"), 50, -1., 4.));
0084 HlistOtherHistos->Add(new TH1F(Form("APVPairCorrectionsTIB2"), Form("APVPairCorrectionsTIB2"), 50, -1., 4.));
0085 HlistOtherHistos->Add(new TH1F(Form("APVPairCorrectionsTOB1"), Form("APVPairCorrectionsTOB1"), 50, -1., 4.));
0086 HlistOtherHistos->Add(new TH1F(Form("APVPairCorrectionsTOB2"), Form("APVPairCorrectionsTOB2"), 50, -1., 4.));
0087 HlistOtherHistos->Add(new TH1F(Form("LocalAngle"), Form("LocalAngle"), 70, -0.1, 3.4));
0088 HlistOtherHistos->Add(new TH1F(Form("LocalAngleAbsoluteCosine"), Form("LocalAngleAbsoluteCosine"), 48, -0.1, 1.1));
0089 HlistOtherHistos->Add(new TH1F(Form("LocalPosition_cm"), Form("LocalPosition_cm"), 100, -5., 5.));
0090 HlistOtherHistos->Add(new TH1F(Form("LocalPosition_normalized"), Form("LocalPosition_normalized"), 100, -1.1, 1.1));
0091 TH1F* local_histo = new TH1F(Form("SiStripRecHitType"), Form("SiStripRecHitType"), 2, 0.5, 2.5);
0092 HlistOtherHistos->Add(local_histo);
0093 local_histo->GetXaxis()->SetBinLabel(1, "simple");
0094 local_histo->GetXaxis()->SetBinLabel(2, "matched");
0095
0096
0097 std::vector<uint32_t> activeDets;
0098 activeDets.clear();
0099 SelectedDetIds.clear();
0100 siStripDetCabling_->addActiveDetectorsRawIds(activeDets);
0101
0102
0103 SiStripSubStructure::getTIBDetectors(
0104 activeDets, SelectedDetIds, tTopo_, 0, 0, 0, 0);
0105 SiStripSubStructure::getTOBDetectors(
0106 activeDets, SelectedDetIds, tTopo_, 0, 0, 0);
0107
0108 for (auto det : tkGeom_->dets()) {
0109 if (dynamic_cast<const StripGeomDetUnit*>(det) != nullptr) {
0110 uint32_t detid = det->geographicalId().rawId();
0111
0112 double module_thickness =
0113 det->surface()
0114 .bounds()
0115 .thickness();
0116 thickness_map.insert(std::make_pair(detid, module_thickness));
0117
0118 const bool is_active_detector =
0119 std::end(SelectedDetIds) != std::find(std::begin(SelectedDetIds), std::end(SelectedDetIds), detid);
0120
0121 const bool exclude_this_detid =
0122 std::end(detModulesToBeExcluded) !=
0123 std::find(std::begin(detModulesToBeExcluded), std::end(detModulesToBeExcluded), detid);
0124
0125 if (is_active_detector &&
0126 (!exclude_this_detid)) {
0127 const StripTopology& p = dynamic_cast<const StripGeomDetUnit*>(det)->specificTopology();
0128 unsigned short NAPVPairs = p.nstrips() / 256;
0129 if (NAPVPairs < 2 || NAPVPairs > 3) {
0130 edm::LogError("SiStripGainCosmicCalculator")
0131 << "Problem with Number of strips in detector: " << p.nstrips() << " Exiting program";
0132 exit(1);
0133 }
0134 for (int iapp = 0; iapp < NAPVPairs; iapp++) {
0135 TString hid = Form("ChargeAPVPair_%i_%i", detid, iapp);
0136 HlistAPVPairs->Add(
0137 new TH1F(hid, hid, 45, 0., 1350.));
0138 }
0139 }
0140 }
0141 }
0142 }
0143
0144
0145 void SiStripGainCosmicCalculator::algoAnalyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0146 using namespace edm;
0147 total_nr_of_events++;
0148
0149
0150
0151
0152
0153
0154
0155
0156 Handle<reco::TrackCollection> trackCollection;
0157 iEvent.getByLabel(TrackProducer, TrackLabel, trackCollection);
0158 const reco::TrackCollection* tracks = trackCollection.product();
0159
0160
0161
0162
0163
0164
0165 for (reco::TrackCollection::const_iterator itr = tracks->begin(); itr != tracks->end();
0166 itr++) {
0167
0168
0169
0170 std::vector<std::pair<const TrackingRecHit*, float> >
0171 hitangle;
0172
0173 for (std::vector<std::pair<const TrackingRecHit*, float> >::const_iterator hitangle_iter = hitangle.begin();
0174 hitangle_iter != hitangle.end();
0175 hitangle_iter++) {
0176 const TrackingRecHit* trechit = hitangle_iter->first;
0177 float local_angle = hitangle_iter->second;
0178 LocalPoint local_position = trechit->localPosition();
0179 const SiStripRecHit2D* sistripsimplehit = dynamic_cast<const SiStripRecHit2D*>(trechit);
0180 const SiStripMatchedRecHit2D* sistripmatchedhit = dynamic_cast<const SiStripMatchedRecHit2D*>(trechit);
0181
0182 ((TH1F*)HlistOtherHistos->FindObject("LocalAngle"))->Fill(local_angle);
0183 ((TH1F*)HlistOtherHistos->FindObject("LocalAngleAbsoluteCosine"))->Fill(fabs(cos(local_angle)));
0184 if (sistripsimplehit) {
0185 ((TH1F*)HlistOtherHistos->FindObject("SiStripRecHitType"))->Fill(1.);
0186 const SiStripRecHit2D::ClusterRef& cluster = sistripsimplehit->cluster();
0187 const auto& ampls = cluster->amplitudes();
0188 uint32_t thedetid = 0;
0189 double module_width = moduleWidth(thedetid);
0190 ((TH1F*)HlistOtherHistos->FindObject("LocalPosition_cm"))->Fill(local_position.x());
0191 ((TH1F*)HlistOtherHistos->FindObject("LocalPosition_normalized"))->Fill(local_position.x() / module_width);
0192 double module_thickness = moduleThickness(thedetid);
0193 int ifirststrip = cluster->firstStrip();
0194 int theapvpairid = int(float(ifirststrip) / 256.);
0195 TH1F* histopointer = (TH1F*)HlistAPVPairs->FindObject(Form("ChargeAPVPair_%i_%i", thedetid, theapvpairid));
0196 if (histopointer) {
0197 short cCharge = 0;
0198 for (unsigned int iampl = 0; iampl < ampls.size(); iampl++) {
0199 cCharge += ampls[iampl];
0200 }
0201 double cluster_charge_over_path = ((double)cCharge) * fabs(cos(local_angle)) / (10. * module_thickness);
0202 histopointer->Fill(cluster_charge_over_path);
0203 }
0204 } else {
0205 if (sistripmatchedhit)
0206 ((TH1F*)HlistOtherHistos->FindObject("SiStripRecHitType"))->Fill(2.);
0207 }
0208 }
0209 }
0210 }
0211
0212
0213 std::pair<double, double> SiStripGainCosmicCalculator::getPeakOfLandau(
0214 TH1F* inputHisto) {
0215
0216 double adcs = -0.5;
0217 double error = 0.;
0218 double nr_of_entries = inputHisto->GetEntries();
0219 if (nr_of_entries < MinNrEntries) {
0220 return std::make_pair(adcs, error);
0221 }
0222
0223
0224
0225
0226
0227
0228
0229
0230
0231
0232 std::unique_ptr<TF1> fitfunction(new TF1("landaufit", "landau"));
0233 inputHisto->Fit(fitfunction.get(), "0Q");
0234 adcs = fitfunction->GetParameter("MPV");
0235 error = fitfunction->GetParError(1);
0236 double chi2 = fitfunction->GetChisquare();
0237 double ndf = fitfunction->GetNDF();
0238 double chi2overndf = chi2 / ndf;
0239
0240 if (adcs < 2. || (error / adcs) > 1.8) {
0241 inputHisto->Fit(fitfunction.get(), "0Q", nullptr, 0., 400.);
0242 std::cout << "refitting landau for histogram " << inputHisto->GetTitle() << std::endl;
0243 std::cout << "initial error/adcs =" << error << " / " << adcs << std::endl;
0244 std::cout << "new error/adcs =" << fitfunction->GetParError(1) << " / " << fitfunction->GetParameter("MPV")
0245 << std::endl;
0246 adcs = fitfunction->GetParameter("MPV");
0247 error = fitfunction->GetParError(1);
0248 chi2 = fitfunction->GetChisquare();
0249 ndf = fitfunction->GetNDF();
0250 chi2overndf = chi2 / ndf;
0251 }
0252
0253 if (adcs < 2. || chi2overndf > MaxChi2OverNDF) {
0254 adcs = -0.5;
0255 error = 0.;
0256 }
0257 return std::make_pair(adcs, error);
0258 }
0259
0260
0261 double SiStripGainCosmicCalculator::moduleWidth(const uint32_t detid)
0262 {
0263 double module_width = 0.;
0264 const GeomDetUnit* it = tkGeom_->idToDetUnit(DetId(detid));
0265 if (dynamic_cast<const StripGeomDetUnit*>(it) == nullptr && dynamic_cast<const PixelGeomDetUnit*>(it) == nullptr) {
0266 std::cout << "this detID doesn't seem to belong to the Tracker" << std::endl;
0267 } else {
0268 module_width = it->surface().bounds().width();
0269 }
0270 return module_width;
0271 }
0272
0273
0274 double SiStripGainCosmicCalculator::moduleThickness(const uint32_t detid)
0275 {
0276 double module_thickness = 0.;
0277 const GeomDetUnit* it = tkGeom_->idToDetUnit(DetId(detid));
0278 if (dynamic_cast<const StripGeomDetUnit*>(it) == nullptr && dynamic_cast<const PixelGeomDetUnit*>(it) == nullptr) {
0279 std::cout << "this detID doesn't seem to belong to the Tracker" << std::endl;
0280 } else {
0281 module_thickness = it->surface().bounds().thickness();
0282 }
0283 return module_thickness;
0284 }
0285
0286
0287 std::unique_ptr<SiStripApvGain> SiStripGainCosmicCalculator::getNewObject() {
0288 std::cout << "SiStripGainCosmicCalculator::getNewObject called" << std::endl;
0289
0290 std::cout << "total_nr_of_events=" << total_nr_of_events << std::endl;
0291
0292 TH1F* ChargeOfEachAPVPair = new TH1F("ChargeOfEachAPVPair", "ChargeOfEachAPVPair", 1, 0, 1);
0293 ChargeOfEachAPVPair->SetCanExtend(TH1::kAllAxes);
0294 TH1F* EntriesApvPairs = new TH1F("EntriesApvPairs", "EntriesApvPairs", 1, 0, 1);
0295 EntriesApvPairs->SetCanExtend(TH1::kAllAxes);
0296 TH1F* NrOfEntries =
0297 new TH1F("NrOfEntries", "NrOfEntries", 351, -0.5, 350.5);
0298 TH1F* ModuleThickness = new TH1F("ModuleThickness", "ModuleThickness", 2, 0.5, 2.5);
0299 HlistOtherHistos->Add(ModuleThickness);
0300 ModuleThickness->GetXaxis()->SetBinLabel(1, "320mu");
0301 ModuleThickness->GetXaxis()->SetBinLabel(2, "500mu");
0302 ModuleThickness->SetYTitle("Nr APVPairs");
0303 TH1F* ModuleWidth = new TH1F("ModuleWidth", "ModuleWidth", 5, 0.5, 5.5);
0304 HlistOtherHistos->Add(ModuleWidth);
0305 ModuleWidth->GetXaxis()->SetBinLabel(1, "6.144cm");
0306 ModuleWidth->GetXaxis()->SetBinLabel(2, "7.14cm");
0307 ModuleWidth->GetXaxis()->SetBinLabel(3, "9.3696cm");
0308 ModuleWidth->GetXaxis()->SetBinLabel(4, "10.49cm");
0309 ModuleWidth->GetXaxis()->SetBinLabel(5, "12.03cm");
0310 ModuleWidth->SetYTitle("Nr APVPairs");
0311
0312 HlistAPVPairs->Sort();
0313 TIter hiterator(HlistAPVPairs);
0314 double MeanCharge = 0.;
0315 double NrOfApvPairs = 0.;
0316 TH1F* MyHisto = (TH1F*)hiterator();
0317 while (MyHisto) {
0318 TString histo_title = MyHisto->GetTitle();
0319 if (histo_title.Contains("ChargeAPVPair_")) {
0320 std::pair<double, double> two_values = getPeakOfLandau(MyHisto);
0321 double local_nrofadcs = two_values.first;
0322 double local_sigma = two_values.second;
0323 ChargeOfEachAPVPair->Fill(histo_title, local_nrofadcs);
0324 int ichbin = ChargeOfEachAPVPair->GetXaxis()->FindBin(histo_title.Data());
0325 ChargeOfEachAPVPair->SetBinError(ichbin, local_sigma);
0326 EntriesApvPairs->Fill(histo_title, MyHisto->GetEntries());
0327 NrOfEntries->Fill(MyHisto->GetEntries());
0328 if (local_nrofadcs > 0) {
0329 MeanCharge += local_nrofadcs;
0330 NrOfApvPairs += 1.;
0331 }
0332 }
0333 MyHisto = (TH1F*)hiterator();
0334 }
0335 ChargeOfEachAPVPair->LabelsDeflate("X");
0336 EntriesApvPairs->LabelsDeflate("X");
0337 HlistOtherHistos->Add(ChargeOfEachAPVPair);
0338 HlistOtherHistos->Add(EntriesApvPairs);
0339 HlistOtherHistos->Add(NrOfEntries);
0340 MeanCharge = MeanCharge / NrOfApvPairs;
0341
0342 TH1F* CorrectionOfEachAPVPair = (TH1F*)ChargeOfEachAPVPair->Clone("CorrectionOfEachAPVPair");
0343 TH1F* ChargeOfEachAPVPairControlView =
0344 new TH1F("ChargeOfEachAPVPairControlView", "ChargeOfEachAPVPairControlView", 1, 0, 1);
0345 ChargeOfEachAPVPairControlView->SetCanExtend(TH1::kAllAxes);
0346 TH1F* CorrectionOfEachAPVPairControlView =
0347 new TH1F("CorrectionOfEachAPVPairControlView", "CorrectionOfEachAPVPairControlView", 1, 0, 1);
0348 CorrectionOfEachAPVPairControlView->SetCanExtend(TH1::kAllAxes);
0349 std::ofstream APVPairTextOutput("apvpair_corrections.txt");
0350 APVPairTextOutput << "# MeanCharge = " << MeanCharge << std::endl;
0351 APVPairTextOutput << "# Nr. of APVPairs = " << NrOfApvPairs << std::endl;
0352 for (int ibin = 1; ibin <= ChargeOfEachAPVPair->GetNbinsX(); ibin++) {
0353 TString local_bin_label = ChargeOfEachAPVPair->GetXaxis()->GetBinLabel(ibin);
0354 double local_charge_over_path = ChargeOfEachAPVPair->GetBinContent(ibin);
0355 if (local_bin_label.Contains("ChargeAPVPair_") &&
0356 local_charge_over_path > 0.0000001) {
0357 uint32_t extracted_detid;
0358 std::istringstream read_label((local_bin_label(14, 9)).Data());
0359 read_label >> extracted_detid;
0360 unsigned short extracted_apvpairid;
0361 std::istringstream read_apvpair((local_bin_label(24, 1)).Data());
0362 read_apvpair >> extracted_apvpairid;
0363 double local_error_of_charge = ChargeOfEachAPVPair->GetBinError(ibin);
0364 double local_correction = -0.5;
0365 double local_error_correction = 0.;
0366 local_correction =
0367 MeanCharge / local_charge_over_path;
0368 local_error_correction = local_correction * local_error_of_charge / local_charge_over_path;
0369 if (local_error_correction > 1.8) {
0370 std::cout << "too large error " << local_error_correction << " for histogram " << local_bin_label << std::endl;
0371 }
0372 double nr_of_entries = EntriesApvPairs->GetBinContent(ibin);
0373 APVPairTextOutput << local_bin_label << " " << local_correction << " " << local_charge_over_path << " "
0374 << nr_of_entries << std::endl;
0375 CorrectionOfEachAPVPair->SetBinContent(ibin, local_correction);
0376 CorrectionOfEachAPVPair->SetBinError(ibin, local_error_correction);
0377 ((TH1F*)HlistOtherHistos->FindObject("APVPairCorrections"))->Fill(local_correction);
0378 DetId thedetId = DetId(extracted_detid);
0379 unsigned int generalized_layer = 0;
0380
0381 if (thedetId.subdetId() == StripSubdetector::TIB) {
0382 generalized_layer =
0383 10 * thedetId.subdetId() + tTopo_->tibLayer(thedetId.rawId()) + tTopo_->tibStereo(thedetId.rawId());
0384 if (tTopo_->tibLayer(thedetId.rawId()) == 2) {
0385 generalized_layer++;
0386 if (tTopo_->tibGlued(thedetId.rawId()))
0387 edm::LogError("ClusterMTCCFilter") << "WRONGGGG" << std::endl;
0388 }
0389 } else {
0390 generalized_layer = 10 * thedetId.subdetId();
0391 if (thedetId.subdetId() == StripSubdetector::TOB) {
0392 generalized_layer += tTopo_->tobLayer(thedetId.rawId());
0393 }
0394 }
0395 if (generalized_layer == 31) {
0396 ((TH1F*)HlistOtherHistos->FindObject("APVPairCorrectionsTIB1mono"))->Fill(local_correction);
0397 }
0398 if (generalized_layer == 32) {
0399 ((TH1F*)HlistOtherHistos->FindObject("APVPairCorrectionsTIB1stereo"))->Fill(local_correction);
0400 }
0401 if (generalized_layer == 33) {
0402 ((TH1F*)HlistOtherHistos->FindObject("APVPairCorrectionsTIB2"))->Fill(local_correction);
0403 }
0404 if (generalized_layer == 51) {
0405 ((TH1F*)HlistOtherHistos->FindObject("APVPairCorrectionsTOB1"))->Fill(local_correction);
0406 }
0407 if (generalized_layer == 52) {
0408 ((TH1F*)HlistOtherHistos->FindObject("APVPairCorrectionsTOB2"))->Fill(local_correction);
0409 }
0410
0411 const FedChannelConnection& fedchannelconnection =
0412 siStripDetCabling_->getConnection(extracted_detid, extracted_apvpairid);
0413 std::ostringstream local_key;
0414
0415 local_key << "fecCrate" << fedchannelconnection.fecCrate() << "_fecSlot" << fedchannelconnection.fecSlot()
0416 << "_fecRing" << fedchannelconnection.fecRing() << "_ccuAddr" << fedchannelconnection.ccuAddr()
0417 << "_ccuChan" << fedchannelconnection.ccuChan() << "_apvPair" << extracted_apvpairid;
0418 TString control_key = local_key.str();
0419 ChargeOfEachAPVPairControlView->Fill(control_key, local_charge_over_path);
0420 int ibin1 = ChargeOfEachAPVPairControlView->GetXaxis()->FindBin(control_key);
0421 ChargeOfEachAPVPairControlView->SetBinError(ibin1, local_error_of_charge);
0422 CorrectionOfEachAPVPairControlView->Fill(control_key, local_correction);
0423 int ibin2 = CorrectionOfEachAPVPairControlView->GetXaxis()->FindBin(control_key);
0424 CorrectionOfEachAPVPairControlView->SetBinError(ibin2, local_error_correction);
0425
0426 double module_thickness = moduleThickness(extracted_detid);
0427 if (fabs(module_thickness - 0.032) < 0.001)
0428 ModuleThickness->Fill(1);
0429 if (fabs(module_thickness - 0.05) < 0.001)
0430 ModuleThickness->Fill(2);
0431
0432 double module_width = moduleWidth(extracted_detid);
0433 if (fabs(module_width - 6.144) < 0.01)
0434 ModuleWidth->Fill(1);
0435 if (fabs(module_width - 7.14) < 0.01)
0436 ModuleWidth->Fill(2);
0437 if (fabs(module_width - 9.3696) < 0.01)
0438 ModuleWidth->Fill(3);
0439 if (fabs(module_width - 10.49) < 0.01)
0440 ModuleWidth->Fill(4);
0441 if (fabs(module_width - 12.03) < 0.01)
0442 ModuleWidth->Fill(5);
0443 }
0444 }
0445 HlistOtherHistos->Add(CorrectionOfEachAPVPair);
0446 ChargeOfEachAPVPairControlView->LabelsDeflate("X");
0447 CorrectionOfEachAPVPairControlView->LabelsDeflate("X");
0448 HlistOtherHistos->Add(ChargeOfEachAPVPairControlView);
0449 HlistOtherHistos->Add(CorrectionOfEachAPVPairControlView);
0450
0451
0452 if (outputHistogramsInRootFile) {
0453 TFile* outputfile = new TFile(outputFileName, "RECREATE");
0454 HlistAPVPairs->Write();
0455 HlistOtherHistos->Write();
0456 outputfile->Close();
0457 }
0458
0459 auto obj = std::make_unique<SiStripApvGain>();
0460
0461
0462
0463
0464
0465
0466
0467
0468
0469
0470
0471
0472
0473
0474
0475
0476
0477
0478
0479
0480
0481
0482
0483
0484
0485 return obj;
0486 }