BeamSpotProducer/interface/BeamFitter.h

0001 #ifndef BeamFitter_H
0002 #define BeamFitter_H
0003
0004 /**_________________________________________________________________
0005    class:   BeamFitter.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 "DataFormats/TrackReco/interface/Track.h"
0023 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0024 #include "DataFormats/TrackCandidate/interface/TrackCandidate.h"
0025 #include "DataFormats/TrackCandidate/interface/TrackCandidateCollection.h"
0026 #include "RecoVertex/BeamSpotProducer/interface/BSTrkParameters.h"
0027 #include "RecoVertex/BeamSpotProducer/interface/BSFitter.h"
0028 #include "RecoVertex/BeamSpotProducer/interface/PVFitter.h"
0029
0030 // ROOT
0031 #include "TFile.h"
0032 #include "TTree.h"
0033 #include "TH1F.h"
0034
0035 #include <fstream>
0036
0037 namespace edm {
0038   class ConsumesCollector;
0039 }
0040
0041 class BeamFitter {
0042 public:
0043   BeamFitter() {}
0044   BeamFitter(const edm::ParameterSet &iConfig, edm::ConsumesCollector &&iColl);
0045   virtual ~BeamFitter();
0046
0047   static void fillDescription(edm::ParameterSetDescription &);
0048
0049   void readEvent(const edm::Event &iEvent);
0050
0051   bool runFitter();
0052   bool runBeamWidthFitter();
0053   bool runPVandTrkFitter();
0054   bool runFitterNoTxt();
0055
0056   reco::BeamSpot getBeamWidth() { return fbeamWidthFit; }
0057   void runAllFitter();
0058   void resetTrkVector() { fBSvector.clear(); }
0059   void resetTotTrk() { ftotal_tracks = 0; }
0060   void resetLSRange() { fbeginLumiOfFit = fendLumiOfFit = -1; }
0061   void resetRefTime() { freftime[0] = freftime[1] = 0; }
0062   void setRefTime(time_t t0, time_t t1) {
0063     freftime[0] = t0;
0064     freftime[1] = t1;
0065     // Make sure the string representation of the time
0066     // is up-to-date
0067     updateBTime();
0068   }
0069
0070   std::pair<time_t, time_t> getRefTime() { return std::make_pair(freftime[0], freftime[1]); }
0071
0072   void resetPVFitter() { MyPVFitter->resetAll(); }
0073
0074   //---these are added to fasciliate BeamMonitor stuff for DIP
0075   std::size_t getPVvectorSize() { return (MyPVFitter->getpvStore()).size(); }
0076   //sc
0077   void resizeBSvector(unsigned int nsize) { fBSvector.erase(fBSvector.begin(), fBSvector.begin() + nsize); }
0078
0079   //ssc
0080   void resizePVvector(unsigned int npvsize) { MyPVFitter->resizepvStore(npvsize); }
0081
0082   //ssc
0083   void SetPVInfo(const std::vector<float> &v1_) {
0084     ForDIPPV_.clear();
0085     ForDIPPV_.assign(v1_.begin(), v1_.end());
0086   }
0087
0088   //----------------
0089
0090   void dumpTxtFile(std::string &, bool);
0091   void dumpBWTxtFile(std::string &);
0092   void write2DB();
0093   reco::BeamSpot getBeamSpot() { return fbeamspot; }
0094   std::map<int, reco::BeamSpot> getBeamSpotMap() { return fbspotPVMap; }
0095   std::vector<BSTrkParameters> getBSvector() { return fBSvector; }
0096   TH1F *getCutFlow() { return h1cutFlow; }
0097   void subtractFromCutFlow(const TH1F *toSubtract) {
0098     h1cutFlow->Add(toSubtract, -1.0);
0099     for (unsigned int i = 0; i < sizeof(countPass) / sizeof(countPass[0]); i++) {
0100       countPass[i] = h1cutFlow->GetBinContent(i + 1);
0101     }
0102   }
0103
0104   void resetCutFlow() {
0105     h1cutFlow->Reset();
0106     ftotal_tracks = 0;
0107     for (unsigned int i = 0; i < sizeof(countPass) / sizeof(countPass[0]); i++)
0108       countPass[i] = 0;
0109   }
0110
0111   //ssc
0112   int getRunNumber() { return frun; }
0113
0114   std::pair<int, int> getFitLSRange() { return std::make_pair(fbeginLumiOfFit, fendLumiOfFit); }
0115   void setFitLSRange(int ls0, int ls1) {
0116     fbeginLumiOfFit = ls0;
0117     fendLumiOfFit = ls1;
0118   }
0119   void setRun(int run) { frun = run; }
0120
0121   int getNTracks() { return fBSvector.size(); }
0122   int getNPVs() { return MyPVFitter->getNPVs(); }
0123   const std::map<int, int> &getNPVsperBX() { return MyPVFitter->getNPVsperBX(); }
0124
0125 private:
0126   // Update the fbeginTimeOfFit etc from the refTime
0127   void updateBTime();
0128   std::vector<BSTrkParameters> fBSvector;
0129   reco::BeamSpot fbeamspot;
0130   reco::BeamSpot fbeamWidthFit;
0131   std::map<int, reco::BeamSpot> fbspotPVMap;
0132   //  BSFitter *fmyalgo;
0133   std::ofstream fasciiFile;
0134   std::ofstream fasciiDIP;
0135
0136   bool debug_;
0137   bool appendRunTxt_;
0138   edm::EDGetTokenT<reco::TrackCollection> tracksToken_;
0139   edm::EDGetTokenT<edm::View<reco::Vertex> > vertexToken_;
0140   edm::EDGetTokenT<reco::BeamSpot> beamSpotToken_;  //offlineBeamSpot
0141   bool writeTxt_;
0142   bool writeDIPTxt_;
0143   bool writeDIPBadFit_;
0144   std::string outputTxt_;
0145   std::string outputDIPTxt_;
0146   double trk_MinpT_;
0147   double trk_MaxZ_;
0148   double trk_MaxEta_;
0149   double trk_MaxIP_;
0150   int trk_MinNTotLayers_;
0151   int trk_MinNPixLayers_;
0152   double trk_MaxNormChi2_;
0153   std::vector<std::string> trk_Algorithm_;
0154   std::vector<std::string> trk_Quality_;
0155   std::vector<reco::TrackBase::TrackQuality> quality_;
0156   std::vector<reco::TrackBase::TrackAlgorithm> algorithm_;
0157   double inputBeamWidth_;
0158   double convergence_;
0159   int ftotal_tracks;
0160   int min_Ntrks_;
0161   bool isMuon_;
0162   bool fitted_;
0163   bool ffilename_changed;
0164
0165   //ssc
0166   std::vector<float> ForDIPPV_;
0167
0168   // ntuple
0169   TH1F *h1z;
0170   bool saveNtuple_;
0171   bool saveBeamFit_;
0172   bool savePVVertices_;
0173   std::string outputfilename_;
0174   TFile *file_;
0175   TTree *ftree_;
0176   double ftheta;
0177   double fpt;
0178   double feta;
0179   int fcharge;
0180   double fnormchi2;
0181   double fphi0;
0182   double fd0;
0183   double fd0bs;
0184   double fsigmad0;
0185   double fz0;
0186   double fsigmaz0;
0187   int fnTotLayerMeas;
0188   int fnPixelLayerMeas;
0189   int fnStripLayerMeas;
0190   int fnTIBLayerMeas;
0191   int fnTIDLayerMeas;
0192   int fnTOBLayerMeas;
0193   int fnTECLayerMeas;
0194   int fnPXBLayerMeas;
0195   int fnPXFLayerMeas;
0196   //  double fd0phi_chi2;
0197   //  double fd0phi_d0;
0198   double fcov[7][7];
0199   double fvx;
0200   double fvy;
0201   int frun;
0202   int flumi;
0203   bool fquality;
0204   bool falgo;
0205   bool fpvValid;
0206   double fpvx, fpvy, fpvz;
0207   std::time_t freftime[2];
0208
0209   //beam fit results
0210   TTree *ftreeFit_;
0211   int frunFit;
0212   int fbeginLumiOfFit;
0213   int fendLumiOfFit;
0214   char fbeginTimeOfFit[32];
0215   char fendTimeOfFit[32];
0216   double fx;
0217   double fy;
0218   double fz;
0219   double fsigmaZ;
0220   double fdxdz;
0221   double fdydz;
0222   double fxErr;
0223   double fyErr;
0224   double fzErr;
0225   double fsigmaZErr;
0226   double fdxdzErr;
0227   double fdydzErr;
0228   double fwidthX;
0229   double fwidthY;
0230   double fwidthXErr;
0231   double fwidthYErr;
0232
0233   TH1F *h1ntrks;
0234   TH1F *h1vz_event;
0235   TH1F *h1cutFlow;
0236   int countPass[9];
0237
0238   PVFitter *MyPVFitter;
0239   TTree *fPVTree_;
0240 };
0241
0242 #endif