File indexing completed on 2024-04-06 12:29:03
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011 #include <string>
0012
0013
0014 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0015 #include "FWCore/Framework/interface/ConsumesCollector.h"
0016 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0017 #include "FWCore/Framework/interface/Event.h"
0018 #include "FWCore/Framework/interface/Frameworkfwd.h"
0019 #include "FWCore/Framework/interface/LuminosityBlock.h"
0020 #include "FWCore/Framework/interface/MakerMacros.h"
0021 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "RecoVertex/BeamSpotProducer/interface/BeamFitter.h"
0024
0025 #include "TMath.h"
0026
0027 class BeamSpotAnalyzer : public edm::one::EDAnalyzer<edm::one::WatchLuminosityBlocks> {
0028 public:
0029 explicit BeamSpotAnalyzer(const edm::ParameterSet&);
0030 ~BeamSpotAnalyzer() override;
0031 static void fillDescriptions(edm::ConfigurationDescriptions&);
0032
0033 private:
0034 void analyze(const edm::Event&, const edm::EventSetup&) override;
0035 void endJob() override;
0036 void beginLuminosityBlock(const edm::LuminosityBlock& lumiSeg, const edm::EventSetup& context) override;
0037 void endLuminosityBlock(const edm::LuminosityBlock& lumiSeg, const edm::EventSetup& c) override;
0038
0039 int ftotalevents;
0040 int fitNLumi_;
0041 int resetFitNLumi_;
0042 int countLumi_;
0043 int org_resetFitNLumi_;
0044 int previousLumi_;
0045 int previousRun_;
0046 int ftmprun0, ftmprun;
0047 int beginLumiOfBSFit_;
0048 int endLumiOfBSFit_;
0049 std::time_t refBStime[2];
0050
0051 bool write2DB_;
0052 bool runbeamwidthfit_;
0053 bool runallfitters_;
0054
0055 BeamFitter* theBeamFitter;
0056 };
0057
0058 BeamSpotAnalyzer::BeamSpotAnalyzer(const edm::ParameterSet& iConfig) {
0059
0060 write2DB_ = iConfig.getParameter<edm::ParameterSet>("BSAnalyzerParameters").getParameter<bool>("WriteToDB");
0061 runallfitters_ = iConfig.getParameter<edm::ParameterSet>("BSAnalyzerParameters").getParameter<bool>("RunAllFitters");
0062 fitNLumi_ =
0063 iConfig.getParameter<edm::ParameterSet>("BSAnalyzerParameters").getUntrackedParameter<int>("fitEveryNLumi", -1);
0064 resetFitNLumi_ =
0065 iConfig.getParameter<edm::ParameterSet>("BSAnalyzerParameters").getUntrackedParameter<int>("resetEveryNLumi", -1);
0066 runbeamwidthfit_ =
0067 iConfig.getParameter<edm::ParameterSet>("BSAnalyzerParameters").getParameter<bool>("RunBeamWidthFit");
0068
0069 theBeamFitter = new BeamFitter(iConfig, consumesCollector());
0070 theBeamFitter->resetTrkVector();
0071 theBeamFitter->resetLSRange();
0072 theBeamFitter->resetCutFlow();
0073 theBeamFitter->resetRefTime();
0074 theBeamFitter->resetPVFitter();
0075
0076 ftotalevents = 0;
0077 ftmprun0 = ftmprun = -1;
0078 countLumi_ = 0;
0079 beginLumiOfBSFit_ = endLumiOfBSFit_ = -1;
0080 previousLumi_ = previousRun_ = 0;
0081 org_resetFitNLumi_ = resetFitNLumi_;
0082 }
0083
0084 BeamSpotAnalyzer::~BeamSpotAnalyzer() { delete theBeamFitter; }
0085
0086 void BeamSpotAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0087 ftotalevents++;
0088 theBeamFitter->readEvent(iEvent);
0089 ftmprun = iEvent.id().run();
0090 }
0091
0092
0093 void BeamSpotAnalyzer::beginLuminosityBlock(const edm::LuminosityBlock& lumiSeg, const edm::EventSetup& context) {
0094 const edm::TimeValue_t fbegintimestamp = lumiSeg.beginTime().value();
0095 const std::time_t ftmptime = fbegintimestamp >> 32;
0096
0097 if (countLumi_ == 0 || (resetFitNLumi_ > 0 && countLumi_ % resetFitNLumi_ == 0)) {
0098 ftmprun0 = lumiSeg.run();
0099 ftmprun = ftmprun0;
0100 beginLumiOfBSFit_ = lumiSeg.luminosityBlock();
0101 refBStime[0] = ftmptime;
0102 }
0103
0104 countLumi_++;
0105 if (ftmprun == previousRun_) {
0106 if ((previousLumi_ + 1) != int(lumiSeg.luminosityBlock()))
0107 edm::LogWarning("BeamSpotAnalyzer") << "LUMI SECTIONS ARE NOT SORTED!";
0108 }
0109 }
0110
0111
0112 void BeamSpotAnalyzer::endLuminosityBlock(const edm::LuminosityBlock& lumiSeg, const edm::EventSetup& iSetup) {
0113 edm::LogPrint("BeamSpotAnalyzer") << "for lumis " << beginLumiOfBSFit_ << " - " << endLumiOfBSFit_ << std::endl
0114 << "number of selected tracks = " << theBeamFitter->getNTracks();
0115 edm::LogPrint("BeamSpotAnalyzer") << "number of selected PVs = " << theBeamFitter->getNPVs();
0116
0117
0118 const edm::TimeValue_t fendtimestamp = lumiSeg.endTime().value();
0119 const std::time_t fendtime = fendtimestamp >> 32;
0120 refBStime[1] = fendtime;
0121
0122 endLumiOfBSFit_ = lumiSeg.luminosityBlock();
0123 previousLumi_ = endLumiOfBSFit_;
0124
0125 if (fitNLumi_ == -1 && resetFitNLumi_ == -1)
0126 return;
0127
0128 if (fitNLumi_ > 0 && countLumi_ % fitNLumi_ != 0)
0129 return;
0130
0131 theBeamFitter->setFitLSRange(beginLumiOfBSFit_, endLumiOfBSFit_);
0132 theBeamFitter->setRefTime(refBStime[0], refBStime[1]);
0133 theBeamFitter->setRun(ftmprun0);
0134
0135 std::pair<int, int> LSRange = theBeamFitter->getFitLSRange();
0136
0137 if (theBeamFitter->runPVandTrkFitter()) {
0138 reco::BeamSpot bs = theBeamFitter->getBeamSpot();
0139 edm::LogPrint("BeamSpotAnalyzer") << "\n RESULTS OF DEFAULT FIT ";
0140 edm::LogPrint("BeamSpotAnalyzer") << " for runs: " << ftmprun0 << " - " << ftmprun;
0141 edm::LogPrint("BeamSpotAnalyzer") << " for lumi blocks : " << LSRange.first << " - " << LSRange.second;
0142 edm::LogPrint("BeamSpotAnalyzer") << " lumi counter # " << countLumi_;
0143 edm::LogPrint("BeamSpotAnalyzer") << bs;
0144 edm::LogPrint("BeamSpotAnalyzer") << "[BeamFitter] fit done. \n";
0145 } else {
0146 reco::BeamSpot bs;
0147 bs.setType(reco::BeamSpot::Fake);
0148 edm::LogPrint("BeamSpotAnalyzer") << "\n Empty Beam spot fit";
0149 edm::LogPrint("BeamSpotAnalyzer") << " for runs: " << ftmprun0 << " - " << ftmprun;
0150 edm::LogPrint("BeamSpotAnalyzer") << " for lumi blocks : " << LSRange.first << " - " << LSRange.second;
0151 edm::LogPrint("BeamSpotAnalyzer") << " lumi counter # " << countLumi_;
0152 edm::LogPrint("BeamSpotAnalyzer") << bs;
0153 edm::LogPrint("BeamSpotAnalyzer") << "[BeamFitter] fit failed \n";
0154
0155
0156
0157
0158 }
0159
0160 if (resetFitNLumi_ > 0 && countLumi_ % resetFitNLumi_ == 0) {
0161 std::vector<BSTrkParameters> theBSvector = theBeamFitter->getBSvector();
0162 edm::LogPrint("BeamSpotAnalyzer") << "Total number of tracks accumulated = " << theBSvector.size();
0163 edm::LogPrint("BeamSpotAnalyzer") << "Reset track collection for beam fit";
0164 theBeamFitter->resetTrkVector();
0165 theBeamFitter->resetLSRange();
0166 theBeamFitter->resetCutFlow();
0167 theBeamFitter->resetRefTime();
0168 theBeamFitter->resetPVFitter();
0169 countLumi_ = 0;
0170
0171 resetFitNLumi_ = org_resetFitNLumi_;
0172 }
0173 }
0174
0175 void BeamSpotAnalyzer::endJob() {
0176 edm::LogPrint("BeamSpotAnalyzer") << "\n-------------------------------------\n";
0177 edm::LogPrint("BeamSpotAnalyzer") << "\n Total number of events processed: " << ftotalevents;
0178 edm::LogPrint("BeamSpotAnalyzer") << "\n-------------------------------------\n\n";
0179
0180 if (fitNLumi_ == -1 && resetFitNLumi_ == -1) {
0181 if (theBeamFitter->runPVandTrkFitter()) {
0182 reco::BeamSpot beam_default = theBeamFitter->getBeamSpot();
0183 std::pair<int, int> LSRange = theBeamFitter->getFitLSRange();
0184
0185 edm::LogPrint("BeamSpotAnalyzer") << "\n RESULTS OF DEFAULT FIT:";
0186 edm::LogPrint("BeamSpotAnalyzer") << " for runs: " << ftmprun0 << " - " << ftmprun;
0187 edm::LogPrint("BeamSpotAnalyzer") << " for lumi blocks : " << LSRange.first << " - " << LSRange.second;
0188 edm::LogPrint("BeamSpotAnalyzer") << " lumi counter # " << countLumi_;
0189 edm::LogPrint("BeamSpotAnalyzer") << beam_default;
0190
0191 if (write2DB_) {
0192 edm::LogPrint("BeamSpotAnalyzer") << "\n-------------------------------------\n\n";
0193 edm::LogPrint("BeamSpotAnalyzer") << " write results to DB...";
0194 theBeamFitter->write2DB();
0195 }
0196
0197 if (runallfitters_) {
0198 theBeamFitter->runAllFitter();
0199 }
0200 }
0201 if ((runbeamwidthfit_)) {
0202 theBeamFitter->runBeamWidthFitter();
0203 reco::BeamSpot beam_width = theBeamFitter->getBeamWidth();
0204 edm::LogPrint("BeamSpotAnalyzer") << beam_width;
0205 } else {
0206 edm::LogPrint("BeamSpotAnalyzer") << "[BeamSpotAnalyzer] beamfit fails !!!";
0207 }
0208 }
0209
0210 edm::LogPrint("BeamSpotAnalyzer") << "[BeamSpotAnalyzer] endJob done \n";
0211 }
0212
0213 void BeamSpotAnalyzer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0214 edm::ParameterSetDescription desc;
0215 desc.setComment("Analyzer of BeamSpot Objects");
0216
0217 edm::ParameterSetDescription bsAnalyzerParamsDesc;
0218 bsAnalyzerParamsDesc.add("WriteToDB", false);
0219 bsAnalyzerParamsDesc.add("RunAllFitters", false);
0220 bsAnalyzerParamsDesc.addUntracked("fitEveryNLumi", -1);
0221 bsAnalyzerParamsDesc.addUntracked("resetEveryNLumi", -1);
0222 bsAnalyzerParamsDesc.add("RunBeamWidthFit", false);
0223 desc.add<edm::ParameterSetDescription>("BSAnalyzerParameters", bsAnalyzerParamsDesc);
0224
0225 BeamFitter::fillDescription(desc);
0226 PVFitter::fillDescription(desc);
0227
0228 descriptions.addWithDefaultLabel(desc);
0229 }
0230
0231
0232 DEFINE_FWK_MODULE(BeamSpotAnalyzer);