OfflineValidation/plugins/PixelBaryCentreAnalyzer.cc

0001 /** \class PixelBaryCentreAnalyzer
0002  *  The analyer works as the following :
0003  *  - Read global tracker position from global tag
0004  *  - Read tracker alignment constants from different ESsource with different labels
0005  *  - Calculate barycentres for different pixel substructures using global tracker position and alignment constants and store them in trees, one for each ESsource label.
0006  *
0007  *  Python script plotBaryCentre_VS_BeamSpot.py under script dir is used to plot barycentres from alignment constants used in Prompt-Reco, End-of-Year Rereco and so-called Run-2 (Ultra)Legacy Rereco. Options of the plotting script can be found from the helper in the script.
0008  *
0009  *  $Date: 2021/01/05 $
0010  *  $Last Modified: 2024/09/23 $
0011  *  $Revision: 1.1 $
0012  *  \author Tongguang Cheng - Beihang University <tongguang.cheng@cern.ch>
0013  *  \author Alberto Mecca - Torino University <alberto.mecca@cern.ch>
0014  *
0015 */
0016
0017 // Framework
0018 #include "FWCore/Framework/interface/Frameworkfwd.h"
0019 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0020 #include "FWCore/Framework/interface/Event.h"
0021 #include "FWCore/Framework/interface/MakerMacros.h"
0022 #include "FWCore/Framework/interface/ESWatcher.h"
0023 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0024 #include "FWCore/Utilities/interface/ESGetToken.h"
0025
0026 // Phase-1 Pixel
0027 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0028 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0029 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0030 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0031 #include "DataFormats/TrackerCommon/interface/PixelBarrelName.h"
0032 #include "DataFormats/TrackerCommon/interface/PixelEndcapName.h"
0033
0034 // pixel quality
0035 #include "CondFormats/SiPixelObjects/interface/SiPixelQuality.h"
0036 #include "CondFormats/DataRecord/interface/SiPixelQualityFromDbRcd.h"
0037 // global postion
0038 #include "CondFormats/Alignment/interface/DetectorGlobalPosition.h"
0039 #include "CondFormats/AlignmentRecord/interface/GlobalPositionRcd.h"
0040 // tracker alignment
0041 #include "CondFormats/AlignmentRecord/interface/TrackerAlignmentRcd.h"
0042 #include "CondFormats/Alignment/interface/Alignments.h"
0043 // beamspot
0044 #include "CondFormats/DataRecord/interface/BeamSpotObjectsRcd.h"
0045 #include "CondFormats/BeamSpotObjects/interface/BeamSpotObjects.h"
0046
0047 // Point and Vector
0048 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0049 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
0050
0051 // TFileService
0052 #include "FWCore/ServiceRegistry/interface/Service.h"
0053 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0054
0055 // ROOT
0056 #include "TTree.h"
0057 #include "TString.h"
0058
0059 //
0060 // class declaration
0061 //
0062
0063 class PixelBaryCentreAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0064 public:
0065   explicit PixelBaryCentreAnalyzer(const edm::ParameterSet&);
0066   ~PixelBaryCentreAnalyzer() override = default;
0067
0068   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0069
0070   struct SimplePoint {
0071     float x, y, z;
0072     SimplePoint(const GlobalPoint& p) : x(p.x()), y(p.y()), z(p.z()) {}
0073     SimplePoint() : x(0), y(0), z(0) {}
0074   };
0075   static const unsigned int nPixelLayers = 4;
0076   static const unsigned int nPixelDiscs = 3;
0077
0078 private:
0079   void beginJob() override;
0080   void analyze(const edm::Event&, const edm::EventSetup&) override;
0081   void endJob() override;
0082
0083   void initBC();
0084   void initBS();
0085
0086   const bool usePixelQuality_;
0087   int phase_;
0088
0089   // ----------member data ---------------------------
0090   edm::ESWatcher<BeamSpotObjectsRcd> watcherBS_;
0091   edm::ESWatcher<TrackerAlignmentRcd> watcherTkAlign_;
0092
0093   // labels of TkAlign tags
0094   std::vector<std::string> bcLabels_;
0095   // labels of beamspot tags
0096   std::vector<std::string> bsLabels_;
0097
0098   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> trackerGeometryToken_;
0099   const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> trackerTopologyToken_;
0100   const edm::ESGetToken<SiPixelQuality, SiPixelQualityFromDbRcd> siPixelQualityToken_;
0101
0102   const edm::ESGetToken<Alignments, GlobalPositionRcd> gprToken_;
0103   std::map<std::string, edm::ESGetToken<Alignments, TrackerAlignmentRcd>> tkAlignTokens_;
0104   std::map<std::string, edm::ESGetToken<BeamSpotObjects, BeamSpotObjectsRcd>> bsTokens_;
0105
0106   // tree content
0107   int run_;
0108   int ls_;
0109
0110   GlobalPoint BS_;
0111
0112   GlobalPoint PIX_, BPIX_, FPIX_;
0113   GlobalPoint BPIX_Flipped_, BPIX_NonFlipped_, BPIX_DiffFlippedNonFlipped_;
0114
0115   GlobalPoint BPIXLayer_[nPixelLayers];
0116   GlobalPoint BPIXLayer_Flipped_[nPixelLayers];
0117   GlobalPoint BPIXLayer_NonFlipped_[nPixelLayers];
0118   GlobalPoint BPIXLayer_DiffFlippedNonFlipped_[nPixelLayers];
0119
0120   GlobalPoint FPIX_plus_, FPIX_minus_;
0121   GlobalPoint FPIXDisks_plus_[nPixelDiscs];
0122   GlobalPoint FPIXDisks_minus_[nPixelDiscs];
0123
0124   SimplePoint vBS_;
0125
0126   SimplePoint vPIX_, vBPIX_, vFPIX_;
0127   SimplePoint vBPIX_Flipped_, vBPIX_NonFlipped_, vBPIX_DiffFlippedNonFlipped_;
0128
0129   SimplePoint vBPIXLayer_[nPixelLayers];
0130   SimplePoint vBPIXLayer_Flipped_[nPixelLayers];
0131   SimplePoint vBPIXLayer_NonFlipped_[nPixelLayers];
0132   SimplePoint vBPIXLayer_DiffFlippedNonFlipped_[nPixelLayers];
0133
0134   SimplePoint vFPIX_plus_, vFPIX_minus_;
0135   SimplePoint vFPIXDisks_plus_[nPixelDiscs];
0136   SimplePoint vFPIXDisks_minus_[nPixelDiscs];
0137
0138   edm::Service<TFileService> tFileService;
0139   std::map<std::string, TTree*> bcTrees_;
0140   std::map<std::string, TTree*> bsTrees_;
0141 };
0142
0143 //
0144 // constructors and destructor
0145 //
0146 PixelBaryCentreAnalyzer::PixelBaryCentreAnalyzer(const edm::ParameterSet& iConfig)
0147     : usePixelQuality_(iConfig.getUntrackedParameter<bool>("usePixelQuality")),
0148       bcLabels_(iConfig.getUntrackedParameter<std::vector<std::string>>("tkAlignLabels")),
0149       bsLabels_(iConfig.getUntrackedParameter<std::vector<std::string>>("beamSpotLabels")),
0150       trackerGeometryToken_(esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>()),
0151       trackerTopologyToken_(esConsumes<TrackerTopology, TrackerTopologyRcd>()),
0152       siPixelQualityToken_(esConsumes<SiPixelQuality, SiPixelQualityFromDbRcd>()),
0153       gprToken_(esConsumes<Alignments, GlobalPositionRcd>()) {
0154   for (const auto& label : bcLabels_) {
0155     bcTrees_[label] = nullptr;
0156     tkAlignTokens_[label] = esConsumes<Alignments, TrackerAlignmentRcd>(edm::ESInputTag{"", label});
0157   }
0158
0159   for (const auto& label : bsLabels_) {
0160     bsTrees_[label] = nullptr;
0161     bsTokens_[label] = esConsumes<BeamSpotObjects, BeamSpotObjectsRcd>(edm::ESInputTag{"", label});
0162   }
0163
0164   usesResource("TFileService");
0165 }
0166
0167 //
0168 // member functions
0169 //
0170
0171 void PixelBaryCentreAnalyzer::initBS() {
0172   double dummy_float = 999999.0;
0173
0174   BS_ = GlobalPoint(dummy_float, dummy_float, dummy_float);
0175   vBS_ = SimplePoint(BS_);
0176 }
0177
0178 void PixelBaryCentreAnalyzer::initBC() {
0179   // init to large number (unreasonable number) not zero
0180   double dummy_float = 999999.0;
0181
0182   PIX_ = GlobalPoint(dummy_float, dummy_float, dummy_float);
0183   BPIX_ = GlobalPoint(dummy_float, dummy_float, dummy_float);
0184   FPIX_ = GlobalPoint(dummy_float, dummy_float, dummy_float);
0185
0186   BPIX_Flipped_ = GlobalPoint(dummy_float, dummy_float, dummy_float);
0187   BPIX_NonFlipped_ = GlobalPoint(dummy_float, dummy_float, dummy_float);
0188   BPIX_DiffFlippedNonFlipped_ = GlobalPoint(dummy_float, dummy_float, dummy_float);
0189
0190   FPIX_plus_ = GlobalPoint(dummy_float, dummy_float, dummy_float);
0191   FPIX_minus_ = GlobalPoint(dummy_float, dummy_float, dummy_float);
0192
0193   for (unsigned int i = 0; i < nPixelLayers; i++) {
0194     BPIXLayer_[i] = GlobalPoint(dummy_float, dummy_float, dummy_float);
0195     BPIXLayer_Flipped_[i] = GlobalPoint(dummy_float, dummy_float, dummy_float);
0196     BPIXLayer_NonFlipped_[i] = GlobalPoint(dummy_float, dummy_float, dummy_float);
0197     BPIXLayer_DiffFlippedNonFlipped_[i] = GlobalPoint(dummy_float, dummy_float, dummy_float);
0198   }
0199
0200   for (unsigned int i = 0; i < nPixelDiscs; i++) {
0201     FPIXDisks_plus_[i] = GlobalPoint(dummy_float, dummy_float, dummy_float);
0202     FPIXDisks_minus_[i] = GlobalPoint(dummy_float, dummy_float, dummy_float);
0203   }
0204
0205   vPIX_ = SimplePoint(PIX_);
0206   vBPIX_ = SimplePoint(BPIX_);
0207   vFPIX_ = SimplePoint(FPIX_);
0208
0209   vBPIX_Flipped_ = SimplePoint(BPIX_Flipped_);
0210   vBPIX_NonFlipped_ = SimplePoint(BPIX_NonFlipped_);
0211   vBPIX_DiffFlippedNonFlipped_ = SimplePoint(BPIX_DiffFlippedNonFlipped_);
0212
0213   vFPIX_plus_ = SimplePoint(FPIX_plus_);
0214   vFPIX_minus_ = SimplePoint(FPIX_minus_);
0215
0216   for (unsigned int i = 0; i < nPixelLayers; i++) {
0217     vBPIXLayer_[i] = SimplePoint(BPIXLayer_[i]);
0218     vBPIXLayer_Flipped_[i] = SimplePoint(BPIXLayer_Flipped_[i]);
0219     vBPIXLayer_NonFlipped_[i] = SimplePoint(BPIXLayer_NonFlipped_[i]);
0220     vBPIXLayer_DiffFlippedNonFlipped_[i] = SimplePoint(BPIXLayer_DiffFlippedNonFlipped_[i]);
0221   }
0222
0223   for (unsigned int i = 0; i < nPixelDiscs; i++) {
0224     vFPIXDisks_plus_[i] = SimplePoint(FPIXDisks_plus_[i]);
0225     vFPIXDisks_minus_[i] = SimplePoint(FPIXDisks_minus_[i]);
0226   }
0227 }
0228
0229 // ------------ method called for each event  ------------
0230 void PixelBaryCentreAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0231   bool prepareTkAlign = false;
0232   bool prepareBS = false;
0233
0234   // ES watcher can noly run once in the same event,
0235   // otherwise it will turn false whatsoever because the condition doesn't change in the second time call.
0236   if (watcherTkAlign_.check(iSetup))
0237     prepareTkAlign = true;
0238   if (watcherBS_.check(iSetup))
0239     prepareBS = true;
0240
0241   if (!prepareTkAlign && !prepareBS)
0242     return;
0243
0244   run_ = iEvent.id().run();
0245   ls_ = iEvent.id().luminosityBlock();
0246
0247   if (prepareTkAlign) {  // check for new IOV for TKAlign
0248
0249     phase_ = -1;
0250
0251     const TrackerGeometry& tkGeom = iSetup.getData(trackerGeometryToken_);
0252     const TrackerTopology& tkTopo = iSetup.getData(trackerTopologyToken_);
0253
0254     if (tkGeom.isThere(GeomDetEnumerators::PixelBarrel) && tkGeom.isThere(GeomDetEnumerators::PixelEndcap)) {
0255       phase_ = 0;
0256     } else if (tkGeom.isThere(GeomDetEnumerators::P1PXB) && tkGeom.isThere(GeomDetEnumerators::P1PXEC)) {
0257       phase_ = 1;
0258     }
0259
0260     // pixel quality
0261     const SiPixelQuality& badPixelInfo = iSetup.getData(siPixelQualityToken_);
0262
0263     // Tracker global position
0264     const Alignments& globalAlignments = iSetup.getData(gprToken_);
0265     const AlignTransform globalCoordinates = align::DetectorGlobalPosition(globalAlignments, DetId(DetId::Tracker));
0266     GlobalVector globalTkPosition(
0267         globalCoordinates.translation().x(), globalCoordinates.translation().y(), globalCoordinates.translation().z());
0268
0269     // loop over bclabels
0270     for (const auto& label : bcLabels_) {
0271       // init tree content
0272       PixelBaryCentreAnalyzer::initBC();
0273
0274       // Get TkAlign from EventSetup:
0275       const Alignments& alignments = iSetup.getData(tkAlignTokens_[label]);
0276       std::vector<AlignTransform> tkAlignments = alignments.m_align;
0277
0278       // PIX
0279       GlobalVector barycentre_PIX(0.0, 0.0, 0.0);
0280       // BPIX
0281       GlobalVector barycentre_BPIX(0.0, 0.0, 0.0);
0282       float nmodules_BPIX(0.);
0283       // FPIX
0284       GlobalVector barycentre_FPIX(0.0, 0.0, 0.0);
0285       float nmodules_FPIX(0.);
0286
0287       // Per-layer/ladder barycentre for BPIX
0288       std::map<int, std::map<int, float>> nmodules_bpix;           // layer-ladder
0289       std::map<int, std::map<int, GlobalVector>> barycentre_bpix;  // layer-ladder
0290
0291       // Per-disk/ring barycentre for FPIX
0292       std::map<int, std::map<int, float>> nmodules_fpix;           // disk-ring
0293       std::map<int, std::map<int, GlobalVector>> barycentre_fpix;  // disk-ring
0294
0295       // Loop over tracker module
0296       for (const auto& ali : tkAlignments) {
0297         //DetId
0298         const DetId& detId = DetId(ali.rawId());
0299         // remove bad module
0300         if (usePixelQuality_ && badPixelInfo.IsModuleBad(detId))
0301           continue;
0302
0303         // alignment for a given module
0304         GlobalVector ali_translation(ali.translation().x(), ali.translation().y(), ali.translation().z());
0305
0306         int subid = DetId(detId).subdetId();
0307         // BPIX
0308         if (subid == PixelSubdetector::PixelBarrel) {
0309           nmodules_BPIX += 1;
0310           barycentre_BPIX += ali_translation;
0311           barycentre_PIX += ali_translation;
0312
0313           int layer = tkTopo.pxbLayer(detId);
0314           int ladder = tkTopo.pxbLadder(detId);
0315           nmodules_bpix[layer][ladder] += 1;
0316           barycentre_bpix[layer][ladder] += ali_translation;
0317
0318         }  // BPIX
0319
0320         // FPIX
0321         if (subid == PixelSubdetector::PixelEndcap) {
0322           nmodules_FPIX += 1;
0323           barycentre_FPIX += ali_translation;
0324           barycentre_PIX += ali_translation;
0325
0326           int disk = tkTopo.pxfDisk(detId);
0327           int quadrant = PixelEndcapName(detId, &tkTopo, phase_).halfCylinder();
0328           if (quadrant < 3)
0329             disk *= -1;
0330
0331           int ring = -9999;
0332           if (phase_ == 0) {
0333             ring = 1 + (tkTopo.pxfPanel(detId) + tkTopo.pxfModule(detId.rawId()) > 3);
0334           } else if (phase_ == 1) {
0335             ring = PixelEndcapName(detId, &tkTopo, phase_).ringName();
0336           }
0337
0338           nmodules_fpix[disk][ring] += 1;
0339           barycentre_fpix[disk][ring] += ali_translation;
0340
0341         }  // FPIX
0342
0343       }  // loop over tracker module
0344
0345       //PIX
0346       float nmodules_PIX = nmodules_BPIX + nmodules_FPIX;
0347       barycentre_PIX *= (1.0 / nmodules_PIX);
0348       barycentre_PIX += globalTkPosition;
0349       PIX_ = GlobalPoint(barycentre_PIX.x(), barycentre_PIX.y(), barycentre_PIX.z());
0350       vPIX_ = SimplePoint(PIX_);
0351
0352       //BPIX
0353       barycentre_BPIX *= (1.0 / nmodules_BPIX);
0354       barycentre_BPIX += globalTkPosition;
0355       BPIX_ = GlobalPoint(barycentre_BPIX.x(), barycentre_BPIX.y(), barycentre_BPIX.z());
0356       vBPIX_ = SimplePoint(BPIX_);
0357       //FPIX
0358       barycentre_FPIX *= (1.0 / nmodules_FPIX);
0359       barycentre_FPIX += globalTkPosition;
0360       FPIX_ = GlobalPoint(barycentre_FPIX.x(), barycentre_FPIX.y(), barycentre_FPIX.z());
0361       vFPIX_ = SimplePoint(FPIX_);
0362       // Pixel substructures
0363
0364       // BPix barycentre per-layer/per-ladder
0365       // assuming each ladder has the same number of modules in the same layer
0366       // inner =  flipped; outer = non-flipped
0367       //
0368       // Phase 0: Outer ladders are odd for layer 1,3 and even for layer 2
0369       // Phase 1: Outer ladders are odd for layer 4 and even for layer 1,2,3
0370       //
0371
0372       int nmodules_BPIX_Flipped = 0;
0373       int nmodules_BPIX_NonFlipped = 0;
0374       GlobalVector BPIX_Flipped(0.0, 0.0, 0.0);
0375       GlobalVector BPIX_NonFlipped(0.0, 0.0, 0.0);
0376
0377       // loop over layers
0378       for (const auto& il : barycentre_bpix) {
0379         int layer = il.first;
0380
0381         int nmodulesLayer = 0;
0382         int nmodulesLayer_Flipped = 0;
0383         int nmodulesLayer_NonFlipped = 0;
0384         GlobalVector BPIXLayer(0.0, 0.0, 0.0);
0385         GlobalVector BPIXLayer_Flipped(0.0, 0.0, 0.0);
0386         GlobalVector BPIXLayer_NonFlipped(0.0, 0.0, 0.0);
0387
0388         // loop over ladder
0389         std::map<int, GlobalVector> barycentreLayer = barycentre_bpix[layer];
0390         for (const auto& it : barycentreLayer) {
0391           int ladder = it.first;
0392           //BPIXLayerLadder_[layer][ladder] = (1.0/nmodules[layer][ladder])*barycentreLayer[ladder] + globalTkPosition;
0393
0394           nmodulesLayer += nmodules_bpix[layer][ladder];
0395           BPIXLayer += barycentreLayer[ladder];
0396
0397           // Phase-1
0398           //
0399           // Phase 1: Outer ladders are odd for layer 4 and even for layer 1,2,3
0400           if (phase_ == 1) {
0401             if (layer != 4) {  // layer 1-3
0402
0403               if (ladder % 2 != 0) {  // odd ladder = outer ladder = unflipped
0404                 nmodulesLayer_NonFlipped += nmodules_bpix[layer][ladder];
0405                 BPIXLayer_NonFlipped += barycentreLayer[ladder];
0406               } else {  // even ladder = inner ladder = flipped
0407                 nmodulesLayer_Flipped += nmodules_bpix[layer][ladder];
0408                 BPIXLayer_Flipped += barycentreLayer[ladder];
0409               }
0410             } else {  // layer-4
0411
0412               if (ladder % 2 != 0) {  // odd ladder = inner = flipped
0413                 nmodulesLayer_Flipped += nmodules_bpix[layer][ladder];
0414                 BPIXLayer_Flipped += barycentreLayer[ladder];
0415               } else {  //even ladder = outer ladder  = unflipped
0416                 nmodulesLayer_NonFlipped += nmodules_bpix[layer][ladder];
0417                 BPIXLayer_NonFlipped += barycentreLayer[ladder];
0418               }
0419             }
0420
0421           }  // phase-1
0422
0423           // Phase-0
0424           //
0425           // Phase 0: Outer ladders are odd for layer 1,3 and even for layer 2
0426           if (phase_ == 0) {
0427             if (layer == 2) {  // layer-2
0428
0429               if (ladder % 2 != 0) {  // odd ladder = inner = flipped
0430                 nmodulesLayer_Flipped += nmodules_bpix[layer][ladder];
0431                 BPIXLayer_Flipped += barycentreLayer[ladder];
0432               } else {
0433                 nmodulesLayer_NonFlipped += nmodules_bpix[layer][ladder];
0434                 BPIXLayer_NonFlipped += barycentreLayer[ladder];
0435               }
0436             } else {  // layer-1,3
0437
0438               if (ladder % 2 == 0) {  // even ladder = inner = flipped
0439                 nmodulesLayer_Flipped += nmodules_bpix[layer][ladder];
0440                 BPIXLayer_Flipped += barycentreLayer[ladder];
0441               } else {  // odd ladder = outer = non-flipped
0442                 nmodulesLayer_NonFlipped += nmodules_bpix[layer][ladder];
0443                 BPIXLayer_NonFlipped += barycentreLayer[ladder];
0444               }
0445             }
0446
0447           }  // phase-0
0448
0449         }  //loop over ladders
0450
0451         // total BPIX flipped/non-flipped
0452         BPIX_Flipped += BPIXLayer_Flipped;
0453         BPIX_NonFlipped += BPIXLayer_NonFlipped;
0454         nmodules_BPIX_Flipped += nmodulesLayer_Flipped;
0455         nmodules_BPIX_NonFlipped += nmodulesLayer_NonFlipped;
0456
0457         //BPIX per-layer
0458         BPIXLayer *= (1.0 / nmodulesLayer);
0459         BPIXLayer += globalTkPosition;
0460         BPIXLayer_Flipped *= (1.0 / nmodulesLayer_Flipped);
0461         BPIXLayer_Flipped += globalTkPosition;
0462         BPIXLayer_NonFlipped *= (1.0 / nmodulesLayer_NonFlipped);
0463         BPIXLayer_NonFlipped += globalTkPosition;
0464
0465         BPIXLayer_[layer - 1] = GlobalPoint(BPIXLayer.x(), BPIXLayer.y(), BPIXLayer.z());
0466         vBPIXLayer_[layer - 1] = SimplePoint(BPIXLayer_[layer - 1]);
0467         BPIXLayer_Flipped_[layer - 1] =
0468             GlobalPoint(BPIXLayer_Flipped.x(), BPIXLayer_Flipped.y(), BPIXLayer_Flipped.z());
0469         vBPIXLayer_Flipped_[layer - 1] = SimplePoint(BPIXLayer_Flipped_[layer - 1]);
0470         BPIXLayer_NonFlipped_[layer - 1] =
0471             GlobalPoint(BPIXLayer_NonFlipped.x(), BPIXLayer_NonFlipped.y(), BPIXLayer_NonFlipped.z());
0472         vBPIXLayer_NonFlipped_[layer - 1] = SimplePoint(BPIXLayer_NonFlipped_[layer - 1]);
0473         BPIXLayer_DiffFlippedNonFlipped_[layer - 1] = GlobalPoint(BPIXLayer_Flipped.x() - BPIXLayer_NonFlipped.x(),
0474                                                                   BPIXLayer_Flipped.y() - BPIXLayer_NonFlipped.y(),
0475                                                                   BPIXLayer_Flipped.z() - BPIXLayer_NonFlipped.z());
0476         vBPIXLayer_DiffFlippedNonFlipped_[layer - 1] = SimplePoint(BPIXLayer_DiffFlippedNonFlipped_[layer - 1]);
0477
0478       }  // loop over layers
0479
0480       BPIX_Flipped *= (1.0 / nmodules_BPIX_Flipped);
0481       BPIX_Flipped += globalTkPosition;
0482       BPIX_Flipped_ = GlobalPoint(BPIX_Flipped.x(), BPIX_Flipped.y(), BPIX_Flipped.z());
0483       vBPIX_Flipped_ = SimplePoint(BPIX_Flipped_);
0484       BPIX_NonFlipped *= (1.0 / nmodules_BPIX_NonFlipped);
0485       BPIX_NonFlipped += globalTkPosition;
0486       BPIX_NonFlipped_ = GlobalPoint(BPIX_NonFlipped.x(), BPIX_NonFlipped.y(), BPIX_NonFlipped.z());
0487       vBPIX_NonFlipped_ = SimplePoint(BPIX_NonFlipped_);
0488       BPIX_DiffFlippedNonFlipped_ = GlobalPoint(BPIX_Flipped.x() - BPIX_NonFlipped.x(),
0489                                                 BPIX_Flipped.y() - BPIX_NonFlipped.y(),
0490                                                 BPIX_Flipped.z() - BPIX_NonFlipped.z());
0491       vBPIX_DiffFlippedNonFlipped_ = SimplePoint(BPIX_DiffFlippedNonFlipped_);
0492
0493       // FPIX substructures per-(signed)disk/per-ring
0494       int nmodules_FPIX_plus = 0;
0495       int nmodules_FPIX_minus = 0;
0496       GlobalVector FPIX_plus(0.0, 0.0, 0.0);
0497       GlobalVector FPIX_minus(0.0, 0.0, 0.0);
0498       // loop over disks
0499
0500       for (const auto& id : barycentre_fpix) {
0501         int disk = id.first;
0502
0503         int nmodulesDisk = 0;
0504         GlobalVector FPIXDisk(0.0, 0.0, 0.0);
0505
0506         std::map<int, GlobalVector> baryCentreDisk = id.second;
0507         for (const auto& ir : baryCentreDisk) {
0508           int ring = ir.first;
0509           nmodulesDisk += nmodules_fpix[disk][ring];
0510           FPIXDisk += ir.second;
0511           if (disk > 0) {
0512             nmodules_FPIX_plus += nmodules_fpix[disk][ring];
0513             FPIX_plus += ir.second;
0514           }
0515           if (disk < 0) {
0516             nmodules_FPIX_minus += nmodules_fpix[disk][ring];
0517             FPIX_minus += ir.second;
0518           }
0519
0520         }  // loop over rings
0521
0522         FPIXDisk *= (1.0 / nmodulesDisk);
0523         FPIXDisk += globalTkPosition;
0524
0525         if (disk > 0) {
0526           FPIXDisks_plus_[disk - 1] = GlobalPoint(FPIXDisk.x(), FPIXDisk.y(), FPIXDisk.z());
0527           vFPIXDisks_plus_[disk - 1] = SimplePoint(FPIXDisks_plus_[disk - 1]);
0528         }
0529         if (disk < 0) {
0530           FPIXDisks_minus_[-disk - 1] = GlobalPoint(FPIXDisk.x(), FPIXDisk.y(), FPIXDisk.z());
0531           vFPIXDisks_minus_[-disk - 1] = SimplePoint(FPIXDisks_minus_[-disk - 1]);
0532         }
0533       }  // loop over disks
0534
0535       FPIX_plus *= (1.0 / nmodules_FPIX_plus);
0536       FPIX_plus += globalTkPosition;
0537       FPIX_plus_ = GlobalPoint(FPIX_plus.x(), FPIX_plus.y(), FPIX_plus.z());
0538       vFPIX_plus_ = SimplePoint(FPIX_plus_);
0539       FPIX_minus *= (1.0 / nmodules_FPIX_minus);
0540       FPIX_minus += globalTkPosition;
0541       FPIX_minus_ = GlobalPoint(FPIX_minus.x(), FPIX_minus.y(), FPIX_minus.z());
0542       vFPIX_minus_ = SimplePoint(FPIX_minus_);
0543
0544       bcTrees_[label]->Fill();
0545
0546     }  // bcLabels_
0547
0548   }  // check for new IOV for TKAlign
0549
0550   // beamspot
0551   if (prepareBS) {
0552     // loop over bsLabels_
0553     for (const auto& label : bsLabels_) {
0554       // init bstree content
0555       PixelBaryCentreAnalyzer::initBS();
0556
0557       // Get BeamSpot from EventSetup
0558       const BeamSpotObjects* mybeamspot = &iSetup.getData(bsTokens_[label]);
0559
0560       BS_ = GlobalPoint(mybeamspot->x(), mybeamspot->y(), mybeamspot->z());
0561       vBS_ = SimplePoint(BS_);
0562
0563       bsTrees_[label]->Fill();
0564     }  // bsLabels_
0565
0566   }  // check for new IOV for BS
0567 }
0568
0569 // ------------ method called once each job just before starting event loop  ------------
0570 void PixelBaryCentreAnalyzer::beginJob() {
0571   // init bc bs trees
0572   for (const auto& label : bsLabels_) {
0573     std::string treeName = "BeamSpot";
0574     if (!label.empty())
0575       treeName = "BeamSpot_";
0576     treeName += label;
0577
0578     bsTrees_[label] = tFileService->make<TTree>(TString(treeName), "PixelBarycentre analyzer ntuple");
0579
0580     bsTrees_[label]->Branch("run", &run_, "run/I");
0581     bsTrees_[label]->Branch("ls", &ls_, "ls/I");
0582
0583     bsTrees_[label]->Branch("BS", &vBS_, "x/F:y/F:z/F");
0584
0585   }  // bsLabels_
0586
0587   for (const auto& label : bcLabels_) {
0588     std::string treeName = "PixelBarycentre";
0589     if (!label.empty())
0590       treeName = "PixelBarycentre_";
0591     treeName += label;
0592     bcTrees_[label] = tFileService->make<TTree>(TString(treeName), "PixelBarycentre analyzer ntuple");
0593
0594     bcTrees_[label]->Branch("run", &run_, "run/I");
0595     bcTrees_[label]->Branch("ls", &ls_, "ls/I");
0596
0597     bcTrees_[label]->Branch("PIX", &vPIX_, "x/F:y/F:z/F");
0598
0599     bcTrees_[label]->Branch("BPIX", &vBPIX_, "x/F:y/F:z/F");
0600     bcTrees_[label]->Branch("BPIX_Flipped", &vBPIX_Flipped_, "x/F:y/F:z/F");
0601     bcTrees_[label]->Branch("BPIX_NonFlipped", &vBPIX_NonFlipped_, "x/F:y/F:z/F");
0602     bcTrees_[label]->Branch("BPIX_DiffFlippedNonFlipped", &vBPIX_DiffFlippedNonFlipped_, "x/F:y/F:z/F");
0603
0604     bcTrees_[label]->Branch("FPIX", &vFPIX_, "x/F:y/F:z/F");
0605     bcTrees_[label]->Branch("FPIX_plus", &vFPIX_plus_, "x/F:y/F:z/F");
0606     bcTrees_[label]->Branch("FPIX_minus", &vFPIX_minus_, "x/F:y/F:z/F");
0607
0608     //per-layer
0609     for (unsigned int i = 0; i < nPixelLayers; i++) {
0610       TString structure = "BPIXLYR";
0611       int layer = i + 1;
0612       structure += layer;
0613
0614       bcTrees_[label]->Branch(structure, &vBPIXLayer_[i], "x/F:y/F:z/F");
0615       bcTrees_[label]->Branch(structure + "_Flipped", &vBPIXLayer_Flipped_[i], "x/F:y/F:z/F");
0616       bcTrees_[label]->Branch(structure + "_NonFlipped", &vBPIXLayer_NonFlipped_[i], "x/F:y/F:z/F");
0617       bcTrees_[label]->Branch(
0618           structure + "_DiffFlippedNonFlipped", &vBPIXLayer_DiffFlippedNonFlipped_[i], "x/F:y/F:z/F");
0619     }
0620
0621     //per-disk/ring
0622     for (unsigned int i = 0; i < nPixelDiscs; i++) {
0623       TString structure = "FPIXDisk_plus";
0624       int disk = i + 1;
0625       structure += disk;
0626       bcTrees_[label]->Branch(structure, &vFPIXDisks_plus_[i], "x/F:y/F:z/F");
0627
0628       structure = "FPIXDisk_minus";
0629       structure += disk;
0630       bcTrees_[label]->Branch(structure, &vFPIXDisks_minus_[i], "x/F:y/F:z/F");
0631     }
0632
0633   }  // bcLabels_
0634 }
0635
0636 // ------------ method called once each job just after ending the event loop  ------------
0637 void PixelBaryCentreAnalyzer::endJob() {
0638   bcLabels_.clear();
0639   bsLabels_.clear();
0640
0641   bcTrees_.clear();
0642   bsTrees_.clear();
0643 }
0644
0645 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0646 void PixelBaryCentreAnalyzer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0647   edm::ParameterSetDescription desc;
0648   desc.setComment("Validates alignment payloads by providing the position of the pixel barycenter positions");
0649   desc.addUntracked<bool>("usePixelQuality", false);
0650   desc.addUntracked<std::vector<std::string>>("tkAlignLabels", {});
0651   desc.addUntracked<std::vector<std::string>>("beamSpotLabels", {});
0652   descriptions.addWithDefaultLabel(desc);
0653 }
0654
0655 //define this as a plug-in
0656 DEFINE_FWK_MODULE(PixelBaryCentreAnalyzer);