BeamSpotProducer/interface/PVFitter.h

0001 #ifndef PVFitter_H
0002 #define PVFitter_H
0003
0004 /**_________________________________________________________________
0005    class:   PVFitter.h
0006    package: RecoVertex/BeamSpotProducer
0007
0008
0009
0010  author: Francisco Yumiceva, Fermilab (yumiceva@fnal.gov)
0011          Geng-Yuan Jeng, UC Riverside (Geng-Yuan.Jeng@cern.ch)
0012
0013
0014  ________________________________________________________________**/
0015
0016 #include "FWCore/Framework/interface/Event.h"
0017 #include "DataFormats/Provenance/interface/Timestamp.h"
0018 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0019 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0020 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0021 #include "DataFormats/TrackReco/interface/TrackBase.h"
0022 #include "RecoVertex/BeamSpotProducer/interface/BSTrkParameters.h"
0023 #include "RecoVertex/BeamSpotProducer/interface/BSFitter.h"
0024
0025 #include "RecoVertex/BeamSpotProducer/interface/BeamSpotTreeData.h"
0026 #include "RecoVertex/BeamSpotProducer/interface/BeamSpotFitPVData.h"
0027
0028 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0029
0030 // ROOT
0031 #include "TFile.h"
0032 #include "TTree.h"
0033 #include "TH2F.h"
0034
0035 #include <fstream>
0036
0037 namespace edm {
0038   class ConsumesCollector;
0039 }
0040
0041 namespace reco {
0042   class Vertex;
0043 }
0044
0045 class PVFitter {
0046 public:
0047   PVFitter() {}
0048   PVFitter(const edm::ParameterSet& iConfig, edm::ConsumesCollector&& iColl);
0049   PVFitter(const edm::ParameterSet& iConfig, edm::ConsumesCollector& iColl);
0050   virtual ~PVFitter();
0051
0052   static void fillDescription(edm::ParameterSetDescription&);
0053
0054   void initialize(const edm::ParameterSet& iConfig, edm::ConsumesCollector& iColl);
0055   void readEvent(const edm::Event& iEvent);
0056   void setTree(TTree* tree);
0057
0058   double getWidthX() { return fwidthX; }
0059   double getWidthY() { return fwidthY; }
0060   double getWidthZ() { return fwidthZ; }
0061   double getWidthXerr() { return fwidthXerr; }
0062   double getWidthYerr() { return fwidthYerr; }
0063   double getWidthZerr() { return fwidthZerr; }
0064   //ssc
0065   std::vector<BeamSpotFitPVData> getpvStore() { return pvStore_; }
0066
0067   void FitPerBunchCrossing() { fFitPerBunchCrossing = true; }
0068   bool runBXFitter();
0069   bool runFitter();
0070   void resetLSRange() { fbeginLumiOfFit = fendLumiOfFit = -1; }
0071   void resetRefTime() { freftime[0] = freftime[1] = 0; }
0072   //ssc
0073   void setRefTime(std::time_t t0, std::time_t t1) {
0074     freftime[0] = t0;
0075     freftime[1] = t1;
0076   }
0077   void setFitLSRange(int ls0, int ls1) {
0078     fbeginLumiOfFit = ls0;
0079     fendLumiOfFit = ls1;
0080   }
0081
0082   void dumpTxtFile();
0083   void resetAll() {
0084     resetLSRange();
0085     resetRefTime();
0086     pvStore_.clear();
0087     bxMap_.clear();
0088     dynamicQualityCut_ = 1.e30;
0089     hPVx->Reset();
0090     hPVy->Reset();
0091     fbeamspot = reco::BeamSpot();
0092     fbspotMap.clear();
0093   };
0094   reco::BeamSpot getBeamSpot() { return fbeamspot; }
0095   std::map<int, reco::BeamSpot> getBeamSpotMap() { return fbspotMap; }
0096   bool IsFitPerBunchCrossing() { return fFitPerBunchCrossing; }
0097   int* getFitLSRange() {
0098     int* tmp = new int[2];
0099     tmp[0] = fbeginLumiOfFit;
0100     tmp[1] = fendLumiOfFit;
0101     return tmp;
0102   }
0103   //ssc
0104   time_t* getRefTime() {
0105     time_t* tmptime = new time_t[2];
0106     tmptime[0] = freftime[0];
0107     tmptime[1] = freftime[1];
0108     return tmptime;
0109   }
0110
0111   //ssc
0112   void resizepvStore(unsigned int rmSize) { pvStore_.erase(pvStore_.begin(), pvStore_.begin() + rmSize); }
0113
0114   /// reduce size of primary vertex cache by increasing quality limit
0115   void compressStore();
0116   /// vertex quality measure
0117   double pvQuality(const reco::Vertex& pv) const;
0118   /// vertex quality measure
0119   double pvQuality(const BeamSpotFitPVData& pv) const;
0120   int getNPVs() { return pvStore_.size(); }
0121
0122   const std::map<int, int>& getNPVsperBX() {
0123     npvsmap_.clear();
0124
0125     for (std::map<int, std::vector<BeamSpotFitPVData> >::const_iterator pvStore = bxMap_.begin();
0126          pvStore != bxMap_.end();
0127          ++pvStore) {
0128       npvsmap_[pvStore->first] = (pvStore->second).size();
0129     }
0130     return npvsmap_;
0131   }
0132
0133 private:
0134   std::map<int, int> npvsmap_;
0135   reco::BeamSpot fbeamspot;
0136   std::map<int, reco::BeamSpot> fbspotMap;
0137   bool fFitPerBunchCrossing;
0138   bool useOnlyFirstPV_;
0139
0140   std::ofstream fasciiFile;
0141
0142   bool debug_;
0143   bool do3DFit_;
0144   edm::EDGetTokenT<reco::VertexCollection> vertexToken_;
0145   std::string outputTxt_;
0146
0147   unsigned int maxNrVertices_;
0148   unsigned int minNrVertices_;
0149   double minVtxNdf_;
0150   double maxVtxNormChi2_;
0151   unsigned int minVtxTracks_;
0152   double minVtxWgt_;
0153   double maxVtxR_;
0154   double maxVtxZ_;
0155   double errorScale_;
0156   double sigmaCut_;
0157   double minSumPt_;
0158
0159   std::time_t freftime[2];
0160
0161   std::unique_ptr<TH2F> hPVx;
0162   std::unique_ptr<TH2F> hPVy;
0163
0164   TTree* ftree_;
0165
0166   //beam fit results
0167   int fbeginLumiOfFit;
0168   int fendLumiOfFit;
0169   double fwidthX;
0170   double fwidthY;
0171   double fwidthZ;
0172   double fwidthXerr;
0173   double fwidthYerr;
0174   double fwidthZerr;
0175
0176   std::vector<BeamSpotFitPVData> pvStore_;                //< cache for PV data
0177   std::map<int, std::vector<BeamSpotFitPVData> > bxMap_;  // store PV data as a function of bunch crossings
0178   double dynamicQualityCut_;                              //< quality cut for vertices (dynamic adjustment)
0179   std::vector<double> pvQualities_;                       //< working space for PV store adjustement
0180
0181   BeamSpotTreeData theBeamSpotTreeData_;
0182 };
0183
0184 #endif