Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-12-14 11:44:40

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