File indexing completed on 2023-03-17 11:21:19
0001 #ifndef RecoParticleFlow_PFTracking_PFDisplacedVertexHelper_h
0002 #define RecoParticleFlow_PFTracking_PFDisplacedVertexHelper_h
0003
0004 #include "FWCore/Framework/interface/ESHandle.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0007 #include "DataFormats/ParticleFlowReco/interface/PFDisplacedVertex.h"
0008 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0009 #include "DataFormats/TrackReco/interface/Track.h"
0010
0011
0012
0013
0014
0015
0016
0017 class PFDisplacedVertexHelper {
0018 public:
0019 PFDisplacedVertexHelper();
0020 ~PFDisplacedVertexHelper();
0021
0022
0023 void setTracksSelector(const edm::ParameterSet& ps) { tracksSelector_ = TracksSelector(ps); };
0024
0025
0026 void setVertexIdentifier(const edm::ParameterSet& ps) { vertexIdentifier_ = VertexIdentifier(ps); };
0027
0028
0029 void setPrimaryVertex(edm::Handle<reco::VertexCollection> mainVertexHandle,
0030 edm::Handle<reco::BeamSpot> beamSpotHandle);
0031
0032
0033 bool isTrackSelected(const reco::Track& trk, const reco::PFDisplacedVertex::VertexTrackType vertexTrackType) const;
0034
0035
0036 reco::PFDisplacedVertex::VertexType identifyVertex(const reco::PFDisplacedVertex& v) const;
0037
0038
0039 math::XYZPoint primaryVertex() const { return pvtx_; }
0040
0041 void Dump(std::ostream& out = std::cout) const;
0042
0043 private:
0044
0045 int lambdaCP(const reco::PFDisplacedVertex& v) const;
0046 bool isKaonMass(const reco::PFDisplacedVertex& v) const;
0047
0048
0049 struct TracksSelector {
0050 TracksSelector()
0051 : bSelectTracks_(false),
0052 nChi2_min_(0),
0053 nChi2_max_(100),
0054 pt_min_(0),
0055 dxy_min_(0),
0056 nHits_min_(3),
0057 nOuterHits_max_(100),
0058 quality_("loose") {}
0059
0060 TracksSelector(const edm::ParameterSet& ps) {
0061 bSelectTracks_ = ps.getParameter<bool>("bSelectTracks");
0062 nChi2_min_ = ps.getParameter<double>("nChi2_min");
0063 nChi2_max_ = ps.getParameter<double>("nChi2_max");
0064 pt_min_ = ps.getParameter<double>("pt_min");
0065 dxy_min_ = ps.getParameter<double>("dxy_min");
0066 nHits_min_ = ps.getParameter<int>("nHits_min");
0067 nOuterHits_max_ = ps.getParameter<int>("nOuterHits_max");
0068 std::string quality_ = ps.getParameter<std::string>("quality");
0069 }
0070
0071 bool selectTracks() const { return bSelectTracks_; }
0072 double nChi2_min() const { return nChi2_min_; }
0073 double nChi2_max() const { return nChi2_max_; }
0074 double pt_min() const { return pt_min_; }
0075 double dxy_min() const { return dxy_min_; }
0076 int nHits_min() const { return nHits_min_; }
0077 int nOuterHits_max() const { return nOuterHits_max_; }
0078 std::string quality() const { return quality_; }
0079 double dxy(const reco::Track& trk) const { return trk.dxy(pvtx_); }
0080
0081 bool bSelectTracks_;
0082 double nChi2_min_;
0083 double nChi2_max_;
0084 double pt_min_;
0085 double dxy_min_;
0086 int nHits_min_;
0087 int nOuterHits_max_;
0088 math::XYZPoint pvtx_;
0089 std::string quality_;
0090
0091 void Dump(std::ostream& out = std::cout) const {
0092 if (!out)
0093 return;
0094 std::string s = bSelectTracks_ ? "On" : "Off";
0095
0096 out << "" << std::endl;
0097 out << " ==== The TrackerSelector is " << s.data() << " ==== " << std::endl;
0098
0099 out << " nChi2_min_ = " << nChi2_min_ << " nChi2_max_ = " << nChi2_max_ << std::endl
0100 << " pt_min_ = " << pt_min_ << " dxy_min_ = " << dxy_min_ << std::endl
0101 << " nHits_min_ = " << nHits_min_ << " nOuterHits_max_ = " << nOuterHits_max_ << std::endl
0102 << " quality = " << quality_ << std::endl;
0103 }
0104 };
0105
0106
0107 struct VertexIdentifier {
0108 VertexIdentifier()
0109 : bIdentifyVertices_(false), pt_min_(0.2), pt_kink_min_(1.4), looper_eta_max_(0.1), logPrimSec_min_(0.2) {
0110 double m[] = {0.050, 0.470, 0.525, 0.470, 0.525, 1.107, 1.125, 0.200};
0111 std::vector<double> masses(m, m + 8);
0112 masses_ = masses;
0113
0114 double a[] = {60, 40};
0115 std::vector<double> angles(a, a + 1);
0116 angles_ = angles;
0117 };
0118
0119 VertexIdentifier(const edm::ParameterSet& ps) {
0120 bIdentifyVertices_ = ps.getParameter<bool>("bIdentifyVertices");
0121 angles_ = ps.getParameter<std::vector<double> >("angles");
0122 masses_ = ps.getParameter<std::vector<double> >("masses");
0123 pt_min_ = ps.getParameter<double>("pt_min");
0124 pt_kink_min_ = ps.getParameter<double>("pt_kink_min");
0125 looper_eta_max_ = ps.getParameter<double>("looper_eta_max");
0126 logPrimSec_min_ = ps.getParameter<double>("logPrimSec_min");
0127 }
0128
0129 bool identifyVertices() const { return bIdentifyVertices_; }
0130
0131 double angle_max() const { return angles_[0]; }
0132 double angle_V0Conv_max() const { return angles_[1]; }
0133
0134 double pt_min() const { return pt_min_; }
0135 double pt_kink_min() const { return pt_kink_min_; }
0136
0137 double mConv_max() const { return masses_[0]; }
0138 double mK0_min() const { return masses_[1]; }
0139 double mK0_max() const { return masses_[2]; }
0140 double mK_min() const { return masses_[3]; }
0141 double mK_max() const { return masses_[4]; }
0142 double mLambda_min() const { return masses_[5]; }
0143 double mLambda_max() const { return masses_[6]; }
0144 double mNucl_min() const { return masses_[7]; }
0145
0146 double looper_eta_max() const { return looper_eta_max_; }
0147 double logPrimSec_min() const { return logPrimSec_min_; }
0148
0149 bool bIdentifyVertices_;
0150 std::vector<double> angles_;
0151 std::vector<double> masses_;
0152 double pt_min_;
0153 double pt_kink_min_;
0154 double looper_eta_max_;
0155 double logPrimSec_min_;
0156
0157 void Dump(std::ostream& out = std::cout) const {
0158 if (!out)
0159 return;
0160 std::string s = bIdentifyVertices_ ? "On" : "Off";
0161 out << "" << std::endl;
0162 out << " ==== The Vertex Identifier is " << s.data() << " ==== " << std::endl;
0163
0164 out << " pt_min_ = " << pt_min_ << " pt_kink_min_ = " << pt_kink_min_ << std::endl
0165 << " looper_eta_max_ = " << looper_eta_max_ << " log10(P_Prim/P_Sec)_min " << logPrimSec_min_ << std::endl
0166 << " Mass_conv > " << mConv_max() << std::endl
0167 << " " << mK0_min() << " < Mass_K0 < " << mK0_max() << std::endl
0168 << " " << mK_min() << " < Mass_K+- < " << mK_max() << std::endl
0169 << " " << mLambda_min() << " < Mass_Lambda < " << mLambda_max() << std::endl
0170 << " Mass_Nucl_ee > " << mNucl_min() << std::endl
0171 << " angle_max = " << angle_max() << " angle_V0Conv_max = " << angle_V0Conv_max() << std::endl;
0172 }
0173 };
0174
0175 TracksSelector tracksSelector_;
0176 VertexIdentifier vertexIdentifier_;
0177
0178 math::XYZPoint pvtx_;
0179
0180
0181 static const double pion_mass2;
0182 static const double muon_mass2;
0183 static const double proton_mass2;
0184 };
0185
0186 #endif