File indexing completed on 2024-04-06 11:59:40
0001
0002
0003
0004 #include <memory>
0005 #include <iostream>
0006
0007 #include "FWCore/Framework/interface/Frameworkfwd.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/Framework/interface/MakerMacros.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012 #include "FWCore/Utilities/interface/Exception.h"
0013 #include "FWCore/Utilities/interface/EDGetToken.h"
0014
0015 #include "CondFormats/RunInfo/interface/RunInfo.h"
0016
0017 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0018 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0019 #include "Geometry/CommonTopologies/interface/StripTopology.h"
0020 #include "DataFormats/GeometrySurface/interface/TrapezoidalPlaneBounds.h"
0021 #include "DataFormats/GeometrySurface/interface/RectangularPlaneBounds.h"
0022
0023 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0024 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0025 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0026 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
0027 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0028 #include "Geometry/TrackerNumberingBuilder/interface/GeometricDet.h"
0029 #include "Geometry/CommonDetUnit/interface/TrackingGeometry.h"
0030
0031 #include "DataFormats/SiStripCluster/interface/SiStripClusterCollection.h"
0032
0033 #include "CalibFormats/SiStripObjects/interface/SiStripDetCabling.h"
0034 #include "CalibTracker/Records/interface/SiStripDetCablingRcd.h"
0035
0036 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
0037 #include "DataFormats/TrackReco/interface/Track.h"
0038 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0039 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h"
0040 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit1D.h"
0041 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
0042 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
0043 #include "DataFormats/DetId/interface/DetId.h"
0044 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0045 #include "DataFormats/TrackReco/interface/DeDxHit.h"
0046 #include "DataFormats/TrackReco/interface/TrackDeDxHits.h"
0047 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0048
0049 #include "CommonTools/ConditionDBWriter/interface/ConditionDBWriter.h"
0050 #include "CondFormats/SiStripObjects/interface/SiStripApvGain.h"
0051
0052 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0053 #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h"
0054
0055 #include "DQMServices/Core/interface/DQMStore.h"
0056 #include "FWCore/ServiceRegistry/interface/Service.h"
0057 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0058
0059 #include "CalibFormats/SiStripObjects/interface/SiStripGain.h"
0060 #include "CalibTracker/Records/interface/SiStripGainRcd.h"
0061
0062 #include "CalibFormats/SiStripObjects/interface/SiStripQuality.h"
0063 #include "CalibTracker/Records/interface/SiStripQualityRcd.h"
0064
0065 #include "TFile.h"
0066 #include "TObjString.h"
0067 #include "TString.h"
0068 #include "TH1F.h"
0069 #include "TH2S.h"
0070 #include "TProfile.h"
0071 #include "TF1.h"
0072 #include "TROOT.h"
0073 #include "TTree.h"
0074 #include "TChain.h"
0075
0076
0077 #include "CalibTracker/SiStripChannelGain/interface/APVGainStruct.h"
0078 #include "CalibTracker/SiStripChannelGain/interface/APVGainHelpers.h"
0079
0080 #include <unordered_map>
0081 #include <array>
0082
0083 using namespace edm;
0084 using namespace reco;
0085 using namespace std;
0086 using namespace APVGain;
0087
0088 class SiStripGainFromCalibTree : public ConditionDBWriter<SiStripApvGain> {
0089 public:
0090 typedef dqm::legacy::MonitorElement MonitorElement;
0091 typedef dqm::legacy::DQMStore DQMStore;
0092 explicit SiStripGainFromCalibTree(const edm::ParameterSet&);
0093 ~SiStripGainFromCalibTree() override;
0094
0095 private:
0096 void algoBeginRun(const edm::Run& run, const edm::EventSetup& iSetup) override;
0097 void algoEndRun(const edm::Run& run, const edm::EventSetup& iSetup) override;
0098 void algoBeginJob(const edm::EventSetup& iSetup) override;
0099 void algoEndJob() override;
0100 void algoAnalyze(const edm::Event&, const edm::EventSetup&) override;
0101
0102 int statCollectionFromMode(const char* tag) const;
0103 void bookDQMHistos(const char* dqm_dir, const char* tag);
0104
0105 bool isBFieldConsistentWithMode(const edm::EventSetup& iSetup) const;
0106 void swapBFieldMode(void);
0107
0108 void merge(TH2* A, TH2* B);
0109 void algoAnalyzeTheTree();
0110 void algoComputeMPVandGain();
0111 void processEvent();
0112
0113 void getPeakOfLandau(TH1* InputHisto, double* FitResults, double LowRange = 50, double HighRange = 5400);
0114 bool IsGoodLandauFit(double* FitResults);
0115 void storeOnTree(TFileService* tfs);
0116 void qualityMonitor();
0117 void MakeCalibrationMap();
0118 bool produceTagFilter();
0119
0120 template <typename T>
0121 inline edm::Handle<T> connect(const T*& ptr, edm::EDGetTokenT<T> token, const edm::Event& evt) {
0122 edm::Handle<T> handle;
0123 evt.getByToken(token, handle);
0124 ptr = handle.product();
0125 return handle;
0126 }
0127
0128 std::unique_ptr<SiStripApvGain> getNewObject() override;
0129
0130 TFileService* tfs;
0131 DQMStore* dbe;
0132 double MagFieldCurrentTh;
0133 double MinNrEntries;
0134 double MaxMPVError;
0135 double MaxChi2OverNDF;
0136 double MinTrackMomentum;
0137 double MaxTrackMomentum;
0138 double MinTrackEta;
0139 double MaxTrackEta;
0140 unsigned int MaxNrStrips;
0141 unsigned int MinTrackHits;
0142 double MaxTrackChiOverNdf;
0143 int MaxTrackingIteration;
0144 bool AllowSaturation;
0145 bool FirstSetOfConstants;
0146 bool Validation;
0147 bool OldGainRemoving;
0148 int CalibrationLevel;
0149
0150 bool saveSummary;
0151 bool useCalibration;
0152 bool m_harvestingMode;
0153 bool m_splitDQMstat;
0154 bool doChargeMonitorPerPlane;
0155 std::string m_calibrationMode;
0156 std::string m_calibrationPath;
0157 std::string m_DQMdir;
0158 std::vector<std::string> VChargeHisto;
0159
0160 double tagCondition_NClusters;
0161 double tagCondition_GoodFrac;
0162
0163 std::string AlgoMode;
0164 std::string OutputGains;
0165 vector<string> VInputFiles;
0166
0167
0168
0169 std::vector<string> dqm_tag_;
0170 std::string booked_dir_;
0171
0172 std::vector<MonitorElement*> Charge_Vs_Index;
0173 std::array<std::vector<APVGain::APVmon>, 7> Charge_1;
0174 std::array<std::vector<APVGain::APVmon>, 7> Charge_2;
0175 std::array<std::vector<APVGain::APVmon>, 7> Charge_3;
0176 std::array<std::vector<APVGain::APVmon>, 7> Charge_4;
0177
0178 std::vector<MonitorElement*> Charge_Vs_PathlengthTIB;
0179 std::vector<MonitorElement*> Charge_Vs_PathlengthTOB;
0180 std::vector<MonitorElement*> Charge_Vs_PathlengthTIDP;
0181 std::vector<MonitorElement*> Charge_Vs_PathlengthTIDM;
0182 std::vector<MonitorElement*> Charge_Vs_PathlengthTECP1;
0183 std::vector<MonitorElement*> Charge_Vs_PathlengthTECP2;
0184 std::vector<MonitorElement*> Charge_Vs_PathlengthTECM1;
0185 std::vector<MonitorElement*> Charge_Vs_PathlengthTECM2;
0186
0187
0188
0189
0190 MonitorElement* MPV_Vs_EtaTIB;
0191 MonitorElement* MPV_Vs_EtaTID;
0192 MonitorElement* MPV_Vs_EtaTOB;
0193 MonitorElement* MPV_Vs_EtaTEC;
0194 MonitorElement* MPV_Vs_EtaTECthin;
0195 MonitorElement* MPV_Vs_EtaTECthick;
0196
0197 MonitorElement* MPV_Vs_PhiTIB;
0198 MonitorElement* MPV_Vs_PhiTID;
0199 MonitorElement* MPV_Vs_PhiTOB;
0200 MonitorElement* MPV_Vs_PhiTEC;
0201 MonitorElement* MPV_Vs_PhiTECthin;
0202 MonitorElement* MPV_Vs_PhiTECthick;
0203
0204 MonitorElement* NoMPVmasked;
0205 MonitorElement* NoMPVfit;
0206
0207 MonitorElement* Gains;
0208 MonitorElement* MPVs;
0209 MonitorElement* MPVs320;
0210 MonitorElement* MPVs500;
0211 MonitorElement* MPVsTIB;
0212 MonitorElement* MPVsTID;
0213 MonitorElement* MPVsTIDP;
0214 MonitorElement* MPVsTIDM;
0215 MonitorElement* MPVsTOB;
0216 MonitorElement* MPVsTEC;
0217 MonitorElement* MPVsTECP;
0218 MonitorElement* MPVsTECM;
0219 MonitorElement* MPVsTECthin;
0220 MonitorElement* MPVsTECthick;
0221 MonitorElement* MPVsTECP1;
0222 MonitorElement* MPVsTECP2;
0223 MonitorElement* MPVsTECM1;
0224 MonitorElement* MPVsTECM2;
0225
0226 MonitorElement* MPVError;
0227 MonitorElement* MPVErrorVsMPV;
0228 MonitorElement* MPVErrorVsEta;
0229 MonitorElement* MPVErrorVsPhi;
0230 MonitorElement* MPVErrorVsN;
0231
0232 MonitorElement* DiffWRTPrevGainTIB;
0233 MonitorElement* DiffWRTPrevGainTID;
0234 MonitorElement* DiffWRTPrevGainTOB;
0235 MonitorElement* DiffWRTPrevGainTEC;
0236
0237 MonitorElement* GainVsPrevGainTIB;
0238 MonitorElement* GainVsPrevGainTID;
0239 MonitorElement* GainVsPrevGainTOB;
0240 MonitorElement* GainVsPrevGainTEC;
0241
0242 std::vector<APVGain::APVmon> newCharge;
0243
0244 unsigned int NEvent;
0245 unsigned int NTrack;
0246 unsigned int NClusterStrip;
0247 unsigned int NClusterPixel;
0248 int NStripAPVs;
0249 int NPixelDets;
0250 unsigned int SRun;
0251 unsigned int ERun;
0252 unsigned int GOOD;
0253 unsigned int BAD;
0254 unsigned int MASKED;
0255
0256
0257
0258
0259 unsigned int eventnumber = 0;
0260 unsigned int runnumber = 0;
0261 const std::vector<bool>* TrigTech = nullptr;
0262 edm::EDGetTokenT<std::vector<bool>> TrigTech_token_;
0263
0264
0265 const std::vector<double>* trackchi2ndof = nullptr;
0266 edm::EDGetTokenT<std::vector<double>> trackchi2ndof_token_;
0267 const std::vector<float>* trackp = nullptr;
0268 edm::EDGetTokenT<std::vector<float>> trackp_token_;
0269 const std::vector<float>* trackpt = nullptr;
0270 edm::EDGetTokenT<std::vector<float>> trackpt_token_;
0271 const std::vector<double>* tracketa = nullptr;
0272 edm::EDGetTokenT<std::vector<double>> tracketa_token_;
0273 const std::vector<double>* trackphi = nullptr;
0274 edm::EDGetTokenT<std::vector<double>> trackphi_token_;
0275 const std::vector<unsigned int>* trackhitsvalid = nullptr;
0276 edm::EDGetTokenT<std::vector<unsigned int>> trackhitsvalid_token_;
0277 const std::vector<int>* trackalgo = nullptr;
0278 edm::EDGetTokenT<std::vector<int>> trackalgo_token_;
0279
0280
0281 const std::vector<int>* trackindex = nullptr;
0282 edm::EDGetTokenT<std::vector<int>> trackindex_token_;
0283 const std::vector<unsigned int>* rawid = nullptr;
0284 edm::EDGetTokenT<std::vector<unsigned int>> rawid_token_;
0285 const std::vector<double>* localdirx = nullptr;
0286 edm::EDGetTokenT<std::vector<double>> localdirx_token_;
0287 const std::vector<double>* localdiry = nullptr;
0288 edm::EDGetTokenT<std::vector<double>> localdiry_token_;
0289 const std::vector<double>* localdirz = nullptr;
0290 edm::EDGetTokenT<std::vector<double>> localdirz_token_;
0291 const std::vector<unsigned short>* firststrip = nullptr;
0292 edm::EDGetTokenT<std::vector<unsigned short>> firststrip_token_;
0293 const std::vector<unsigned short>* nstrips = nullptr;
0294 edm::EDGetTokenT<std::vector<unsigned short>> nstrips_token_;
0295 const std::vector<bool>* saturation = nullptr;
0296 edm::EDGetTokenT<std::vector<bool>> saturation_token_;
0297 const std::vector<bool>* overlapping = nullptr;
0298 edm::EDGetTokenT<std::vector<bool>> overlapping_token_;
0299 const std::vector<bool>* farfromedge = nullptr;
0300 edm::EDGetTokenT<std::vector<bool>> farfromedge_token_;
0301 const std::vector<unsigned int>* charge = nullptr;
0302 edm::EDGetTokenT<std::vector<unsigned int>> charge_token_;
0303 const std::vector<double>* path = nullptr;
0304 edm::EDGetTokenT<std::vector<double>> path_token_;
0305 const std::vector<double>* chargeoverpath = nullptr;
0306 edm::EDGetTokenT<std::vector<double>> chargeoverpath_token_;
0307 const std::vector<unsigned char>* amplitude = nullptr;
0308 edm::EDGetTokenT<std::vector<unsigned char>> amplitude_token_;
0309 const std::vector<double>* gainused = nullptr;
0310 edm::EDGetTokenT<std::vector<double>> gainused_token_;
0311 const std::vector<double>* gainusedTick = nullptr;
0312 edm::EDGetTokenT<std::vector<double>> gainusedTick_token_;
0313
0314 string EventPrefix_;
0315 string EventSuffix_;
0316 string TrackPrefix_;
0317 string TrackSuffix_;
0318 string CalibPrefix_;
0319 string CalibSuffix_;
0320
0321 private:
0322 std::vector<stAPVGain*> APVsCollOrdered;
0323 std::unordered_map<unsigned int, stAPVGain*> APVsColl;
0324
0325 edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> tTopoToken_;
0326 edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> tkGeomToken_;
0327 edm::ESGetToken<RunInfo, RunInfoRcd> runInfoToken_;
0328 edm::ESGetToken<SiStripGain, SiStripGainRcd> gainToken_;
0329 edm::ESGetToken<SiStripQuality, SiStripQualityRcd> qualityToken_;
0330 const TrackerTopology* tTopo_ = nullptr;
0331 };
0332
0333 inline int SiStripGainFromCalibTree::statCollectionFromMode(const char* tag) const {
0334 std::vector<string>::const_iterator it = dqm_tag_.begin();
0335 while (it != dqm_tag_.end()) {
0336 if (*it == std::string(tag))
0337 return it - dqm_tag_.begin();
0338 it++;
0339 }
0340
0341 if (std::string(tag).empty())
0342 return 0;
0343
0344 return None;
0345 }
0346
0347 void SiStripGainFromCalibTree::merge(TH2* A, TH2* B) {
0348 if (A->GetNbinsX() == B->GetNbinsX()) {
0349 A->Add(B);
0350 } else {
0351 for (int x = 0; x <= B->GetNbinsX() + 1; x++) {
0352 for (int y = 0; y <= B->GetNbinsY() + 1; y++) {
0353 A->SetBinContent(x, y, A->GetBinContent(x, y) + B->GetBinContent(x, y));
0354 }
0355 }
0356 }
0357 }
0358
0359 SiStripGainFromCalibTree::SiStripGainFromCalibTree(const edm::ParameterSet& iConfig)
0360 : ConditionDBWriter<SiStripApvGain>(iConfig) {
0361 usesResource(TFileService::kSharedResource);
0362 OutputGains = iConfig.getParameter<std::string>("OutputGains");
0363 AlgoMode = iConfig.getUntrackedParameter<std::string>("AlgoMode", "CalibTree");
0364 MagFieldCurrentTh = iConfig.getUntrackedParameter<double>("MagFieldCurrentTh", 2000.);
0365 MinNrEntries = iConfig.getUntrackedParameter<double>("minNrEntries", 20);
0366 MaxMPVError = iConfig.getUntrackedParameter<double>("maxMPVError", 500.0);
0367 MaxChi2OverNDF = iConfig.getUntrackedParameter<double>("maxChi2OverNDF", 5.0);
0368 MinTrackMomentum = iConfig.getUntrackedParameter<double>("minTrackMomentum", 3.0);
0369 MaxTrackMomentum = iConfig.getUntrackedParameter<double>("maxTrackMomentum", 99999.0);
0370 MinTrackEta = iConfig.getUntrackedParameter<double>("minTrackEta", -5.0);
0371 MaxTrackEta = iConfig.getUntrackedParameter<double>("maxTrackEta", 5.0);
0372 MaxNrStrips = iConfig.getUntrackedParameter<unsigned>("maxNrStrips", 2);
0373 MinTrackHits = iConfig.getUntrackedParameter<unsigned>("MinTrackHits", 8);
0374 MaxTrackChiOverNdf = iConfig.getUntrackedParameter<double>("MaxTrackChiOverNdf", 3);
0375 MaxTrackingIteration = iConfig.getUntrackedParameter<int>("MaxTrackingIteration", 7);
0376 AllowSaturation = iConfig.getUntrackedParameter<bool>("AllowSaturation", false);
0377 FirstSetOfConstants = iConfig.getUntrackedParameter<bool>("FirstSetOfConstants", true);
0378 Validation = iConfig.getUntrackedParameter<bool>("Validation", false);
0379 OldGainRemoving = iConfig.getUntrackedParameter<bool>("OldGainRemoving", false);
0380
0381 CalibrationLevel = iConfig.getUntrackedParameter<int>("CalibrationLevel", 0);
0382 VInputFiles = iConfig.getUntrackedParameter<vector<string>>("InputFiles");
0383 VChargeHisto = iConfig.getUntrackedParameter<vector<string>>("ChargeHisto");
0384
0385 useCalibration = iConfig.getUntrackedParameter<bool>("UseCalibration", false);
0386 m_harvestingMode = iConfig.getUntrackedParameter<bool>("harvestingMode", false);
0387 m_splitDQMstat = iConfig.getUntrackedParameter<bool>("splitDQMstat", false);
0388 m_calibrationMode = iConfig.getUntrackedParameter<string>("calibrationMode", "StdBunch");
0389 m_calibrationPath = iConfig.getUntrackedParameter<string>("calibrationPath");
0390 m_DQMdir = iConfig.getUntrackedParameter<string>("DQMdir", "AlCaReco/SiStripGains");
0391
0392 tagCondition_NClusters = iConfig.getUntrackedParameter<double>("NClustersForTagProd", 2E8);
0393 tagCondition_GoodFrac = iConfig.getUntrackedParameter<double>("GoodFracForTagProd", 0.95);
0394
0395 saveSummary = iConfig.getUntrackedParameter<bool>("saveSummary", false);
0396
0397 doChargeMonitorPerPlane = iConfig.getUntrackedParameter<bool>("doChargeMonitorPerPlane", false);
0398
0399
0400 dbe = edm::Service<DQMStore>().operator->();
0401
0402
0403 dqm_tag_.reserve(7);
0404 dqm_tag_.clear();
0405 dqm_tag_.push_back("StdBunch");
0406 dqm_tag_.push_back("StdBunch0T");
0407 dqm_tag_.push_back("AagBunch");
0408 dqm_tag_.push_back("AagBunch0T");
0409 dqm_tag_.push_back("IsoMuon");
0410 dqm_tag_.push_back("IsoMuon0T");
0411 dqm_tag_.push_back("Harvest");
0412
0413 Charge_Vs_Index.insert(Charge_Vs_Index.begin(), dqm_tag_.size(), nullptr);
0414
0415 Charge_Vs_PathlengthTIB.insert(Charge_Vs_PathlengthTIB.begin(), dqm_tag_.size(), nullptr);
0416 Charge_Vs_PathlengthTOB.insert(Charge_Vs_PathlengthTOB.begin(), dqm_tag_.size(), nullptr);
0417 Charge_Vs_PathlengthTIDP.insert(Charge_Vs_PathlengthTIDP.begin(), dqm_tag_.size(), nullptr);
0418 Charge_Vs_PathlengthTIDM.insert(Charge_Vs_PathlengthTIDM.begin(), dqm_tag_.size(), nullptr);
0419 Charge_Vs_PathlengthTECP1.insert(Charge_Vs_PathlengthTECP1.begin(), dqm_tag_.size(), nullptr);
0420 Charge_Vs_PathlengthTECP2.insert(Charge_Vs_PathlengthTECP2.begin(), dqm_tag_.size(), nullptr);
0421 Charge_Vs_PathlengthTECM1.insert(Charge_Vs_PathlengthTECM1.begin(), dqm_tag_.size(), nullptr);
0422 Charge_Vs_PathlengthTECM2.insert(Charge_Vs_PathlengthTECM2.begin(), dqm_tag_.size(), nullptr);
0423
0424
0425 edm::ParameterSet swhallowgain_pset = iConfig.getUntrackedParameter<edm::ParameterSet>("gain");
0426
0427 string label = swhallowgain_pset.getUntrackedParameter<string>("label");
0428 CalibPrefix_ = swhallowgain_pset.getUntrackedParameter<string>("prefix");
0429 CalibSuffix_ = swhallowgain_pset.getUntrackedParameter<string>("suffix");
0430
0431 trackindex_token_ = consumes<std::vector<int>>(edm::InputTag(label, CalibPrefix_ + "trackindex" + CalibSuffix_));
0432 rawid_token_ = consumes<std::vector<unsigned int>>(edm::InputTag(label, CalibPrefix_ + "rawid" + CalibSuffix_));
0433 localdirx_token_ = consumes<std::vector<double>>(edm::InputTag(label, CalibPrefix_ + "localdirx" + CalibSuffix_));
0434 localdiry_token_ = consumes<std::vector<double>>(edm::InputTag(label, CalibPrefix_ + "localdiry" + CalibSuffix_));
0435 localdirz_token_ = consumes<std::vector<double>>(edm::InputTag(label, CalibPrefix_ + "localdirz" + CalibSuffix_));
0436 firststrip_token_ =
0437 consumes<std::vector<unsigned short>>(edm::InputTag(label, CalibPrefix_ + "firststrip" + CalibSuffix_));
0438 nstrips_token_ = consumes<std::vector<unsigned short>>(edm::InputTag(label, CalibPrefix_ + "nstrips" + CalibSuffix_));
0439 saturation_token_ = consumes<std::vector<bool>>(edm::InputTag(label, CalibPrefix_ + "saturation" + CalibSuffix_));
0440 overlapping_token_ = consumes<std::vector<bool>>(edm::InputTag(label, CalibPrefix_ + "overlapping" + CalibSuffix_));
0441 farfromedge_token_ = consumes<std::vector<bool>>(edm::InputTag(label, CalibPrefix_ + "farfromedge" + CalibSuffix_));
0442 charge_token_ = consumes<std::vector<unsigned int>>(edm::InputTag(label, CalibPrefix_ + "charge" + CalibSuffix_));
0443 path_token_ = consumes<std::vector<double>>(edm::InputTag(label, CalibPrefix_ + "path" + CalibSuffix_));
0444 chargeoverpath_token_ =
0445 consumes<std::vector<double>>(edm::InputTag(label, CalibPrefix_ + "chargeoverpath" + CalibSuffix_));
0446 amplitude_token_ =
0447 consumes<std::vector<unsigned char>>(edm::InputTag(label, CalibPrefix_ + "amplitude" + CalibSuffix_));
0448 gainused_token_ = consumes<std::vector<double>>(edm::InputTag(label, CalibPrefix_ + "gainused" + CalibSuffix_));
0449 gainusedTick_token_ =
0450 consumes<std::vector<double>>(edm::InputTag(label, CalibPrefix_ + "gainusedTick" + CalibSuffix_));
0451
0452 edm::ParameterSet evtinfo_pset = iConfig.getUntrackedParameter<edm::ParameterSet>("evtinfo");
0453 label = evtinfo_pset.getUntrackedParameter<string>("label");
0454 EventPrefix_ = evtinfo_pset.getUntrackedParameter<string>("prefix");
0455 EventSuffix_ = evtinfo_pset.getUntrackedParameter<string>("suffix");
0456 TrigTech_token_ = consumes<std::vector<bool>>(edm::InputTag(label, EventPrefix_ + "TrigTech" + EventSuffix_));
0457
0458 edm::ParameterSet track_pset = iConfig.getUntrackedParameter<edm::ParameterSet>("tracks");
0459 label = track_pset.getUntrackedParameter<string>("label");
0460 TrackPrefix_ = track_pset.getUntrackedParameter<string>("prefix");
0461 TrackSuffix_ = track_pset.getUntrackedParameter<string>("suffix");
0462
0463 trackchi2ndof_token_ = consumes<std::vector<double>>(edm::InputTag(label, TrackPrefix_ + "chi2ndof" + TrackSuffix_));
0464 trackp_token_ = consumes<std::vector<float>>(edm::InputTag(label, TrackPrefix_ + "momentum" + TrackSuffix_));
0465 trackpt_token_ = consumes<std::vector<float>>(edm::InputTag(label, TrackPrefix_ + "pt" + TrackSuffix_));
0466 tracketa_token_ = consumes<std::vector<double>>(edm::InputTag(label, TrackPrefix_ + "eta" + TrackSuffix_));
0467 trackphi_token_ = consumes<std::vector<double>>(edm::InputTag(label, TrackPrefix_ + "phi" + TrackSuffix_));
0468 trackhitsvalid_token_ =
0469 consumes<std::vector<unsigned int>>(edm::InputTag(label, TrackPrefix_ + "hitsvalid" + TrackSuffix_));
0470 trackalgo_token_ = consumes<std::vector<int>>(edm::InputTag(label, TrackPrefix_ + "algo" + TrackSuffix_));
0471
0472 tTopoToken_ = esConsumes<edm::Transition::BeginRun>();
0473 tkGeomToken_ = esConsumes<edm::Transition::BeginRun>();
0474 runInfoToken_ = esConsumes<edm::Transition::BeginRun>();
0475 gainToken_ = esConsumes<edm::Transition::BeginRun>();
0476 qualityToken_ = esConsumes<edm::Transition::BeginRun>();
0477 }
0478
0479 void SiStripGainFromCalibTree::bookDQMHistos(const char* dqm_dir, const char* tag) {
0480 edm::LogInfo("SiStripGainFromCalibTree")
0481 << "Setting " << dqm_dir << "in DQM and booking histograms for tag " << tag << std::endl;
0482
0483 if (strcmp(booked_dir_.c_str(), dqm_dir) != 0) {
0484 booked_dir_ = dqm_dir;
0485 dbe->setCurrentFolder(dqm_dir);
0486 }
0487
0488 std::string stag(tag);
0489 if (!stag.empty() && stag[0] != '_')
0490 stag.insert(0, 1, '_');
0491
0492 std::string cvi = std::string("Charge_Vs_Index") + stag;
0493
0494 std::string cvpTIB = std::string("Charge_Vs_PathlengthTIB") + stag;
0495 std::string cvpTOB = std::string("Charge_Vs_PathlengthTOB") + stag;
0496 std::string cvpTIDP = std::string("Charge_Vs_PathlengthTIDP") + stag;
0497 std::string cvpTIDM = std::string("Charge_Vs_PathlengthTIDM") + stag;
0498 std::string cvpTECP1 = std::string("Charge_Vs_PathlengthTECP1") + stag;
0499 std::string cvpTECP2 = std::string("Charge_Vs_PathlengthTECP2") + stag;
0500 std::string cvpTECM1 = std::string("Charge_Vs_PathlengthTECM1") + stag;
0501 std::string cvpTECM2 = std::string("Charge_Vs_PathlengthTECM2") + stag;
0502
0503 int elepos = (m_harvestingMode && AlgoMode == "PCL") ? Harvest : statCollectionFromMode(tag);
0504
0505
0506
0507
0508
0509
0510
0511
0512
0513 std::vector<float> binXarray;
0514 binXarray.reserve(NStripAPVs + 1);
0515 for (int a = 0; a <= NStripAPVs; a++) {
0516 binXarray.push_back((float)a);
0517 }
0518
0519 std::array<float, 688> binYarray;
0520 double p0 = 5.445;
0521 double p1 = 0.002113;
0522 double p2 = 69.01576;
0523 double y = 0.;
0524 for (int b = 0; b < 687; b++) {
0525 binYarray[b] = y;
0526 if (y <= 902.)
0527 y = y + 2.;
0528 else
0529 y = (p0 - log(exp(p0 - p1 * y) - p2 * p1)) / p1;
0530 }
0531 binYarray[687] = 4000.;
0532
0533 Charge_Vs_Index[elepos] = dbe->book2S(cvi.c_str(), cvi.c_str(), NStripAPVs, &binXarray[0], 687, binYarray.data());
0534
0535 Charge_Vs_PathlengthTIB[elepos] = dbe->book2S(cvpTIB.c_str(), cvpTIB.c_str(), 20, 0.3, 1.3, 250, 0, 2000);
0536 Charge_Vs_PathlengthTOB[elepos] = dbe->book2S(cvpTOB.c_str(), cvpTOB.c_str(), 20, 0.3, 1.3, 250, 0, 2000);
0537 Charge_Vs_PathlengthTIDP[elepos] = dbe->book2S(cvpTIDP.c_str(), cvpTIDP.c_str(), 20, 0.3, 1.3, 250, 0, 2000);
0538 Charge_Vs_PathlengthTIDM[elepos] = dbe->book2S(cvpTIDM.c_str(), cvpTIDM.c_str(), 20, 0.3, 1.3, 250, 0, 2000);
0539 Charge_Vs_PathlengthTECP1[elepos] = dbe->book2S(cvpTECP1.c_str(), cvpTECP1.c_str(), 20, 0.3, 1.3, 250, 0, 2000);
0540 Charge_Vs_PathlengthTECP2[elepos] = dbe->book2S(cvpTECP2.c_str(), cvpTECP2.c_str(), 20, 0.3, 1.3, 250, 0, 2000);
0541 Charge_Vs_PathlengthTECM1[elepos] = dbe->book2S(cvpTECM1.c_str(), cvpTECM1.c_str(), 20, 0.3, 1.3, 250, 0, 2000);
0542 Charge_Vs_PathlengthTECM2[elepos] = dbe->book2S(cvpTECM2.c_str(), cvpTECM2.c_str(), 20, 0.3, 1.3, 250, 0, 2000);
0543
0544
0545 std::vector<std::pair<std::string, std::string>> hnames =
0546 APVGain::monHnames(VChargeHisto, doChargeMonitorPerPlane, "");
0547 for (unsigned int i = 0; i < hnames.size(); i++) {
0548 std::string htag = (hnames[i]).first + stag;
0549 MonitorElement* monitor = dbe->book1DD(htag.c_str(), (hnames[i]).second.c_str(), 100, 0., 1000.);
0550 int thick = APVGain::thickness((hnames[i]).first);
0551 int id = APVGain::subdetectorId((hnames[i]).first);
0552 int side = APVGain::subdetectorSide((hnames[i]).first);
0553 int plane = APVGain::subdetectorPlane((hnames[i]).first);
0554 Charge_1[elepos].push_back(APVGain::APVmon(thick, id, side, plane, monitor));
0555 }
0556
0557 hnames = APVGain::monHnames(VChargeHisto, doChargeMonitorPerPlane, "woG2");
0558 for (unsigned int i = 0; i < hnames.size(); i++) {
0559 std::string htag = (hnames[i]).first + stag;
0560 MonitorElement* monitor = dbe->book1DD(htag.c_str(), (hnames[i]).second.c_str(), 100, 0., 1000.);
0561 int thick = APVGain::thickness((hnames[i]).first);
0562 int id = APVGain::subdetectorId((hnames[i]).first);
0563 int side = APVGain::subdetectorSide((hnames[i]).first);
0564 int plane = APVGain::subdetectorPlane((hnames[i]).first);
0565 Charge_2[elepos].push_back(APVGain::APVmon(thick, id, side, plane, monitor));
0566 }
0567
0568 hnames = APVGain::monHnames(VChargeHisto, doChargeMonitorPerPlane, "woG1");
0569 for (unsigned int i = 0; i < hnames.size(); i++) {
0570 std::string htag = (hnames[i]).first + stag;
0571 MonitorElement* monitor = dbe->book1DD(htag.c_str(), (hnames[i]).second.c_str(), 100, 0., 1000.);
0572 int thick = APVGain::thickness((hnames[i]).first);
0573 int id = APVGain::subdetectorId((hnames[i]).first);
0574 int side = APVGain::subdetectorSide((hnames[i]).first);
0575 int plane = APVGain::subdetectorPlane((hnames[i]).first);
0576 Charge_3[elepos].push_back(APVGain::APVmon(thick, id, side, plane, monitor));
0577 }
0578
0579 hnames = APVGain::monHnames(VChargeHisto, doChargeMonitorPerPlane, "woG1G2");
0580 for (unsigned int i = 0; i < hnames.size(); i++) {
0581 std::string htag = (hnames[i]).first + stag;
0582 MonitorElement* monitor = dbe->book1DD(htag.c_str(), (hnames[i]).second.c_str(), 100, 0., 1000.);
0583 int thick = APVGain::thickness((hnames[i]).first);
0584 int id = APVGain::subdetectorId((hnames[i]).first);
0585 int side = APVGain::subdetectorSide((hnames[i]).first);
0586 int plane = APVGain::subdetectorPlane((hnames[i]).first);
0587 Charge_4[elepos].push_back(APVGain::APVmon(thick, id, side, plane, monitor));
0588 }
0589
0590
0591 if (m_harvestingMode) {
0592 int MPVbin = 300;
0593 float MPVmin = 0.;
0594 float MPVmax = 600.;
0595
0596 MPV_Vs_EtaTIB = dbe->book2DD("MPV_vs_EtaTIB", "MPV vs Eta TIB", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
0597 MPV_Vs_EtaTID = dbe->book2DD("MPV_vs_EtaTID", "MPV vs Eta TID", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
0598 MPV_Vs_EtaTOB = dbe->book2DD("MPV_vs_EtaTOB", "MPV vs Eta TOB", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
0599 MPV_Vs_EtaTEC = dbe->book2DD("MPV_vs_EtaTEC", "MPV vs Eta TEC", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
0600 MPV_Vs_EtaTECthin = dbe->book2DD("MPV_vs_EtaTEC1", "MPV vs Eta TEC-thin", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
0601 MPV_Vs_EtaTECthick = dbe->book2DD("MPV_vs_EtaTEC2", "MPV vs Eta TEC-thick", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
0602
0603 MPV_Vs_PhiTIB = dbe->book2DD("MPV_vs_PhiTIB", "MPV vs Phi TIB", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
0604 MPV_Vs_PhiTID = dbe->book2DD("MPV_vs_PhiTID", "MPV vs Phi TID", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
0605 MPV_Vs_PhiTOB = dbe->book2DD("MPV_vs_PhiTOB", "MPV vs Phi TOB", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
0606 MPV_Vs_PhiTEC = dbe->book2DD("MPV_vs_PhiTEC", "MPV vs Phi TEC", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
0607 MPV_Vs_PhiTECthin = dbe->book2DD("MPV_vs_PhiTEC1", "MPV vs Phi TEC-thin", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
0608 MPV_Vs_PhiTECthick = dbe->book2DD("MPV_vs_PhiTEC2", "MPV vs Phi TEC-thick", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
0609
0610 NoMPVfit = dbe->book2DD("NoMPVfit", "Modules with bad Landau Fit", 350, -350, 350, 240, 0, 120);
0611 NoMPVmasked = dbe->book2DD("NoMPVmasked", "Masked Modules", 350, -350, 350, 240, 0, 120);
0612
0613 Gains = dbe->book1DD("Gains", "Gains", 300, 0, 2);
0614 MPVs = dbe->book1DD("MPVs", "MPVs", MPVbin, MPVmin, MPVmax);
0615 MPVs320 = dbe->book1DD("MPV_320", "MPV 320 thickness", MPVbin, MPVmin, MPVmax);
0616 MPVs500 = dbe->book1DD("MPV_500", "MPV 500 thickness", MPVbin, MPVmin, MPVmax);
0617 MPVsTIB = dbe->book1DD("MPV_TIB", "MPV TIB", MPVbin, MPVmin, MPVmax);
0618 MPVsTID = dbe->book1DD("MPV_TID", "MPV TID", MPVbin, MPVmin, MPVmax);
0619 MPVsTIDP = dbe->book1DD("MPV_TIDP", "MPV TIDP", MPVbin, MPVmin, MPVmax);
0620 MPVsTIDM = dbe->book1DD("MPV_TIDM", "MPV TIDM", MPVbin, MPVmin, MPVmax);
0621 MPVsTOB = dbe->book1DD("MPV_TOB", "MPV TOB", MPVbin, MPVmin, MPVmax);
0622 MPVsTEC = dbe->book1DD("MPV_TEC", "MPV TEC", MPVbin, MPVmin, MPVmax);
0623 MPVsTECP = dbe->book1DD("MPV_TECP", "MPV TECP", MPVbin, MPVmin, MPVmax);
0624 MPVsTECM = dbe->book1DD("MPV_TECM", "MPV TECM", MPVbin, MPVmin, MPVmax);
0625 MPVsTECthin = dbe->book1DD("MPV_TEC1", "MPV TEC1", MPVbin, MPVmin, MPVmax);
0626 MPVsTECthick = dbe->book1DD("MPV_TEC2", "MPV TEC2", MPVbin, MPVmin, MPVmax);
0627 MPVsTECP1 = dbe->book1DD("MPV_TECP1", "MPV TECP1", MPVbin, MPVmin, MPVmax);
0628 MPVsTECP2 = dbe->book1DD("MPV_TECP2", "MPV TECP2", MPVbin, MPVmin, MPVmax);
0629 MPVsTECM1 = dbe->book1DD("MPV_TECM1", "MPV TECM1", MPVbin, MPVmin, MPVmax);
0630 MPVsTECM2 = dbe->book1DD("MPV_TECM2", "MPV TECM2", MPVbin, MPVmin, MPVmax);
0631
0632 MPVError = dbe->book1DD("MPVError", "MPV Error", 150, 0, 150);
0633 MPVErrorVsMPV = dbe->book2DD("MPVErrorVsMPV", "MPV Error vs MPV", 300, 0, 600, 150, 0, 150);
0634 MPVErrorVsEta = dbe->book2DD("MPVErrorVsEta", "MPV Error vs Eta", 50, -3.0, 3.0, 150, 0, 150);
0635 MPVErrorVsPhi = dbe->book2DD("MPVErrorVsPhi", "MPV Error vs Phi", 50, -3.4, 3.4, 150, 0, 150);
0636 MPVErrorVsN = dbe->book2DD("MPVErrorVsN", "MPV Error vs N", 500, 0, 1000, 150, 0, 150);
0637
0638 DiffWRTPrevGainTIB = dbe->book1DD("DiffWRTPrevGainTIB", "Diff w.r.t. PrevGain TIB", 250, 0, 2);
0639 DiffWRTPrevGainTID = dbe->book1DD("DiffWRTPrevGainTID", "Diff w.r.t. PrevGain TID", 250, 0, 2);
0640 DiffWRTPrevGainTOB = dbe->book1DD("DiffWRTPrevGainTOB", "Diff w.r.t. PrevGain TOB", 250, 0, 2);
0641 DiffWRTPrevGainTEC = dbe->book1DD("DiffWRTPrevGainTEC", "Diff w.r.t. PrevGain TEC", 250, 0, 2);
0642
0643 GainVsPrevGainTIB = dbe->book2DD("GainVsPrevGainTIB", "Gain vs PrevGain TIB", 100, 0, 2, 100, 0, 2);
0644 GainVsPrevGainTID = dbe->book2DD("GainVsPrevGainTID", "Gain vs PrevGain TID", 100, 0, 2, 100, 0, 2);
0645 GainVsPrevGainTOB = dbe->book2DD("GainVsPrevGainTOB", "Gain vs PrevGain TOB", 100, 0, 2, 100, 0, 2);
0646 GainVsPrevGainTEC = dbe->book2DD("GainVsPrevGainTEC", "Gain vs PrevGain TEC", 100, 0, 2, 100, 0, 2);
0647
0648 std::vector<std::pair<std::string, std::string>> hnames =
0649 APVGain::monHnames(VChargeHisto, doChargeMonitorPerPlane, "newG2");
0650 for (unsigned int i = 0; i < hnames.size(); i++) {
0651 MonitorElement* monitor = dbe->book1DD((hnames[i]).first.c_str(), (hnames[i]).second.c_str(), 100, 0., 1000.);
0652 int thick = APVGain::thickness((hnames[i]).first);
0653 int id = APVGain::subdetectorId((hnames[i]).first);
0654 int side = APVGain::subdetectorSide((hnames[i]).first);
0655 int plane = APVGain::subdetectorPlane((hnames[i]).first);
0656 newCharge.push_back(APVGain::APVmon(thick, id, side, plane, monitor));
0657 }
0658 }
0659 }
0660
0661 void SiStripGainFromCalibTree::algoBeginJob(const edm::EventSetup& iSetup) {
0662 edm::LogInfo("SiStripGainFromCalibTree") << "AlgoMode : " << AlgoMode << "\n"
0663 << "CalibrationMode : " << m_calibrationMode << "\n"
0664 << "HarvestingMode : " << m_harvestingMode << std::endl;
0665
0666 if (AlgoMode != "PCL" or m_harvestingMode) {
0667 const char* dqm_dir = "AlCaReco/SiStripGainsHarvesting/";
0668 this->bookDQMHistos(dqm_dir, dqm_tag_[statCollectionFromMode(m_calibrationMode.c_str())].c_str());
0669 } else {
0670
0671 if (!isBFieldConsistentWithMode(iSetup)) {
0672 string prevMode = m_calibrationMode;
0673 swapBFieldMode();
0674 edm::LogInfo("SiStripGainFromCalibTree") << "Switching calibration mode for endorsing BField status: " << prevMode
0675 << " ==> " << m_calibrationMode << std::endl;
0676 }
0677 std::string dqm_dir = m_DQMdir + ((m_splitDQMstat) ? m_calibrationMode : "") + "/";
0678 int elem = statCollectionFromMode(m_calibrationMode.c_str());
0679 this->bookDQMHistos(dqm_dir.c_str(), dqm_tag_[elem].c_str());
0680 }
0681
0682 tTopo_ = &iSetup.getData(tTopoToken_);
0683
0684 auto const& Det = iSetup.getData(tkGeomToken_).dets();
0685
0686 NPixelDets = 0;
0687 NStripAPVs = 0;
0688 unsigned int Index = 0;
0689 for (unsigned int i = 0; i < Det.size(); i++) {
0690 DetId Detid = Det[i]->geographicalId();
0691 int SubDet = Detid.subdetId();
0692
0693 if (SubDet == StripSubdetector::TIB || SubDet == StripSubdetector::TID || SubDet == StripSubdetector::TOB ||
0694 SubDet == StripSubdetector::TEC) {
0695 auto DetUnit = dynamic_cast<const StripGeomDetUnit*>(Det[i]);
0696 if (!DetUnit)
0697 continue;
0698
0699 const StripTopology& Topo = DetUnit->specificTopology();
0700 unsigned int NAPV = Topo.nstrips() / 128;
0701
0702 for (unsigned int j = 0; j < NAPV; j++) {
0703 stAPVGain* APV = new stAPVGain;
0704 APV->Index = Index;
0705 APV->Bin = -1;
0706 APV->DetId = Detid.rawId();
0707 APV->APVId = j;
0708 APV->SubDet = SubDet;
0709 APV->FitMPV = -1;
0710 APV->FitMPVErr = -1;
0711 APV->FitWidth = -1;
0712 APV->FitWidthErr = -1;
0713 APV->FitChi2 = -1;
0714 APV->FitNorm = -1;
0715 APV->Gain = -1;
0716 APV->PreviousGain = 1;
0717 APV->PreviousGainTick = 1;
0718 APV->x = DetUnit->position().basicVector().x();
0719 APV->y = DetUnit->position().basicVector().y();
0720 APV->z = DetUnit->position().basicVector().z();
0721 APV->Eta = DetUnit->position().basicVector().eta();
0722 APV->Phi = DetUnit->position().basicVector().phi();
0723 APV->R = DetUnit->position().basicVector().transverse();
0724 APV->Thickness = DetUnit->surface().bounds().thickness();
0725 APV->NEntries = 0;
0726 APV->isMasked = false;
0727
0728 APVsCollOrdered.push_back(APV);
0729 APVsColl[(APV->DetId << 4) | APV->APVId] = APV;
0730 Index++;
0731 NStripAPVs++;
0732 }
0733 }
0734 }
0735
0736 for (unsigned int i = 0; i < Det.size();
0737 i++) {
0738 DetId Detid = Det[i]->geographicalId();
0739 int SubDet = Detid.subdetId();
0740 if (SubDet == PixelSubdetector::PixelBarrel || SubDet == PixelSubdetector::PixelEndcap) {
0741 auto DetUnit = dynamic_cast<const PixelGeomDetUnit*>(Det[i]);
0742 if (!DetUnit)
0743 continue;
0744
0745 const PixelTopology& Topo = DetUnit->specificTopology();
0746 unsigned int NROCRow = Topo.nrows() / (80.);
0747 unsigned int NROCCol = Topo.ncolumns() / (52.);
0748
0749 for (unsigned int j = 0; j < NROCRow; j++) {
0750 for (unsigned int i = 0; i < NROCCol; i++) {
0751 stAPVGain* APV = new stAPVGain;
0752 APV->Index = Index;
0753 APV->Bin = -1;
0754 APV->DetId = Detid.rawId();
0755 APV->APVId = (j << 3 | i);
0756 APV->SubDet = SubDet;
0757 APV->FitMPV = -1;
0758 APV->FitMPVErr = -1;
0759 APV->FitWidth = -1;
0760 APV->FitWidthErr = -1;
0761 APV->FitChi2 = -1;
0762 APV->Gain = -1;
0763 APV->PreviousGain = 1;
0764 APV->PreviousGainTick = 1;
0765 APV->x = DetUnit->position().basicVector().x();
0766 APV->y = DetUnit->position().basicVector().y();
0767 APV->z = DetUnit->position().basicVector().z();
0768 APV->Eta = DetUnit->position().basicVector().eta();
0769 APV->Phi = DetUnit->position().basicVector().phi();
0770 APV->R = DetUnit->position().basicVector().transverse();
0771 APV->Thickness = DetUnit->surface().bounds().thickness();
0772 APV->isMasked = false;
0773 APV->NEntries = 0;
0774
0775 APVsCollOrdered.push_back(APV);
0776 APVsColl[(APV->DetId << 4) | APV->APVId] = APV;
0777 Index++;
0778 NPixelDets++;
0779 }
0780 }
0781 }
0782 }
0783
0784 MakeCalibrationMap();
0785
0786 NEvent = 0;
0787 NTrack = 0;
0788 NClusterStrip = 0;
0789 NClusterPixel = 0;
0790 SRun = 1 << 31;
0791 ERun = 0;
0792 GOOD = 0;
0793 BAD = 0;
0794 MASKED = 0;
0795 }
0796
0797 bool SiStripGainFromCalibTree::isBFieldConsistentWithMode(const edm::EventSetup& iSetup) const {
0798 const auto& runInfo = iSetup.getData(runInfoToken_);
0799
0800 double average_current = runInfo.m_avg_current;
0801 bool isOn = (average_current > MagFieldCurrentTh);
0802 bool is0T = (m_calibrationMode.substr(m_calibrationMode.length() - 2) == "0T");
0803
0804 return ((isOn && !is0T) || (!isOn && is0T));
0805 }
0806
0807 void SiStripGainFromCalibTree::swapBFieldMode() {
0808 if (m_calibrationMode.substr(m_calibrationMode.length() - 2) == "0T") {
0809 m_calibrationMode.erase(m_calibrationMode.length() - 2, 2);
0810 } else {
0811 m_calibrationMode.append("0T");
0812 }
0813 }
0814
0815 void SiStripGainFromCalibTree::algoBeginRun(const edm::Run& run, const edm::EventSetup& iSetup) {
0816 if (!m_harvestingMode && AlgoMode == "PCL") {
0817
0818 if (!isBFieldConsistentWithMode(iSetup)) {
0819 string prevMode = m_calibrationMode;
0820 swapBFieldMode();
0821 edm::LogInfo("SiStripGainFromCalibTree") << "Switching calibration mode for endorsing BField status: " << prevMode
0822 << " ==> " << m_calibrationMode << std::endl;
0823 }
0824 }
0825
0826 const auto gainHandle = iSetup.getHandle(gainToken_);
0827 if (!gainHandle.isValid()) {
0828 edm::LogError("SiStripGainFromCalibTree") << "gainHandle is not valid\n";
0829 exit(0);
0830 }
0831
0832 const auto& siStripQuality = iSetup.getData(qualityToken_);
0833 for (unsigned int a = 0; a < APVsCollOrdered.size(); a++) {
0834 stAPVGain* APV = APVsCollOrdered[a];
0835
0836
0837 if (APV->SubDet == PixelSubdetector::PixelBarrel || APV->SubDet == PixelSubdetector::PixelEndcap)
0838 continue;
0839
0840 APV->isMasked = siStripQuality.IsApvBad(APV->DetId, APV->APVId);
0841
0842
0843 if (gainHandle->getNumberOfTags() != 2) {
0844 edm::LogError("SiStripGainFromCalibTree") << "NUMBER OF GAIN TAG IS EXPECTED TO BE 2\n";
0845 fflush(stdout);
0846 exit(0);
0847 };
0848 float newPreviousGain = gainHandle->getApvGain(APV->APVId, gainHandle->getRange(APV->DetId, 1), 1);
0849 if (APV->PreviousGain != 1 and newPreviousGain != APV->PreviousGain)
0850 edm::LogWarning("SiStripGainFromCalibTree") << "WARNING: ParticleGain in the global tag changed\n";
0851 APV->PreviousGain = newPreviousGain;
0852
0853 float newPreviousGainTick = gainHandle->getApvGain(APV->APVId, gainHandle->getRange(APV->DetId, 0), 0);
0854 if (APV->PreviousGainTick != 1 and newPreviousGainTick != APV->PreviousGainTick) {
0855 edm::LogWarning("SiStripGainFromCalibTree")
0856 << "WARNING: TickMarkGain in the global tag changed\n"
0857 << std::endl
0858 << " APV->SubDet: " << APV->SubDet << " APV->APVId:" << APV->APVId << std::endl
0859 << " APV->PreviousGainTick: " << APV->PreviousGainTick << " newPreviousGainTick: " << newPreviousGainTick
0860 << std::endl;
0861 }
0862 APV->PreviousGainTick = newPreviousGainTick;
0863
0864
0865
0866 }
0867 }
0868
0869 void SiStripGainFromCalibTree::algoEndRun(const edm::Run& run, const edm::EventSetup& iSetup) {
0870 if (AlgoMode == "PCL" && !m_harvestingMode)
0871 return;
0872
0873 if (AlgoMode == "PCL" and m_harvestingMode) {
0874
0875
0876
0877
0878 edm::LogInfo("SiStripGainFromCalibTree") << "Starting harvesting statistics" << std::endl;
0879
0880
0881 int elepos = statCollectionFromMode(m_calibrationMode.c_str());
0882 if (elepos != Harvest) {
0883
0884 std::string stag = m_calibrationMode;
0885 if (!stag.empty() && stag[0] != '_')
0886 stag.insert(0, 1, '_');
0887
0888 if (elepos == -1) {
0889
0890 elepos = 0;
0891 stag = "";
0892 }
0893
0894 std::string DQM_dir = m_DQMdir + ((m_splitDQMstat) ? m_calibrationMode : "");
0895
0896 std::string cvi = DQM_dir + std::string("/Charge_Vs_Index") + stag;
0897
0898 std::string cvpTIB = DQM_dir + std::string("/Charge_Vs_PathlengthTIB") + stag;
0899 std::string cvpTOB = DQM_dir + std::string("/Charge_Vs_PathlengthTOB") + stag;
0900 std::string cvpTIDP = DQM_dir + std::string("/Charge_Vs_PathlengthTIDP") + stag;
0901 std::string cvpTIDM = DQM_dir + std::string("/Charge_Vs_PathlengthTIDM") + stag;
0902 std::string cvpTECP1 = DQM_dir + std::string("/Charge_Vs_PathlengthTECP1") + stag;
0903 std::string cvpTECP2 = DQM_dir + std::string("/Charge_Vs_PathlengthTECP2") + stag;
0904 std::string cvpTECM1 = DQM_dir + std::string("/Charge_Vs_PathlengthTECM1") + stag;
0905 std::string cvpTECM2 = DQM_dir + std::string("/Charge_Vs_PathlengthTECM2") + stag;
0906
0907 Charge_Vs_Index[elepos] = dbe->get(cvi);
0908
0909 Charge_Vs_PathlengthTIB[elepos] = dbe->get(cvpTIB);
0910 Charge_Vs_PathlengthTOB[elepos] = dbe->get(cvpTOB);
0911 Charge_Vs_PathlengthTIDP[elepos] = dbe->get(cvpTIDP);
0912 Charge_Vs_PathlengthTIDM[elepos] = dbe->get(cvpTIDM);
0913 Charge_Vs_PathlengthTECP1[elepos] = dbe->get(cvpTECP1);
0914 Charge_Vs_PathlengthTECP2[elepos] = dbe->get(cvpTECP2);
0915 Charge_Vs_PathlengthTECM1[elepos] = dbe->get(cvpTECM1);
0916 Charge_Vs_PathlengthTECM2[elepos] = dbe->get(cvpTECM2);
0917
0918 if (Charge_Vs_Index[elepos] == nullptr) {
0919 edm::LogError("SiStripGainFromCalibTree")
0920 << "Harvesting: could not retrieve " << cvi.c_str() << ", statistics will not be summed!" << std::endl;
0921 } else {
0922 merge((Charge_Vs_Index[Harvest])->getTH2S(), (Charge_Vs_Index[elepos])->getTH2S());
0923 edm::LogInfo("SiStripGainFromCalibTree")
0924 << "Harvesting " << (Charge_Vs_Index[elepos])->getTH2S()->GetEntries() << " more clusters" << std::endl;
0925 }
0926
0927
0928
0929
0930
0931
0932 if (Charge_Vs_PathlengthTIB[elepos] == nullptr) {
0933 edm::LogError("SiStripGainFromCalibTree")
0934 << "Harvesting: could not retrieve " << cvpTIB.c_str() << ", statistics will not be summed!" << std::endl;
0935 } else
0936 (Charge_Vs_PathlengthTIB[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTIB[elepos])->getTH2S());
0937
0938 if (Charge_Vs_PathlengthTOB[elepos] == nullptr) {
0939 edm::LogError("SiStripGainFromCalibTree")
0940 << "Harvesting: could not retrieve " << cvpTOB.c_str() << ", statistics will not be summed!" << std::endl;
0941 } else
0942 (Charge_Vs_PathlengthTOB[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTOB[elepos])->getTH2S());
0943
0944 if (Charge_Vs_PathlengthTIDP[elepos] == nullptr) {
0945 edm::LogError("SiStripGainFromCalibTree")
0946 << "Harvesting: could not retrieve " << cvpTIDP.c_str() << ", statistics will not be summed!" << std::endl;
0947 } else
0948 (Charge_Vs_PathlengthTIDP[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTIDP[elepos])->getTH2S());
0949
0950 if (Charge_Vs_PathlengthTIDM[elepos] == nullptr) {
0951 edm::LogError("SiStripGainFromCalibTree")
0952 << "Harvesting: could not retrieve " << cvpTIDM.c_str() << ", statistics will not be summed!" << std::endl;
0953 } else
0954 (Charge_Vs_PathlengthTIDM[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTIDM[elepos])->getTH2S());
0955
0956 if (Charge_Vs_PathlengthTECP1[elepos] == nullptr) {
0957 edm::LogError("SiStripGainFromCalibTree")
0958 << "Harvesting: could not retrieve " << cvpTECP1.c_str() << ", statistics will not be summed!" << std::endl;
0959 } else
0960 (Charge_Vs_PathlengthTECP1[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTECP1[elepos])->getTH2S());
0961
0962 if (Charge_Vs_PathlengthTECP2[elepos] == nullptr) {
0963 edm::LogError("SiStripGainFromCalibTree")
0964 << "Harvesting: could not retrieve " << cvpTECP2.c_str() << ", statistics will not be summed!" << std::endl;
0965 } else
0966 (Charge_Vs_PathlengthTECP2[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTECP2[elepos])->getTH2S());
0967
0968 if (Charge_Vs_PathlengthTECM1[elepos] == nullptr) {
0969 edm::LogError("SiStripGainFromCalibTree")
0970 << "Harvesting: could not retrieve " << cvpTECM1.c_str() << ", statistics will not be summed!" << std::endl;
0971 } else
0972 (Charge_Vs_PathlengthTECM1[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTECM1[elepos])->getTH2S());
0973
0974 if (Charge_Vs_PathlengthTECM2[elepos] == nullptr) {
0975 edm::LogError("SiStripGainFromCalibTree")
0976 << "Harvesting: could not retrieve " << cvpTECM2.c_str() << ", statistics will not be summed!" << std::endl;
0977 } else
0978 (Charge_Vs_PathlengthTECM2[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTECM2[elepos])->getTH2S());
0979
0980
0981 std::vector<std::pair<std::string, std::string>> tags =
0982 APVGain::monHnames(VChargeHisto, doChargeMonitorPerPlane, "");
0983 for (unsigned int i = 0; i < tags.size(); i++) {
0984 std::string tag = DQM_dir + "/" + (tags[i]).first + stag;
0985 Charge_1[elepos].push_back(APVGain::APVmon(0, 0, 0, 0, dbe->get(tag)));
0986 if ((Charge_1[elepos][i]).getMonitor() == nullptr) {
0987 edm::LogError("SiStripGainFromCalibTree")
0988 << "Harvesting: could not retrieve " << tag.c_str() << ", statistics will not be summed!" << std::endl;
0989 } else
0990 (Charge_1[Harvest][i]).getMonitor()->getTH1D()->Add((Charge_1[elepos][i]).getMonitor()->getTH1D());
0991 }
0992
0993 tags = APVGain::monHnames(VChargeHisto, doChargeMonitorPerPlane, "woG2");
0994 for (unsigned int i = 0; i < tags.size(); i++) {
0995 std::string tag = DQM_dir + "/" + (tags[i]).first + stag;
0996 Charge_2[elepos].push_back(APVGain::APVmon(0, 0, 0, 0, dbe->get(tag)));
0997 if ((Charge_2[elepos][i]).getMonitor() == nullptr) {
0998 edm::LogError("SiStripGainFromCalibTree")
0999 << "Harvesting: could not retrieve " << tag.c_str() << ", statistics will not be summed!" << std::endl;
1000 } else
1001 (Charge_2[Harvest][i]).getMonitor()->getTH1D()->Add((Charge_2[elepos][i]).getMonitor()->getTH1D());
1002 }
1003
1004 tags = APVGain::monHnames(VChargeHisto, doChargeMonitorPerPlane, "woG1");
1005 for (unsigned int i = 0; i < tags.size(); i++) {
1006 std::string tag = DQM_dir + "/" + (tags[i]).first + stag;
1007 Charge_3[elepos].push_back(APVGain::APVmon(0, 0, 0, 0, dbe->get(tag)));
1008 if ((Charge_3[elepos][i]).getMonitor() == nullptr) {
1009 edm::LogError("SiStripGainFromCalibTree")
1010 << "Harvesting: could not retrieve " << tag.c_str() << ", statistics will not be summed!" << std::endl;
1011 } else
1012 (Charge_3[Harvest][i]).getMonitor()->getTH1D()->Add((Charge_3[elepos][i]).getMonitor()->getTH1D());
1013 }
1014
1015 tags = APVGain::monHnames(VChargeHisto, doChargeMonitorPerPlane, "woG1G2");
1016 for (unsigned int i = 0; i < tags.size(); i++) {
1017 std::string tag = DQM_dir + "/" + (tags[i]).first + stag;
1018 Charge_4[elepos].push_back(APVGain::APVmon(0, 0, 0, 0, dbe->get(tag)));
1019 if ((Charge_4[elepos][i]).getMonitor() == nullptr) {
1020 edm::LogError("SiStripGainFromCalibTree")
1021 << "Harvesting: could not retrieve " << tag.c_str() << ", statistics will not be summed!" << std::endl;
1022 } else
1023 (Charge_4[Harvest][i]).getMonitor()->getTH1D()->Add((Charge_4[elepos][i]).getMonitor()->getTH1D());
1024 }
1025 }
1026 }
1027 }
1028
1029 void SiStripGainFromCalibTree::algoEndJob() {
1030 if (AlgoMode == "PCL" && !m_harvestingMode)
1031 return;
1032
1033 if (AlgoMode == "CalibTree") {
1034 edm::LogInfo("SiStripGainFromCalibTree") << "Analyzing calibration tree" << std::endl;
1035
1036 algoAnalyzeTheTree();
1037 } else if (m_harvestingMode) {
1038 NClusterStrip = (Charge_Vs_Index[Harvest])->getTH2S()->Integral(0, NStripAPVs + 1, 0, 99999);
1039
1040 }
1041
1042
1043 algoComputeMPVandGain();
1044
1045
1046 qualityMonitor();
1047
1048
1049
1050 storeOnDbNow();
1051
1052 if (AlgoMode != "PCL" or saveSummary) {
1053 edm::LogInfo("SiStripGainFromCalibTree") << "Saving summary into root file" << std::endl;
1054
1055
1056 tfs = edm::Service<TFileService>().operator->();
1057
1058
1059 int elepos = statCollectionFromMode(m_calibrationMode.c_str());
1060
1061 if (Charge_Vs_Index[elepos] != nullptr)
1062 tfs->make<TH2S>(*(Charge_Vs_Index[elepos])->getTH2S());
1063
1064 if (Charge_Vs_PathlengthTIB[elepos] != nullptr)
1065 tfs->make<TH2S>(*(Charge_Vs_PathlengthTIB[elepos])->getTH2S());
1066 if (Charge_Vs_PathlengthTOB[elepos] != nullptr)
1067 tfs->make<TH2S>(*(Charge_Vs_PathlengthTOB[elepos])->getTH2S());
1068 if (Charge_Vs_PathlengthTIDP[elepos] != nullptr)
1069 tfs->make<TH2S>(*(Charge_Vs_PathlengthTIDP[elepos])->getTH2S());
1070 if (Charge_Vs_PathlengthTIDM[elepos] != nullptr)
1071 tfs->make<TH2S>(*(Charge_Vs_PathlengthTIDM[elepos])->getTH2S());
1072 if (Charge_Vs_PathlengthTECP1[elepos] != nullptr)
1073 tfs->make<TH2S>(*(Charge_Vs_PathlengthTECP1[elepos])->getTH2S());
1074 if (Charge_Vs_PathlengthTECP2[elepos] != nullptr)
1075 tfs->make<TH2S>(*(Charge_Vs_PathlengthTECP2[elepos])->getTH2S());
1076 if (Charge_Vs_PathlengthTECM1[elepos] != nullptr)
1077 tfs->make<TH2S>(*(Charge_Vs_PathlengthTECM1[elepos])->getTH2S());
1078 if (Charge_Vs_PathlengthTECM2[elepos] != nullptr)
1079 tfs->make<TH2S>(*(Charge_Vs_PathlengthTECM2[elepos])->getTH2S());
1080
1081 storeOnTree(tfs);
1082 }
1083 }
1084
1085 void SiStripGainFromCalibTree::getPeakOfLandau(TH1* InputHisto, double* FitResults, double LowRange, double HighRange) {
1086 FitResults[0] = -0.5;
1087 FitResults[1] = 0;
1088 FitResults[2] = -0.5;
1089 FitResults[3] = 0;
1090 FitResults[4] = -0.5;
1091 FitResults[5] = 0;
1092
1093 if (InputHisto->GetEntries() < MinNrEntries)
1094 return;
1095
1096
1097 TF1* MyLandau = new TF1("MyLandau", "landau", LowRange, HighRange);
1098 MyLandau->SetParameter(1, 300);
1099 InputHisto->Fit(MyLandau, "0QR WW");
1100
1101
1102 FitResults[0] = MyLandau->GetParameter(1);
1103 FitResults[1] = MyLandau->GetParError(1);
1104 FitResults[2] = MyLandau->GetParameter(2);
1105 FitResults[3] = MyLandau->GetParError(2);
1106 FitResults[4] = MyLandau->GetChisquare() / MyLandau->GetNDF();
1107 FitResults[5] = MyLandau->GetParameter(0);
1108
1109 delete MyLandau;
1110 }
1111
1112 bool SiStripGainFromCalibTree::IsGoodLandauFit(double* FitResults) {
1113 if (FitResults[0] <= 0)
1114 return false;
1115
1116
1117 return true;
1118 }
1119
1120 void SiStripGainFromCalibTree::processEvent() {
1121 edm::LogInfo("SiStripGainFromCalibTree") << "Processing run " << runnumber << " and event " << eventnumber << " for "
1122 << m_calibrationMode << " calibration." << std::endl;
1123
1124 if (runnumber < SRun)
1125 SRun = runnumber;
1126 if (runnumber > ERun)
1127 ERun = runnumber;
1128
1129 NEvent++;
1130 NTrack += (*trackp).size();
1131
1132 int elepos = statCollectionFromMode(m_calibrationMode.c_str());
1133
1134 unsigned int FirstAmplitude = 0;
1135 for (unsigned int i = 0; i < (*chargeoverpath).size(); i++) {
1136 FirstAmplitude += (*nstrips)[i];
1137 int TI = (*trackindex)[i];
1138
1139
1140 if ((*tracketa)[TI] < MinTrackEta)
1141 continue;
1142 if ((*tracketa)[TI] > MaxTrackEta)
1143 continue;
1144 if ((*trackp)[TI] < MinTrackMomentum)
1145 continue;
1146 if ((*trackp)[TI] > MaxTrackMomentum)
1147 continue;
1148 if ((*trackhitsvalid)[TI] < MinTrackHits)
1149 continue;
1150 if ((*trackchi2ndof)[TI] > MaxTrackChiOverNdf)
1151 continue;
1152 if ((*trackalgo)[TI] > MaxTrackingIteration)
1153 continue;
1154
1155 stAPVGain* APV =
1156 APVsColl[((*rawid)[i] << 4) |
1157 ((*firststrip)[i] /
1158 128)];
1159
1160 if (APV->SubDet > 2 && (*farfromedge)[i] == false)
1161 continue;
1162 if (APV->SubDet > 2 && (*overlapping)[i] == true)
1163 continue;
1164 if (APV->SubDet > 2 && (*saturation)[i] && !AllowSaturation)
1165 continue;
1166 if (APV->SubDet > 2 && (*nstrips)[i] > MaxNrStrips)
1167 continue;
1168
1169
1170
1171
1172
1173
1174
1175
1176
1177
1178
1179
1180 if (APV->SubDet > 2) {
1181 NClusterStrip++;
1182 } else {
1183 NClusterPixel++;
1184 }
1185
1186 int Charge = 0;
1187 if (APV->SubDet > 2 && (useCalibration || !FirstSetOfConstants)) {
1188 bool Saturation = false;
1189 for (unsigned int s = 0; s < (*nstrips)[i]; s++) {
1190 int StripCharge = (*amplitude)[FirstAmplitude - (*nstrips)[i] + s];
1191 if (useCalibration && !FirstSetOfConstants) {
1192 StripCharge = (int)(StripCharge * (APV->PreviousGain / APV->CalibGain));
1193 } else if (useCalibration) {
1194 StripCharge = (int)(StripCharge / APV->CalibGain);
1195 } else if (!FirstSetOfConstants) {
1196 StripCharge = (int)(StripCharge * APV->PreviousGain);
1197 }
1198 if (StripCharge > 1024) {
1199 StripCharge = 255;
1200 Saturation = true;
1201 } else if (StripCharge > 254) {
1202 StripCharge = 254;
1203 Saturation = true;
1204 }
1205 Charge += StripCharge;
1206 }
1207 if (Saturation && !AllowSaturation)
1208 continue;
1209 } else if (APV->SubDet > 2) {
1210 Charge = (*charge)[i];
1211 } else {
1212 Charge = (*charge)[i] / 265.0;
1213 }
1214
1215
1216
1217 double ClusterChargeOverPath = ((double)Charge) / (*path)[i];
1218 if (APV->SubDet > 2) {
1219 if (Validation) {
1220 ClusterChargeOverPath /= (*gainused)[i];
1221 }
1222 if (OldGainRemoving) {
1223 ClusterChargeOverPath *= (*gainused)[i];
1224 }
1225 }
1226
1227
1228 if (APV->SubDet <= 2)
1229 continue;
1230
1231 (Charge_Vs_Index[elepos])->Fill(APV->Index, ClusterChargeOverPath);
1232
1233
1234 int mCharge1 = 0;
1235 int mCharge2 = 0;
1236 int mCharge3 = 0;
1237 int mCharge4 = 0;
1238 if (APV->SubDet > 2) {
1239 for (unsigned int s = 0; s < (*nstrips)[i]; s++) {
1240 int StripCharge = (*amplitude)[FirstAmplitude - (*nstrips)[i] + s];
1241 if (StripCharge > 1024)
1242 StripCharge = 255;
1243 else if (StripCharge > 254)
1244 StripCharge = 254;
1245 mCharge1 += StripCharge;
1246 mCharge2 += StripCharge;
1247 mCharge3 += StripCharge;
1248 mCharge4 += StripCharge;
1249 }
1250
1251 mCharge2 *= (*gainused)[i];
1252 mCharge3 *= (*gainusedTick)[i];
1253 mCharge4 *= ((*gainused)[i] * (*gainusedTick)[i]);
1254 }
1255 std::vector<APVGain::APVmon>& v1 = Charge_1[elepos];
1256 std::vector<MonitorElement*> cmon1 = APVGain::FetchMonitor(v1, (*rawid)[i], tTopo_);
1257 for (unsigned int m = 0; m < cmon1.size(); m++)
1258 cmon1[m]->Fill(((double)mCharge1) / (*path)[i]);
1259
1260 std::vector<APVGain::APVmon>& v2 = Charge_2[elepos];
1261 std::vector<MonitorElement*> cmon2 = APVGain::FetchMonitor(v2, (*rawid)[i], tTopo_);
1262 for (unsigned int m = 0; m < cmon2.size(); m++)
1263 cmon2[m]->Fill(((double)mCharge2) / (*path)[i]);
1264
1265 std::vector<APVGain::APVmon>& v3 = Charge_3[elepos];
1266 std::vector<MonitorElement*> cmon3 = APVGain::FetchMonitor(v3, (*rawid)[i], tTopo_);
1267 for (unsigned int m = 0; m < cmon3.size(); m++)
1268 cmon3[m]->Fill(((double)mCharge3) / (*path)[i]);
1269
1270 std::vector<APVGain::APVmon>& v4 = Charge_4[elepos];
1271 std::vector<MonitorElement*> cmon4 = APVGain::FetchMonitor(v4, (*rawid)[i], tTopo_);
1272 for (unsigned int m = 0; m < cmon4.size(); m++)
1273 cmon4[m]->Fill(((double)mCharge4) / (*path)[i]);
1274
1275
1276 if (APV->SubDet == StripSubdetector::TIB) {
1277 (Charge_Vs_PathlengthTIB[elepos])->Fill((*path)[i], Charge);
1278
1279 } else if (APV->SubDet == StripSubdetector::TOB) {
1280 (Charge_Vs_PathlengthTOB[elepos])->Fill((*path)[i], Charge);
1281
1282 } else if (APV->SubDet == StripSubdetector::TID) {
1283 if (APV->Eta < 0) {
1284 (Charge_Vs_PathlengthTIDM[elepos])->Fill((*path)[i], Charge);
1285 }
1286 else if (APV->Eta > 0) {
1287 (Charge_Vs_PathlengthTIDP[elepos])->Fill((*path)[i], Charge);
1288 }
1289
1290 } else if (APV->SubDet == StripSubdetector::TEC) {
1291 if (APV->Eta < 0) {
1292 if (APV->Thickness < 0.04) {
1293 (Charge_Vs_PathlengthTECM1[elepos])->Fill((*path)[i], Charge);
1294 }
1295 else if (APV->Thickness > 0.04) {
1296 (Charge_Vs_PathlengthTECM2[elepos])->Fill((*path)[i], Charge);
1297 }
1298 } else if (APV->Eta > 0) {
1299 if (APV->Thickness < 0.04) {
1300 (Charge_Vs_PathlengthTECP1[elepos])->Fill((*path)[i], Charge);
1301 }
1302 else if (APV->Thickness > 0.04) {
1303 (Charge_Vs_PathlengthTECP2[elepos])->Fill((*path)[i], Charge);
1304 }
1305 }
1306 }
1307
1308 }
1309 }
1310
1311 void SiStripGainFromCalibTree::algoAnalyzeTheTree() {
1312 for (unsigned int i = 0; i < VInputFiles.size(); i++) {
1313 printf("Openning file %3i/%3i --> %s\n", i + 1, (int)VInputFiles.size(), (char*)(VInputFiles[i].c_str()));
1314 fflush(stdout);
1315 TFile* tfile = TFile::Open(VInputFiles[i].c_str());
1316 TString tree_path = TString::Format("gainCalibrationTree%s/tree", m_calibrationMode.c_str());
1317 TTree* tree = dynamic_cast<TTree*>(tfile->Get(tree_path.Data()));
1318
1319 tree->SetBranchAddress((EventPrefix_ + "event" + EventSuffix_).c_str(), &eventnumber, nullptr);
1320 tree->SetBranchAddress((EventPrefix_ + "run" + EventSuffix_).c_str(), &runnumber, nullptr);
1321 tree->SetBranchAddress((EventPrefix_ + "TrigTech" + EventSuffix_).c_str(), &TrigTech, nullptr);
1322
1323 tree->SetBranchAddress((TrackPrefix_ + "chi2ndof" + TrackSuffix_).c_str(), &trackchi2ndof, nullptr);
1324 tree->SetBranchAddress((TrackPrefix_ + "momentum" + TrackSuffix_).c_str(), &trackp, nullptr);
1325 tree->SetBranchAddress((TrackPrefix_ + "pt" + TrackSuffix_).c_str(), &trackpt, nullptr);
1326 tree->SetBranchAddress((TrackPrefix_ + "eta" + TrackSuffix_).c_str(), &tracketa, nullptr);
1327 tree->SetBranchAddress((TrackPrefix_ + "phi" + TrackSuffix_).c_str(), &trackphi, nullptr);
1328 tree->SetBranchAddress((TrackPrefix_ + "hitsvalid" + TrackSuffix_).c_str(), &trackhitsvalid, nullptr);
1329 tree->SetBranchAddress((TrackPrefix_ + "algo" + TrackSuffix_).c_str(), &trackalgo, nullptr);
1330
1331 tree->SetBranchAddress((CalibPrefix_ + "trackindex" + CalibSuffix_).c_str(), &trackindex, nullptr);
1332 tree->SetBranchAddress((CalibPrefix_ + "rawid" + CalibSuffix_).c_str(), &rawid, nullptr);
1333 tree->SetBranchAddress((CalibPrefix_ + "localdirx" + CalibSuffix_).c_str(), &localdirx, nullptr);
1334 tree->SetBranchAddress((CalibPrefix_ + "localdiry" + CalibSuffix_).c_str(), &localdiry, nullptr);
1335 tree->SetBranchAddress((CalibPrefix_ + "localdirz" + CalibSuffix_).c_str(), &localdirz, nullptr);
1336 tree->SetBranchAddress((CalibPrefix_ + "firststrip" + CalibSuffix_).c_str(), &firststrip, nullptr);
1337 tree->SetBranchAddress((CalibPrefix_ + "nstrips" + CalibSuffix_).c_str(), &nstrips, nullptr);
1338 tree->SetBranchAddress((CalibPrefix_ + "saturation" + CalibSuffix_).c_str(), &saturation, nullptr);
1339 tree->SetBranchAddress((CalibPrefix_ + "overlapping" + CalibSuffix_).c_str(), &overlapping, nullptr);
1340 tree->SetBranchAddress((CalibPrefix_ + "farfromedge" + CalibSuffix_).c_str(), &farfromedge, nullptr);
1341 tree->SetBranchAddress((CalibPrefix_ + "charge" + CalibSuffix_).c_str(), &charge, nullptr);
1342 tree->SetBranchAddress((CalibPrefix_ + "path" + CalibSuffix_).c_str(), &path, nullptr);
1343 tree->SetBranchAddress((CalibPrefix_ + "chargeoverpath" + CalibSuffix_).c_str(), &chargeoverpath, nullptr);
1344 tree->SetBranchAddress((CalibPrefix_ + "amplitude" + CalibSuffix_).c_str(), &litude, nullptr);
1345 tree->SetBranchAddress((CalibPrefix_ + "gainused" + CalibSuffix_).c_str(), &gainused, nullptr);
1346 tree->SetBranchAddress((CalibPrefix_ + "gainusedTick" + CalibSuffix_).c_str(), &gainusedTick, nullptr);
1347
1348 unsigned int nentries = tree->GetEntries();
1349 printf("Number of Events = %i + %i = %i\n", NEvent, nentries, (NEvent + nentries));
1350 printf("Progressing Bar :0%% 20%% 40%% 60%% 80%% 100%%\n");
1351 printf("Looping on the Tree :");
1352 int TreeStep = nentries / 50;
1353 if (TreeStep <= 1)
1354 TreeStep = 1;
1355 for (unsigned int ientry = 0; ientry < tree->GetEntries(); ientry++) {
1356 if (ientry % TreeStep == 0) {
1357 printf(".");
1358 fflush(stdout);
1359 }
1360 tree->GetEntry(ientry);
1361 processEvent();
1362 }
1363 printf("\n");
1364 }
1365 }
1366
1367 void SiStripGainFromCalibTree::algoComputeMPVandGain() {
1368 unsigned int I = 0;
1369 TH1F* Proj = nullptr;
1370 double FitResults[6];
1371 double MPVmean = 300;
1372
1373 int elepos = (AlgoMode == "PCL") ? Harvest : statCollectionFromMode(m_calibrationMode.c_str());
1374
1375 if (Charge_Vs_Index[elepos] == nullptr) {
1376 edm::LogError("SiStripGainFromCalibTree")
1377 << "Harvesting: could not execute algoComputeMPVandGain method because " << m_calibrationMode.c_str()
1378 << " statistics cannot be retrieved.\n"
1379 << "Please check if input contains " << m_calibrationMode.c_str() << " data." << std::endl;
1380 return;
1381 }
1382
1383 TH2S* chvsidx = (Charge_Vs_Index[elepos])->getTH2S();
1384
1385 printf("Progressing Bar :0%% 20%% 40%% 60%% 80%% 100%%\n");
1386 printf("Fitting Charge Distribution :");
1387 int TreeStep = APVsColl.size() / 50;
1388 for (auto it = APVsColl.begin(); it != APVsColl.end(); it++, I++) {
1389 if (I % TreeStep == 0) {
1390 printf(".");
1391 fflush(stdout);
1392 }
1393 stAPVGain* APV = it->second;
1394 if (APV->Bin < 0)
1395 APV->Bin = chvsidx->GetXaxis()->FindBin(APV->Index);
1396
1397 if (APV->isMasked) {
1398 APV->Gain = APV->PreviousGain;
1399 MASKED++;
1400 continue;
1401 }
1402
1403 Proj = (TH1F*)(chvsidx->ProjectionY(
1404 "", chvsidx->GetXaxis()->FindBin(APV->Index), chvsidx->GetXaxis()->FindBin(APV->Index), "e"));
1405 if (!Proj)
1406 continue;
1407
1408 if (CalibrationLevel == 0) {
1409 } else if (CalibrationLevel == 1) {
1410 int SecondAPVId = APV->APVId;
1411 if (SecondAPVId % 2 == 0) {
1412 SecondAPVId = SecondAPVId + 1;
1413 } else {
1414 SecondAPVId = SecondAPVId - 1;
1415 }
1416 stAPVGain* APV2 = APVsColl[(APV->DetId << 4) | SecondAPVId];
1417 if (APV2->Bin < 0)
1418 APV2->Bin = chvsidx->GetXaxis()->FindBin(APV2->Index);
1419 TH1F* Proj2 = (TH1F*)(chvsidx->ProjectionY("", APV2->Bin, APV2->Bin, "e"));
1420 if (Proj2) {
1421 Proj->Add(Proj2, 1);
1422 delete Proj2;
1423 }
1424 } else if (CalibrationLevel == 2) {
1425 for (unsigned int i = 0; i < 16; i++) {
1426 auto tmpit = APVsColl.find((APV->DetId << 4) | i);
1427 if (tmpit == APVsColl.end())
1428 continue;
1429 stAPVGain* APV2 = tmpit->second;
1430 if (APV2->DetId != APV->DetId || APV2->APVId == APV->APVId)
1431 continue;
1432 if (APV2->Bin < 0)
1433 APV2->Bin = chvsidx->GetXaxis()->FindBin(APV2->Index);
1434 TH1F* Proj2 = (TH1F*)(chvsidx->ProjectionY("", APV2->Bin, APV2->Bin, "e"));
1435 if (Proj2) {
1436 Proj->Add(Proj2, 1);
1437 delete Proj2;
1438 }
1439 }
1440 } else {
1441 CalibrationLevel = 0;
1442 printf("Unknown Calibration Level, will assume %i\n", CalibrationLevel);
1443 }
1444
1445 getPeakOfLandau(Proj, FitResults);
1446 APV->FitMPV = FitResults[0];
1447 APV->FitMPVErr = FitResults[1];
1448 APV->FitWidth = FitResults[2];
1449 APV->FitWidthErr = FitResults[3];
1450 APV->FitChi2 = FitResults[4];
1451 APV->FitNorm = FitResults[5];
1452 APV->NEntries = Proj->GetEntries();
1453
1454 if (IsGoodLandauFit(FitResults)) {
1455 APV->Gain = APV->FitMPV / MPVmean;
1456 if (APV->SubDet > 2)
1457 GOOD++;
1458 } else {
1459 APV->Gain = APV->PreviousGain;
1460 if (APV->SubDet > 2)
1461 BAD++;
1462 }
1463 if (APV->Gain <= 0)
1464 APV->Gain = 1;
1465
1466
1467 delete Proj;
1468 }
1469 printf("\n");
1470 }
1471
1472 void SiStripGainFromCalibTree::qualityMonitor() {
1473 int elepos = (AlgoMode == "PCL") ? Harvest : statCollectionFromMode(m_calibrationMode.c_str());
1474
1475 for (unsigned int a = 0; a < APVsCollOrdered.size(); a++) {
1476 stAPVGain* APV = APVsCollOrdered[a];
1477 if (APV == nullptr)
1478 continue;
1479
1480 unsigned int Index = APV->Index;
1481 unsigned int SubDet = APV->SubDet;
1482 unsigned int DetId = APV->DetId;
1483 float z = APV->z;
1484 float Eta = APV->Eta;
1485 float R = APV->R;
1486 float Phi = APV->Phi;
1487 float Thickness = APV->Thickness;
1488 double FitMPV = APV->FitMPV;
1489 double FitMPVErr = APV->FitMPVErr;
1490 double Gain = APV->Gain;
1491 double NEntries = APV->NEntries;
1492 double PreviousGain = APV->PreviousGain;
1493
1494 if (SubDet < 3)
1495 continue;
1496
1497 if (Gain != 1.) {
1498 std::vector<MonitorElement*> charge_histos = APVGain::FetchMonitor(newCharge, DetId, tTopo_);
1499 TH2S* chvsidx = (Charge_Vs_Index[elepos])->getTH2S();
1500 int bin = chvsidx->GetXaxis()->FindBin(Index);
1501 TH1D* Proj = chvsidx->ProjectionY("proj", bin, bin);
1502 for (int binId = 0; binId < Proj->GetXaxis()->GetNbins(); binId++) {
1503 double new_charge = Proj->GetXaxis()->GetBinCenter(binId) / Gain;
1504 if (Proj->GetBinContent(binId) != 0.) {
1505 for (unsigned int h = 0; h < charge_histos.size(); h++) {
1506 TH1D* chisto = (charge_histos[h])->getTH1D();
1507 for (int e = 0; e < Proj->GetBinContent(binId); e++)
1508 chisto->Fill(new_charge);
1509 }
1510 }
1511 }
1512 }
1513
1514 if (FitMPV <= 0.) {
1515 if (APV->isMasked)
1516 NoMPVmasked->Fill(z, R);
1517 else
1518 NoMPVfit->Fill(z, R);
1519
1520 } else {
1521 if (FitMPV > 0.)
1522 Gains->Fill(Gain);
1523
1524 MPVs->Fill(FitMPV);
1525 if (Thickness < 0.04)
1526 MPVs320->Fill(FitMPV);
1527 if (Thickness > 0.04)
1528 MPVs500->Fill(FitMPV);
1529
1530 MPVError->Fill(FitMPVErr);
1531 MPVErrorVsMPV->Fill(FitMPV, FitMPVErr);
1532 MPVErrorVsEta->Fill(Eta, FitMPVErr);
1533 MPVErrorVsPhi->Fill(Phi, FitMPVErr);
1534 MPVErrorVsN->Fill(NEntries, FitMPVErr);
1535
1536 if (SubDet == 3) {
1537 MPV_Vs_EtaTIB->Fill(Eta, FitMPV);
1538 MPV_Vs_PhiTIB->Fill(Phi, FitMPV);
1539 MPVsTIB->Fill(FitMPV);
1540
1541 } else if (SubDet == 4) {
1542 MPV_Vs_EtaTID->Fill(Eta, FitMPV);
1543 MPV_Vs_PhiTID->Fill(Phi, FitMPV);
1544 MPVsTID->Fill(FitMPV);
1545 if (Eta < 0.)
1546 MPVsTIDM->Fill(FitMPV);
1547 if (Eta > 0.)
1548 MPVsTIDP->Fill(FitMPV);
1549
1550 } else if (SubDet == 5) {
1551 MPV_Vs_EtaTOB->Fill(Eta, FitMPV);
1552 MPV_Vs_PhiTOB->Fill(Phi, FitMPV);
1553 MPVsTOB->Fill(FitMPV);
1554
1555 } else if (SubDet == 6) {
1556 MPV_Vs_EtaTEC->Fill(Eta, FitMPV);
1557 MPV_Vs_PhiTEC->Fill(Phi, FitMPV);
1558 MPVsTEC->Fill(FitMPV);
1559 if (Eta < 0.)
1560 MPVsTECM->Fill(FitMPV);
1561 if (Eta > 0.)
1562 MPVsTECP->Fill(FitMPV);
1563 if (Thickness < 0.04) {
1564 MPV_Vs_EtaTECthin->Fill(Eta, FitMPV);
1565 MPV_Vs_PhiTECthin->Fill(Phi, FitMPV);
1566 MPVsTECthin->Fill(FitMPV);
1567 if (Eta > 0.)
1568 MPVsTECP1->Fill(FitMPV);
1569 if (Eta < 0.)
1570 MPVsTECM1->Fill(FitMPV);
1571 }
1572 if (Thickness > 0.04) {
1573 MPV_Vs_EtaTECthick->Fill(Eta, FitMPV);
1574 MPV_Vs_PhiTECthick->Fill(Phi, FitMPV);
1575 MPVsTECthick->Fill(FitMPV);
1576 if (Eta > 0.)
1577 MPVsTECP2->Fill(FitMPV);
1578 if (Eta < 0.)
1579 MPVsTECM2->Fill(FitMPV);
1580 }
1581 }
1582 }
1583
1584 if (SubDet == 3 && PreviousGain != 0.)
1585 DiffWRTPrevGainTIB->Fill(Gain / PreviousGain);
1586 else if (SubDet == 4 && PreviousGain != 0.)
1587 DiffWRTPrevGainTID->Fill(Gain / PreviousGain);
1588 else if (SubDet == 5 && PreviousGain != 0.)
1589 DiffWRTPrevGainTOB->Fill(Gain / PreviousGain);
1590 else if (SubDet == 6 && PreviousGain != 0.)
1591 DiffWRTPrevGainTEC->Fill(Gain / PreviousGain);
1592
1593 if (SubDet == 3)
1594 GainVsPrevGainTIB->Fill(PreviousGain, Gain);
1595 else if (SubDet == 4)
1596 GainVsPrevGainTID->Fill(PreviousGain, Gain);
1597 else if (SubDet == 5)
1598 GainVsPrevGainTOB->Fill(PreviousGain, Gain);
1599 else if (SubDet == 6)
1600 GainVsPrevGainTEC->Fill(PreviousGain, Gain);
1601 }
1602 }
1603
1604 void SiStripGainFromCalibTree::storeOnTree(TFileService* tfs) {
1605 unsigned int tree_Index;
1606 unsigned int tree_Bin;
1607 unsigned int tree_DetId;
1608 unsigned char tree_APVId;
1609 unsigned char tree_SubDet;
1610 float tree_x;
1611 float tree_y;
1612 float tree_z;
1613 float tree_Eta;
1614 float tree_R;
1615 float tree_Phi;
1616 float tree_Thickness;
1617 float tree_FitMPV;
1618 float tree_FitMPVErr;
1619 float tree_FitWidth;
1620 float tree_FitWidthErr;
1621 float tree_FitChi2NDF;
1622 float tree_FitNorm;
1623 double tree_Gain;
1624 double tree_PrevGain;
1625 double tree_PrevGainTick;
1626 double tree_NEntries;
1627 bool tree_isMasked;
1628
1629 TTree* MyTree;
1630 MyTree = tfs->make<TTree>("APVGain", "APVGain");
1631 MyTree->Branch("Index", &tree_Index, "Index/i");
1632 MyTree->Branch("Bin", &tree_Bin, "Bin/i");
1633 MyTree->Branch("DetId", &tree_DetId, "DetId/i");
1634 MyTree->Branch("APVId", &tree_APVId, "APVId/b");
1635 MyTree->Branch("SubDet", &tree_SubDet, "SubDet/b");
1636 MyTree->Branch("x", &tree_x, "x/F");
1637 MyTree->Branch("y", &tree_y, "y/F");
1638 MyTree->Branch("z", &tree_z, "z/F");
1639 MyTree->Branch("Eta", &tree_Eta, "Eta/F");
1640 MyTree->Branch("R", &tree_R, "R/F");
1641 MyTree->Branch("Phi", &tree_Phi, "Phi/F");
1642 MyTree->Branch("Thickness", &tree_Thickness, "Thickness/F");
1643 MyTree->Branch("FitMPV", &tree_FitMPV, "FitMPV/F");
1644 MyTree->Branch("FitMPVErr", &tree_FitMPVErr, "FitMPVErr/F");
1645 MyTree->Branch("FitWidth", &tree_FitWidth, "FitWidth/F");
1646 MyTree->Branch("FitWidthErr", &tree_FitWidthErr, "FitWidthErr/F");
1647 MyTree->Branch("FitChi2NDF", &tree_FitChi2NDF, "FitChi2NDF/F");
1648 MyTree->Branch("FitNorm", &tree_FitNorm, "FitNorm/F");
1649 MyTree->Branch("Gain", &tree_Gain, "Gain/D");
1650 MyTree->Branch("PrevGain", &tree_PrevGain, "PrevGain/D");
1651 MyTree->Branch("PrevGainTick", &tree_PrevGainTick, "PrevGainTick/D");
1652 MyTree->Branch("NEntries", &tree_NEntries, "NEntries/D");
1653 MyTree->Branch("isMasked", &tree_isMasked, "isMasked/O");
1654
1655 FILE* Gains = stdout;
1656 fprintf(Gains, "NEvents = %i\n", NEvent);
1657 fprintf(Gains, "NTracks = %i\n", NTrack);
1658
1659 fprintf(Gains, "NClustersStrip = %i\n", NClusterStrip);
1660
1661 fprintf(Gains, "Number of Strip APVs = %lu\n", static_cast<unsigned long>(NStripAPVs));
1662 fprintf(Gains,
1663 "GoodFits = %i BadFits = %i ratio = %f%% (MASKED=%i)\n",
1664 GOOD,
1665 BAD,
1666 (100.0 * GOOD) / (GOOD + BAD),
1667 MASKED);
1668
1669 Gains = fopen(OutputGains.c_str(), "w");
1670 fprintf(Gains, "NEvents = %i\n", NEvent);
1671 fprintf(Gains, "NTracks = %i\n", NTrack);
1672
1673 fprintf(Gains, "NClustersStrip = %i\n", NClusterStrip);
1674 fprintf(Gains, "Number of Strip APVs = %lu\n", static_cast<unsigned long>(NStripAPVs));
1675
1676 fprintf(Gains,
1677 "GoodFits = %i BadFits = %i ratio = %f%% (MASKED=%i)\n",
1678 GOOD,
1679 BAD,
1680 (100.0 * GOOD) / (GOOD + BAD),
1681 MASKED);
1682
1683 int elepos = statCollectionFromMode(m_calibrationMode.c_str());
1684
1685 for (unsigned int a = 0; a < APVsCollOrdered.size(); a++) {
1686 stAPVGain* APV = APVsCollOrdered[a];
1687 if (APV == nullptr)
1688 continue;
1689
1690 fprintf(Gains,
1691 "%i | %i | PreviousGain = %7.5f(tick) x %7.5f(particle) NewGain (particle) = %7.5f (#clusters=%8.0f)\n",
1692 APV->DetId,
1693 APV->APVId,
1694 APV->PreviousGainTick,
1695 APV->PreviousGain,
1696 APV->Gain,
1697 APV->NEntries);
1698
1699 tree_Index = APV->Index;
1700 tree_Bin = (Charge_Vs_Index[elepos])->getTH2S()->GetXaxis()->FindBin(APV->Index);
1701 tree_DetId = APV->DetId;
1702 tree_APVId = APV->APVId;
1703 tree_SubDet = APV->SubDet;
1704 tree_x = APV->x;
1705 tree_y = APV->y;
1706 tree_z = APV->z;
1707 tree_Eta = APV->Eta;
1708 tree_R = APV->R;
1709 tree_Phi = APV->Phi;
1710 tree_Thickness = APV->Thickness;
1711 tree_FitMPV = APV->FitMPV;
1712 tree_FitMPVErr = APV->FitMPVErr;
1713 tree_FitWidth = APV->FitWidth;
1714 tree_FitWidthErr = APV->FitWidthErr;
1715 tree_FitChi2NDF = APV->FitChi2;
1716 tree_FitNorm = APV->FitNorm;
1717 tree_Gain = APV->Gain;
1718 tree_PrevGain = APV->PreviousGain;
1719 tree_PrevGainTick = APV->PreviousGainTick;
1720 tree_NEntries = APV->NEntries;
1721 tree_isMasked = APV->isMasked;
1722
1723 if (tree_DetId == 402673324) {
1724 printf("%i | %i : %f --> %f (%f)\n", tree_DetId, tree_APVId, tree_PrevGain, tree_Gain, tree_NEntries);
1725 }
1726
1727 MyTree->Fill();
1728 }
1729 if (Gains)
1730 fclose(Gains);
1731 }
1732
1733 bool SiStripGainFromCalibTree::produceTagFilter() {
1734
1735
1736 int elepos = (AlgoMode == "PCL") ? Harvest : statCollectionFromMode(m_calibrationMode.c_str());
1737 if (Charge_Vs_Index[elepos] == nullptr) {
1738 edm::LogError("SiStripGainFromCalibTree")
1739 << "produceTagFilter -> Return false: could not retrieve the " << m_calibrationMode.c_str() << " statistics.\n"
1740 << "Please check if input contains " << m_calibrationMode.c_str() << " data." << std::endl;
1741 return false;
1742 }
1743
1744 float integral = (Charge_Vs_Index[elepos])->getTH2S()->Integral();
1745 if ((Charge_Vs_Index[elepos])->getTH2S()->Integral(0, NStripAPVs + 1, 0, 99999) < tagCondition_NClusters) {
1746 edm::LogWarning("SiStripGainFromCalibTree")
1747 << "calibrationMode -> " << m_calibrationMode << "\n"
1748 << "produceTagFilter -> Return false: Statistics is too low : " << integral << endl;
1749 return false;
1750 }
1751 if ((1.0 * GOOD) / (GOOD + BAD) < tagCondition_GoodFrac) {
1752 edm::LogWarning("SiStripGainFromCalibTree")
1753 << "calibrationMode -> " << m_calibrationMode << "\n"
1754 << "produceTagFilter -> Return false: ratio of GOOD/TOTAL is too low: " << (1.0 * GOOD) / (GOOD + BAD) << endl;
1755 return false;
1756 }
1757 return true;
1758 }
1759
1760 std::unique_ptr<SiStripApvGain> SiStripGainFromCalibTree::getNewObject() {
1761 auto obj = std::make_unique<SiStripApvGain>();
1762 if (!m_harvestingMode)
1763 return obj;
1764
1765 if (!produceTagFilter()) {
1766 edm::LogWarning("SiStripGainFromCalibTree")
1767 << "getNewObject -> will not produce a paylaod because produceTagFilter returned false " << endl;
1768 setDoStore(false);
1769 return obj;
1770 }
1771
1772 std::vector<float>* theSiStripVector = nullptr;
1773 unsigned int PreviousDetId = 0;
1774 for (unsigned int a = 0; a < APVsCollOrdered.size(); a++) {
1775 stAPVGain* APV = APVsCollOrdered[a];
1776 if (APV == nullptr) {
1777 printf("Bug\n");
1778 continue;
1779 }
1780 if (APV->SubDet <= 2)
1781 continue;
1782 if (APV->DetId != PreviousDetId) {
1783 if (theSiStripVector != nullptr) {
1784 SiStripApvGain::Range range(theSiStripVector->begin(), theSiStripVector->end());
1785 if (!obj->put(PreviousDetId, range))
1786 printf("Bug to put detId = %i\n", PreviousDetId);
1787 }
1788 theSiStripVector = new std::vector<float>;
1789 PreviousDetId = APV->DetId;
1790 }
1791 if (theSiStripVector != nullptr) {
1792 theSiStripVector->push_back(APV->Gain);
1793 }
1794 }
1795 if (theSiStripVector != nullptr) {
1796 SiStripApvGain::Range range(theSiStripVector->begin(), theSiStripVector->end());
1797 if (!obj->put(PreviousDetId, range))
1798 printf("Bug to put detId = %i\n", PreviousDetId);
1799 }
1800
1801 if (theSiStripVector != nullptr)
1802 delete theSiStripVector;
1803
1804 return obj;
1805 }
1806
1807 SiStripGainFromCalibTree::~SiStripGainFromCalibTree() {
1808 APVsColl.clear();
1809 for (unsigned int a = 0; a < APVsCollOrdered.size(); a++) {
1810 stAPVGain* APV = APVsCollOrdered[a];
1811 if (APV != nullptr)
1812 delete APV;
1813 }
1814 APVsCollOrdered.clear();
1815 }
1816
1817 void SiStripGainFromCalibTree::MakeCalibrationMap() {
1818 if (!useCalibration)
1819 return;
1820
1821 TChain* t1 = new TChain("SiStripCalib/APVGain");
1822 t1->Add(m_calibrationPath.c_str());
1823
1824 unsigned int tree_DetId;
1825 unsigned char tree_APVId;
1826 double tree_Gain;
1827
1828 t1->SetBranchAddress("DetId", &tree_DetId);
1829 t1->SetBranchAddress("APVId", &tree_APVId);
1830 t1->SetBranchAddress("Gain", &tree_Gain);
1831
1832 for (unsigned int ientry = 0; ientry < t1->GetEntries(); ientry++) {
1833 t1->GetEntry(ientry);
1834 stAPVGain* APV = APVsColl[(tree_DetId << 4) | (unsigned int)tree_APVId];
1835 APV->CalibGain = tree_Gain;
1836 }
1837
1838 delete t1;
1839 }
1840
1841 void SiStripGainFromCalibTree::algoAnalyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
1842
1843 if (m_harvestingMode)
1844 return;
1845
1846 if (AlgoMode == "CalibTree")
1847 return;
1848
1849 eventnumber = iEvent.id().event();
1850 runnumber = iEvent.id().run();
1851 auto handle01 = connect(TrigTech, TrigTech_token_, iEvent);
1852 auto handle02 = connect(trackchi2ndof, trackchi2ndof_token_, iEvent);
1853 auto handle03 = connect(trackp, trackp_token_, iEvent);
1854 auto handle04 = connect(trackpt, trackpt_token_, iEvent);
1855 auto handle05 = connect(tracketa, tracketa_token_, iEvent);
1856 auto handle06 = connect(trackphi, trackphi_token_, iEvent);
1857 auto handle07 = connect(trackhitsvalid, trackhitsvalid_token_, iEvent);
1858 auto handle08 = connect(trackindex, trackindex_token_, iEvent);
1859 auto handle09 = connect(rawid, rawid_token_, iEvent);
1860 auto handle11 = connect(localdirx, localdirx_token_, iEvent);
1861 auto handle12 = connect(localdiry, localdiry_token_, iEvent);
1862 auto handle13 = connect(localdirz, localdirz_token_, iEvent);
1863 auto handle14 = connect(firststrip, firststrip_token_, iEvent);
1864 auto handle15 = connect(nstrips, nstrips_token_, iEvent);
1865 auto handle16 = connect(saturation, saturation_token_, iEvent);
1866 auto handle17 = connect(overlapping, overlapping_token_, iEvent);
1867 auto handle18 = connect(farfromedge, farfromedge_token_, iEvent);
1868 auto handle19 = connect(charge, charge_token_, iEvent);
1869 auto handle21 = connect(path, path_token_, iEvent);
1870 auto handle22 = connect(chargeoverpath, chargeoverpath_token_, iEvent);
1871 auto handle23 = connect(amplitude, amplitude_token_, iEvent);
1872 auto handle24 = connect(gainused, gainused_token_, iEvent);
1873 auto handle25 = connect(gainusedTick, gainusedTick_token_, iEvent);
1874 auto handle26 = connect(trackalgo, trackalgo_token_, iEvent);
1875
1876 processEvent();
1877 }
1878
1879 DEFINE_FWK_MODULE(SiStripGainFromCalibTree);