Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-04-06 03:10:46

0001 // -*- C++ -*-
0002 //
0003 // Package:    CalibTracker/SiPixelLorentzAnglePCLWorker
0004 // Class:      SiPixelLorentzAnglePCLWorker
0005 //
0006 /**\class SiPixelLorentzAnglePCLWorker SiPixelLorentzAnglePCLWorker.cc CalibTracker/SiPixelLorentzAnglePCLWorker/src/SiPixelLorentzAnglePCLWorker.cc
0007  Description: generates the intermediate ALCAPROMPT dataset for the measurement of the SiPixel Lorentz Angle in the Prompt Calibration Loop
0008  Implementation:
0009      Books and fills 2D histograms of the drift vs depth in bins of pixel module rings to be fed into the SiPixelLorentzAnglePCLHarvester
0010 */
0011 //
0012 // Original Author:  mmusich
0013 //         Created:  Sat, 29 May 2021 14:46:19 GMT
0014 //
0015 //
0016 
0017 // system includes
0018 #include <string>
0019 #include <fmt/printf.h>
0020 
0021 // user include files
0022 #include "CalibTracker/Records/interface/SiPixelTemplateDBObjectESProducerRcd.h"
0023 #include "CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h"
0024 #include "CondFormats/SiPixelObjects/interface/SiPixelTemplateDBObject.h"
0025 #include "CondFormats/SiPixelTransient/interface/SiPixelTemplate.h"
0026 #include "CondFormats/SiPixelTransient/interface/SiPixelTemplateDefs.h"
0027 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0028 #include "DQMServices/Core/interface/DQMStore.h"
0029 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
0030 #include "DataFormats/GeometryVector/interface/LocalVector.h"
0031 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0032 #include "DataFormats/TrackReco/interface/TrackExtra.h"
0033 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0034 #include "DataFormats/TrackerCommon/interface/PixelBarrelName.h"
0035 #include "DataFormats/TrackerCommon/interface/PixelEndcapName.h"
0036 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0037 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h"
0038 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
0039 #include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h"
0040 #include "FWCore/Framework/interface/Event.h"
0041 #include "FWCore/Framework/interface/ESWatcher.h"
0042 #include "FWCore/Framework/interface/Frameworkfwd.h"
0043 #include "FWCore/Framework/interface/MakerMacros.h"
0044 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0045 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0046 #include "Geometry/CommonTopologies/interface/StripTopology.h"
0047 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0048 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0049 #include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyMap.h"
0050 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0051 #include "Geometry/CommonTopologies/interface/SimplePixelTopology.h"
0052 #include "RecoTracker/TransientTrackingRecHit/interface/TkTransientTrackingRecHitBuilder.h"
0053 #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h"
0054 #include "TrackingTools/Records/interface/TransientRecHitRecord.h"
0055 #include "TrackingTools/TrackFitters/interface/TrajectoryStateCombiner.h"
0056 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0057 
0058 // ROOT includes
0059 #include <TTree.h>
0060 #include <TFile.h>
0061 #include <fstream>
0062 
0063 //
0064 // class declaration
0065 //
0066 
0067 static const int maxpix = 1000;
0068 struct Pixinfo {
0069   int npix;
0070   float row[maxpix];
0071   float col[maxpix];
0072   float adc[maxpix];
0073   float x[maxpix];
0074   float y[maxpix];
0075 };
0076 
0077 struct Hit {
0078   float x;
0079   float y;
0080   double alpha;
0081   double beta;
0082   double gamma;
0083 };
0084 struct Clust {
0085   float x;
0086   float y;
0087   float charge;
0088   int size_x;
0089   int size_y;
0090   int maxPixelCol;
0091   int maxPixelRow;
0092   int minPixelCol;
0093   int minPixelRow;
0094 };
0095 struct Rechit {
0096   float x;
0097   float y;
0098 };
0099 
0100 class SiPixelLorentzAnglePCLWorker : public DQMEDAnalyzer {
0101 public:
0102   explicit SiPixelLorentzAnglePCLWorker(const edm::ParameterSet&);
0103   ~SiPixelLorentzAnglePCLWorker() override = default;
0104 
0105   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0106 
0107 private:
0108   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
0109 
0110   void analyze(edm::Event const&, edm::EventSetup const&) override;
0111 
0112   void dqmBeginRun(edm::Run const&, edm::EventSetup const&) override;
0113 
0114   void dqmEndRun(edm::Run const&, edm::EventSetup const&);
0115 
0116   const Pixinfo fillPix(const SiPixelCluster& LocPix, const PixelTopology* topol) const;
0117   const std::pair<LocalPoint, LocalPoint> surface_deformation(const PixelTopology* topol,
0118                                                               TrajectoryStateOnSurface& tsos,
0119                                                               const SiPixelRecHit* recHitPix) const;
0120   // ------------ member data ------------
0121   SiPixelLorentzAngleCalibrationHistograms iHists;
0122 
0123   // template stuff
0124   edm::ESWatcher<SiPixelTemplateDBObjectESProducerRcd> watchSiPixelTemplateRcd_;
0125   const SiPixelTemplateDBObject* templateDBobject_;
0126   std::vector<SiPixelTemplateStore> thePixelTemp_;
0127 
0128   std::string folder_;
0129   bool notInPCL_;
0130   std::string filename_;
0131   std::vector<std::string> newmodulelist_;
0132 
0133   // tree branches barrel
0134   int run_;
0135   long int event_;
0136   int lumiblock_;
0137   int bx_;
0138   int orbit_;
0139   int module_;
0140   int ladder_;
0141   int layer_;
0142   int isflipped_;
0143   float pt_;
0144   float eta_;
0145   float phi_;
0146   double chi2_;
0147   double ndof_;
0148   Pixinfo pixinfo_;
0149   Hit simhit_, trackhit_;
0150   Clust clust_;
0151   Rechit rechit_;
0152   Rechit rechitCorr_;
0153   float trackhitCorrX_;
0154   float trackhitCorrY_;
0155   float qScale_;
0156   float rQmQt_;
0157 
0158   // tree branches forward
0159   int sideF_;
0160   int diskF_;
0161   int bladeF_;
0162   int panelF_;
0163   int moduleF_;
0164   Pixinfo pixinfoF_;
0165   Hit simhitF_, trackhitF_;
0166   Clust clustF_;
0167   Rechit rechitF_;
0168   Rechit rechitCorrF_;
0169   float trackhitCorrXF_;
0170   float trackhitCorrYF_;
0171   float qScaleF_;
0172   float rQmQtF_;
0173 
0174   // parameters from config file
0175   double ptmin_;
0176   double normChi2Max_;
0177   std::vector<int> clustSizeYMin_;
0178   int clustSizeXMax_;
0179   double residualMax_;
0180   double clustChargeMaxPerLength_;
0181   int hist_depth_;
0182   int hist_drift_;
0183 
0184   std::unique_ptr<TFile> hFile_;
0185   std::unique_ptr<TTree> SiPixelLorentzAngleTreeBarrel_;
0186   std::unique_ptr<TTree> SiPixelLorentzAngleTreeForward_;
0187 
0188   // es consumes
0189   edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geomEsToken_;
0190   edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> topoEsToken_;
0191   edm::ESGetToken<SiPixelTemplateDBObject, SiPixelTemplateDBObjectESProducerRcd> siPixelTemplateEsToken_;
0192   edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> topoPerEventEsToken_;
0193   edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geomPerEventEsToken_;
0194 
0195   // event consumes
0196   edm::EDGetTokenT<TrajTrackAssociationCollection> t_trajTrack;
0197 };
0198 
0199 //
0200 // constructors and destructor
0201 //
0202 SiPixelLorentzAnglePCLWorker::SiPixelLorentzAnglePCLWorker(const edm::ParameterSet& iConfig)
0203     : folder_(iConfig.getParameter<std::string>("folder")),
0204       notInPCL_(iConfig.getParameter<bool>("notInPCL")),
0205       filename_(iConfig.getParameter<std::string>("fileName")),
0206       newmodulelist_(iConfig.getParameter<std::vector<std::string>>("newmodulelist")),
0207       ptmin_(iConfig.getParameter<double>("ptMin")),
0208       normChi2Max_(iConfig.getParameter<double>("normChi2Max")),
0209       clustSizeYMin_(iConfig.getParameter<std::vector<int>>("clustSizeYMin")),
0210       clustSizeXMax_(iConfig.getParameter<int>("clustSizeXMax")),
0211       residualMax_(iConfig.getParameter<double>("residualMax")),
0212       clustChargeMaxPerLength_(iConfig.getParameter<double>("clustChargeMaxPerLength")),
0213       hist_depth_(iConfig.getParameter<int>("binsDepth")),
0214       hist_drift_(iConfig.getParameter<int>("binsDrift")),
0215       geomEsToken_(esConsumes<edm::Transition::BeginRun>()),
0216       topoEsToken_(esConsumes<edm::Transition::BeginRun>()),
0217       siPixelTemplateEsToken_(esConsumes<edm::Transition::BeginRun>()),
0218       topoPerEventEsToken_(esConsumes()),
0219       geomPerEventEsToken_(esConsumes()) {
0220   t_trajTrack = consumes<TrajTrackAssociationCollection>(iConfig.getParameter<edm::InputTag>("src"));
0221 
0222   // now do what ever initialization is needed
0223   int bufsize = 64000;
0224 
0225   //    create tree structure
0226   //    Barrel pixel
0227   if (notInPCL_) {
0228     hFile_ = std::make_unique<TFile>(filename_.c_str(), "RECREATE");
0229     SiPixelLorentzAngleTreeBarrel_ =
0230         std::make_unique<TTree>("SiPixelLorentzAngleTreeBarrel_", "SiPixel LorentzAngle tree barrel", bufsize);
0231     SiPixelLorentzAngleTreeBarrel_->Branch("run", &run_, "run/I", bufsize);
0232     SiPixelLorentzAngleTreeBarrel_->Branch("event", &event_, "event/l", bufsize);
0233     SiPixelLorentzAngleTreeBarrel_->Branch("lumiblock", &lumiblock_, "lumiblock/I", bufsize);
0234     SiPixelLorentzAngleTreeBarrel_->Branch("bx", &bx_, "bx/I", bufsize);
0235     SiPixelLorentzAngleTreeBarrel_->Branch("orbit", &orbit_, "orbit/I", bufsize);
0236     SiPixelLorentzAngleTreeBarrel_->Branch("module", &module_, "module/I", bufsize);
0237     SiPixelLorentzAngleTreeBarrel_->Branch("ladder", &ladder_, "ladder/I", bufsize);
0238     SiPixelLorentzAngleTreeBarrel_->Branch("layer", &layer_, "layer/I", bufsize);
0239     SiPixelLorentzAngleTreeBarrel_->Branch("isflipped", &isflipped_, "isflipped/I", bufsize);
0240     SiPixelLorentzAngleTreeBarrel_->Branch("pt", &pt_, "pt/F", bufsize);
0241     SiPixelLorentzAngleTreeBarrel_->Branch("eta", &eta_, "eta/F", bufsize);
0242     SiPixelLorentzAngleTreeBarrel_->Branch("phi", &phi_, "phi/F", bufsize);
0243     SiPixelLorentzAngleTreeBarrel_->Branch("chi2", &chi2_, "chi2/D", bufsize);
0244     SiPixelLorentzAngleTreeBarrel_->Branch("ndof", &ndof_, "ndof/D", bufsize);
0245     SiPixelLorentzAngleTreeBarrel_->Branch("trackhit", &trackhit_, "x/F:y/F:alpha/D:beta/D:gamma_/D", bufsize);
0246     SiPixelLorentzAngleTreeBarrel_->Branch("npix", &pixinfo_.npix, "npix/I", bufsize);
0247     SiPixelLorentzAngleTreeBarrel_->Branch("rowpix", pixinfo_.row, "row[npix]/F", bufsize);
0248     SiPixelLorentzAngleTreeBarrel_->Branch("colpix", pixinfo_.col, "col[npix]/F", bufsize);
0249     SiPixelLorentzAngleTreeBarrel_->Branch("adc", pixinfo_.adc, "adc[npix]/F", bufsize);
0250     SiPixelLorentzAngleTreeBarrel_->Branch("xpix", pixinfo_.x, "x[npix]/F", bufsize);
0251     SiPixelLorentzAngleTreeBarrel_->Branch("ypix", pixinfo_.y, "y[npix]/F", bufsize);
0252 
0253     SiPixelLorentzAngleTreeBarrel_->Branch(
0254         "clust",
0255         &clust_,
0256         "x/F:y/F:charge/F:size_x/I:size_y/I:maxPixelCol/I:maxPixelRow:minPixelCol/I:minPixelRow/I",
0257         bufsize);
0258     SiPixelLorentzAngleTreeBarrel_->Branch("rechit", &rechit_, "x/F:y/F", bufsize);
0259     SiPixelLorentzAngleTreeBarrel_->Branch("rechit_corr", &rechitCorr_, "x/F:y/F", bufsize);
0260     SiPixelLorentzAngleTreeBarrel_->Branch("trackhitcorr_x", &trackhitCorrX_, "trackhitcorr_x/F", bufsize);
0261     SiPixelLorentzAngleTreeBarrel_->Branch("trackhitcorr_y", &trackhitCorrY_, "trackhitcorr_y/F", bufsize);
0262     SiPixelLorentzAngleTreeBarrel_->Branch("qScale", &qScale_, "qScale/F", bufsize);
0263     SiPixelLorentzAngleTreeBarrel_->Branch("rQmQt", &rQmQt_, "rQmQt/F", bufsize);
0264     //    Forward pixel
0265 
0266     SiPixelLorentzAngleTreeForward_ =
0267         std::make_unique<TTree>("SiPixelLorentzAngleTreeForward_", "SiPixel LorentzAngle tree forward", bufsize);
0268     SiPixelLorentzAngleTreeForward_->Branch("run", &run_, "run/I", bufsize);
0269     SiPixelLorentzAngleTreeForward_->Branch("event", &event_, "event/l", bufsize);
0270     SiPixelLorentzAngleTreeForward_->Branch("lumiblock", &lumiblock_, "lumiblock/I", bufsize);
0271     SiPixelLorentzAngleTreeForward_->Branch("bx", &bx_, "bx/I", bufsize);
0272     SiPixelLorentzAngleTreeForward_->Branch("orbit", &orbit_, "orbit/I", bufsize);
0273     SiPixelLorentzAngleTreeForward_->Branch("side", &sideF_, "side/I", bufsize);
0274     SiPixelLorentzAngleTreeForward_->Branch("disk", &diskF_, "disk/I", bufsize);
0275     SiPixelLorentzAngleTreeForward_->Branch("blade", &bladeF_, "blade/I", bufsize);
0276     SiPixelLorentzAngleTreeForward_->Branch("panel", &panelF_, "panel/I", bufsize);
0277     SiPixelLorentzAngleTreeForward_->Branch("module", &moduleF_, "module/I", bufsize);
0278     SiPixelLorentzAngleTreeForward_->Branch("pt", &pt_, "pt/F", bufsize);
0279     SiPixelLorentzAngleTreeForward_->Branch("eta", &eta_, "eta/F", bufsize);
0280     SiPixelLorentzAngleTreeForward_->Branch("phi", &phi_, "phi/F", bufsize);
0281     SiPixelLorentzAngleTreeForward_->Branch("chi2", &chi2_, "chi2/D", bufsize);
0282     SiPixelLorentzAngleTreeForward_->Branch("ndof", &ndof_, "ndof/D", bufsize);
0283     SiPixelLorentzAngleTreeForward_->Branch("trackhit", &trackhitF_, "x/F:y/F:alpha/D:beta/D:gamma_/D", bufsize);
0284     SiPixelLorentzAngleTreeForward_->Branch("npix", &pixinfoF_.npix, "npix/I", bufsize);
0285     SiPixelLorentzAngleTreeForward_->Branch("rowpix", pixinfoF_.row, "row[npix]/F", bufsize);
0286     SiPixelLorentzAngleTreeForward_->Branch("colpix", pixinfoF_.col, "col[npix]/F", bufsize);
0287     SiPixelLorentzAngleTreeForward_->Branch("adc", pixinfoF_.adc, "adc[npix]/F", bufsize);
0288     SiPixelLorentzAngleTreeForward_->Branch("xpix", pixinfoF_.x, "x[npix]/F", bufsize);
0289     SiPixelLorentzAngleTreeForward_->Branch("ypix", pixinfoF_.y, "y[npix]/F", bufsize);
0290 
0291     SiPixelLorentzAngleTreeForward_->Branch(
0292         "clust",
0293         &clustF_,
0294         "x/F:y/F:charge/F:size_x/I:size_y/I:maxPixelCol/I:maxPixelRow:minPixelCol/I:minPixelRow/I",
0295         bufsize);
0296     SiPixelLorentzAngleTreeForward_->Branch("rechit", &rechitF_, "x/F:y/F", bufsize);
0297     SiPixelLorentzAngleTreeForward_->Branch("rechit_corr", &rechitCorrF_, "x/F:y/F", bufsize);
0298     SiPixelLorentzAngleTreeForward_->Branch("trackhitcorr_x", &trackhitCorrXF_, "trackhitcorr_x/F", bufsize);
0299     SiPixelLorentzAngleTreeForward_->Branch("trackhitcorr_y", &trackhitCorrYF_, "trackhitcorr_y/F", bufsize);
0300     SiPixelLorentzAngleTreeForward_->Branch("qScale", &qScaleF_, "qScale/F", bufsize);
0301     SiPixelLorentzAngleTreeForward_->Branch("rQmQt", &rQmQtF_, "rQmQt/F", bufsize);
0302   }
0303 }
0304 
0305 //
0306 // member functions
0307 //
0308 
0309 // ------------ method called for each event  ------------
0310 
0311 void SiPixelLorentzAnglePCLWorker::analyze(edm::Event const& iEvent, edm::EventSetup const& iSetup) {
0312   static constexpr float cmToum = 10000.;
0313 
0314   // Retrieve tracker topology from geometry
0315   const TrackerTopology* const tTopo = &iSetup.getData(topoPerEventEsToken_);
0316 
0317   // Retrieve track geometry
0318   const TrackerGeometry* tracker = &iSetup.getData(geomPerEventEsToken_);
0319 
0320   // get the association map between tracks and trajectories
0321   edm::Handle<TrajTrackAssociationCollection> trajTrackCollectionHandle;
0322   iEvent.getByToken(t_trajTrack, trajTrackCollectionHandle);
0323 
0324   module_ = -1;
0325   layer_ = -1;
0326   ladder_ = -1;
0327   isflipped_ = -1;
0328   pt_ = -999;
0329   eta_ = 999;
0330   phi_ = 999;
0331   pixinfo_.npix = 0;
0332 
0333   run_ = iEvent.id().run();
0334   event_ = iEvent.id().event();
0335   lumiblock_ = iEvent.luminosityBlock();
0336   bx_ = iEvent.bunchCrossing();
0337   orbit_ = iEvent.orbitNumber();
0338 
0339   if (!trajTrackCollectionHandle->empty()) {
0340     for (TrajTrackAssociationCollection::const_iterator it = trajTrackCollectionHandle->begin();
0341          it != trajTrackCollectionHandle->end();
0342          ++it) {
0343       const reco::Track& track = *it->val;
0344       const Trajectory& traj = *it->key;
0345 
0346       // get the trajectory measurements
0347       std::vector<TrajectoryMeasurement> tmColl = traj.measurements();
0348       pt_ = track.pt();
0349       eta_ = track.eta();
0350       phi_ = track.phi();
0351       chi2_ = traj.chiSquared();
0352       ndof_ = traj.ndof();
0353 
0354       if (pt_ < ptmin_)
0355         continue;
0356 
0357       iHists.h_trackEta_->Fill(eta_);
0358       iHists.h_trackPhi_->Fill(phi_);
0359       iHists.h_trackPt_->Fill(pt_);
0360       iHists.h_trackChi2_->Fill(chi2_ / ndof_);
0361       iHists.h_tracks_->Fill(0);
0362       bool pixeltrack = false;
0363 
0364       // iterate over trajectory measurements
0365       for (const auto& itTraj : tmColl) {
0366         if (!itTraj.updatedState().isValid())
0367           continue;
0368         const TransientTrackingRecHit::ConstRecHitPointer& recHit = itTraj.recHit();
0369         if (!recHit->isValid() || recHit->geographicalId().det() != DetId::Tracker)
0370           continue;
0371         unsigned int subDetID = (recHit->geographicalId().subdetId());
0372         if (subDetID == PixelSubdetector::PixelBarrel || subDetID == PixelSubdetector::PixelEndcap) {
0373           if (!pixeltrack) {
0374             iHists.h_tracks_->Fill(1);
0375           }
0376           pixeltrack = true;
0377         }
0378 
0379         if (subDetID == PixelSubdetector::PixelBarrel) {
0380           DetId detIdObj = recHit->geographicalId();
0381           const PixelGeomDetUnit* theGeomDet = dynamic_cast<const PixelGeomDetUnit*>(tracker->idToDet(detIdObj));
0382           if (!theGeomDet)
0383             continue;
0384 
0385           const PixelTopology* topol = &(theGeomDet->specificTopology());
0386 
0387           float ypitch_ = topol->pitch().second;
0388           float width_ = theGeomDet->surface().bounds().thickness();
0389 
0390           if (!topol)
0391             continue;
0392 
0393           layer_ = tTopo->pxbLayer(detIdObj);
0394           ladder_ = tTopo->pxbLadder(detIdObj);
0395           module_ = tTopo->pxbModule(detIdObj);
0396 
0397           float tmp1 = theGeomDet->surface().toGlobal(Local3DPoint(0., 0., 0.)).perp();
0398           float tmp2 = theGeomDet->surface().toGlobal(Local3DPoint(0., 0., 1.)).perp();
0399 
0400           isflipped_ = (tmp2 < tmp1) ? 1 : 0;
0401 
0402           const SiPixelRecHit* recHitPix = dynamic_cast<const SiPixelRecHit*>((*recHit).hit());
0403           if (!recHitPix)
0404             continue;
0405           rechit_.x = recHitPix->localPosition().x();
0406           rechit_.y = recHitPix->localPosition().y();
0407           SiPixelRecHit::ClusterRef const& cluster = recHitPix->cluster();
0408 
0409           pixinfo_ = fillPix(*cluster, topol);
0410 
0411           // fill entries in clust_
0412 
0413           clust_.x = (cluster)->x();
0414           clust_.y = (cluster)->y();
0415           clust_.charge = (cluster->charge()) / 1000.;  // clust_.charge: in the unit of 1000e
0416           clust_.size_x = cluster->sizeX();
0417           clust_.size_y = cluster->sizeY();
0418           clust_.maxPixelCol = cluster->maxPixelCol();
0419           clust_.maxPixelRow = cluster->maxPixelRow();
0420           clust_.minPixelCol = cluster->minPixelCol();
0421           clust_.minPixelRow = cluster->minPixelRow();
0422 
0423           // fill the trackhit info
0424           TrajectoryStateOnSurface tsos = itTraj.updatedState();
0425           if (!tsos.isValid()) {
0426             edm::LogWarning("SiPixelLorentzAnglePCLWorker") << "tsos not valid";
0427             continue;
0428           }
0429           LocalVector trackdirection = tsos.localDirection();
0430           LocalPoint trackposition = tsos.localPosition();
0431 
0432           if (trackdirection.z() == 0)
0433             continue;
0434           // the local position and direction
0435           trackhit_.alpha = atan2(trackdirection.z(), trackdirection.x());
0436           trackhit_.beta = atan2(trackdirection.z(), trackdirection.y());
0437           trackhit_.gamma = atan2(trackdirection.x(), trackdirection.y());
0438           trackhit_.x = trackposition.x();
0439           trackhit_.y = trackposition.y();
0440 
0441           // get qScale_ = templ.qscale() and  templ.r_qMeas_qTrue();
0442           float cotalpha = trackdirection.x() / trackdirection.z();
0443           float cotbeta = trackdirection.y() / trackdirection.z();
0444           float cotbeta_min = clustSizeYMin_[layer_ - 1] * ypitch_ / width_;
0445           if (fabs(cotbeta) <= cotbeta_min)
0446             continue;
0447           double drdz = sqrt(1. + cotalpha * cotalpha + cotbeta * cotbeta);
0448           double clusterCharge_cut = clustChargeMaxPerLength_ * drdz;
0449 
0450           auto detId = detIdObj.rawId();
0451           int DetId_index = -1;
0452 
0453           const auto& newModIt = (std::find(iHists.BPixnewDetIds_.begin(), iHists.BPixnewDetIds_.end(), detId));
0454           bool isNewMod = (newModIt != iHists.BPixnewDetIds_.end());
0455           if (isNewMod) {
0456             DetId_index = std::distance(iHists.BPixnewDetIds_.begin(), newModIt);
0457           }
0458 
0459           if (notInPCL_) {
0460             // fill the template from the store (from dqmBeginRun)
0461             SiPixelTemplate theTemplate(thePixelTemp_);
0462 
0463             float locBx = (cotbeta < 0.) ? -1 : 1.;
0464             float locBz = (cotalpha < 0.) ? -locBx : locBx;
0465 
0466             int TemplID = templateDBobject_->getTemplateID(detId);
0467             theTemplate.interpolate(TemplID, cotalpha, cotbeta, locBz, locBx);
0468             qScale_ = theTemplate.qscale();
0469             rQmQt_ = theTemplate.r_qMeas_qTrue();
0470           }
0471 
0472           // Surface deformation
0473           const auto& lp_pair = surface_deformation(topol, tsos, recHitPix);
0474 
0475           LocalPoint lp_track = lp_pair.first;
0476           LocalPoint lp_rechit = lp_pair.second;
0477 
0478           rechitCorr_.x = lp_rechit.x();
0479           rechitCorr_.y = lp_rechit.y();
0480           trackhitCorrX_ = lp_track.x();
0481           trackhitCorrY_ = lp_track.y();
0482 
0483           if (notInPCL_) {
0484             SiPixelLorentzAngleTreeBarrel_->Fill();
0485           }
0486 
0487           // is one pixel in cluster a large pixel ? (hit will be excluded)
0488           bool large_pix = false;
0489           for (int j = 0; j < pixinfo_.npix; j++) {
0490             int colpos = static_cast<int>(pixinfo_.col[j]);
0491             if (pixinfo_.row[j] == 0 || pixinfo_.row[j] == 79 || pixinfo_.row[j] == 80 || pixinfo_.row[j] == 159 ||
0492                 colpos % 52 == 0 || colpos % 52 == 51) {
0493               large_pix = true;
0494             }
0495           }
0496 
0497           double residualsq = (trackhitCorrX_ - rechitCorr_.x) * (trackhitCorrX_ - rechitCorr_.x) +
0498                               (trackhitCorrY_ - rechitCorr_.y) * (trackhitCorrY_ - rechitCorr_.y);
0499 
0500           double xlim1 = trackhitCorrX_ - width_ * cotalpha / 2.;
0501           double hypitch_ = ypitch_ / 2.;
0502           double ylim1 = trackhitCorrY_ - width_ * cotbeta / 2.;
0503           double ylim2 = trackhitCorrY_ + width_ * cotbeta / 2.;
0504 
0505           int clustSizeY_cut = clustSizeYMin_[layer_ - 1];
0506 
0507           if (!large_pix && (chi2_ / ndof_) < normChi2Max_ && cluster->sizeY() >= clustSizeY_cut &&
0508               residualsq < residualMax_ * residualMax_ && cluster->charge() < clusterCharge_cut &&
0509               cluster->sizeX() < clustSizeXMax_) {
0510             // iterate over pixels in hit
0511             for (int j = 0; j < pixinfo_.npix; j++) {
0512               // use trackhits and include bowing correction
0513               float ypixlow = pixinfo_.y[j] - hypitch_;
0514               float ypixhigh = pixinfo_.y[j] + hypitch_;
0515               if (cotbeta > 0.) {
0516                 if (ylim1 > ypixlow)
0517                   ypixlow = ylim1;
0518                 if (ylim2 < ypixhigh)
0519                   ypixhigh = ylim2;
0520               } else {
0521                 if (ylim2 > ypixlow)
0522                   ypixlow = ylim2;
0523                 if (ylim1 < ypixhigh)
0524                   ypixhigh = ylim1;
0525               }
0526               float ypixavg = 0.5f * (ypixlow + ypixhigh);
0527 
0528               float dx = (pixinfo_.x[j] - xlim1) * cmToum;  // dx: in the unit of micrometer
0529               float dy = (ypixavg - ylim1) * cmToum;        // dy: in the unit of micrometer
0530               float depth = dy * tan(trackhit_.beta);
0531               float drift = dx - dy * tan(trackhit_.gamma);
0532 
0533               if (isNewMod == false) {
0534                 int i_index = module_ + (layer_ - 1) * iHists.nModules_[layer_ - 1];
0535                 iHists.h_drift_depth_adc_[i_index]->Fill(drift, depth, pixinfo_.adc[j]);
0536                 iHists.h_drift_depth_adc2_[i_index]->Fill(drift, depth, pixinfo_.adc[j] * pixinfo_.adc[j]);
0537                 iHists.h_drift_depth_noadc_[i_index]->Fill(drift, depth, 1.);
0538                 iHists.h_bySectOccupancy_->Fill(i_index - 1);  // histogram starts at 0
0539 
0540                 if (tracker->getDetectorType(subDetID) == TrackerGeometry::ModuleType::Ph1PXB) {
0541                   if ((module_ == 3 || module_ == 5) && (layer_ == 3 || layer_ == 4)) {
0542                     int i_index_merge = i_index + 1;
0543                     iHists.h_drift_depth_adc_[i_index_merge]->Fill(drift, depth, pixinfo_.adc[j]);
0544                     iHists.h_drift_depth_adc2_[i_index_merge]->Fill(drift, depth, pixinfo_.adc[j] * pixinfo_.adc[j]);
0545                     iHists.h_drift_depth_noadc_[i_index_merge]->Fill(drift, depth, 1.);
0546                     iHists.h_bySectOccupancy_->Fill(i_index_merge - 1);
0547                   }
0548                   if ((module_ == 4 || module_ == 6) && (layer_ == 3 || layer_ == 4)) {
0549                     int i_index_merge = i_index - 1;
0550                     iHists.h_drift_depth_adc_[i_index_merge]->Fill(drift, depth, pixinfo_.adc[j]);
0551                     iHists.h_drift_depth_adc2_[i_index_merge]->Fill(drift, depth, pixinfo_.adc[j] * pixinfo_.adc[j]);
0552                     iHists.h_drift_depth_noadc_[i_index_merge]->Fill(drift, depth, 1.);
0553                     iHists.h_bySectOccupancy_->Fill(i_index_merge - 1);
0554                   }
0555                 }
0556 
0557               } else {
0558                 int new_index = iHists.nModules_[iHists.nlay - 1] +
0559                                 (iHists.nlay - 1) * iHists.nModules_[iHists.nlay - 1] + 1 + DetId_index;
0560 
0561                 iHists.h_drift_depth_adc_[new_index]->Fill(drift, depth, pixinfo_.adc[j]);
0562                 iHists.h_drift_depth_adc2_[new_index]->Fill(drift, depth, pixinfo_.adc[j] * pixinfo_.adc[j]);
0563                 iHists.h_drift_depth_noadc_[new_index]->Fill(drift, depth, 1.);
0564                 iHists.h_bySectOccupancy_->Fill(new_index - 1);  // histogram starts at 0
0565               }
0566             }
0567           }
0568         } else if (subDetID == PixelSubdetector::PixelEndcap) {
0569           DetId detIdObj = recHit->geographicalId();
0570           const PixelGeomDetUnit* theGeomDet = dynamic_cast<const PixelGeomDetUnit*>(tracker->idToDet(detIdObj));
0571           if (!theGeomDet)
0572             continue;
0573 
0574           const PixelTopology* topol = &(theGeomDet->specificTopology());
0575 
0576           if (!topol)
0577             continue;
0578 
0579           sideF_ = tTopo->pxfSide(detIdObj);
0580           diskF_ = tTopo->pxfDisk(detIdObj);
0581           bladeF_ = tTopo->pxfBlade(detIdObj);
0582           panelF_ = tTopo->pxfPanel(detIdObj);
0583           moduleF_ = tTopo->pxfModule(detIdObj);
0584 
0585           const SiPixelRecHit* recHitPix = dynamic_cast<const SiPixelRecHit*>((*recHit).hit());
0586           if (!recHitPix)
0587             continue;
0588           rechitF_.x = recHitPix->localPosition().x();
0589           rechitF_.y = recHitPix->localPosition().y();
0590           SiPixelRecHit::ClusterRef const& cluster = recHitPix->cluster();
0591 
0592           pixinfoF_ = fillPix(*cluster, topol);
0593 
0594           // fill entries in clust_
0595 
0596           clustF_.x = (cluster)->x();
0597           clustF_.y = (cluster)->y();
0598           clustF_.charge = (cluster->charge()) / 1000.;  // clustF_.charge: in the unit of 1000e
0599           clustF_.size_x = cluster->sizeX();
0600           clustF_.size_y = cluster->sizeY();
0601           clustF_.maxPixelCol = cluster->maxPixelCol();
0602           clustF_.maxPixelRow = cluster->maxPixelRow();
0603           clustF_.minPixelCol = cluster->minPixelCol();
0604           clustF_.minPixelRow = cluster->minPixelRow();
0605 
0606           // fill the trackhit info
0607           TrajectoryStateOnSurface tsos = itTraj.updatedState();
0608           if (!tsos.isValid()) {
0609             edm::LogWarning("SiPixelLorentzAnglePCLWorker") << "tsos not valid";
0610             continue;
0611           }
0612           LocalVector trackdirection = tsos.localDirection();
0613           LocalPoint trackposition = tsos.localPosition();
0614 
0615           if (trackdirection.z() == 0)
0616             continue;
0617           // the local position and direction
0618           trackhitF_.alpha = atan2(trackdirection.z(), trackdirection.x());
0619           trackhitF_.beta = atan2(trackdirection.z(), trackdirection.y());
0620           trackhitF_.gamma = atan2(trackdirection.x(), trackdirection.y());
0621           trackhitF_.x = trackposition.x();
0622           trackhitF_.y = trackposition.y();
0623 
0624           float cotalpha = trackdirection.x() / trackdirection.z();
0625           float cotbeta = trackdirection.y() / trackdirection.z();
0626 
0627           auto detId = detIdObj.rawId();
0628 
0629           if (notInPCL_) {
0630             // fill the template from the store (from dqmBeginRun)
0631             SiPixelTemplate theTemplate(thePixelTemp_);
0632 
0633             float locBx = (cotbeta < 0.) ? -1 : 1.;
0634             float locBz = (cotalpha < 0.) ? -locBx : locBx;
0635 
0636             int TemplID = templateDBobject_->getTemplateID(detId);
0637             theTemplate.interpolate(TemplID, cotalpha, cotbeta, locBz, locBx);
0638             qScaleF_ = theTemplate.qscale();
0639             rQmQtF_ = theTemplate.r_qMeas_qTrue();
0640           }
0641 
0642           // Surface deformation
0643           const auto& lp_pair = surface_deformation(topol, tsos, recHitPix);
0644 
0645           LocalPoint lp_track = lp_pair.first;
0646           LocalPoint lp_rechit = lp_pair.second;
0647 
0648           rechitCorrF_.x = lp_rechit.x();
0649           rechitCorrF_.y = lp_rechit.y();
0650           trackhitCorrXF_ = lp_track.x();
0651           trackhitCorrYF_ = lp_track.y();
0652           if (notInPCL_) {
0653             SiPixelLorentzAngleTreeForward_->Fill();
0654           }
0655         }
0656       }  //end iteration over trajectory measurements
0657     }    //end iteration over trajectories
0658   }
0659 }
0660 
0661 void SiPixelLorentzAnglePCLWorker::dqmBeginRun(edm::Run const& run, edm::EventSetup const& iSetup) {
0662   // geometry
0663   const TrackerGeometry* geom = &iSetup.getData(geomEsToken_);
0664   const TrackerTopology* tTopo = &iSetup.getData(topoEsToken_);
0665 
0666   if (notInPCL_) {
0667     // Initialize 1D templates
0668     if (watchSiPixelTemplateRcd_.check(iSetup)) {
0669       templateDBobject_ = &iSetup.getData(siPixelTemplateEsToken_);
0670       if (!SiPixelTemplate::pushfile(*templateDBobject_, thePixelTemp_)) {
0671         edm::LogError("SiPixelLorentzAnglePCLWorker")
0672             << "Templates not filled correctly. Check the sqlite file. Using SiPixelTemplateDBObject version "
0673             << (*templateDBobject_).version();
0674       }
0675     }
0676   }
0677 
0678   PixelTopologyMap map = PixelTopologyMap(geom, tTopo);
0679   iHists.nlay = geom->numberOfLayers(PixelSubdetector::PixelBarrel);
0680   iHists.nModules_.resize(iHists.nlay);
0681   for (int i = 0; i < iHists.nlay; i++) {
0682     iHists.nModules_[i] = map.getPXBModules(i + 1);
0683   }
0684 
0685   // list of modules already filled, then return (we already entered here)
0686   if (!iHists.BPixnewDetIds_.empty() || !iHists.FPixnewDetIds_.empty())
0687     return;
0688 
0689   if (!newmodulelist_.empty()) {
0690     for (auto const& modulename : newmodulelist_) {
0691       if (modulename.find("BPix_") != std::string::npos) {
0692         PixelBarrelName bn(modulename, true);
0693         const auto& detId = bn.getDetId(tTopo);
0694         iHists.BPixnewmodulename_.push_back(modulename);
0695         iHists.BPixnewDetIds_.push_back(detId.rawId());
0696         iHists.BPixnewModule_.push_back(bn.moduleName());
0697         iHists.BPixnewLayer_.push_back(bn.layerName());
0698       } else if (modulename.find("FPix_") != std::string::npos) {
0699         PixelEndcapName en(modulename, true);
0700         const auto& detId = en.getDetId(tTopo);
0701         iHists.FPixnewmodulename_.push_back(modulename);
0702         iHists.FPixnewDetIds_.push_back(detId.rawId());
0703         iHists.FPixnewDisk_.push_back(en.diskName());
0704         iHists.FPixnewBlade_.push_back(en.bladeName());
0705       }
0706     }
0707   }
0708 }
0709 
0710 void SiPixelLorentzAnglePCLWorker::bookHistograms(DQMStore::IBooker& iBooker,
0711                                                   edm::Run const& run,
0712                                                   edm::EventSetup const& iSetup) {
0713   // book the by partition monitoring
0714   const auto maxSect = iHists.nlay * iHists.nModules_[iHists.nlay - 1] + (int)iHists.BPixnewDetIds_.size();
0715 
0716   iBooker.setCurrentFolder(fmt::sprintf("%s/SectorMonitoring", folder_.data()));
0717   iHists.h_bySectOccupancy_ = iBooker.book1D(
0718       "h_bySectorOccupancy", "hit occupancy by sector;pixel sector;hits on track", maxSect, -0.5, maxSect - 0.5);
0719 
0720   iBooker.setCurrentFolder(folder_);
0721   static constexpr double min_depth_ = -100.;
0722   static constexpr double max_depth_ = 400.;
0723   static constexpr double min_drift_ = -500.;
0724   static constexpr double max_drift_ = 500.;
0725 
0726   // book the mean values projections and set the bin names of the by sector monitoring
0727   char name[128];
0728   char title[256];
0729   for (int i_layer = 1; i_layer <= iHists.nlay; i_layer++) {
0730     for (int i_module = 1; i_module <= iHists.nModules_[i_layer - 1]; i_module++) {
0731       unsigned int i_index = i_module + (i_layer - 1) * iHists.nModules_[i_layer - 1];
0732       std::string binName = fmt::sprintf("BPix Lay%i Mod%i", i_layer, i_module);
0733       LogDebug("SiPixelLorentzAnglePCLWorker") << " i_index: " << i_index << " bin name: " << binName
0734                                                << " (i_layer: " << i_layer << " i_module:" << i_module << ")";
0735 
0736       iHists.h_bySectOccupancy_->setBinLabel(i_index, binName);
0737 
0738       sprintf(name, "h_mean_layer%i_module%i", i_layer, i_module);
0739       sprintf(title,
0740               "average drift vs depth layer%i module%i; production depth [#mum]; #LTdrift#GT [#mum]",
0741               i_layer,
0742               i_module);
0743       iHists.h_mean_[i_index] = iBooker.book1D(name, title, hist_depth_, min_depth_, max_depth_);
0744     }
0745   }
0746   for (int i = 0; i < (int)iHists.BPixnewDetIds_.size(); i++) {
0747     sprintf(name, "h_BPixnew_mean_%s", iHists.BPixnewmodulename_[i].c_str());
0748     sprintf(title,
0749             "average drift vs depth %s; production depth [#mum]; #LTdrift#GT [#mum]",
0750             iHists.BPixnewmodulename_[i].c_str());
0751     int new_index = iHists.nModules_[iHists.nlay - 1] + (iHists.nlay - 1) * iHists.nModules_[iHists.nlay - 1] + 1 + i;
0752     iHists.h_mean_[new_index] = iBooker.book1D(name, title, hist_depth_, min_depth_, max_depth_);
0753 
0754     LogDebug("SiPixelLorentzAnglePCLWorker") << "i_index" << new_index << " bin name: " << iHists.BPixnewmodulename_[i];
0755 
0756     iHists.h_bySectOccupancy_->setBinLabel(new_index, iHists.BPixnewmodulename_[i]);
0757   }
0758 
0759   //book the 2D histograms
0760   for (int i_layer = 1; i_layer <= iHists.nlay; i_layer++) {
0761     iBooker.setCurrentFolder(fmt::sprintf("%s/BPix/BPixLayer%i", folder_.data(), i_layer));
0762     for (int i_module = 1; i_module <= iHists.nModules_[i_layer - 1]; i_module++) {
0763       unsigned int i_index = i_module + (i_layer - 1) * iHists.nModules_[i_layer - 1];
0764 
0765       sprintf(name, "h_drift_depth_adc_layer%i_module%i", i_layer, i_module);
0766       sprintf(title, "depth vs drift (ADC) layer%i module%i; drift [#mum]; production depth [#mum]", i_layer, i_module);
0767       iHists.h_drift_depth_adc_[i_index] =
0768           iBooker.book2D(name, title, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
0769 
0770       sprintf(name, "h_drift_depth_adc2_layer%i_module%i", i_layer, i_module);
0771       sprintf(
0772           title, "depth vs drift (ADC^{2}) layer%i module%i; drift [#mum]; production depth [#mum]", i_layer, i_module);
0773       iHists.h_drift_depth_adc2_[i_index] =
0774           iBooker.book2D(name, title, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
0775 
0776       sprintf(name, "h_drift_depth_noadc_layer%i_module%i", i_layer, i_module);
0777       sprintf(
0778           title, "depth vs drift (no ADC) layer%i module%i; drift [#mum]; production depth [#mum]", i_layer, i_module);
0779       iHists.h_drift_depth_noadc_[i_index] =
0780           iBooker.book2D(name, title, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
0781 
0782       sprintf(name, "h_drift_depth_layer%i_module%i", i_layer, i_module);
0783       sprintf(title, "depth vs drift layer%i module%i; drift [#mum]; production depth [#mum]", i_layer, i_module);
0784       iHists.h_drift_depth_[i_index] =
0785           iBooker.book2D(name, title, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
0786     }
0787   }
0788 
0789   // book the "new" modules
0790   iBooker.setCurrentFolder(fmt::sprintf("%s/BPix/NewModules", folder_.data()));
0791   for (int i = 0; i < (int)iHists.BPixnewDetIds_.size(); i++) {
0792     int new_index = iHists.nModules_[iHists.nlay - 1] + (iHists.nlay - 1) * iHists.nModules_[iHists.nlay - 1] + 1 + i;
0793 
0794     sprintf(name, "h_BPixnew_drift_depth_adc_%s", iHists.BPixnewmodulename_[i].c_str());
0795     sprintf(
0796         title, "depth vs drift (ADC) %s; drift [#mum]; production depth [#mum]", iHists.BPixnewmodulename_[i].c_str());
0797     iHists.h_drift_depth_adc_[new_index] =
0798         iBooker.book2D(name, title, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
0799 
0800     sprintf(name, "h_BPixnew_drift_depth_adc2_%s", iHists.BPixnewmodulename_[i].c_str());
0801     sprintf(title,
0802             "depth vs drift (ADC^{2}) %s; drift [#mum]; production depth [#mum]",
0803             iHists.BPixnewmodulename_[i].c_str());
0804     iHists.h_drift_depth_adc2_[new_index] =
0805         iBooker.book2D(name, title, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
0806 
0807     sprintf(name, "h_BPixnew_drift_depth_noadc_%s", iHists.BPixnewmodulename_[i].c_str());
0808     sprintf(title,
0809             "depth vs drift (no ADC)%s; drift [#mum]; production depth [#mum]",
0810             iHists.BPixnewmodulename_[i].c_str());
0811     iHists.h_drift_depth_noadc_[new_index] =
0812         iBooker.book2D(name, title, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
0813 
0814     sprintf(name, "h_BPixnew_drift_depth_%s", iHists.BPixnewmodulename_[i].c_str());
0815     sprintf(title, "depth vs drift %s; drift [#mum]; production depth [#mum]", iHists.BPixnewmodulename_[i].c_str());
0816     iHists.h_drift_depth_[new_index] =
0817         iBooker.book2D(name, title, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
0818   }
0819 
0820   // book the track monitoring plots
0821   iBooker.setCurrentFolder(fmt::sprintf("%s/TrackMonitoring", folder_.data()));
0822   iHists.h_tracks_ = iBooker.book1D("h_tracks", ";tracker volume;tracks", 2, -0.5, 1.5);
0823   iHists.h_tracks_->setBinLabel(1, "all tracks", 1);
0824   iHists.h_tracks_->setBinLabel(2, "has pixel hits", 1);
0825   iHists.h_trackEta_ = iBooker.book1D("h_trackEta", ";track #eta; #tracks", 30, -3., 3.);
0826   iHists.h_trackPhi_ = iBooker.book1D("h_trackPhi", ";track #phi; #tracks", 48, -M_PI, M_PI);
0827   iHists.h_trackPt_ = iBooker.book1D("h_trackPt", ";track p_{T} [GeV]; #tracks", 100, 0., 100.);
0828   iHists.h_trackChi2_ = iBooker.book1D("h_trackChi2ndof", ";track #chi^{2}/ndof; #tracks", 100, 0., 10.);
0829 }
0830 
0831 void SiPixelLorentzAnglePCLWorker::dqmEndRun(edm::Run const& run, edm::EventSetup const& iSetup) {
0832   if (notInPCL_) {
0833     hFile_->cd();
0834     hFile_->Write();
0835     hFile_->Close();
0836   }
0837 }
0838 
0839 // method used to fill per pixel info
0840 const Pixinfo SiPixelLorentzAnglePCLWorker::fillPix(const SiPixelCluster& LocPix, const PixelTopology* topol) const {
0841   Pixinfo pixinfo;
0842   const std::vector<SiPixelCluster::Pixel>& pixvector = LocPix.pixels();
0843   pixinfo.npix = 0;
0844   for (std::vector<SiPixelCluster::Pixel>::const_iterator itPix = pixvector.begin(); itPix != pixvector.end();
0845        itPix++) {
0846     pixinfo.row[pixinfo.npix] = itPix->x;
0847     pixinfo.col[pixinfo.npix] = itPix->y;
0848     pixinfo.adc[pixinfo.npix] = itPix->adc;
0849     LocalPoint lp = topol->localPosition(MeasurementPoint(itPix->x + 0.5, itPix->y + 0.5));
0850     pixinfo.x[pixinfo.npix] = lp.x();
0851     pixinfo.y[pixinfo.npix] = lp.y();
0852     pixinfo.npix++;
0853   }
0854   return pixinfo;
0855 }
0856 
0857 // method used to correct for the surface deformation
0858 const std::pair<LocalPoint, LocalPoint> SiPixelLorentzAnglePCLWorker::surface_deformation(
0859     const PixelTopology* topol, TrajectoryStateOnSurface& tsos, const SiPixelRecHit* recHitPix) const {
0860   LocalPoint trackposition = tsos.localPosition();
0861   const LocalTrajectoryParameters& ltp = tsos.localParameters();
0862   const Topology::LocalTrackAngles localTrackAngles(ltp.dxdz(), ltp.dydz());
0863 
0864   std::pair<float, float> pixels_track = topol->pixel(trackposition, localTrackAngles);
0865   std::pair<float, float> pixels_rechit = topol->pixel(recHitPix->localPosition(), localTrackAngles);
0866 
0867   LocalPoint lp_track = topol->localPosition(MeasurementPoint(pixels_track.first, pixels_track.second));
0868 
0869   LocalPoint lp_rechit = topol->localPosition(MeasurementPoint(pixels_rechit.first, pixels_rechit.second));
0870 
0871   std::pair<LocalPoint, LocalPoint> lps = std::make_pair(lp_track, lp_rechit);
0872   return lps;
0873 }
0874 
0875 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0876 void SiPixelLorentzAnglePCLWorker::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0877   edm::ParameterSetDescription desc;
0878   desc.setComment("Worker module of the SiPixel Lorentz Angle PCL monitoring workflow");
0879   desc.add<std::string>("folder", "AlCaReco/SiPixelLorentzAngle")->setComment("directory of PCL Worker output");
0880   desc.add<bool>("notInPCL", false)->setComment("create TTree (true) or not (false)");
0881   desc.add<std::string>("fileName", "testrun.root")->setComment("name of the TTree file if notInPCL = true");
0882   desc.add<std::vector<std::string>>("newmodulelist", {})->setComment("the list of DetIds for new sensors");
0883   desc.add<edm::InputTag>("src", edm::InputTag("TrackRefitter"))->setComment("input track collections");
0884   desc.add<double>("ptMin", 3.)->setComment("minimum pt on tracks");
0885   desc.add<double>("normChi2Max", 2.)->setComment("maximum reduced chi squared");
0886   desc.add<std::vector<int>>("clustSizeYMin", {4, 3, 3, 2})
0887       ->setComment("minimum cluster size on Y axis for all Barrel Layers");
0888   desc.add<int>("clustSizeXMax", 5)->setComment("maximum cluster size on X axis");
0889   desc.add<double>("residualMax", 0.005)->setComment("maximum residual");
0890   desc.add<double>("clustChargeMaxPerLength", 50000)
0891       ->setComment("maximum cluster charge per unit length of pixel depth (z)");
0892   desc.add<int>("binsDepth", 50)->setComment("# bins for electron production depth axis");
0893   desc.add<int>("binsDrift", 100)->setComment("# bins for electron drift axis");
0894   descriptions.addWithDefaultLabel(desc);
0895 }
0896 
0897 // define this as a plug-in
0898 DEFINE_FWK_MODULE(SiPixelLorentzAnglePCLWorker);