File indexing completed on 2024-04-06 11:58:33
0001
0002
0003
0004
0005
0006
0007 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010 #include "FWCore/Framework/interface/EventSetup.h"
0011 #include "FWCore/Framework/interface/MakerMacros.h"
0012 #include "FWCore/Framework/interface/ESWatcher.h"
0013 #include "FWCore/Utilities/interface/ESGetToken.h"
0014 #include "FWCore/Utilities/interface/InputTag.h"
0015
0016 #include "CalibPPS/AlignmentRelative/interface/StraightTrackAlignment.h"
0017
0018 #include "CondFormats/AlignmentRecord/interface/RPRealAlignmentRecord.h"
0019 #include "Geometry/Records/interface/VeryForwardMisalignedGeometryRecord.h"
0020 #include "Geometry/Records/interface/VeryForwardRealGeometryRecord.h"
0021 #include "Geometry/VeryForwardGeometryBuilder/interface/CTPPSGeometry.h"
0022
0023 #include "DataFormats/Common/interface/DetSetVector.h"
0024 #include "DataFormats/CTPPSReco/interface/TotemRPRecHit.h"
0025 #include "DataFormats/CTPPSReco/interface/TotemRPUVPattern.h"
0026 #include "DataFormats/CTPPSReco/interface/CTPPSDiamondRecHit.h"
0027 #include "DataFormats/CTPPSReco/interface/CTPPSPixelRecHit.h"
0028 #include "DataFormats/CTPPSReco/interface/CTPPSPixelLocalTrack.h"
0029
0030
0031
0032
0033 class PPSStraightTrackAligner : public edm::one::EDAnalyzer<> {
0034 public:
0035 PPSStraightTrackAligner(const edm::ParameterSet &ps);
0036 ~PPSStraightTrackAligner() override {}
0037
0038 private:
0039 unsigned int verbosity_;
0040
0041 edm::InputTag tagUVPatternsStrip_;
0042 edm::InputTag tagDiamondHits_;
0043 edm::InputTag tagPixelHits_;
0044 edm::InputTag tagPixelLocalTracks_;
0045
0046 edm::EDGetTokenT<edm::DetSetVector<TotemRPUVPattern>> tokenUVPatternsStrip_;
0047 edm::EDGetTokenT<edm::DetSetVector<CTPPSDiamondRecHit>> tokenDiamondHits_;
0048 edm::EDGetTokenT<edm::DetSetVector<CTPPSPixelRecHit>> tokenPixelHits_;
0049 edm::EDGetTokenT<edm::DetSetVector<CTPPSPixelLocalTrack>> tokenPixelLocalTracks_;
0050
0051 edm::ESWatcher<VeryForwardRealGeometryRecord> geometryWatcher_;
0052
0053 edm::ESGetToken<CTPPSGeometry, VeryForwardRealGeometryRecord> esTokenRealGeometry_;
0054 edm::ESGetToken<CTPPSGeometry, VeryForwardMisalignedGeometryRecord> esTokenMisalignedGeometry_;
0055 edm::ESGetToken<CTPPSRPAlignmentCorrectionsData, RPRealAlignmentRecord> esTokenRealAlignment_;
0056
0057 bool worker_initialized_;
0058 StraightTrackAlignment worker_;
0059
0060 void analyze(const edm::Event &e, const edm::EventSetup &es) override;
0061
0062 void endJob() override;
0063 };
0064
0065 using namespace std;
0066 using namespace edm;
0067
0068
0069
0070 PPSStraightTrackAligner::PPSStraightTrackAligner(const ParameterSet &ps)
0071 : verbosity_(ps.getUntrackedParameter<unsigned int>("verbosity", 0)),
0072
0073 tagUVPatternsStrip_(ps.getParameter<edm::InputTag>("tagUVPatternsStrip")),
0074 tagDiamondHits_(ps.getParameter<edm::InputTag>("tagDiamondHits")),
0075 tagPixelHits_(ps.getParameter<edm::InputTag>("tagPixelHits")),
0076 tagPixelLocalTracks_(ps.getParameter<edm::InputTag>("tagPixelLocalTracks")),
0077
0078 tokenUVPatternsStrip_(consumes<DetSetVector<TotemRPUVPattern>>(tagUVPatternsStrip_)),
0079 tokenDiamondHits_(consumes<DetSetVector<CTPPSDiamondRecHit>>(tagDiamondHits_)),
0080 tokenPixelHits_(consumes<DetSetVector<CTPPSPixelRecHit>>(tagPixelHits_)),
0081 tokenPixelLocalTracks_(consumes<DetSetVector<CTPPSPixelLocalTrack>>(tagPixelLocalTracks_)),
0082
0083 esTokenRealGeometry_(esConsumes()),
0084 esTokenMisalignedGeometry_(esConsumes()),
0085 esTokenRealAlignment_(esConsumes()),
0086
0087 worker_initialized_(false),
0088 worker_(ps) {
0089 if (!tagPixelHits_.label().empty() && !tagPixelLocalTracks_.label().empty())
0090 LogWarning("PPS")
0091 << "Both tagPixelHits and tagPixelLocalTracks are not empty - most likely this not what you want.";
0092 }
0093
0094
0095
0096 void PPSStraightTrackAligner::analyze(const edm::Event &event, const edm::EventSetup &es) {
0097
0098 if (geometryWatcher_.check(es)) {
0099 if (worker_initialized_)
0100 throw cms::Exception("PPS") << "PPSStraightTrackAligner can't cope with changing geometry - change in event "
0101 << event.id() << endl;
0102 }
0103
0104
0105 if (!worker_initialized_) {
0106 auto hRealGeometry = es.getHandle(esTokenRealGeometry_);
0107 auto hMisalignedGeometry = es.getHandle(esTokenMisalignedGeometry_);
0108 auto hRealAlignment = es.getHandle(esTokenRealAlignment_);
0109
0110 worker_.begin(hRealAlignment, hRealGeometry, hMisalignedGeometry);
0111 worker_initialized_ = true;
0112 }
0113
0114
0115 DetSetVector<TotemRPUVPattern> defaultStripUVPatterns;
0116 const DetSetVector<TotemRPUVPattern> *pStripUVPatterns = &defaultStripUVPatterns;
0117 Handle<DetSetVector<TotemRPUVPattern>> inputStripUVPatterns;
0118 if (!tagUVPatternsStrip_.label().empty()) {
0119 event.getByToken(tokenUVPatternsStrip_, inputStripUVPatterns);
0120 pStripUVPatterns = &(*inputStripUVPatterns);
0121 }
0122
0123 DetSetVector<CTPPSDiamondRecHit> defaultDiamondHits;
0124 const DetSetVector<CTPPSDiamondRecHit> *pDiamondHits = &defaultDiamondHits;
0125 Handle<DetSetVector<CTPPSDiamondRecHit>> inputDiamondHits;
0126 if (!tagDiamondHits_.label().empty()) {
0127 event.getByToken(tokenDiamondHits_, inputDiamondHits);
0128 pDiamondHits = &(*inputDiamondHits);
0129 }
0130
0131 DetSetVector<CTPPSPixelRecHit> defaultPixelHits;
0132 const DetSetVector<CTPPSPixelRecHit> *pPixelHits = &defaultPixelHits;
0133 Handle<DetSetVector<CTPPSPixelRecHit>> inputPixelHits;
0134 if (!tagPixelHits_.label().empty()) {
0135 event.getByToken(tokenPixelHits_, inputPixelHits);
0136 pPixelHits = &(*inputPixelHits);
0137 }
0138
0139 DetSetVector<CTPPSPixelLocalTrack> defaultPixelLocalTracks;
0140 const DetSetVector<CTPPSPixelLocalTrack> *pPixelLocalTracks = &defaultPixelLocalTracks;
0141 Handle<DetSetVector<CTPPSPixelLocalTrack>> inputPixelLocalTracks;
0142 if (!tagPixelLocalTracks_.label().empty()) {
0143 event.getByToken(tokenPixelLocalTracks_, inputPixelLocalTracks);
0144 pPixelLocalTracks = &(*inputPixelLocalTracks);
0145 }
0146
0147
0148 worker_.processEvent(event.id(), *pStripUVPatterns, *pDiamondHits, *pPixelHits, *pPixelLocalTracks);
0149 }
0150
0151
0152
0153 void PPSStraightTrackAligner::endJob() {
0154 if (worker_initialized_)
0155 worker_.finish();
0156 else
0157 throw cms::Exception("PPS") << "worker not initialized." << endl;
0158 }
0159
0160 DEFINE_FWK_MODULE(PPSStraightTrackAligner);