File indexing completed on 2021-07-03 01:57:38
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014 #include "Alignment/CommonAlignmentMonitor/interface/AlignmentMonitorPluginFactory.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "FWCore/Utilities/interface/InputTag.h"
0017 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h"
0018 #include "DataFormats/GeometrySurface/interface/Surface.h"
0019 #include "TH1.h"
0020 #include "TObject.h"
0021 #include "Alignment/CommonAlignmentMonitor/interface/AlignmentMonitorBase.h"
0022
0023 #include "RecoMuon/TrackingTools/interface/MuonServiceProxy.h"
0024 #include "RecoMuon/TrackingTools/interface/MuonUpdatorAtVertex.h"
0025 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0026
0027 #include <fstream>
0028
0029
0030
0031
0032
0033
0034
0035 class AlignmentMonitorTracksFromTrajectories : public AlignmentMonitorBase {
0036 public:
0037 AlignmentMonitorTracksFromTrajectories(const edm::ParameterSet &cfg, edm::ConsumesCollector iC);
0038 ~AlignmentMonitorTracksFromTrajectories(){};
0039
0040 void book();
0041 void event(const edm::Event &iEvent, const edm::EventSetup &iSetup, const ConstTrajTrackPairCollection &iTrajTracks);
0042 void afterAlignment();
0043
0044 private:
0045 MuonServiceProxy *theMuonServiceProxy;
0046 MuonUpdatorAtVertex *theMuonUpdatorAtVertex;
0047 bool m_vertexConstraint;
0048 edm::InputTag m_beamSpot;
0049
0050 TH1F *m_diMuon_Z;
0051 TH1F *m_diMuon_Zforward;
0052 TH1F *m_diMuon_Zbarrel;
0053 TH1F *m_diMuon_Zbackward;
0054 TH1F *m_diMuon_Ups;
0055 TH1F *m_diMuon_Jpsi;
0056 TH1F *m_diMuon_log;
0057 TH1F *m_chi2_100;
0058 TH1F *m_chi2_10000;
0059 TH1F *m_chi2_1000000;
0060 TH1F *m_chi2_log;
0061 TH1F *m_chi2DOF_5;
0062 TH1F *m_chi2DOF_100;
0063 TH1F *m_chi2DOF_1000;
0064 TH1F *m_chi2DOF_log;
0065 TH1F *m_chi2_improvement;
0066 TH1F *m_chi2DOF_improvement;
0067 TH1F *m_pt[36];
0068 };
0069
0070
0071
0072
0073
0074
0075
0076
0077
0078
0079
0080
0081
0082
0083
0084
0085
0086
0087
0088
0089
0090
0091
0092
0093
0094
0095
0096
0097
0098
0099 AlignmentMonitorTracksFromTrajectories::AlignmentMonitorTracksFromTrajectories(const edm::ParameterSet &cfg,
0100 edm::ConsumesCollector iC)
0101 : AlignmentMonitorBase(cfg, iC, "AlignmentMonitorTracksFromTrajectories"),
0102 m_vertexConstraint(cfg.getParameter<bool>("vertexConstraint")),
0103 m_beamSpot(cfg.getParameter<edm::InputTag>("beamSpot")) {
0104 theMuonServiceProxy = new MuonServiceProxy(cfg.getParameter<edm::ParameterSet>("ServiceParameters"));
0105 theMuonUpdatorAtVertex = new MuonUpdatorAtVertex(cfg.getParameter<edm::ParameterSet>("MuonUpdatorAtVertexParameters"),
0106 theMuonServiceProxy);
0107 }
0108
0109
0110
0111
0112
0113
0114
0115
0116
0117 void AlignmentMonitorTracksFromTrajectories::book() {
0118 m_diMuon_Z = book1D("/iterN/", "diMuon_Z", "Di-muon mass (GeV)", 200, 90. - 50., 90. + 50.);
0119 m_diMuon_Zforward = book1D("/iterN/", "diMuon_Zforward", "Di-muon mass (GeV) eta > 1.4", 200, 90. - 50., 90. + 50.);
0120 m_diMuon_Zbarrel =
0121 book1D("/iterN/", "diMuon_Zbarrel", "Di-muon mass (GeV) -1.4 < eta < 1.4", 200, 90. - 50., 90. + 50.);
0122 m_diMuon_Zbackward =
0123 book1D("/iterN/", "diMuon_Zbackward", "Di-muon mass (GeV) eta < -1.4", 200, 90. - 50., 90. + 50.);
0124 m_diMuon_Ups = book1D("/iterN/", "diMuon_Ups", "Di-muon mass (GeV)", 200, 9. - 3., 9. + 3.);
0125 m_diMuon_Jpsi = book1D("/iterN/", "diMuon_Jpsi", "Di-muon mass (GeV)", 200, 3. - 1., 3. + 1.);
0126 m_diMuon_log = book1D("/iterN/", "diMuon_log", "Di-muon mass (log GeV)", 300, 0, 3);
0127 m_chi2_100 = book1D("/iterN/", "m_chi2_100", "Track chi^2", 100, 0, 100);
0128 m_chi2_10000 = book1D("/iterN/", "m_chi2_10000", "Track chi^2", 100, 0, 10000);
0129 m_chi2_1000000 = book1D("/iterN/", "m_chi2_1000000", "Track chi^2", 100, 1, 1000000);
0130 m_chi2_log = book1D("/iterN/", "m_chi2_log", "Log track chi^2", 100, -3, 7);
0131 m_chi2DOF_5 = book1D("/iterN/", "m_chi2DOF_5", "Track chi^2/nDOF", 100, 0, 5);
0132 m_chi2DOF_100 = book1D("/iterN/", "m_chi2DOF_100", "Track chi^2/nDOF", 100, 0, 100);
0133 m_chi2DOF_1000 = book1D("/iterN/", "m_chi2DOF_1000", "Track chi^2/nDOF", 100, 0, 1000);
0134 m_chi2DOF_log = book1D("/iterN/", "m_chi2DOF_log", "Log track chi^2/nDOF", 100, -3, 7);
0135 m_chi2_improvement = book1D("/iterN/", "m_chi2_improvement", "Track-by-track chi^2/original chi^2", 100, 0., 10.);
0136 m_chi2DOF_improvement = book1D(
0137 "/iterN/", "m_chi2DOF_improvement", "Track-by-track (chi^2/DOF)/(original chi^2/original DOF)", 100, 0., 10.);
0138 for (int i = 0; i < 36; i++) {
0139 char name[100], title[100];
0140 snprintf(name, sizeof(name), "m_pt_phi%d", i);
0141 snprintf(title, sizeof(title), "Track pt (GeV) in phi bin %d/36", i);
0142 m_pt[i] = book1D("/iterN/", name, title, 100, 0, 100);
0143 }
0144 }
0145
0146
0147
0148
0149
0150 void AlignmentMonitorTracksFromTrajectories::event(const edm::Event &iEvent,
0151 const edm::EventSetup &iSetup,
0152 const ConstTrajTrackPairCollection &tracks) {
0153 theMuonServiceProxy->update(iSetup);
0154
0155 edm::Handle<reco::BeamSpot> beamSpot;
0156 iEvent.getByLabel(m_beamSpot, beamSpot);
0157
0158 GlobalVector p1, p2;
0159 double e1 = 0.;
0160 double e2 = 0.;
0161
0162 for (ConstTrajTrackPairCollection::const_iterator it = tracks.begin(); it != tracks.end(); ++it) {
0163 const Trajectory *traj = it->first;
0164 const reco::Track *track = it->second;
0165
0166 int nDOF = 0;
0167 TrajectoryStateOnSurface closestTSOS;
0168 double closest = 10000.;
0169
0170 std::vector<TrajectoryMeasurement> measurements = traj->measurements();
0171 for (std::vector<TrajectoryMeasurement>::const_iterator im = measurements.begin(); im != measurements.end(); ++im) {
0172 const TrajectoryMeasurement meas = *im;
0173 const TransientTrackingRecHit *hit = &(*meas.recHit());
0174
0175
0176 nDOF += hit->dimension();
0177
0178 TrajectoryStateOnSurface TSOS = meas.backwardPredictedState();
0179 GlobalPoint where = TSOS.surface().toGlobal(LocalPoint(0, 0, 0));
0180
0181 if (where.mag() < closest) {
0182 closest = where.mag();
0183 closestTSOS = TSOS;
0184 }
0185
0186 }
0187 nDOF -= 5;
0188
0189 if (closest != 10000.) {
0190 std::pair<bool, FreeTrajectoryState> state;
0191
0192 if (m_vertexConstraint) {
0193 state = theMuonUpdatorAtVertex->propagateWithUpdate(closestTSOS, *beamSpot);
0194
0195 } else {
0196 state = theMuonUpdatorAtVertex->propagate(closestTSOS, *beamSpot);
0197 }
0198
0199 if (state.first) {
0200 double chi2 = traj->chiSquared();
0201 double chi2DOF = chi2 / double(nDOF);
0202
0203 m_chi2_100->Fill(chi2);
0204 m_chi2_10000->Fill(chi2);
0205 m_chi2_1000000->Fill(chi2);
0206 m_chi2_log->Fill(log10(chi2));
0207
0208 m_chi2DOF_5->Fill(chi2DOF);
0209 m_chi2DOF_100->Fill(chi2DOF);
0210 m_chi2DOF_1000->Fill(chi2DOF);
0211 m_chi2DOF_log->Fill(log10(chi2DOF));
0212 m_chi2_improvement->Fill(chi2 / track->chi2());
0213 m_chi2DOF_improvement->Fill(chi2DOF / track->normalizedChi2());
0214
0215 GlobalVector momentum = state.second.momentum();
0216 double energy = momentum.mag();
0217
0218 if (energy > e1) {
0219 e2 = e1;
0220 e1 = energy;
0221 p2 = p1;
0222 p1 = momentum;
0223 } else if (energy > e2) {
0224 e2 = energy;
0225 p2 = momentum;
0226 }
0227
0228 double pt = momentum.perp();
0229 double phi = momentum.phi();
0230 while (phi < -M_PI)
0231 phi += 2. * M_PI;
0232 while (phi > M_PI)
0233 phi -= 2. * M_PI;
0234
0235 double phibin = (phi + M_PI) / (2. * M_PI) * 36.;
0236
0237 for (int i = 0; i < 36; i++) {
0238 if (phibin < i + 1) {
0239 m_pt[i]->Fill(pt);
0240 break;
0241 }
0242 }
0243
0244 }
0245 }
0246 }
0247
0248 if (e1 > 0. && e2 > 0.) {
0249 double energy_tot = e1 + e2;
0250 GlobalVector momentum_tot = p1 + p2;
0251 double mass = sqrt(energy_tot * energy_tot - momentum_tot.mag2());
0252 double eta = momentum_tot.eta();
0253
0254 m_diMuon_Z->Fill(mass);
0255 if (eta > 1.4)
0256 m_diMuon_Zforward->Fill(mass);
0257 else if (eta > -1.4)
0258 m_diMuon_Zbarrel->Fill(mass);
0259 else
0260 m_diMuon_Zbackward->Fill(mass);
0261
0262 m_diMuon_Ups->Fill(mass);
0263 m_diMuon_Jpsi->Fill(mass);
0264 m_diMuon_log->Fill(log10(mass));
0265 }
0266 }
0267
0268
0269
0270
0271
0272 void AlignmentMonitorTracksFromTrajectories::afterAlignment() {}
0273
0274
0275
0276
0277
0278
0279
0280
0281
0282
0283
0284
0285
0286 DEFINE_EDM_PLUGIN(AlignmentMonitorPluginFactory,
0287 AlignmentMonitorTracksFromTrajectories,
0288 "AlignmentMonitorTracksFromTrajectories");