From da6a318d9d54ea31e4755a1f620c51017f1e2779 Mon Sep 17 00:00:00 2001 From: mmusich Date: Sat, 29 May 2021 18:31:11 +0200 Subject: [PATCH 01/12] add stub of SiPixelLorentzAnglePCLWorker --- .../SiPixelLorentzAngle/BuildFile.xml | 3 + .../SiPixelLorentzAngleCalibrationStruct.h | 26 + .../SiPixelLorentzAnglePCLWorker_cfi.py | 15 + .../src/SiPixelLorentzAnglePCLWorker.cc | 476 ++++++++++++++++++ ...CalibProdSiPixelLorentzAngle_Output_cff.py | 15 + ...OPromptCalibProdSiPixelLorentzAngle_cff.py | 70 +++ Configuration/AlCa/python/autoAlca.py | 2 +- .../EventContent/python/AlCaRecoOutput_cff.py | 2 +- .../python/relval_steps.py | 4 +- .../python/AlCaRecoStreams_cff.py | 12 +- 10 files changed, 620 insertions(+), 5 deletions(-) create mode 100644 CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h create mode 100644 CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py create mode 100644 CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc create mode 100644 Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_Output_cff.py create mode 100644 Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py diff --git a/CalibTracker/SiPixelLorentzAngle/BuildFile.xml b/CalibTracker/SiPixelLorentzAngle/BuildFile.xml index 28ad08af08f7b..68ed8d0fb831d 100644 --- a/CalibTracker/SiPixelLorentzAngle/BuildFile.xml +++ b/CalibTracker/SiPixelLorentzAngle/BuildFile.xml @@ -8,5 +8,8 @@ + + + diff --git a/CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h b/CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h new file mode 100644 index 0000000000000..c68bf7885c136 --- /dev/null +++ b/CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h @@ -0,0 +1,26 @@ +#ifndef CalibTracker_SiPixelLorentzAngle_SiPixelLorentzAngleCalibrationStruct_h +#define CalibTracker_SiPixelLorentzAngle_SiPixelLorentzAngleCalibrationStruct_h + +#include "DQMServices/Core/interface/DQMStore.h" +#include + +struct SiPixelLorentzAngleCalibrationHistograms { +public: + SiPixelLorentzAngleCalibrationHistograms() = default; + + using MonitorMap = std::unordered_map; + + // hardcode 4 BPix layers + int nlay; + int nModules_[4]; + + MonitorMap h_drift_depth_adc_; + MonitorMap h_drift_depth_adc2_; + MonitorMap h_drift_depth_noadc_; + MonitorMap h_drift_depth_; + MonitorMap h_mean_; + + dqm::reco::MonitorElement* h_tracks_; +}; + +#endif diff --git a/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py new file mode 100644 index 0000000000000..04e139b5aba43 --- /dev/null +++ b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py @@ -0,0 +1,15 @@ +import FWCore.ParameterSet.Config as cms + +from DQMServices.Core.DQMEDAnalyzer import DQMEDAnalyzer +SiPixelLorentzAnglePCLWorker = DQMEDAnalyzer( + "SiPixelLorentzAnglePCLWorker", + folder = cms.string('AlCaReco/SiPixelLorentzAngle'), + src = cms.InputTag("TrackRefitter"), + binsDepth = cms.int32(50), + binsDrift = cms.int32(200), + ptMin = cms.double(3), + normChi2Max = cms.double(2), + clustSizeYMin = cms.int32(4), + residualMax = cms.double(0.005), + clustChargeMax = cms.double(120000) +) diff --git a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc new file mode 100644 index 0000000000000..55daa1a914d60 --- /dev/null +++ b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc @@ -0,0 +1,476 @@ +// -*- C++ -*- +// +// Package: CalibTracker/SiPixelLorentzAnglePCLWorker +// Class: SiPixelLorentzAnglePCLWorker +// +/**\class SiPixelLorentzAnglePCLWorker SiPixelLorentzAnglePCLWorker.cc CalibTracker/SiPixelLorentzAnglePCLWorker/plugins/SiPixelLorentzAnglePCLWorker.cc + + Description: [one line class summary] + + Implementation: + [Notes on implementation] +*/ +// +// Original Author: mmusich +// Created: Sat, 29 May 2021 14:46:19 GMT +// +// + +#include + +// user include files +#include "DQMServices/Core/interface/DQMStore.h" +#include "FWCore/Framework/interface/Frameworkfwd.h" +#include "DQMServices/Core/interface/DQMGlobalEDAnalyzer.h" + +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +#include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h" +#include "FWCore/Framework/interface/ESHandle.h" +#include "DataFormats/GeometryVector/interface/GlobalVector.h" +#include "DataFormats/GeometryVector/interface/LocalVector.h" +#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" +#include "Geometry/CommonDetUnit/interface/GeomDetType.h" +#include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h" +#include "DataFormats/TrackReco/interface/TrackExtra.h" +#include "DataFormats/SiStripDetId/interface/StripSubdetector.h" +#include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h" +#include "Geometry/CommonTopologies/interface/StripTopology.h" +#include "DataFormats/TrackerCommon/interface/TrackerTopology.h" +#include "Geometry/Records/interface/IdealGeometryRecord.h" +#include "DataFormats/TrackReco/interface/TrackFwd.h" + +#include "TrackingTools/Records/interface/TransientRecHitRecord.h" +#include "TrackingTools/TrackFitters/interface/TrajectoryStateCombiner.h" +#include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h" +#include "RecoTracker/TransientTrackingRecHit/interface/TkTransientTrackingRecHitBuilder.h" +#include "TrackingTools/TransientTrack/interface/TransientTrack.h" +#include "CalibTracker/Records/interface/SiPixelTemplateDBObjectESProducerRcd.h" +#include "CondFormats/SiPixelObjects/interface/SiPixelTemplateDBObject.h" +#include "CondFormats/SiPixelTransient/interface/SiPixelTemplateDefs.h" +#include "CondFormats/SiPixelTransient/interface/SiPixelTemplate.h" +#include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyMap.h" +#include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" + +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "../interface/SiPixelLorentzAngleCalibrationStruct.h" + +// +// class declaration +// + +static const int maxpix = 1000; +struct Pixinfo { + int npix; + float row[maxpix]; + float col[maxpix]; + float adc[maxpix]; + float x[maxpix]; + float y[maxpix]; +}; + +class SiPixelLorentzAnglePCLWorker : public DQMGlobalEDAnalyzer { +public: + explicit SiPixelLorentzAnglePCLWorker(const edm::ParameterSet&); + ~SiPixelLorentzAnglePCLWorker() override; + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +private: + void bookHistograms(DQMStore::IBooker&, + edm::Run const&, + edm::EventSetup const&, + SiPixelLorentzAngleCalibrationHistograms&) const override; + + void dqmAnalyze(edm::Event const&, + edm::EventSetup const&, + SiPixelLorentzAngleCalibrationHistograms const&) const override; + + void dqmBeginRun(edm::Run const&, edm::EventSetup const&, SiPixelLorentzAngleCalibrationHistograms&) const override; + + const Pixinfo fillPix(const SiPixelCluster& LocPix, const PixelTopology* topol) const; + const std::pair surface_deformation(const PixelTopology* topol, + TrajectoryStateOnSurface& tsos, + const SiPixelRecHit* recHitPix) const; + // ------------ member data ------------ + std::string folder_; + + // histogram etc + int hist_x_; + int hist_y_; + double min_x_; + double max_x_; + double min_y_; + double max_y_; + double width_; + double min_depth_; + double max_depth_; + double min_drift_; + double max_drift_; + + // parameters from config file + double ptmin_; + double normChi2Max_; + int clustSizeYMin_; + double residualMax_; + double clustChargeMax_; + int hist_depth_; + int hist_drift_; + + // es consumes + edm::ESGetToken geomEsToken_; + edm::ESGetToken topoEsToken_; + edm::ESGetToken siPixelTemplateEsToken_; + edm::ESGetToken topoPerEventEsToken_; + edm::ESGetToken geomPerEventEsToken_; + + // event consumes + edm::EDGetTokenT t_trajTrack; +}; + +// +// constants, enums and typedefs +// + +// +// static data member definitions +// + +// +// constructors and destructor +// +SiPixelLorentzAnglePCLWorker::SiPixelLorentzAnglePCLWorker(const edm::ParameterSet& iConfig) + : folder_(iConfig.getParameter("folder")), + ptmin_(iConfig.getParameter("ptMin")), + normChi2Max_(iConfig.getParameter("normChi2Max")), + clustSizeYMin_(iConfig.getParameter("clustSizeYMin")), + residualMax_(iConfig.getParameter("residualMax")), + clustChargeMax_(iConfig.getParameter("clustChargeMax")), + hist_depth_(iConfig.getParameter("binsDepth")), + hist_drift_(iConfig.getParameter("binsDrift")), + geomEsToken_(esConsumes()), + topoEsToken_(esConsumes()), + siPixelTemplateEsToken_(esConsumes()), + topoPerEventEsToken_(esConsumes()), + geomPerEventEsToken_(esConsumes()) { + t_trajTrack = consumes(iConfig.getParameter("src")); + + // now do what ever initialization is needed + hist_x_ = 50; + hist_y_ = 100; + min_x_ = -500.; + max_x_ = 500.; + min_y_ = -1500.; + max_y_ = 500.; + width_ = 0.0285; + min_depth_ = -100.; + max_depth_ = 400.; + min_drift_ = -1000.; //-200.;(iConfig.getParameter("residualMax")) + max_drift_ = 1000.; //400.; +} + +SiPixelLorentzAnglePCLWorker::~SiPixelLorentzAnglePCLWorker() { + // do anything here that needs to be done at desctruction time + // (e.g. close files, deallocate resources etc.) +} + +// +// member functions +// + +// ------------ method called for each event ------------ + +void SiPixelLorentzAnglePCLWorker::dqmAnalyze(edm::Event const& iEvent, + edm::EventSetup const& iSetup, + SiPixelLorentzAngleCalibrationHistograms const& iHists) const { + // Retrieve template stuff + const SiPixelTemplateDBObject* templateDBobject_ = &iSetup.getData(siPixelTemplateEsToken_); + std::vector thePixelTemp_; + SiPixelTemplate templ(thePixelTemp_); + if (!SiPixelTemplate::pushfile(*templateDBobject_, thePixelTemp_)) { + edm::LogError("SiPixelLorentzAnglePCLWorker") + << "\nERROR: Templates not filled correctly. Check the sqlite file." + << "Using SiPixelTemplateDBObject version " << (*templateDBobject_).version() << "\n\n"; + } + + // Retrieve tracker topology from geometry + const TrackerTopology* const tTopo = &iSetup.getData(topoPerEventEsToken_); + + // Retrieve track geometry + const TrackerGeometry* tracker = &iSetup.getData(geomPerEventEsToken_); + + // get the association map between tracks and trajectories + edm::Handle trajTrackCollectionHandle; + iEvent.getByToken(t_trajTrack, trajTrackCollectionHandle); + if (!trajTrackCollectionHandle->empty()) { + for (TrajTrackAssociationCollection::const_iterator it = trajTrackCollectionHandle->begin(); + it != trajTrackCollectionHandle->end(); + ++it) { + const reco::Track& track = *it->val; + const Trajectory& traj = *it->key; + + // get the trajectory measurements + std::vector tmColl = traj.measurements(); + float pt_ = track.pt(); + //float eta_ = track.eta(); + //float phi_ = track.phi(); + float chi2_ = traj.chiSquared(); + float ndof_ = traj.ndof(); + + if (pt_ < ptmin_) + continue; + // iterate over trajectory measurements + iHists.h_tracks_->Fill(0); + bool pixeltrack = false; + for (std::vector::const_iterator itTraj = tmColl.begin(); itTraj != tmColl.end(); + itTraj++) { + if (!itTraj->updatedState().isValid()) + continue; + TransientTrackingRecHit::ConstRecHitPointer recHit = itTraj->recHit(); + if (!recHit->isValid() || recHit->geographicalId().det() != DetId::Tracker) + continue; + unsigned int subDetID = (recHit->geographicalId().subdetId()); + if (subDetID == PixelSubdetector::PixelBarrel || subDetID == PixelSubdetector::PixelEndcap) { + if (!pixeltrack) { + iHists.h_tracks_->Fill(1); + } + pixeltrack = true; + } + + if (subDetID == PixelSubdetector::PixelBarrel) { + DetId detIdObj = recHit->geographicalId(); + const PixelGeomDetUnit* theGeomDet = dynamic_cast(tracker->idToDet(detIdObj)); + if (!theGeomDet) + continue; + + const PixelTopology* topol = &(theGeomDet->specificTopology()); + + if (!topol) + continue; + + const auto& layer_ = tTopo->pxbLayer(detIdObj); + //const auto& ladder_ = tTopo->pxbLadder(detIdObj); + const auto& module_ = tTopo->pxbModule(detIdObj); + + /* + float tmp1 = theGeomDet->surface().toGlobal(Local3DPoint(0.,0.,0.)).perp(); + float tmp2 = theGeomDet->surface().toGlobal(Local3DPoint(0.,0.,1.)).perp(); + + int isflipped_; + if ( tmp2((*recHit).hit()); + if (!recHitPix) + continue; + const auto& rechit_x = recHitPix->localPosition().x(); + const auto& rechit_y = recHitPix->localPosition().y(); + SiPixelRecHit::ClusterRef const& cluster = recHitPix->cluster(); + + Pixinfo pixinfo_ = fillPix(*cluster, topol); + + // fill entries in clust_ + /* + const auto& clust_x = (cluster)->x(); + const auto& clust_y = (cluster)->y(); + const auto& clust_charge = (cluster->charge())/1000.; + const auto& clust_size_x = cluster->sizeX(); + const auto& clust_size_y = cluster->sizeY(); + const auto& clust_maxPixelCol = cluster->maxPixelCol(); + const auto& clust_maxPixelRow = cluster->maxPixelRow(); + const auto& clust_minPixelCol = cluster->minPixelCol(); + const auto& clust_minPixelRow = cluster->minPixelRow(); + */ + + // fill the trackhit info + TrajectoryStateOnSurface tsos = itTraj->updatedState(); + if (!tsos.isValid()) { + edm::LogWarning("SiPixelLorentzAnglePCLWorker") << "tsos not valid" << std::endl; + continue; + } + LocalVector trackdirection = tsos.localDirection(); + LocalPoint trackposition = tsos.localPosition(); + + if (trackdirection.z() == 0) + continue; + // the local position and direction + const auto& trackhit_alpha = atan2(trackdirection.z(), trackdirection.x()); + const auto& trackhit_beta = atan2(trackdirection.z(), trackdirection.y()); + const auto& trackhit_gamma = atan2(trackdirection.x(), trackdirection.y()); + const auto& trackhit_x = trackposition.x(); + const auto& trackhit_y = trackposition.y(); + + // get qScale_ = templ.qscale() and templ.r_qMeas_qTrue(); + + float cotalpha = 1. / TMath::Tan(trackhit_alpha); + float cotbeta = 1. / TMath::Tan(trackhit_beta); + float locBx = 1.; + if (cotbeta < 0.) + locBx = -1.; + float locBz = locBx; + if (cotalpha < 0.) + locBz = -locBx; + + auto detId = detIdObj.rawId(); + int TemplID = -9999; + TemplID = templateDBobject_->getTemplateID(detId); + templ.interpolate(TemplID, cotalpha, cotbeta, locBz, locBx); + //const auto& qScale_ = templ.qscale(); + //const auto& rQmQt_ = templ.r_qMeas_qTrue(); + + // Surface deformation + /* + const auto& lp_pair = surface_deformation(topol, tsos, recHitPix); + + LocalPoint lp_track = lp_pair.first; + LocalPoint lp_rechit = lp_pair.second; + + const auto& rechitCorr_x = lp_rechit.x(); + const auto& rechitCorr_y = lp_rechit.y(); + const auto& trackhitCorrX_ = lp_track.x(); + const auto& trackhitCorrY_ = lp_track.y(); + */ + + // is one pixel in cluster a large pixel ? (hit will be excluded) + bool large_pix = false; + for (int j = 0; j < pixinfo_.npix; j++) { + int colpos = static_cast(pixinfo_.col[j]); + if (pixinfo_.row[j] == 0 || pixinfo_.row[j] == 79 || pixinfo_.row[j] == 80 || pixinfo_.row[j] == 159 || + colpos % 52 == 0 || colpos % 52 == 51) { + large_pix = true; + } + } + + double residual = TMath::Sqrt((trackhit_x - rechit_x) * (trackhit_x - rechit_x) + + (trackhit_y - rechit_y) * (trackhit_y - rechit_y)); + + if (!large_pix && (chi2_ / ndof_) < normChi2Max_ && cluster->sizeY() >= clustSizeYMin_ && + residual < residualMax_ && (cluster->charge() < clustChargeMax_)) { + // iterate over pixels in hit + for (int j = 0; j < pixinfo_.npix; j++) { + // use trackhits + float dx = (pixinfo_.x[j] - (trackhit_x - width_ / 2. / TMath::Tan(trackhit_alpha))) * 10000.; + float dy = (pixinfo_.y[j] - (trackhit_y - width_ / 2. / TMath::Tan(trackhit_beta))) * 10000.; + float depth = dy * tan(trackhit_beta); + float drift = dx - dy * tan(trackhit_gamma); + + int i_index = module_ + (layer_ - 1) * iHists.nModules_[layer_ - 1]; + iHists.h_drift_depth_adc_.at(i_index)->Fill(drift, depth, pixinfo_.adc[j]); + iHists.h_drift_depth_adc2_.at(i_index)->Fill(drift, depth, pixinfo_.adc[j] * pixinfo_.adc[j]); + iHists.h_drift_depth_noadc_.at(i_index)->Fill(drift, depth); + } + } + } else if (subDetID == PixelSubdetector::PixelEndcap) { + // what to do here? + edm::LogInfo("SiPixelLorentzAnglePCLWorker") << "do Nothing here?" << std::endl; + } + } //end iteration over trajectory measurements + } //end iteration over trajectories + } +} + +void SiPixelLorentzAnglePCLWorker::dqmBeginRun(edm::Run const& run, + edm::EventSetup const& iSetup, + SiPixelLorentzAngleCalibrationHistograms& iHists) const { + // geometry + const TrackerGeometry* geom = &iSetup.getData(geomEsToken_); + const TrackerTopology* tTopo = &iSetup.getData(topoEsToken_); + + PixelTopologyMap map = PixelTopologyMap(geom, tTopo); + iHists.nlay = geom->numberOfLayers(PixelSubdetector::PixelBarrel); + + for (int i = 0; i < iHists.nlay; i++) { + iHists.nModules_[i] = map.getPXBModules(i + 1); + } +} + +void SiPixelLorentzAnglePCLWorker::bookHistograms(DQMStore::IBooker& iBooker, + edm::Run const& run, + edm::EventSetup const& iSetup, + SiPixelLorentzAngleCalibrationHistograms& iHists) const { + iBooker.setCurrentFolder(folder_); + iHists.h_tracks_ = iBooker.book1D("h_tracks", "h_tracks", 2, 0., 2.); + + //book histograms + char name[128]; + for (int i_layer = 1; i_layer <= iHists.nlay; i_layer++) { + for (int i_module = 1; i_module <= iHists.nModules_[i_layer - 1]; i_module++) { + unsigned int i_index = i_module + (i_layer - 1) * iHists.nModules_[i_layer - 1]; + + sprintf(name, "h_drift_depth_adc_layer%i_module%i", i_layer, i_module); + iHists.h_drift_depth_adc_[i_index] = + iBooker.book2D(name, name, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_); + + sprintf(name, "h_drift_depth_adc2_layer%i_module%i", i_layer, i_module); + iHists.h_drift_depth_adc2_[i_index] = + iBooker.book2D(name, name, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_); + + sprintf(name, "h_drift_depth_noadc_layer%i_module%i", i_layer, i_module); + iHists.h_drift_depth_noadc_[i_index] = + iBooker.book2D(name, name, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_); + + sprintf(name, "h_drift_depth_layer%i_module%i", i_layer, i_module); + iHists.h_drift_depth_[i_index] = + iBooker.book2D(name, name, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_); + + sprintf(name, "h_mean_layer%i_module%i", i_layer, i_module); + iHists.h_mean_[i_index] = iBooker.book1D(name, name, hist_depth_, min_depth_, max_depth_); + } + } +} + +// method used to fill per pixel info +const Pixinfo SiPixelLorentzAnglePCLWorker::fillPix(const SiPixelCluster& LocPix, const PixelTopology* topol) const { + static Pixinfo pixinfo; + const std::vector& pixvector = LocPix.pixels(); + pixinfo.npix = 0; + for (std::vector::const_iterator itPix = pixvector.begin(); itPix != pixvector.end(); + itPix++) { + pixinfo.row[pixinfo.npix] = itPix->x; + pixinfo.col[pixinfo.npix] = itPix->y; + pixinfo.adc[pixinfo.npix] = itPix->adc; + LocalPoint lp = topol->localPosition(MeasurementPoint(itPix->x + 0.5, itPix->y + 0.5)); + pixinfo.x[pixinfo.npix] = lp.x(); + pixinfo.y[pixinfo.npix] = lp.y(); + pixinfo.npix++; + } + return pixinfo; +} + +// method used to correct for the surface deformation +const std::pair SiPixelLorentzAnglePCLWorker::surface_deformation( + const PixelTopology* topol, TrajectoryStateOnSurface& tsos, const SiPixelRecHit* recHitPix) const { + LocalPoint trackposition = tsos.localPosition(); + const LocalTrajectoryParameters& ltp = tsos.localParameters(); + const Topology::LocalTrackAngles localTrackAngles(ltp.dxdz(), ltp.dydz()); + + std::pair pixels_track = topol->pixel(trackposition, localTrackAngles); + std::pair pixels_rechit = topol->pixel(recHitPix->localPosition(), localTrackAngles); + + LocalPoint lp_track = topol->localPosition(MeasurementPoint(pixels_track.first, pixels_track.second)); + + LocalPoint lp_rechit = topol->localPosition(MeasurementPoint(pixels_rechit.first, pixels_rechit.second)); + + static std::pair lps = std::make_pair(lp_track, lp_rechit); + return lps; +} + +// ------------ method fills 'descriptions' with the allowed parameters for the module ------------ +void SiPixelLorentzAnglePCLWorker::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + desc.add("folder", "AlCaReco/SiPixelLorentzAngle"); + desc.add("src", edm::InputTag("TrackRefitter")); + desc.add("ptMin", 3.); + desc.add("normChi2Max", 2.); + desc.add("clustSizeYMin", 4); + desc.add("residualMax", 0.005); + desc.add("clustChargeMax", 120000); + desc.add("binsDepth", 50); + desc.add("binsDrift", 200); + descriptions.addWithDefaultLabel(desc); +} + +// define this as a plug-in +DEFINE_FWK_MODULE(SiPixelLorentzAnglePCLWorker); diff --git a/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_Output_cff.py b/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_Output_cff.py new file mode 100644 index 0000000000000..af2ec0ab2aa4b --- /dev/null +++ b/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_Output_cff.py @@ -0,0 +1,15 @@ +import FWCore.ParameterSet.Config as cms + +OutALCARECOPromptCalibProdSiPixelLorentzAngle_noDrop = cms.PSet( + SelectEvents = cms.untracked.PSet( + SelectEvents = cms.vstring('pathALCARECOPromptCalibProdSiPixelLorentzAngle') + ), + outputCommands = cms.untracked.vstring( + 'keep *_alcaBeamSpotProducer_*_*', + 'keep *_MEtoEDMConvertSiPixelLorentzAngle_*_*', + ) +) + +import copy +OutALCARECOPromptCalibProdSiPixelLorentzAngle=copy.deepcopy(OutALCARECOPromptCalibProdSiPixelLorentzAngle_noDrop) +OutALCARECOPromptCalibProdSiPixelLorentzAngle.outputCommands.insert(0, "drop *") diff --git a/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py b/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py new file mode 100644 index 0000000000000..91fe152ab9b67 --- /dev/null +++ b/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py @@ -0,0 +1,70 @@ +import FWCore.ParameterSet.Config as cms + +# ------------------------------------------------------------------------------ +# configure a filter to run only on the events selected by SiPixelCalSingleMuonLoose AlcaReco +import copy +from HLTrigger.HLTfilters.hltHighLevel_cfi import * +ALCARECOCalSignleMuonFilterForSiPixelLorentzAngle = copy.deepcopy(hltHighLevel) +ALCARECOCalSignleMuonFilterForSiPixelLorentzAngle.HLTPaths = ['pathALCARECOSiPixelCalSingleMuon'] +ALCARECOCalSignleMuonFilterForSiPixelLorentzAngle.throw = True ## dont throw on unknown path names +ALCARECOCalSignleMuonFilterForSiPixelLorentzAngle.TriggerResultsTag = cms.InputTag("TriggerResults","","RECO") + +# ------------------------------------------------------------------------------ +# This is the sequence for track refitting of the track saved by SiPixelCalSingleMuonLoose +# to have access to transient objects produced during RECO step and not saved + +from Alignment.CommonAlignmentProducer.AlignmentTrackSelector_cfi import * +ALCARECOPixelLACalibrationTracks = AlignmentTrackSelector.clone( + src = 'ALCARECOSiPixelCalSingleMuon', + filter = True, + applyBasicCuts = True, + ptMin = 3. +) + +# FIXME: the beam-spot should be kept in the AlCaReco (if not already there) and dropped from here +from RecoVertex.BeamSpotProducer.BeamSpot_cff import * + +from RecoTracker.IterativeTracking.InitialStep_cff import * +from RecoTracker.Configuration.RecoTrackerP5_cff import * +from RecoTracker.TrackProducer.TrackRefitter_cfi import * +from TrackingTools.TransientTrack.TransientTrackBuilder_cfi import * + +ALCARECOPixelLACalibrationTracksRefit = TrackRefitter.clone(src = cms.InputTag("ALCARECOCalibrationTracks"), + TrajectoryInEvent = cms.bool(True), + NavigationSchool = cms.string("") + ) + +# refit and BS can be dropped if done together with RECO. +# track filter can be moved in acalreco if no otehr users +ALCARECOPixelLATrackFilterRefit = cms.Sequence(ALCARECOPixelLACalibrationTracks + + offlineBeamSpot + + ALCARECOPixelLACalibrationTracksRefit) + +# ------------------------------------------------------------------------------ +# This is the module actually doing the calibration +from CalibTracker.SiPixelLorentzAngle.SiPixelLorentzAnglePCLWorker_cfi import SiPixelLorentzAnglePCLWorker +ALCARECOSiPixelLACalib = SiPixelLorentzAnglePCLWorker.clone( + folder = cms.string('AlCaReco/SiPixelLorentzAngle'), + src = cms.InputTag('ALCARECOPixelLACalibrationTracksRefit') +) +# ---------------------------------------------------------------------------- + +# **************************************************************************** +# ** Conversion for the SiPixelLorentzAngle DQM dir ** +# **************************************************************************** +MEtoEDMConvertSiPixelLorentzAngle = cms.EDProducer("MEtoEDMConverter", + Name = cms.untracked.string('MEtoEDMConverter'), + Verbosity = cms.untracked.int32(0), # 0 provides no output + # 1 provides basic output + # 2 provide more detailed output + Frequency = cms.untracked.int32(50), + MEPathToSave = cms.untracked.string('AlCaReco/SiPixelLorentzAngle'), + ) + +# The actual sequence +seqALCARECOPromptCalibProdSiPixelLorentzAngle = cms.Sequence( + ALCARECOCalSignleMuonFilterForSiPixelLorentzAngle * + ALCARECOPixelLATrackFilterRefit * + ALCARECOSiPixelLACalib * + MEtoEDMConvertSiPixelLorentzAngle + ) diff --git a/Configuration/AlCa/python/autoAlca.py b/Configuration/AlCa/python/autoAlca.py index 51f166826ffe7..0fafdd26dacb8 100644 --- a/Configuration/AlCa/python/autoAlca.py +++ b/Configuration/AlCa/python/autoAlca.py @@ -35,7 +35,7 @@ "DoubleMuParked" : "MuAlCalIsolatedMu+MuAlOverlaps+TkAlZMuMu", "MuOniaParked" : "TkAlJpsiMuMu+TkAlUpsilonMuMu", "DoubleElectron" : "EcalCalZElectron+EcalUncalZElectron+HcalCalIsoTrkFilter", - "StreamExpress" : "SiStripCalZeroBias+TkAlMinBias+SiStripPCLHistos+SiStripCalMinBias+SiStripCalMinBiasAAG+Hotline+LumiPixelsMinBias+SiPixelCalZeroBias", + "StreamExpress" : "SiStripCalZeroBias+TkAlMinBias+SiStripPCLHistos+SiStripCalMinBias+SiStripCalMinBiasAAG+Hotline+LumiPixelsMinBias+SiPixelCalZeroBias+SiPixelCalSingleMuon", "StreamExpressHI" : "SiStripCalZeroBias+TkAlMinBiasHI+SiStripPCLHistos+SiStripCalMinBias+SiStripCalMinBiasAAG+SiPixelCalZeroBias" } diff --git a/Configuration/EventContent/python/AlCaRecoOutput_cff.py b/Configuration/EventContent/python/AlCaRecoOutput_cff.py index fa6443a21035a..7b5e0e36c1af5 100644 --- a/Configuration/EventContent/python/AlCaRecoOutput_cff.py +++ b/Configuration/EventContent/python/AlCaRecoOutput_cff.py @@ -158,7 +158,7 @@ from Calibration.TkAlCaRecoProducers.ALCARECOPromptCalibProdSiStrip_Output_cff import * from Calibration.TkAlCaRecoProducers.ALCARECOPromptCalibProdSiStripGains_Output_cff import * from Calibration.TkAlCaRecoProducers.ALCARECOPromptCalibProdSiStripGainsAAG_Output_cff import * - +from Calibration.TkAlCaRecoProducers.ALCARECOPromptCalibProdSiPixelLorentzAngle_Output_cff import * from Calibration.TkAlCaRecoProducers.ALCARECOSiStripPCLHistos_Output_cff import * from Alignment.CommonAlignmentProducer.ALCARECOPromptCalibProdSiPixelAli_Output_cff import * diff --git a/Configuration/PyReleaseValidation/python/relval_steps.py b/Configuration/PyReleaseValidation/python/relval_steps.py index 5586e1657abfd..7923acfdac2fa 100644 --- a/Configuration/PyReleaseValidation/python/relval_steps.py +++ b/Configuration/PyReleaseValidation/python/relval_steps.py @@ -2608,13 +2608,13 @@ def gen2021HiMix(fragment,howMuch): '--conditions':'auto:run1_data', '--datatier':'ALCARECO', '--eventcontent':'ALCARECO'} -steps['ALCAEXP']={'-s':'ALCAOUTPUT:SiStripCalZeroBias+TkAlMinBias+Hotline+LumiPixelsMinBias+AlCaPCCZeroBiasFromRECO+AlCaPCCRandomFromRECO,ALCA:PromptCalibProd+PromptCalibProdSiStrip+PromptCalibProdSiStripGains+PromptCalibProdSiStripGainsAAG+PromptCalibProdSiPixelAli+PromptCalibProdSiPixel', +steps['ALCAEXP']={'-s':'ALCAOUTPUT:SiStripCalZeroBias+TkAlMinBias+Hotline+LumiPixelsMinBias+AlCaPCCZeroBiasFromRECO+AlCaPCCRandomFromRECO+SiPixelCalSingleMuon,ALCA:PromptCalibProd+PromptCalibProdSiStrip+PromptCalibProdSiStripGains+PromptCalibProdSiStripGainsAAG+PromptCalibProdSiPixelAli+PromptCalibProdSiPixel+PromptCalibProdSiPixelLorentzAngle', '--conditions':'auto:run1_data', '--datatier':'ALCARECO', '--eventcontent':'ALCARECO', '--triggerResultsProcess': 'RECO'} -steps['ALCAEXPRUN2']={'-s':'ALCAOUTPUT:SiStripCalZeroBias+TkAlMinBias+LumiPixelsMinBias+AlCaPCCZeroBiasFromRECO+AlCaPCCRandomFromRECO+SiPixelCalZeroBias,ALCA:PromptCalibProd+PromptCalibProdSiStrip+PromptCalibProdSiStripGains+PromptCalibProdSiStripGainsAAG+PromptCalibProdSiPixelAli+PromptCalibProdSiPixel', +steps['ALCAEXPRUN2']={'-s':'ALCAOUTPUT:SiStripCalZeroBias+TkAlMinBias+LumiPixelsMinBias+AlCaPCCZeroBiasFromRECO+AlCaPCCRandomFromRECO+SiPixelCalZeroBias+SiPixelCalSingleMuon,ALCA:PromptCalibProd+PromptCalibProdSiStrip+PromptCalibProdSiStripGains+PromptCalibProdSiStripGainsAAG+PromptCalibProdSiPixelAli+PromptCalibProdSiPixel+PromptCalibProdSiPixelLorentzAngle', '--conditions':'auto:run2_data', '--datatier':'ALCARECO', '--eventcontent':'ALCARECO', diff --git a/Configuration/StandardSequences/python/AlCaRecoStreams_cff.py b/Configuration/StandardSequences/python/AlCaRecoStreams_cff.py index c999d0a17b723..b7a9757a2714a 100644 --- a/Configuration/StandardSequences/python/AlCaRecoStreams_cff.py +++ b/Configuration/StandardSequences/python/AlCaRecoStreams_cff.py @@ -152,6 +152,8 @@ from Calibration.TkAlCaRecoProducers.ALCARECOPromptCalibProdSiStripGains_cff import * from Calibration.TkAlCaRecoProducers.ALCARECOPromptCalibProdSiStripGainsAAG_cff import * +from Calibration.TkAlCaRecoProducers.ALCARECOPromptCalibProdSiPixelLorentzAngle_cff import * + from Calibration.TkAlCaRecoProducers.ALCARECOSiStripPCLHistos_cff import * from Alignment.CommonAlignmentProducer.ALCARECOPromptCalibProdSiPixelAli_cff import * @@ -265,6 +267,7 @@ pathALCARECOPromptCalibProdSiStrip = cms.Path(seqALCARECOPromptCalibProdSiStrip) pathALCARECOPromptCalibProdSiStripGains = cms.Path(seqALCARECOPromptCalibProdSiStripGains) pathALCARECOPromptCalibProdSiStripGainsAAG = cms.Path(seqALCARECOPromptCalibProdSiStripGainsAAG) +pathALCARECOPromptCalibProdSiPixelLorentzAngle = cms.Path(seqALCARECOPromptCalibProdSiPixelLorentzAngle) pathALCARECOPromptCalibProdSiPixelAli = cms.Path(seqALCARECOPromptCalibProdSiPixelAli) pathALCARECOPromptCalibProdSiPixel = cms.Path(seqALCARECOPromptCalibProdSiPixel) pathALCARECOPromptCalibProdEcalPedestals = cms.Path(seqALCARECOPromptCalibProdEcalPedestals) @@ -953,7 +956,14 @@ dataTier = cms.untracked.string('ALCARECO') ) - +ALCARECOStreamPromptCalibProdSiPixelLorentzAngle = cms.FilteredStream( + responsible = 'Marco Musich', + name = 'PromptCalibProdSiPixelLorentzAngle', + paths = (pathALCARECOPromptCalibProdSiPixelLorentzAngle), + content = OutALCARECOPromptCalibProdSiPixelLorentzAngle.outputCommands, + selectEvents = OutALCARECOPromptCalibProdSiPixelLorentzAngle.SelectEvents, + dataTier = cms.untracked.string('ALCARECO') + ) ALCARECOStreamPromptCalibProdSiPixelAli = cms.FilteredStream( responsible = 'Gianluca Cerminara', From 7fcb87b2e52e22c559dee686a11d9aba5716daee Mon Sep 17 00:00:00 2001 From: mmusich Date: Thu, 3 Jun 2021 15:51:21 +0200 Subject: [PATCH 02/12] add stub of SiPixelLorentzAnglePCLHarvester --- .../SiPixelLorentzAngle/BuildFile.xml | 1 + .../SiPixelLorentzAnglePCLHarvester_cfi.py | 7 + .../src/SiPixelLorentzAnglePCLHarvester.cc | 227 ++++++++++++++++++ .../AlcaSiPixelLorentzAngleHarvester_cff.py | 17 ++ .../AlcaSiPixelLorentzAngleHarvester_cfi.py | 4 + Configuration/AlCa/python/autoPCL.py | 4 +- .../python/relval_production.py | 2 +- .../python/relval_steps.py | 6 + .../python/AlCaHarvesting_cff.py | 11 + 9 files changed, 277 insertions(+), 2 deletions(-) create mode 100644 CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py create mode 100644 CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc create mode 100644 Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cff.py create mode 100644 Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cfi.py diff --git a/CalibTracker/SiPixelLorentzAngle/BuildFile.xml b/CalibTracker/SiPixelLorentzAngle/BuildFile.xml index 68ed8d0fb831d..d01ba32f1258e 100644 --- a/CalibTracker/SiPixelLorentzAngle/BuildFile.xml +++ b/CalibTracker/SiPixelLorentzAngle/BuildFile.xml @@ -11,5 +11,6 @@ + diff --git a/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py new file mode 100644 index 0000000000000..04dd577d5edf9 --- /dev/null +++ b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py @@ -0,0 +1,7 @@ +import FWCore.ParameterSet.Config as cms +from DQMServices.Core.DQMEDHarvester import DQMEDHarvester + +SiPixelLorentzAnglePCLHarvester = DQMEDHarvester( + "SiPixelLorentzAnglePCLHarvester", + dqmDir = cms.string('AlCaReco/SiPixelLorentzAngle') +) diff --git a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc new file mode 100644 index 0000000000000..9df5c76f33fd8 --- /dev/null +++ b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc @@ -0,0 +1,227 @@ +#include "DQMServices/Core/interface/DQMEDHarvester.h" + +#include "FWCore/Framework/interface/Frameworkfwd.h" +#include "FWCore/Framework/interface/MakerMacros.h" +#include "FWCore/Framework/interface/ESHandle.h" +#include "FWCore/Framework/interface/EventSetup.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ServiceRegistry/interface/Service.h" + +#include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyMap.h" +#include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" + +#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" +#include "../interface/SiPixelLorentzAngleCalibrationStruct.h" +#include "CondCore/DBOutputService/interface/PoolDBOutputService.h" + +#include +#include + +//------------------------------------------------------------------------------ + +class SiPixelLorentzAnglePCLHarvester : public DQMEDHarvester { +public: + SiPixelLorentzAnglePCLHarvester(const edm::ParameterSet&); + void beginRun(const edm::Run&, const edm::EventSetup&) override; + + static void fillDescriptions(edm::ConfigurationDescriptions&); + +private: + void dqmEndJob(DQMStore::IBooker&, DQMStore::IGetter&) override; + void findMean(MonitorElement* h_drift_depth_adc_slice_, int i, int i_ring); + + // es tokens + edm::ESGetToken geomEsToken_; + edm::ESGetToken topoEsToken_; + + const std::string dqmDir_; + SiPixelLorentzAngleCalibrationHistograms hists; +}; + +//------------------------------------------------------------------------------ + +SiPixelLorentzAnglePCLHarvester::SiPixelLorentzAnglePCLHarvester(const edm::ParameterSet& iConfig) + : geomEsToken_(esConsumes()), + topoEsToken_(esConsumes()), + dqmDir_(iConfig.getParameter("dqmDir")) { + // first ensure DB output service is available + edm::Service poolDbService; + if (!poolDbService.isAvailable()) + throw cms::Exception("SiPixelLorentzAnglePCLHarvester") << "PoolDBService required"; +} + +//------------------------------------------------------------------------------ + +void SiPixelLorentzAnglePCLHarvester::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) { + // geometry + const TrackerGeometry* geom = &iSetup.getData(geomEsToken_); + const TrackerTopology* tTopo = &iSetup.getData(topoEsToken_); + + PixelTopologyMap map = PixelTopologyMap(geom, tTopo); + hists.nlay = geom->numberOfLayers(PixelSubdetector::PixelBarrel); + + for (int i = 0; i < hists.nlay; i++) { + hists.nModules_[i] = map.getPXBModules(i + 1); + } +} + +//------------------------------------------------------------------------------ + +void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter) { + // go in the right directory + iGetter.cd(); + iGetter.setCurrentFolder(dqmDir_); + + /* + const auto listOfHistos = iGetter.getMEs(); + for(const auto& hname : listOfHistos){ + const auto& histo = iGetter.get(dqmDir_+"/"+hname); + std::cout << hname << " name: " << histo->getName() << std::endl; + } + */ + + for (int i_layer = 1; i_layer <= hists.nlay; i_layer++) { + for (int i_module = 1; i_module <= hists.nModules_[i_layer - 1]; i_module++) { + int i_index = i_module + (i_layer - 1) * hists.nModules_[i_layer - 1]; + + hists.h_drift_depth_[i_index] = + iGetter.get(fmt::format("{}/h_drift_depth_layer{}_module{}", dqmDir_, i_layer, i_module)); + + hists.h_drift_depth_adc_[i_index] = + iGetter.get(fmt::format("{}/h_drift_depth_adc_layer{}_module{}", dqmDir_, i_layer, i_module)); + + hists.h_drift_depth_adc2_[i_index] = + iGetter.get(fmt::format("{}/h_drift_depth_adc2_layer{}_module{}", dqmDir_, i_layer, i_module)); + + hists.h_drift_depth_noadc_[i_index] = + iGetter.get(fmt::format("{}/h_drift_depth_noadc_layer{}_module{}", dqmDir_, i_layer, i_module)); + + hists.h_mean_[i_index] = iGetter.get(fmt::format("{}/h_mean_layer{}_module{}", dqmDir_, i_layer, i_module)); + + hists.h_drift_depth_[i_index]->divide( + hists.h_drift_depth_adc_[i_index], hists.h_drift_depth_noadc_[i_index], 1., 1., ""); + } + } + + /* + for(const auto& [index,histo] : hists.h_drift_depth_adc_){ + std::cout << index << " => " << histo->getName(); + } + */ + + int hist_drift_ = hists.h_drift_depth_adc_[1]->getNbinsX(); + int hist_depth_ = hists.h_drift_depth_adc_[1]->getNbinsY(); + double min_drift_ = hists.h_drift_depth_adc_[1]->getAxisMin(1); + double max_drift_ = hists.h_drift_depth_adc_[1]->getAxisMax(1); + + iBooker.setCurrentFolder("AlCaReco/SiPixelLorentzAngleHarvesting/"); + MonitorElement* h_drift_depth_adc_slice_ = + iBooker.book1D("h_drift_depth_adc_slice", "slice of adc histogram", hist_drift_, min_drift_, max_drift_); + + TF1* f1 = new TF1("f1", "[0] + [1]*x", 50., 235.); + f1->SetParName(0, "p0"); + f1->SetParName(1, "p1"); + f1->SetParameter(0, 0); + f1->SetParameter(1, 0.4); + std::cout << "module" + << "\t" + << "layer" + << "\t" + << "offset" + << "\t" + << "error" + << "\t" + << "slope" + << "\t" + << "error" + << "\t" + "rel.err" + << "\t" + "pull" + << "\t" + << "chi2" + << "\t" + << "prob" << std::endl; + //loop over modlues and layers to fit the lorentz angle + for (int i_layer = 1; i_layer <= hists.nlay; i_layer++) { + for (int i_module = 1; i_module <= hists.nModules_[i_layer - 1]; i_module++) { + int i_index = i_module + (i_layer - 1) * hists.nModules_[i_layer - 1]; + //loop over bins in depth (z-local-coordinate) (in order to fit slices) + for (int i = 1; i <= hist_depth_; i++) { + //std::cout << i_layer << " " << i_module << " " << i << std::endl; + + findMean(h_drift_depth_adc_slice_, i, i_index); + } // end loop over bins in depth + hists.h_mean_[i_index]->getTH1()->Fit(f1, "ERQ"); + double p0 = f1->GetParameter(0); + double e0 = f1->GetParError(0); + double p1 = f1->GetParameter(1); + double e1 = f1->GetParError(1); + double chi2 = f1->GetChisquare(); + double prob = f1->GetProb(); + std::cout << std::setprecision(4) << i_module << "\t" << i_layer << "\t" << p0 << "\t" << e0 << "\t" << p1 + << std::setprecision(3) << "\t" << e1 << "\t" << e1 / p1 * 100. << "\t" << (p1 - 0.424) / e1 << "\t" + << chi2 << "\t" << prob << std::endl; + } + } // end loop over modules and layers + + /* + // fill the DB object record + + // write the object + edm::Service poolDbService; + poolDbService->writeOne(& , poolDbService->currentTime(), "SiPixelLorentzAngleRcd"); + */ +} + +void SiPixelLorentzAnglePCLHarvester::findMean(MonitorElement* h_drift_depth_adc_slice_, int i, int i_ring) { + double nentries = 0; + h_drift_depth_adc_slice_->Reset(); + int hist_drift_ = h_drift_depth_adc_slice_->getNbinsX(); + + // determine sigma and sigma^2 of the adc counts and average adc counts + //loop over bins in drift width + for (int j = 1; j <= hist_drift_; j++) { + if (hists.h_drift_depth_noadc_[i_ring]->getBinContent(j, i) >= 1) { + /* + std::cout << hists.h_drift_depth_adc_[i_ring]->getBinContent(j, i) << std::endl; + std::cout << hists.h_drift_depth_noadc_[i_ring]->getBinContent(j, i) << std::endl; + std::cout << hists.h_drift_depth_adc2_[i_ring]->getBinContent(j, i) << std::endl; + */ + + double adc_error2 = (hists.h_drift_depth_adc2_[i_ring]->getBinContent(j, i) - + hists.h_drift_depth_adc_[i_ring]->getBinContent(j, i) * + hists.h_drift_depth_adc_[i_ring]->getBinContent(j, i) / + hists.h_drift_depth_noadc_[i_ring]->getBinContent(j, i)) / + hists.h_drift_depth_noadc_[i_ring]->getBinContent(j, i); + + hists.h_drift_depth_adc_[i_ring]->setBinError(j, i, sqrt(adc_error2)); + double error2 = adc_error2 / (hists.h_drift_depth_noadc_[i_ring]->getBinContent(j, i) - 1.); + hists.h_drift_depth_[i_ring]->setBinError(j, i, sqrt(error2)); + } else { + hists.h_drift_depth_[i_ring]->setBinError(j, i, 0); + hists.h_drift_depth_adc_[i_ring]->setBinError(j, i, 0); + } + h_drift_depth_adc_slice_->setBinContent(j, hists.h_drift_depth_adc_[i_ring]->getBinContent(j, i)); + h_drift_depth_adc_slice_->setBinError(j, hists.h_drift_depth_adc_[i_ring]->getBinError(j, i)); + nentries += hists.h_drift_depth_noadc_[i_ring]->getBinContent(j, i); + } // end loop over bins in drift width + + double mean = h_drift_depth_adc_slice_->getMean(1); + double error = 0; + if (nentries != 0) { + error = h_drift_depth_adc_slice_->getRMS(1) / std::sqrt(nentries); + } + + hists.h_mean_[i_ring]->setBinContent(i, mean); + hists.h_mean_[i_ring]->setBinError(i, error); +} + +//------------------------------------------------------------------------------ +void SiPixelLorentzAnglePCLHarvester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + desc.add("dqmDir", "AlCaReco/SiPixelLorentzAngle"); + descriptions.addWithDefaultLabel(desc); +} + +DEFINE_FWK_MODULE(SiPixelLorentzAnglePCLHarvester); diff --git a/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cff.py b/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cff.py new file mode 100644 index 0000000000000..d1920180bfb91 --- /dev/null +++ b/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cff.py @@ -0,0 +1,17 @@ +import FWCore.ParameterSet.Config as cms + +from Calibration.TkAlCaRecoProducers.AlcaSiPixelLorentzAngleHarvester_cfi import * +from DQMServices.Components.EDMtoMEConverter_cfi import * + +EDMtoMEConvertSiPixelLorentzAngle = EDMtoMEConverter.clone() +EDMtoMEConvertSiPixelLorentzAngle.lumiInputTag = cms.InputTag("MEtoEDMConvertSiPixelLorentzAngle","MEtoEDMConverterLumi") +EDMtoMEConvertSiPixelLorentzAngle.runInputTag = cms.InputTag("MEtoEDMConvertSiPixelLorentzAngle","MEtoEDMConverterRun") + +DQMStore = cms.Service("DQMStore") + +from DQMServices.Core.DQMEDHarvester import DQMEDHarvester +dqmEnvSiPixelLorentzAngle = DQMEDHarvester('DQMHarvestingMetadata', + subSystemFolder = cms.untracked.string('AlCaReco'), + ) + +ALCAHARVESTSiPixelLorentzAngle = cms.Sequence( EDMtoMEConvertSiPixelLorentzAngle + alcaSiPixelLorentzAngleHarvester + dqmEnvSiPixelLorentzAngle ) diff --git a/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cfi.py b/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cfi.py new file mode 100644 index 0000000000000..085169f4ac0f0 --- /dev/null +++ b/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cfi.py @@ -0,0 +1,4 @@ +import FWCore.ParameterSet.Config as cms + +from CalibTracker.SiPixelLorentzAngle.SiPixelLorentzAnglePCLHarvester_cfi import SiPixelLorentzAnglePCLHarvester +alcaSiPixelLorentzAngleHarvester = SiPixelLorentzAnglePCLHarvester.clone() diff --git a/Configuration/AlCa/python/autoPCL.py b/Configuration/AlCa/python/autoPCL.py index 5d5a06f7ecd14..e0763522ec0bb 100644 --- a/Configuration/AlCa/python/autoPCL.py +++ b/Configuration/AlCa/python/autoPCL.py @@ -3,9 +3,11 @@ 'PromptCalibProdBeamSpotHPLowPU' : 'BeamSpotHPLowPUByRun+BeamSpotHPLowPUByLumi', 'PromptCalibProdSiStrip' : 'SiStripQuality', 'PromptCalibProdSiStripGains' : 'SiStripGains', + 'PromptCalibProdSiStripGainsAAG' : 'SiStripGainsAAG', 'PromptCalibProdSiPixelAli' : 'SiPixelAli', + 'PromptCalibProdSiPixel' : 'SiPixelQuality', + 'PromptCalibProdSiPixelLA' : 'SiPixelLA', 'PromptCalibProdEcalPedestals': 'EcalPedestals', - 'PromptCalibProdSiStripGainsAAG' : 'SiStripGainsAAG', 'PromptCalibProdLumiPCC': 'LumiPCC', 'PromptCalibProdSiPixel' : 'SiPixelQuality', 'PromptCalibProdPPS' : 'PPSTimingCalibration', diff --git a/Configuration/PyReleaseValidation/python/relval_production.py b/Configuration/PyReleaseValidation/python/relval_production.py index a7bbc199faf7d..f0f540b6cbc1a 100644 --- a/Configuration/PyReleaseValidation/python/relval_production.py +++ b/Configuration/PyReleaseValidation/python/relval_production.py @@ -11,7 +11,7 @@ ## data production test workflows[1000] = [ '',['RunMinBias2011A','TIER0','SKIMD','HARVESTDfst2','ALCASPLIT']] -workflows[1001] = [ '',['RunMinBias2011A','TIER0EXP','ALCAEXP','ALCAHARVDSIPIXELCALRUN1','ALCAHARVD1','ALCAHARVD2','ALCAHARVD3','ALCAHARVD4','ALCAHARVD5']] +workflows[1001] = [ '',['RunMinBias2011A','TIER0EXP','ALCAEXP','ALCAHARVDSIPIXELCALRUN1','ALCAHARVD1','ALCAHARVD2','ALCAHARVD3','ALCAHARVD4','ALCAHARVD5','ALCAHARVD7']] workflows[1001.2] = [ '',['RunZeroBias2017F','TIER0EXPRUN2','ALCAEXPRUN2','ALCAHARVDSIPIXELCAL']] workflows[1002]=['RRD',['RunMinBias2011A','RECODR1','COPYPASTE']] diff --git a/Configuration/PyReleaseValidation/python/relval_steps.py b/Configuration/PyReleaseValidation/python/relval_steps.py index 7923acfdac2fa..af8689b01b987 100644 --- a/Configuration/PyReleaseValidation/python/relval_steps.py +++ b/Configuration/PyReleaseValidation/python/relval_steps.py @@ -2712,6 +2712,12 @@ def gen2021HiMix(fragment,howMuch): '--data':'', '--filein':'file:PromptCalibProdSiPixel.root'} +steps['ALCAHARVD7']={'-s':'ALCAHARVEST:%s'%(autoPCL['PromptCalibProdSiPixelLA']), + '--conditions':'auto:run1_data', + '--scenario':'pp', + '--data':'', + '--filein':'file:PromptCalibProdSiPixelLorentzAngle.root'} + steps['ALCAHARVD5HI']=merge([{'--scenario':'HeavyIons'},steps['ALCAHARVD5']]) steps['ALCAHARVDTE']={'-s':'ALCAHARVEST:%s'%(autoPCL['PromptCalibProdEcalPedestals']), '--conditions':'auto:run2_data', diff --git a/Configuration/StandardSequences/python/AlCaHarvesting_cff.py b/Configuration/StandardSequences/python/AlCaHarvesting_cff.py index 3d4ced64d0675..583a05feb7453 100644 --- a/Configuration/StandardSequences/python/AlCaHarvesting_cff.py +++ b/Configuration/StandardSequences/python/AlCaHarvesting_cff.py @@ -5,6 +5,7 @@ from Calibration.TkAlCaRecoProducers.AlcaSiStripQualityHarvester_cff import * from Calibration.TkAlCaRecoProducers.AlcaSiStripGainsHarvester_cff import * from Calibration.TkAlCaRecoProducers.AlcaSiStripGainsAAGHarvester_cff import * +from Calibration.TkAlCaRecoProducers.AlcaSiPixelLorentzAngleHarvester_cff import * from Alignment.CommonAlignmentProducer.AlcaSiPixelAliHarvester_cff import * from Calibration.EcalCalibAlgos.AlcaEcalPedestalsHarvester_cff import * from Calibration.LumiAlCaRecoProducers.AlcaLumiPCCHarvester_cff import * @@ -155,6 +156,15 @@ timetype = cms.untracked.string('runnumber') ) +# -------------------------------------------------------------------------------------- +# SiPixel Lorentz Angle +ALCAHARVESTSiPixelLA_metadata = cms.PSet(record = cms.untracked.string('SiPixelLorentzAngleRcd')) + +ALCAHARVESTSiPixelLA_dbOutput = cms.PSet(record = cms.string('SiPixelLorentzAngleRcd'), + tag = cms.string('SiPixelLA_pcl'), + timetype = cms.untracked.string('runnumber') + ) + # -------------------------------------------------------------------------------------- # ECAL Pedestals ALCAHARVESTEcalPedestals_metadata = cms.PSet(record = cms.untracked.string('EcalPedestalsRcd')) @@ -239,6 +249,7 @@ SiStripQuality = cms.Path(ALCAHARVESTSiStripQuality) SiStripGains = cms.Path(ALCAHARVESTSiStripGains) SiPixelAli = cms.Path(ALCAHARVESTSiPixelAli) +SiPixelLA = cms.Path(ALCAHARVESTSiPixelLorentzAngle) EcalPedestals = cms.Path(ALCAHARVESTEcalPedestals) SiStripGainsAAG = cms.Path(ALCAHARVESTSiStripGainsAAG) LumiPCC = cms.Path(ALCAHARVESTLumiPCC) From f230a91e66e856041a191e26b1d0385f7915f1b7 Mon Sep 17 00:00:00 2001 From: mmusich Date: Fri, 4 Jun 2021 18:35:21 +0200 Subject: [PATCH 03/12] fixed collection of tracks passed to track refitter --- .../python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py b/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py index 91fe152ab9b67..2c723527933a3 100644 --- a/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py +++ b/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py @@ -29,7 +29,7 @@ from RecoTracker.TrackProducer.TrackRefitter_cfi import * from TrackingTools.TransientTrack.TransientTrackBuilder_cfi import * -ALCARECOPixelLACalibrationTracksRefit = TrackRefitter.clone(src = cms.InputTag("ALCARECOCalibrationTracks"), +ALCARECOPixelLACalibrationTracksRefit = TrackRefitter.clone(src = cms.InputTag("ALCARECOPixelLACalibrationTracks"), TrajectoryInEvent = cms.bool(True), NavigationSchool = cms.string("") ) From 00184e9d6be4dad08c646bb14d869084fa3f3032 Mon Sep 17 00:00:00 2001 From: mmusich Date: Mon, 7 Jun 2021 16:12:07 +0200 Subject: [PATCH 04/12] push code to write out payload, improve functionality, follow-up to review --- .../SiPixelLorentzAngleCalibrationStruct.h | 13 +- .../SiPixelLorentzAnglePCLHarvester_cfi.py | 13 +- .../SiPixelLorentzAnglePCLWorker_cfi.py | 23 +- .../src/SiPixelLorentzAnglePCLHarvester.cc | 408 +++++++++++-- .../src/SiPixelLorentzAnglePCLWorker.cc | 568 ++++++++++++++---- .../test/testPCLAlCaHarvesting.py | 7 +- Configuration/AlCa/python/autoPCL.py | 1 - 7 files changed, 854 insertions(+), 179 deletions(-) diff --git a/CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h b/CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h index c68bf7885c136..9f9be5aca4604 100644 --- a/CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h +++ b/CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h @@ -10,9 +10,18 @@ struct SiPixelLorentzAngleCalibrationHistograms { using MonitorMap = std::unordered_map; - // hardcode 4 BPix layers int nlay; - int nModules_[4]; + std::vector nModules_; + std::vector BPixnewmodulename_; + std::vector BPixnewDetIds_; + std::vector BPixnewModule_; + std::vector BPixnewLayer_; + + std::vector FPixnewmodulename_; + std::vector FPixnewDetIds_; + std::vector FPixnewDisk_; + std::vector FPixnewBlade_; + std::unordered_map > detIdsList; MonitorMap h_drift_depth_adc_; MonitorMap h_drift_depth_adc2_; diff --git a/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py index 04dd577d5edf9..92766f0917920 100644 --- a/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py +++ b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py @@ -3,5 +3,16 @@ SiPixelLorentzAnglePCLHarvester = DQMEDHarvester( "SiPixelLorentzAnglePCLHarvester", - dqmDir = cms.string('AlCaReco/SiPixelLorentzAngle') + newmodulelist = cms.vstring( + "BPix_BmI_SEC7_LYR2_LDR12F_MOD1", + "BPix_BmI_SEC8_LYR2_LDR14F_MOD1", + "BPix_BmO_SEC3_LYR2_LDR5F_MOD1", + "BPix_BmO_SEC3_LYR2_LDR5F_MOD2", + "BPix_BmO_SEC3_LYR2_LDR5F_MOD3", + "BPix_BpO_SEC1_LYR2_LDR1F_MOD1", + "BPix_BpO_SEC1_LYR2_LDR1F_MOD2", + "BPix_BpO_SEC1_LYR2_LDR1F_MOD3"), + dqmDir = cms.string('AlCaReco/SiPixelLorentzAngle'), + record = cms.string("SiPixelLorentzAngleRcd"), + fitProbCut = cms.double(0.5) ) diff --git a/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py index 04e139b5aba43..0ae9692e6b247 100644 --- a/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py +++ b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py @@ -1,15 +1,28 @@ import FWCore.ParameterSet.Config as cms from DQMServices.Core.DQMEDAnalyzer import DQMEDAnalyzer -SiPixelLorentzAnglePCLWorker = DQMEDAnalyzer( +SiPixelLorentzAnglePCLWorker = DQMEDAnalyzer( "SiPixelLorentzAnglePCLWorker", folder = cms.string('AlCaReco/SiPixelLorentzAngle'), - src = cms.InputTag("TrackRefitter"), - binsDepth = cms.int32(50), - binsDrift = cms.int32(200), + notInPCL = cms.bool(False), + fileName = cms.string('testrun.root'), + newmodulelist = cms.vstring( + "BPix_BmI_SEC7_LYR2_LDR12F_MOD1", + "BPix_BmI_SEC8_LYR2_LDR14F_MOD1", + "BPix_BmO_SEC3_LYR2_LDR5F_MOD1", + "BPix_BmO_SEC3_LYR2_LDR5F_MOD2", + "BPix_BmO_SEC3_LYR2_LDR5F_MOD3", + "BPix_BpO_SEC1_LYR2_LDR1F_MOD1", + "BPix_BpO_SEC1_LYR2_LDR1F_MOD2", + "BPix_BpO_SEC1_LYR2_LDR1F_MOD3"), + src = cms.InputTag("TrackRefitter"), + binsDepth = cms.int32(50), + binsDrift = cms.int32(200), ptMin = cms.double(3), normChi2Max = cms.double(2), clustSizeYMin = cms.int32(4), + clustSizeYMinL4 = cms.int32(3), + clustSizeXMax = cms.int32(5), residualMax = cms.double(0.005), - clustChargeMax = cms.double(120000) + clustChargeMaxPerLength = cms.double(50000) ) diff --git a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc index 9df5c76f33fd8..97c372bf6c293 100644 --- a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc +++ b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc @@ -6,16 +6,23 @@ #include "FWCore/Framework/interface/EventSetup.h" #include "FWCore/ParameterSet/interface/ParameterSet.h" #include "FWCore/ServiceRegistry/interface/Service.h" - +#include "DataFormats/TrackerCommon/interface/PixelBarrelName.h" +#include "DataFormats/TrackerCommon/interface/PixelEndcapName.h" #include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyMap.h" #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" -#include "../interface/SiPixelLorentzAngleCalibrationStruct.h" +#include "CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h" #include "CondCore/DBOutputService/interface/PoolDBOutputService.h" +#include "CondFormats/SiPixelObjects/interface/SiPixelLorentzAngle.h" +#include "CondFormats/DataRecord/interface/SiPixelLorentzAngleRcd.h" + +#include "MagneticField/Engine/interface/MagneticField.h" +#include "MagneticField/Records/interface/IdealMagneticFieldRecord.h" #include #include +#include //------------------------------------------------------------------------------ @@ -33,9 +40,18 @@ class SiPixelLorentzAnglePCLHarvester : public DQMEDHarvester { // es tokens edm::ESGetToken geomEsToken_; edm::ESGetToken topoEsToken_; + edm::ESGetToken siPixelLAEsToken_; + edm::ESGetToken magneticFieldToken_; + std::vector newmodulelist_; const std::string dqmDir_; + const double fitProbCut_; + const std::string recordName_; + std::unique_ptr f1; + SiPixelLorentzAngleCalibrationHistograms hists; + const SiPixelLorentzAngle* currentLorentzAngle; + const MagneticField* magField; }; //------------------------------------------------------------------------------ @@ -43,7 +59,12 @@ class SiPixelLorentzAnglePCLHarvester : public DQMEDHarvester { SiPixelLorentzAnglePCLHarvester::SiPixelLorentzAnglePCLHarvester(const edm::ParameterSet& iConfig) : geomEsToken_(esConsumes()), topoEsToken_(esConsumes()), - dqmDir_(iConfig.getParameter("dqmDir")) { + siPixelLAEsToken_(esConsumes()), + magneticFieldToken_(esConsumes()), + newmodulelist_(iConfig.getParameter>("newmodulelist")), + dqmDir_(iConfig.getParameter("dqmDir")), + fitProbCut_(iConfig.getParameter("fitProbCut")), + recordName_(iConfig.getParameter("record")) { // first ensure DB output service is available edm::Service poolDbService; if (!poolDbService.isAvailable()) @@ -57,12 +78,53 @@ void SiPixelLorentzAnglePCLHarvester::beginRun(const edm::Run& iRun, const edm:: const TrackerGeometry* geom = &iSetup.getData(geomEsToken_); const TrackerTopology* tTopo = &iSetup.getData(topoEsToken_); + magField = &iSetup.getData(magneticFieldToken_); + currentLorentzAngle = &iSetup.getData(siPixelLAEsToken_); + PixelTopologyMap map = PixelTopologyMap(geom, tTopo); hists.nlay = geom->numberOfLayers(PixelSubdetector::PixelBarrel); - + hists.nModules_.resize(hists.nlay); for (int i = 0; i < hists.nlay; i++) { hists.nModules_[i] = map.getPXBModules(i + 1); } + + if (!newmodulelist_.empty()) { + for (auto const& modulename : newmodulelist_) { + if (modulename.find("BPix_") != std::string::npos) { + PixelBarrelName bn(modulename, true); + const auto& detId = bn.getDetId(tTopo); + hists.BPixnewmodulename_.push_back(modulename); + hists.BPixnewDetIds_.push_back(detId.rawId()); + hists.BPixnewModule_.push_back(bn.moduleName()); + hists.BPixnewLayer_.push_back(bn.layerName()); + } else if (modulename.find("FPix_") != std::string::npos) { + PixelEndcapName en(modulename, true); + const auto& detId = en.getDetId(tTopo); + hists.FPixnewmodulename_.push_back(modulename); + hists.FPixnewDetIds_.push_back(detId.rawId()); + hists.FPixnewDisk_.push_back(en.diskName()); + hists.FPixnewBlade_.push_back(en.bladeName()); + } + } + } + + std::vector treatedIndices; + + for (auto det : geom->detsPXB()) { + const PixelGeomDetUnit* pixelDet = dynamic_cast(det); + const auto& layer = tTopo->pxbLayer(pixelDet->geographicalId()); + const auto& module = tTopo->pxbModule(pixelDet->geographicalId()); + int i_index = module + (layer - 1) * hists.nModules_[layer - 1]; + + uint32_t rawId = pixelDet->geographicalId().rawId(); + + if (std::find(treatedIndices.begin(), treatedIndices.end(), i_index) != treatedIndices.end()) { + hists.detIdsList.at(i_index).push_back(rawId); + } else { + hists.detIdsList.insert(std::pair>(i_index, {rawId})); + treatedIndices.push_back(i_index); + } + } } //------------------------------------------------------------------------------ @@ -72,14 +134,6 @@ void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMS iGetter.cd(); iGetter.setCurrentFolder(dqmDir_); - /* - const auto listOfHistos = iGetter.getMEs(); - for(const auto& hname : listOfHistos){ - const auto& histo = iGetter.get(dqmDir_+"/"+hname); - std::cout << hname << " name: " << histo->getName() << std::endl; - } - */ - for (int i_layer = 1; i_layer <= hists.nlay; i_layer++) { for (int i_module = 1; i_module <= hists.nModules_[i_layer - 1]; i_module++) { int i_index = i_module + (i_layer - 1) * hists.nModules_[i_layer - 1]; @@ -87,6 +141,12 @@ void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMS hists.h_drift_depth_[i_index] = iGetter.get(fmt::format("{}/h_drift_depth_layer{}_module{}", dqmDir_, i_layer, i_module)); + if (hists.h_drift_depth_[i_index] == nullptr) { + edm::LogError("SiPixelLorentzAnglePCLHarvester::dqmEndJob") + << "Failed to retrieve electron drift over depth for layer " << i_layer << ", module " << i_module << "."; + continue; + } + hists.h_drift_depth_adc_[i_index] = iGetter.get(fmt::format("{}/h_drift_depth_adc_layer{}_module{}", dqmDir_, i_layer, i_module)); @@ -103,75 +163,297 @@ void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMS } } - /* - for(const auto& [index,histo] : hists.h_drift_depth_adc_){ - std::cout << index << " => " << histo->getName(); + for (int i = 0; i < (int)hists.BPixnewDetIds_.size(); i++) { + int new_index = i + 1 + hists.nModules_[hists.nlay - 1] + (hists.nlay - 1) * hists.nModules_[hists.nlay - 1]; + + hists.h_drift_depth_adc_[new_index] = + iGetter.get(fmt::format("{}/h_BPixnew_drift_depth_{}", dqmDir_, hists.BPixnewmodulename_[i])); + + if (hists.h_drift_depth_adc_[new_index] == nullptr) { + edm::LogError("SiPixelLorentzAnglePCLHarvester::dqmEndJob") + << "Failed to retrieve electron drift over depth for new module " << hists.BPixnewmodulename_[i] << "."; + continue; + } + + hists.h_drift_depth_adc2_[new_index] = + iGetter.get(fmt::format("{}/h_BPixnew_drift_depth_adc_{}", dqmDir_, hists.BPixnewmodulename_[i])); + + hists.h_drift_depth_noadc_[new_index] = + iGetter.get(fmt::format("{}/h_BPixnew_drift_depth_adc2_{}", dqmDir_, hists.BPixnewmodulename_[i])); + + hists.h_drift_depth_[new_index] = + iGetter.get(fmt::format("{}/h_BPixnew_drift_depth_noadc_{}", dqmDir_, hists.BPixnewmodulename_[i])); + + hists.h_mean_[new_index] = iGetter.get(fmt::format("{}/h_BPixnew_mean_{}", dqmDir_, hists.BPixnewmodulename_[i])); + + hists.h_drift_depth_[new_index]->divide( + hists.h_drift_depth_adc_[new_index], hists.h_drift_depth_noadc_[new_index], 1., 1., ""); } - */ - int hist_drift_ = hists.h_drift_depth_adc_[1]->getNbinsX(); - int hist_depth_ = hists.h_drift_depth_adc_[1]->getNbinsY(); - double min_drift_ = hists.h_drift_depth_adc_[1]->getAxisMin(1); - double max_drift_ = hists.h_drift_depth_adc_[1]->getAxisMax(1); + int hist_drift_; + int hist_depth_; + double min_drift_; + double max_drift_; + + if (hists.h_drift_depth_adc_[1] != nullptr) { + hist_drift_ = hists.h_drift_depth_adc_[1]->getNbinsX(); + hist_depth_ = hists.h_drift_depth_adc_[1]->getNbinsY(); + min_drift_ = hists.h_drift_depth_adc_[1]->getAxisMin(1); + max_drift_ = hists.h_drift_depth_adc_[1]->getAxisMax(1); + } else { + hist_drift_ = 200; + hist_depth_ = 50; + min_drift_ = -1000.; + max_drift_ = 1000.; + } iBooker.setCurrentFolder("AlCaReco/SiPixelLorentzAngleHarvesting/"); MonitorElement* h_drift_depth_adc_slice_ = iBooker.book1D("h_drift_depth_adc_slice", "slice of adc histogram", hist_drift_, min_drift_, max_drift_); - TF1* f1 = new TF1("f1", "[0] + [1]*x", 50., 235.); - f1->SetParName(0, "p0"); - f1->SetParName(1, "p1"); - f1->SetParameter(0, 0); - f1->SetParameter(1, 0.4); - std::cout << "module" - << "\t" - << "layer" - << "\t" - << "offset" - << "\t" - << "error" - << "\t" - << "slope" - << "\t" - << "error" - << "\t" - "rel.err" - << "\t" - "pull" - << "\t" - << "chi2" - << "\t" - << "prob" << std::endl; + edm::LogPrint("LorentzAngle") << "module" + << "\t" + << "layer" + << "\t" + << "offset" + << "\t" + << "e0" + << "\t" + << "slope" + << "\t" + << "e1" + << "\t" + "rel.err" + << "\t" + "pull" + << "\t" + << "p2" + << "\t" + << "e2" + << "\t" + << "p3" + << "\t" + << "e3" + << "\t" + << "p4" + << "\t" + << "e4" + << "\t" + << "p5" + << "\t" + << "e5" + << "\t" + << "chi2" + << "\t" + << "prob" + << "\t" + << "newDetId" << std::endl; + + SiPixelLorentzAngle* LorentzAngle = new SiPixelLorentzAngle(); + + f1 = std::make_unique("f1", "[0] + [1]*x + [2]*x*x + [3]*x*x*x + [4]*x*x*x*x + [5]*x*x*x*x*x", 5., 280.); + f1->SetParName(0, "offset"); + f1->SetParName(1, "tan#theta_{LA}"); + f1->SetParName(2, "quad term"); + f1->SetParName(3, "cubic term"); + f1->SetParName(4, "quartic term"); + f1->SetParName(5, "quintic term"); + + double p1_simul_newmodule = 0.294044; + + for (int j = 0; j < (int)hists.BPixnewDetIds_.size(); j++) { + int new_index = j + 1 + hists.nModules_[hists.nlay - 1] + (hists.nlay - 1) * hists.nModules_[hists.nlay - 1]; + if (hists.h_drift_depth_adc_[new_index] == nullptr) + continue; + for (int i = 1; i <= hist_depth_; i++) { + findMean(h_drift_depth_adc_slice_, i, new_index); + } + + f1->SetParameter(0, 0); + f1->SetParError(0, 0); + f1->SetParameter(1, 0.4); + f1->SetParError(1, 0); + f1->SetParameter(2, 0.0); + f1->SetParError(2, 0); + f1->SetParameter(3, 0.0); + f1->SetParError(3, 0); + f1->SetParameter(4, 0.0); + f1->SetParError(4, 0); + f1->SetParameter(5, 0.0); + f1->SetParError(5, 0); + f1->SetChisquare(0); + + hists.h_mean_[new_index]->getTH1()->Fit(f1.get(), "ERQ"); + + double p0 = f1->GetParameter(0); + double e0 = f1->GetParError(0); + double p1 = f1->GetParameter(1); + double e1 = f1->GetParError(1); + double p2 = f1->GetParameter(2); + double e2 = f1->GetParError(2); + double p3 = f1->GetParameter(3); + double e3 = f1->GetParError(3); + double p4 = f1->GetParameter(4); + double e4 = f1->GetParError(4); + double p5 = f1->GetParameter(5); + double e5 = f1->GetParError(5); + double chi2 = f1->GetChisquare(); + double prob = f1->GetProb(); + + edm::LogPrint("LorentzAngle") << std::setprecision(4) << hists.BPixnewModule_[j] << "\t" << hists.BPixnewLayer_[j] + << "\t" << p0 << "\t" << e0 << "\t" << p1 << std::setprecision(3) << "\t" << e1 + << "\t" << e1 / p1 * 100. << "\t" << (p1 - p1_simul_newmodule) / e1 << "\t" << p2 + << "\t" << e2 << "\t" << p3 << "\t" << e3 << "\t" << p4 << "\t" << e4 << "\t" << p5 + << "\t" << e5 << "\t" << chi2 << "\t" << prob << "\t" << hists.BPixnewDetIds_[j] + << std::endl; + } + + double p1_simul[hists.nlay][hists.nModules_[hists.nlay - 1]]; + for (int i_layer = 1; i_layer <= hists.nlay; i_layer++) { + for (int i_module = 1; i_module <= hists.nModules_[i_layer - 1]; i_module++) { + if (i_layer == 1) + p1_simul[i_layer - 1][i_module - 1] = 0.436848; + if (i_layer == 2) + p1_simul[i_layer - 1][i_module - 1] = 0.25802; + if (i_layer == 3 && i_module <= 4) + p1_simul[i_layer - 1][i_module - 1] = 0.29374; + if (i_layer == 3 && i_module >= 5) + p1_simul[i_layer - 1][i_module - 1] = 0.31084; + if (i_layer == 4 && i_module <= 4) + p1_simul[i_layer - 1][i_module - 1] = 0.29944; + if (i_layer == 4 && i_module >= 5) + p1_simul[i_layer - 1][i_module - 1] = 0.31426; + } + } + //loop over modlues and layers to fit the lorentz angle for (int i_layer = 1; i_layer <= hists.nlay; i_layer++) { for (int i_module = 1; i_module <= hists.nModules_[i_layer - 1]; i_module++) { int i_index = i_module + (i_layer - 1) * hists.nModules_[i_layer - 1]; + if (hists.h_drift_depth_adc_[i_index] == nullptr) + continue; //loop over bins in depth (z-local-coordinate) (in order to fit slices) for (int i = 1; i <= hist_depth_; i++) { - //std::cout << i_layer << " " << i_module << " " << i << std::endl; - findMean(h_drift_depth_adc_slice_, i, i_index); } // end loop over bins in depth - hists.h_mean_[i_index]->getTH1()->Fit(f1, "ERQ"); + + f1->SetParameter(0, 0); + f1->SetParError(0, 0); + f1->SetParameter(1, 0.4); + f1->SetParError(1, 0); + f1->SetParameter(2, 0.0); + f1->SetParError(2, 0); + f1->SetParameter(3, 0.0); + f1->SetParError(3, 0); + f1->SetParameter(4, 0.0); + f1->SetParError(4, 0); + f1->SetParameter(5, 0.0); + f1->SetParError(5, 0); + f1->SetChisquare(0); + + hists.h_mean_[i_index]->getTH1()->Fit(f1.get(), "ERQ"); double p0 = f1->GetParameter(0); double e0 = f1->GetParError(0); double p1 = f1->GetParameter(1); double e1 = f1->GetParError(1); + double p2 = f1->GetParameter(2); + double e2 = f1->GetParError(2); + double p3 = f1->GetParameter(3); + double e3 = f1->GetParError(3); + double p4 = f1->GetParameter(4); + double e4 = f1->GetParError(4); + double p5 = f1->GetParameter(5); + double e5 = f1->GetParError(5); double chi2 = f1->GetChisquare(); double prob = f1->GetProb(); - std::cout << std::setprecision(4) << i_module << "\t" << i_layer << "\t" << p0 << "\t" << e0 << "\t" << p1 - << std::setprecision(3) << "\t" << e1 << "\t" << e1 / p1 * 100. << "\t" << (p1 - 0.424) / e1 << "\t" - << chi2 << "\t" << prob << std::endl; + + edm::LogPrint("LorentzAngle") << std::setprecision(4) << i_module << "\t" << i_layer << "\t" << p0 << "\t" << e0 + << "\t" << p1 << std::setprecision(3) << "\t" << e1 << "\t" << e1 / p1 * 100. + << "\t" << (p1 - p1_simul[i_layer - 1][i_module - 1]) / e1 << "\t" << p2 << "\t" + << e2 << "\t" << p3 << "\t" << e3 << "\t" << p4 << "\t" << e4 << "\t" << p5 << "\t" + << e5 << "\t" << chi2 << "\t" << prob << "\t" + << "null" << std::endl; + + const auto& detIdsToFill = hists.detIdsList.at(i_index); + + GlobalPoint center(0.0, 0.0, 0.0); + float theMagField = magField->inTesla(center).mag(); + + float bPixLorentzAnglePerTesla_; + // if the fit quality is OK + if (prob > fitProbCut_) { + for (const auto& id : detIdsToFill) { + bPixLorentzAnglePerTesla_ = p1 / theMagField; + if (!LorentzAngle->putLorentzAngle(id, bPixLorentzAnglePerTesla_)) { + edm::LogError("SiPixelLorentzAnglePCLHarvester") + << "[SiPixelLorentzAnglePCLHarvester::dqmEndRun] detid already exists" << std::endl; + } + } + } else { + // just copy the values from the existing payload + for (const auto& id : detIdsToFill) { + bPixLorentzAnglePerTesla_ = currentLorentzAngle->getLorentzAngle(id); + if (!LorentzAngle->putLorentzAngle(id, bPixLorentzAnglePerTesla_)) { + edm::LogError("SiPixelLorentzAnglePCLHarvester") + << "[SiPixelLorentzAnglePCLHarvester::dqmEndRun] detid already exists" << std::endl; + } + } + } } } // end loop over modules and layers - /* + // fill the rest of DetIds not filled above (for the moment FPix) + const auto& currentLAMap = currentLorentzAngle->getLorentzAngles(); + const auto& newLAMap = LorentzAngle->getLorentzAngles(); + std::vector currentLADets; + std::vector newLADets; + + std::transform(currentLAMap.begin(), + currentLAMap.end(), + std::back_inserter(currentLADets), + [](const std::map::value_type& pair) { return pair.first; }); + + std::transform(newLAMap.begin(), + newLAMap.end(), + std::back_inserter(newLADets), + [](const std::map::value_type& pair) { return pair.first; }); + + std::vector notCommon; + std::set_symmetric_difference( + currentLADets.begin(), currentLADets.end(), newLADets.begin(), newLADets.end(), std::back_inserter(notCommon)); + + for (const auto& id : notCommon) { + float fPixLorentzAnglePerTesla_ = currentLorentzAngle->getLorentzAngle(id); + if (!LorentzAngle->putLorentzAngle(id, fPixLorentzAnglePerTesla_)) { + edm::LogError("SiPixelLorentzAnglePCLHarvester") + << "[SiPixelLorentzAnglePCLHarvester::dqmEndRun] detid already exists" << std::endl; + } + } + + // book histogram of differences + MonitorElement* h_diffLA = iBooker.book1D( + "h_diffLA", "difference in #mu_{H}; #Delta #mu_{H}/#mu_{H} (old-new)/old [%];n. modules", 100, -10, 10); + + for (const auto& id : newLADets) { + float deltaMuHoverMuH = (currentLorentzAngle->getLorentzAngle(id) - LorentzAngle->getLorentzAngle(id)) / + currentLorentzAngle->getLorentzAngle(id); + h_diffLA->Fill(deltaMuHoverMuH); + } + // fill the DB object record - - // write the object - edm::Service poolDbService; - poolDbService->writeOne(& , poolDbService->currentTime(), "SiPixelLorentzAngleRcd"); - */ + edm::Service mydbservice; + if (mydbservice.isAvailable()) { + try { + mydbservice->writeOne(LorentzAngle, mydbservice->currentTime(), recordName_); + } catch (const cond::Exception& er) { + edm::LogError("SiPixelLorentzAngleDB") << er.what() << std::endl; + } catch (const std::exception& er) { + edm::LogError("SiPixelLorentzAngleDB") << "caught std::exception " << er.what() << std::endl; + } + } else { + edm::LogError("SiPixelLorentzAngleDB") << "Service is unavailable" << std::endl; + } + delete LorentzAngle; } void SiPixelLorentzAnglePCLHarvester::findMean(MonitorElement* h_drift_depth_adc_slice_, int i, int i_ring) { @@ -183,12 +465,6 @@ void SiPixelLorentzAnglePCLHarvester::findMean(MonitorElement* h_drift_depth_adc //loop over bins in drift width for (int j = 1; j <= hist_drift_; j++) { if (hists.h_drift_depth_noadc_[i_ring]->getBinContent(j, i) >= 1) { - /* - std::cout << hists.h_drift_depth_adc_[i_ring]->getBinContent(j, i) << std::endl; - std::cout << hists.h_drift_depth_noadc_[i_ring]->getBinContent(j, i) << std::endl; - std::cout << hists.h_drift_depth_adc2_[i_ring]->getBinContent(j, i) << std::endl; - */ - double adc_error2 = (hists.h_drift_depth_adc2_[i_ring]->getBinContent(j, i) - hists.h_drift_depth_adc_[i_ring]->getBinContent(j, i) * hists.h_drift_depth_adc_[i_ring]->getBinContent(j, i) / @@ -212,7 +488,6 @@ void SiPixelLorentzAnglePCLHarvester::findMean(MonitorElement* h_drift_depth_adc if (nentries != 0) { error = h_drift_depth_adc_slice_->getRMS(1) / std::sqrt(nentries); } - hists.h_mean_[i_ring]->setBinContent(i, mean); hists.h_mean_[i_ring]->setBinError(i, error); } @@ -220,7 +495,12 @@ void SiPixelLorentzAnglePCLHarvester::findMean(MonitorElement* h_drift_depth_adc //------------------------------------------------------------------------------ void SiPixelLorentzAnglePCLHarvester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { edm::ParameterSetDescription desc; - desc.add("dqmDir", "AlCaReco/SiPixelLorentzAngle"); + desc.add>( + "newmodulelist", + {"BPix_BmO_SEC3_LYR2_LDR5F_MOD2", "BPix_BpO_SEC1_LYR2_LDR1F_MOD1"}); // the layer 2 new module names + desc.add("dqmDir", "AlCaReco/SiPixelLorentzAngle"); // the directory of PCL Worker output + desc.add("fitProbCut", 0.5); + desc.add("record", "SiPixelLorentzAngleRcd"); descriptions.addWithDefaultLabel(desc); } diff --git a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc index 55daa1a914d60..c4c8da7ccb9c5 100644 --- a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc +++ b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc @@ -4,9 +4,7 @@ // Class: SiPixelLorentzAnglePCLWorker // /**\class SiPixelLorentzAnglePCLWorker SiPixelLorentzAnglePCLWorker.cc CalibTracker/SiPixelLorentzAnglePCLWorker/plugins/SiPixelLorentzAnglePCLWorker.cc - Description: [one line class summary] - Implementation: [Notes on implementation] */ @@ -21,7 +19,7 @@ // user include files #include "DQMServices/Core/interface/DQMStore.h" #include "FWCore/Framework/interface/Frameworkfwd.h" -#include "DQMServices/Core/interface/DQMGlobalEDAnalyzer.h" +#include "DQMServices/Core/interface/DQMEDAnalyzer.h" #include "FWCore/Framework/interface/Event.h" #include "FWCore/Framework/interface/MakerMacros.h" @@ -40,7 +38,8 @@ #include "DataFormats/TrackerCommon/interface/TrackerTopology.h" #include "Geometry/Records/interface/IdealGeometryRecord.h" #include "DataFormats/TrackReco/interface/TrackFwd.h" - +#include "DataFormats/TrackerCommon/interface/PixelBarrelName.h" +#include "DataFormats/TrackerCommon/interface/PixelEndcapName.h" #include "TrackingTools/Records/interface/TransientRecHitRecord.h" #include "TrackingTools/TrackFitters/interface/TrajectoryStateCombiner.h" #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h" @@ -51,11 +50,14 @@ #include "CondFormats/SiPixelTransient/interface/SiPixelTemplateDefs.h" #include "CondFormats/SiPixelTransient/interface/SiPixelTemplate.h" #include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyMap.h" +#include "Geometry/TrackerGeometryBuilder/interface/phase1PixelTopology.h" #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" #include "FWCore/ParameterSet/interface/ParameterSet.h" -#include "../interface/SiPixelLorentzAngleCalibrationStruct.h" - +#include "CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h" +#include +#include +#include // // class declaration // @@ -70,31 +72,97 @@ struct Pixinfo { float y[maxpix]; }; -class SiPixelLorentzAnglePCLWorker : public DQMGlobalEDAnalyzer { +struct Hit { + float x; + float y; + double alpha; + double beta; + double gamma; +}; +struct Clust { + float x; + float y; + float charge; + int size_x; + int size_y; + int maxPixelCol; + int maxPixelRow; + int minPixelCol; + int minPixelRow; +}; +struct Rechit { + float x; + float y; +}; + +class SiPixelLorentzAnglePCLWorker : public DQMEDAnalyzer { public: explicit SiPixelLorentzAnglePCLWorker(const edm::ParameterSet&); - ~SiPixelLorentzAnglePCLWorker() override; + ~SiPixelLorentzAnglePCLWorker() override = default; static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); private: - void bookHistograms(DQMStore::IBooker&, - edm::Run const&, - edm::EventSetup const&, - SiPixelLorentzAngleCalibrationHistograms&) const override; + void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override; + + void analyze(edm::Event const&, edm::EventSetup const&) override; - void dqmAnalyze(edm::Event const&, - edm::EventSetup const&, - SiPixelLorentzAngleCalibrationHistograms const&) const override; + void dqmBeginRun(edm::Run const&, edm::EventSetup const&) override; - void dqmBeginRun(edm::Run const&, edm::EventSetup const&, SiPixelLorentzAngleCalibrationHistograms&) const override; + void dqmEndRun(edm::Run const&, edm::EventSetup const&); const Pixinfo fillPix(const SiPixelCluster& LocPix, const PixelTopology* topol) const; const std::pair surface_deformation(const PixelTopology* topol, TrajectoryStateOnSurface& tsos, const SiPixelRecHit* recHitPix) const; // ------------ member data ------------ + SiPixelLorentzAngleCalibrationHistograms iHists; + std::string folder_; + bool notInPCL_; + std::string filename_; + std::vector newmodulelist_; + + // tree branches barrel + int run_; + long int event_; + int lumiblock_; + int bx_; + int orbit_; + int module_; + int ladder_; + int layer_; + int isflipped_; + float pt_; + float eta_; + float phi_; + double chi2_; + double ndof_; + Pixinfo pixinfo_; + Hit simhit_, trackhit_; + Clust clust_; + Rechit rechit_; + Rechit rechitCorr_; + float trackhitCorrX_; + float trackhitCorrY_; + float qScale_; + float rQmQt_; + + // tree branches forward + int sideF_; + int diskF_; + int bladeF_; + int panelF_; + int moduleF_; + Pixinfo pixinfoF_; + Hit simhitF_, trackhitF_; + Clust clustF_; + Rechit rechitF_; + Rechit rechitCorrF_; + float trackhitCorrXF_; + float trackhitCorrYF_; + float qScaleF_; + float rQmQtF_; // histogram etc int hist_x_; @@ -108,16 +176,23 @@ class SiPixelLorentzAnglePCLWorker : public DQMGlobalEDAnalyzer hFile_; + std::unique_ptr SiPixelLorentzAngleTreeBarrel_; + std::unique_ptr SiPixelLorentzAngleTreeForward_; + // es consumes edm::ESGetToken geomEsToken_; edm::ESGetToken topoEsToken_; @@ -142,11 +217,16 @@ class SiPixelLorentzAnglePCLWorker : public DQMGlobalEDAnalyzer("folder")), + notInPCL_(iConfig.getParameter("notInPCL")), + filename_(iConfig.getParameter("fileName")), + newmodulelist_(iConfig.getParameter>("newmodulelist")), ptmin_(iConfig.getParameter("ptMin")), normChi2Max_(iConfig.getParameter("normChi2Max")), clustSizeYMin_(iConfig.getParameter("clustSizeYMin")), + clustSizeYMinL4_(iConfig.getParameter("clustSizeYMinL4")), + clustSizeXMax_(iConfig.getParameter("clustSizeXMax")), residualMax_(iConfig.getParameter("residualMax")), - clustChargeMax_(iConfig.getParameter("clustChargeMax")), + clustChargeMaxPerLength_(iConfig.getParameter("clustChargeMaxPerLength")), hist_depth_(iConfig.getParameter("binsDepth")), hist_drift_(iConfig.getParameter("binsDrift")), geomEsToken_(esConsumes()), @@ -168,11 +248,87 @@ SiPixelLorentzAnglePCLWorker::SiPixelLorentzAnglePCLWorker(const edm::ParameterS max_depth_ = 400.; min_drift_ = -1000.; //-200.;(iConfig.getParameter("residualMax")) max_drift_ = 1000.; //400.; -} -SiPixelLorentzAnglePCLWorker::~SiPixelLorentzAnglePCLWorker() { - // do anything here that needs to be done at desctruction time - // (e.g. close files, deallocate resources etc.) + bufsize = 64000; + + // create tree structure + // Barrel pixel + if (notInPCL_) { + hFile_ = std::make_unique(filename_.c_str(), "RECREATE"); + SiPixelLorentzAngleTreeBarrel_ = + std::make_unique("SiPixelLorentzAngleTreeBarrel_", "SiPixel LorentzAngle tree barrel", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("run", &run_, "run/I", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("event", &event_, "event/l", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("lumiblock", &lumiblock_, "lumiblock/I", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("bx", &bx_, "bx/I", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("orbit", &orbit_, "orbit/I", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("module", &module_, "module/I", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("ladder", &ladder_, "ladder/I", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("layer", &layer_, "layer/I", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("isflipped", &isflipped_, "isflipped/I", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("pt", &pt_, "pt/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("eta", &eta_, "eta/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("phi", &phi_, "phi/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("chi2", &chi2_, "chi2/D", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("ndof", &ndof_, "ndof/D", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("trackhit", &trackhit_, "x/F:y/F:alpha/D:beta/D:gamma_/D", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("npix", &pixinfo_.npix, "npix/I", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("rowpix", pixinfo_.row, "row[npix]/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("colpix", pixinfo_.col, "col[npix]/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("adc", pixinfo_.adc, "adc[npix]/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("xpix", pixinfo_.x, "x[npix]/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("ypix", pixinfo_.y, "y[npix]/F", bufsize); + + SiPixelLorentzAngleTreeBarrel_->Branch( + "clust", + &clust_, + "x/F:y/F:charge/F:size_x/I:size_y/I:maxPixelCol/I:maxPixelRow:minPixelCol/I:minPixelRow/I", + bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("rechit", &rechit_, "x/F:y/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("rechit_corr", &rechitCorr_, "x/F:y/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("trackhitcorr_x", &trackhitCorrX_, "trackhitcorr_x/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("trackhitcorr_y", &trackhitCorrY_, "trackhitcorr_y/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("qScale", &qScale_, "qScale/F", bufsize); + SiPixelLorentzAngleTreeBarrel_->Branch("rQmQt", &rQmQt_, "rQmQt/F", bufsize); + // Forward pixel + + SiPixelLorentzAngleTreeForward_ = + std::make_unique("SiPixelLorentzAngleTreeForward_", "SiPixel LorentzAngle tree forward", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("run", &run_, "run/I", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("event", &event_, "event/l", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("lumiblock", &lumiblock_, "lumiblock/I", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("bx", &bx_, "bx/I", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("orbit", &orbit_, "orbit/I", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("side", &sideF_, "side/I", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("disk", &diskF_, "disk/I", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("blade", &bladeF_, "blade/I", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("panel", &panelF_, "panel/I", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("module", &moduleF_, "module/I", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("pt", &pt_, "pt/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("eta", &eta_, "eta/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("phi", &phi_, "phi/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("chi2", &chi2_, "chi2/D", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("ndof", &ndof_, "ndof/D", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("trackhit", &trackhitF_, "x/F:y/F:alpha/D:beta/D:gamma_/D", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("npix", &pixinfoF_.npix, "npix/I", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("rowpix", pixinfoF_.row, "row[npix]/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("colpix", pixinfoF_.col, "col[npix]/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("adc", pixinfoF_.adc, "adc[npix]/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("xpix", pixinfoF_.x, "x[npix]/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("ypix", pixinfoF_.y, "y[npix]/F", bufsize); + + SiPixelLorentzAngleTreeForward_->Branch( + "clust", + &clustF_, + "x/F:y/F:charge/F:size_x/I:size_y/I:maxPixelCol/I:maxPixelRow:minPixelCol/I:minPixelRow/I", + bufsize); + SiPixelLorentzAngleTreeForward_->Branch("rechit", &rechitF_, "x/F:y/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("rechit_corr", &rechitCorrF_, "x/F:y/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("trackhitcorr_x", &trackhitCorrXF_, "trackhitcorr_x/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("trackhitcorr_y", &trackhitCorrYF_, "trackhitcorr_y/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("qScale", &qScaleF_, "qScale/F", bufsize); + SiPixelLorentzAngleTreeForward_->Branch("rQmQt", &rQmQtF_, "rQmQt/F", bufsize); + } } // @@ -181,9 +337,7 @@ SiPixelLorentzAnglePCLWorker::~SiPixelLorentzAnglePCLWorker() { // ------------ method called for each event ------------ -void SiPixelLorentzAnglePCLWorker::dqmAnalyze(edm::Event const& iEvent, - edm::EventSetup const& iSetup, - SiPixelLorentzAngleCalibrationHistograms const& iHists) const { +void SiPixelLorentzAnglePCLWorker::analyze(edm::Event const& iEvent, edm::EventSetup const& iSetup) { // Retrieve template stuff const SiPixelTemplateDBObject* templateDBobject_ = &iSetup.getData(siPixelTemplateEsToken_); std::vector thePixelTemp_; @@ -203,6 +357,22 @@ void SiPixelLorentzAnglePCLWorker::dqmAnalyze(edm::Event const& iEvent, // get the association map between tracks and trajectories edm::Handle trajTrackCollectionHandle; iEvent.getByToken(t_trajTrack, trajTrackCollectionHandle); + + module_ = -1; + layer_ = -1; + ladder_ = -1; + isflipped_ = -1; + pt_ = -999; + eta_ = 999; + phi_ = 999; + pixinfo_.npix = 0; + + run_ = iEvent.id().run(); + event_ = iEvent.id().event(); + lumiblock_ = iEvent.luminosityBlock(); + bx_ = iEvent.bunchCrossing(); + orbit_ = iEvent.orbitNumber(); + if (!trajTrackCollectionHandle->empty()) { for (TrajTrackAssociationCollection::const_iterator it = trajTrackCollectionHandle->begin(); it != trajTrackCollectionHandle->end(); @@ -212,11 +382,11 @@ void SiPixelLorentzAnglePCLWorker::dqmAnalyze(edm::Event const& iEvent, // get the trajectory measurements std::vector tmColl = traj.measurements(); - float pt_ = track.pt(); - //float eta_ = track.eta(); - //float phi_ = track.phi(); - float chi2_ = traj.chiSquared(); - float ndof_ = traj.ndof(); + pt_ = track.pt(); + eta_ = track.eta(); + phi_ = track.phi(); + chi2_ = traj.chiSquared(); + ndof_ = traj.ndof(); if (pt_ < ptmin_) continue; @@ -246,43 +416,40 @@ void SiPixelLorentzAnglePCLWorker::dqmAnalyze(edm::Event const& iEvent, const PixelTopology* topol = &(theGeomDet->specificTopology()); + float ypitch_ = topol->pitch().second; + if (!topol) continue; - const auto& layer_ = tTopo->pxbLayer(detIdObj); - //const auto& ladder_ = tTopo->pxbLadder(detIdObj); - const auto& module_ = tTopo->pxbModule(detIdObj); + layer_ = tTopo->pxbLayer(detIdObj); + ladder_ = tTopo->pxbLadder(detIdObj); + module_ = tTopo->pxbModule(detIdObj); - /* - float tmp1 = theGeomDet->surface().toGlobal(Local3DPoint(0.,0.,0.)).perp(); - float tmp2 = theGeomDet->surface().toGlobal(Local3DPoint(0.,0.,1.)).perp(); + float tmp1 = theGeomDet->surface().toGlobal(Local3DPoint(0., 0., 0.)).perp(); + float tmp2 = theGeomDet->surface().toGlobal(Local3DPoint(0., 0., 1.)).perp(); - int isflipped_; - if ( tmp2((*recHit).hit()); if (!recHitPix) continue; - const auto& rechit_x = recHitPix->localPosition().x(); - const auto& rechit_y = recHitPix->localPosition().y(); + rechit_.x = recHitPix->localPosition().x(); + rechit_.y = recHitPix->localPosition().y(); SiPixelRecHit::ClusterRef const& cluster = recHitPix->cluster(); - Pixinfo pixinfo_ = fillPix(*cluster, topol); + pixinfo_ = fillPix(*cluster, topol); // fill entries in clust_ - /* - const auto& clust_x = (cluster)->x(); - const auto& clust_y = (cluster)->y(); - const auto& clust_charge = (cluster->charge())/1000.; - const auto& clust_size_x = cluster->sizeX(); - const auto& clust_size_y = cluster->sizeY(); - const auto& clust_maxPixelCol = cluster->maxPixelCol(); - const auto& clust_maxPixelRow = cluster->maxPixelRow(); - const auto& clust_minPixelCol = cluster->minPixelCol(); - const auto& clust_minPixelRow = cluster->minPixelRow(); - */ + + clust_.x = (cluster)->x(); + clust_.y = (cluster)->y(); + clust_.charge = (cluster->charge()) / 1000.; // clust_.charge: in the unit of 1000e + clust_.size_x = cluster->sizeX(); + clust_.size_y = cluster->sizeY(); + clust_.maxPixelCol = cluster->maxPixelCol(); + clust_.maxPixelRow = cluster->maxPixelRow(); + clust_.minPixelCol = cluster->minPixelCol(); + clust_.minPixelRow = cluster->minPixelRow(); // fill the trackhit info TrajectoryStateOnSurface tsos = itTraj->updatedState(); @@ -296,16 +463,22 @@ void SiPixelLorentzAnglePCLWorker::dqmAnalyze(edm::Event const& iEvent, if (trackdirection.z() == 0) continue; // the local position and direction - const auto& trackhit_alpha = atan2(trackdirection.z(), trackdirection.x()); - const auto& trackhit_beta = atan2(trackdirection.z(), trackdirection.y()); - const auto& trackhit_gamma = atan2(trackdirection.x(), trackdirection.y()); - const auto& trackhit_x = trackposition.x(); - const auto& trackhit_y = trackposition.y(); + trackhit_.alpha = atan2(trackdirection.z(), trackdirection.x()); + trackhit_.beta = atan2(trackdirection.z(), trackdirection.y()); + trackhit_.gamma = atan2(trackdirection.x(), trackdirection.y()); + trackhit_.x = trackposition.x(); + trackhit_.y = trackposition.y(); // get qScale_ = templ.qscale() and templ.r_qMeas_qTrue(); - float cotalpha = 1. / TMath::Tan(trackhit_alpha); - float cotbeta = 1. / TMath::Tan(trackhit_beta); + float cotalpha = trackdirection.x() / trackdirection.z(); + float cotbeta = trackdirection.y() / trackdirection.z(); + float cotbeta_min = clustSizeYMin_ * ypitch_ / width_; + if (fabs(cotbeta) <= cotbeta_min) + continue; + double drdz = sqrt(1. + cotalpha * cotalpha + cotbeta * cotbeta); + double clusterCharge_cut = clustChargeMaxPerLength_ * drdz; + float locBx = 1.; if (cotbeta < 0.) locBx = -1.; @@ -314,24 +487,34 @@ void SiPixelLorentzAnglePCLWorker::dqmAnalyze(edm::Event const& iEvent, locBz = -locBx; auto detId = detIdObj.rawId(); - int TemplID = -9999; - TemplID = templateDBobject_->getTemplateID(detId); + bool Isnew = false; + int DetId_index = -1; + for (int i = 0; i < (int)iHists.BPixnewDetIds_.size(); i++) { + if (detId == (unsigned int)iHists.BPixnewDetIds_[i]) + Isnew = true; + DetId_index = i; + } + + int TemplID = templateDBobject_->getTemplateID(detId); templ.interpolate(TemplID, cotalpha, cotbeta, locBz, locBx); - //const auto& qScale_ = templ.qscale(); - //const auto& rQmQt_ = templ.r_qMeas_qTrue(); + qScale_ = templ.qscale(); + rQmQt_ = templ.r_qMeas_qTrue(); // Surface deformation - /* - const auto& lp_pair = surface_deformation(topol, tsos, recHitPix); - LocalPoint lp_track = lp_pair.first; + const auto& lp_pair = surface_deformation(topol, tsos, recHitPix); + + LocalPoint lp_track = lp_pair.first; LocalPoint lp_rechit = lp_pair.second; - const auto& rechitCorr_x = lp_rechit.x(); - const auto& rechitCorr_y = lp_rechit.y(); - const auto& trackhitCorrX_ = lp_track.x(); - const auto& trackhitCorrY_ = lp_track.y(); - */ + rechitCorr_.x = lp_rechit.x(); + rechitCorr_.y = lp_rechit.y(); + trackhitCorrX_ = lp_track.x(); + trackhitCorrY_ = lp_track.y(); + + if (notInPCL_) { + SiPixelLorentzAngleTreeBarrel_->Fill(); + } // is one pixel in cluster a large pixel ? (hit will be excluded) bool large_pix = false; @@ -343,53 +526,190 @@ void SiPixelLorentzAnglePCLWorker::dqmAnalyze(edm::Event const& iEvent, } } - double residual = TMath::Sqrt((trackhit_x - rechit_x) * (trackhit_x - rechit_x) + - (trackhit_y - rechit_y) * (trackhit_y - rechit_y)); + double residual = TMath::Sqrt((trackhitCorrX_ - rechitCorr_.x) * (trackhitCorrX_ - rechitCorr_.x) + + (trackhitCorrY_ - rechitCorr_.y) * (trackhitCorrY_ - rechitCorr_.y)); + + double xlim1 = trackhitCorrX_ - width_ * cotalpha / 2.; + double hypitch_ = ypitch_ / 2.; + double ylim1 = trackhitCorrY_ - width_ * cotbeta / 2.; + double ylim2 = trackhitCorrY_ + width_ * cotbeta / 2.; + + int clustSizeY_cut; + if (layer_ < 4) { + clustSizeY_cut = clustSizeYMin_; + } else { + clustSizeY_cut = clustSizeYMinL4_; + } - if (!large_pix && (chi2_ / ndof_) < normChi2Max_ && cluster->sizeY() >= clustSizeYMin_ && - residual < residualMax_ && (cluster->charge() < clustChargeMax_)) { + if (!large_pix && (chi2_ / ndof_) < normChi2Max_ && cluster->sizeY() >= clustSizeY_cut && + residual < residualMax_ && cluster->charge() < clusterCharge_cut && cluster->sizeX() < clustSizeXMax_) { // iterate over pixels in hit for (int j = 0; j < pixinfo_.npix; j++) { - // use trackhits - float dx = (pixinfo_.x[j] - (trackhit_x - width_ / 2. / TMath::Tan(trackhit_alpha))) * 10000.; - float dy = (pixinfo_.y[j] - (trackhit_y - width_ / 2. / TMath::Tan(trackhit_beta))) * 10000.; - float depth = dy * tan(trackhit_beta); - float drift = dx - dy * tan(trackhit_gamma); - - int i_index = module_ + (layer_ - 1) * iHists.nModules_[layer_ - 1]; - iHists.h_drift_depth_adc_.at(i_index)->Fill(drift, depth, pixinfo_.adc[j]); - iHists.h_drift_depth_adc2_.at(i_index)->Fill(drift, depth, pixinfo_.adc[j] * pixinfo_.adc[j]); - iHists.h_drift_depth_noadc_.at(i_index)->Fill(drift, depth); + // use trackhits and include bowing correction + float ypixlow = pixinfo_.y[j] - hypitch_; + float ypixhigh = pixinfo_.y[j] + hypitch_; + if (cotbeta > 0.) { + if (ylim1 > ypixlow) + ypixlow = ylim1; + if (ylim2 < ypixhigh) + ypixhigh = ylim2; + } else { + if (ylim2 > ypixlow) + ypixlow = ylim2; + if (ylim1 < ypixhigh) + ypixhigh = ylim1; + } + float ypixavg = 0.5 * (ypixlow + ypixhigh); + + float changeUnit = 10000.; + float dx = (pixinfo_.x[j] - xlim1) * changeUnit; // dx: in the unit of micrometer + float dy = (ypixavg - ylim1) * changeUnit; // dy: in the unit of micrometer + float depth = dy * tan(trackhit_.beta); + float drift = dx - dy * tan(trackhit_.gamma); + + if (Isnew == false) { + int i_index = module_ + (layer_ - 1) * iHists.nModules_[layer_ - 1]; + iHists.h_drift_depth_adc_.at(i_index)->Fill(drift, depth, pixinfo_.adc[j]); + iHists.h_drift_depth_adc2_.at(i_index)->Fill(drift, depth, pixinfo_.adc[j] * pixinfo_.adc[j]); + iHists.h_drift_depth_noadc_.at(i_index)->Fill(drift, depth); + } else { + int new_index = iHists.nModules_[iHists.nlay - 1] + + (iHists.nlay - 1) * iHists.nModules_[iHists.nlay - 1] + 1 + DetId_index; + iHists.h_drift_depth_adc_.at(new_index)->Fill(drift, depth, pixinfo_.adc[j]); + iHists.h_drift_depth_adc2_.at(new_index)->Fill(drift, depth, pixinfo_.adc[j] * pixinfo_.adc[j]); + iHists.h_drift_depth_noadc_.at(new_index)->Fill(drift, depth); + } } } } else if (subDetID == PixelSubdetector::PixelEndcap) { - // what to do here? - edm::LogInfo("SiPixelLorentzAnglePCLWorker") << "do Nothing here?" << std::endl; + DetId detIdObj = recHit->geographicalId(); + const PixelGeomDetUnit* theGeomDet = dynamic_cast(tracker->idToDet(detIdObj)); + if (!theGeomDet) + continue; + + const PixelTopology* topol = &(theGeomDet->specificTopology()); + + if (!topol) + continue; + + sideF_ = tTopo->pxfSide(detIdObj); + diskF_ = tTopo->pxfDisk(detIdObj); + bladeF_ = tTopo->pxfBlade(detIdObj); + panelF_ = tTopo->pxfPanel(detIdObj); + moduleF_ = tTopo->pxfModule(detIdObj); + + const SiPixelRecHit* recHitPix = dynamic_cast((*recHit).hit()); + if (!recHitPix) + continue; + rechitF_.x = recHitPix->localPosition().x(); + rechitF_.y = recHitPix->localPosition().y(); + SiPixelRecHit::ClusterRef const& cluster = recHitPix->cluster(); + + pixinfoF_ = fillPix(*cluster, topol); + + // fill entries in clust_ + + clustF_.x = (cluster)->x(); + clustF_.y = (cluster)->y(); + clustF_.charge = (cluster->charge()) / 1000.; // clustF_.charge: in the unit of 1000e + clustF_.size_x = cluster->sizeX(); + clustF_.size_y = cluster->sizeY(); + clustF_.maxPixelCol = cluster->maxPixelCol(); + clustF_.maxPixelRow = cluster->maxPixelRow(); + clustF_.minPixelCol = cluster->minPixelCol(); + clustF_.minPixelRow = cluster->minPixelRow(); + + // fill the trackhit info + TrajectoryStateOnSurface tsos = itTraj->updatedState(); + if (!tsos.isValid()) { + edm::LogWarning("SiPixelLorentzAnglePCLWorker") << "tsos not valid" << std::endl; + continue; + } + LocalVector trackdirection = tsos.localDirection(); + LocalPoint trackposition = tsos.localPosition(); + + if (trackdirection.z() == 0) + continue; + // the local position and direction + trackhitF_.alpha = atan2(trackdirection.z(), trackdirection.x()); + trackhitF_.beta = atan2(trackdirection.z(), trackdirection.y()); + trackhitF_.gamma = atan2(trackdirection.x(), trackdirection.y()); + trackhitF_.x = trackposition.x(); + trackhitF_.y = trackposition.y(); + + float cotalpha = trackdirection.x() / trackdirection.z(); + float cotbeta = trackdirection.y() / trackdirection.z(); + + float locBx = 1.; + if (cotbeta < 0.) + locBx = -1.; + float locBz = locBx; + if (cotalpha < 0.) + locBz = -locBx; + + auto detId = detIdObj.rawId(); + + int TemplID = templateDBobject_->getTemplateID(detId); + templ.interpolate(TemplID, cotalpha, cotbeta, locBz, locBx); + qScaleF_ = templ.qscale(); + rQmQtF_ = templ.r_qMeas_qTrue(); + + // Surface deformation + + const auto& lp_pair = surface_deformation(topol, tsos, recHitPix); + + LocalPoint lp_track = lp_pair.first; + LocalPoint lp_rechit = lp_pair.second; + + rechitCorrF_.x = lp_rechit.x(); + rechitCorrF_.y = lp_rechit.y(); + trackhitCorrXF_ = lp_track.x(); + trackhitCorrYF_ = lp_track.y(); + if (notInPCL_) { + SiPixelLorentzAngleTreeForward_->Fill(); + } } } //end iteration over trajectory measurements } //end iteration over trajectories } } -void SiPixelLorentzAnglePCLWorker::dqmBeginRun(edm::Run const& run, - edm::EventSetup const& iSetup, - SiPixelLorentzAngleCalibrationHistograms& iHists) const { +void SiPixelLorentzAnglePCLWorker::dqmBeginRun(edm::Run const& run, edm::EventSetup const& iSetup) { // geometry const TrackerGeometry* geom = &iSetup.getData(geomEsToken_); const TrackerTopology* tTopo = &iSetup.getData(topoEsToken_); PixelTopologyMap map = PixelTopologyMap(geom, tTopo); iHists.nlay = geom->numberOfLayers(PixelSubdetector::PixelBarrel); - + iHists.nModules_.resize(iHists.nlay); for (int i = 0; i < iHists.nlay; i++) { iHists.nModules_[i] = map.getPXBModules(i + 1); } + + if (!newmodulelist_.empty()) { + for (auto const& modulename : newmodulelist_) { + if (modulename.find("BPix_") != std::string::npos) { + PixelBarrelName bn(modulename, true); + const auto& detId = bn.getDetId(tTopo); + iHists.BPixnewmodulename_.push_back(modulename); + iHists.BPixnewDetIds_.push_back(detId.rawId()); + iHists.BPixnewModule_.push_back(bn.moduleName()); + iHists.BPixnewLayer_.push_back(bn.layerName()); + } else if (modulename.find("FPix_") != std::string::npos) { + PixelEndcapName en(modulename, true); + const auto& detId = en.getDetId(tTopo); + iHists.FPixnewmodulename_.push_back(modulename); + iHists.FPixnewDetIds_.push_back(detId.rawId()); + iHists.FPixnewDisk_.push_back(en.diskName()); + iHists.FPixnewBlade_.push_back(en.bladeName()); + } + } + } } void SiPixelLorentzAnglePCLWorker::bookHistograms(DQMStore::IBooker& iBooker, edm::Run const& run, - edm::EventSetup const& iSetup, - SiPixelLorentzAngleCalibrationHistograms& iHists) const { + edm::EventSetup const& iSetup) { iBooker.setCurrentFolder(folder_); iHists.h_tracks_ = iBooker.book1D("h_tracks", "h_tracks", 2, 0., 2.); @@ -419,11 +739,42 @@ void SiPixelLorentzAnglePCLWorker::bookHistograms(DQMStore::IBooker& iBooker, iHists.h_mean_[i_index] = iBooker.book1D(name, name, hist_depth_, min_depth_, max_depth_); } } + + for (int i = 0; i < (int)iHists.BPixnewDetIds_.size(); i++) { + int new_index = iHists.nModules_[iHists.nlay - 1] + (iHists.nlay - 1) * iHists.nModules_[iHists.nlay - 1] + 1 + i; + + sprintf(name, "h_BPixnew_drift_depth_adc_%s", iHists.BPixnewmodulename_[i].c_str()); + iHists.h_drift_depth_adc_[new_index] = + iBooker.book2D(name, name, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_); + + sprintf(name, "h_BPixnew_drift_depth_adc2_%s", iHists.BPixnewmodulename_[i].c_str()); + iHists.h_drift_depth_adc2_[new_index] = + iBooker.book2D(name, name, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_); + + sprintf(name, "h_BPixnew_drift_depth_noadc_%s", iHists.BPixnewmodulename_[i].c_str()); + iHists.h_drift_depth_noadc_[new_index] = + iBooker.book2D(name, name, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_); + + sprintf(name, "h_BPixnew_drift_depth_%s", iHists.BPixnewmodulename_[i].c_str()); + iHists.h_drift_depth_[new_index] = + iBooker.book2D(name, name, hist_drift_, min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_); + + sprintf(name, "h_BPixnew_mean_%s", iHists.BPixnewmodulename_[i].c_str()); + iHists.h_mean_[new_index] = iBooker.book1D(name, name, hist_depth_, min_depth_, max_depth_); + } +} + +void SiPixelLorentzAnglePCLWorker::dqmEndRun(edm::Run const& run, edm::EventSetup const& iSetup) { + if (notInPCL_) { + hFile_->cd(); + hFile_->Write(); + hFile_->Close(); + } } // method used to fill per pixel info const Pixinfo SiPixelLorentzAnglePCLWorker::fillPix(const SiPixelCluster& LocPix, const PixelTopology* topol) const { - static Pixinfo pixinfo; + Pixinfo pixinfo; const std::vector& pixvector = LocPix.pixels(); pixinfo.npix = 0; for (std::vector::const_iterator itPix = pixvector.begin(); itPix != pixvector.end(); @@ -453,22 +804,29 @@ const std::pair SiPixelLorentzAnglePCLWorker::surface_de LocalPoint lp_rechit = topol->localPosition(MeasurementPoint(pixels_rechit.first, pixels_rechit.second)); - static std::pair lps = std::make_pair(lp_track, lp_rechit); + std::pair lps = std::make_pair(lp_track, lp_rechit); return lps; } // ------------ method fills 'descriptions' with the allowed parameters for the module ------------ void SiPixelLorentzAnglePCLWorker::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { edm::ParameterSetDescription desc; - desc.add("folder", "AlCaReco/SiPixelLorentzAngle"); + desc.add("folder", "AlCaReco/SiPixelLorentzAngle"); // directory of PCL Worker output + desc.add("notInPCL", false); // create TTree (true) or not (false) + desc.add("fileName", "testrun.root"); // name of the TTree file if notInPCL = true + desc.add>( + "newmodulelist", + {"BPix_BmO_SEC3_LYR2_LDR5F_MOD2", "BPix_BpO_SEC1_LYR2_LDR1F_MOD1"}); // list of layer 2 new module names desc.add("src", edm::InputTag("TrackRefitter")); - desc.add("ptMin", 3.); - desc.add("normChi2Max", 2.); - desc.add("clustSizeYMin", 4); - desc.add("residualMax", 0.005); - desc.add("clustChargeMax", 120000); - desc.add("binsDepth", 50); - desc.add("binsDrift", 200); + desc.add("ptMin", 3.); // minimum pt on tracks + desc.add("normChi2Max", 2.); // maximum reduced chi squared + desc.add("clustSizeYMin", 4); // minimum cluster size on Y axis for Layer 1-3 + desc.add("clustSizeYMinL4", 3); // minimum cluster size on Y axis for Layer 4 + desc.add("clustSizeXMax", 5); // maximum cluster size on X axis + desc.add("residualMax", 0.005); // maximum residual + desc.add("clustChargeMaxPerLength", 50000); // maximum cluster charge per unit length of pixel depth (z) + desc.add("binsDepth", 50); // # bins for electron production depth axis + desc.add("binsDrift", 200); // # bins for electron drift axis descriptions.addWithDefaultLabel(desc); } diff --git a/Calibration/TkAlCaRecoProducers/test/testPCLAlCaHarvesting.py b/Calibration/TkAlCaRecoProducers/test/testPCLAlCaHarvesting.py index 25c4d6f6aff00..d72de7123bc39 100644 --- a/Calibration/TkAlCaRecoProducers/test/testPCLAlCaHarvesting.py +++ b/Calibration/TkAlCaRecoProducers/test/testPCLAlCaHarvesting.py @@ -81,6 +81,7 @@ def findRunStopTime(run_number): process.PoolDBOutputService.toPut.append(process.ALCAHARVESTSiStripGains_dbOutput) process.PoolDBOutputService.toPut.append(process.ALCAHARVESTSiStripGainsAAG_dbOutput ) process.PoolDBOutputService.toPut.append(process.ALCAHARVESTSiPixelAli_dbOutput) +process.PoolDBOutputService.toPut.append(process.ALCAHARVESTSiPixelLA_dbOutput) process.PoolDBOutputService.toPut.extend(process.ALCAHARVESTSiPixelQuality_dbOutput) process.PoolDBOutputService.toPut.append(process.ALCAHARVESTBeamSpotByRun_dbOutput) process.PoolDBOutputService.toPut.append(process.ALCAHARVESTBeamSpotByLumi_dbOutput) @@ -94,6 +95,7 @@ def findRunStopTime(run_number): process.pclMetadataWriter.recordsToMap.append(process.ALCAHARVESTSiStripGains_metadata ) process.pclMetadataWriter.recordsToMap.append(process.ALCAHARVESTSiStripGainsAAG_metadata) process.pclMetadataWriter.recordsToMap.append(process.ALCAHARVESTSiPixelAli_metadata) +process.pclMetadataWriter.recordsToMap.append(process.ALCAHARVESTSiPixelLA_metadata) process.pclMetadataWriter.recordsToMap.extend(process.ALCAHARVESTSiPixelQuality_metadata) process.pclMetadataWriter.recordsToMap.append(process.ALCAHARVESTBeamSpotByRun_metadata) process.pclMetadataWriter.recordsToMap.append(process.ALCAHARVESTBeamSpotByLumi_metadata) @@ -124,6 +126,8 @@ def findRunStopTime(run_number): process.SiPixelAliMilleFileExtractor.outputBinaryFile = cms.string('') process.SiPixelAliPedeAlignmentProducer.algoConfig.mergeBinaryFiles=[] +process.SiPixelLA = cms.Path(process.ALCAHARVESTSiPixelLorentzAngle) + process.SiPixelQuality = cms.Path(process.ALCAHARVESTSiPixelQuality) process.ALCAHARVESTDQMSaveAndMetadataWriter = cms.Path(process.dqmSaver+process.pclMetadataWriter) @@ -136,7 +140,8 @@ def findRunStopTime(run_number): process.schedule = cms.Schedule(process.SiStripQuality, process.SiStripGains, process.SiStripGainsAAG, - process.SiPixelAli, + process.SiPixelAli, + process.SiPixelLA, process.SiPixelQuality, process.BeamSpotByRun, process.BeamSpotByLumi, diff --git a/Configuration/AlCa/python/autoPCL.py b/Configuration/AlCa/python/autoPCL.py index e0763522ec0bb..0876ddef0fc4c 100644 --- a/Configuration/AlCa/python/autoPCL.py +++ b/Configuration/AlCa/python/autoPCL.py @@ -9,7 +9,6 @@ 'PromptCalibProdSiPixelLA' : 'SiPixelLA', 'PromptCalibProdEcalPedestals': 'EcalPedestals', 'PromptCalibProdLumiPCC': 'LumiPCC', - 'PromptCalibProdSiPixel' : 'SiPixelQuality', 'PromptCalibProdPPS' : 'PPSTimingCalibration', 'PromptCalibProdPPSDiamondSampic' : 'PPSDiamondSampicTimingCalibration' } From 1cd279b7126febb062ba4ace351c5ea66f1ae414 Mon Sep 17 00:00:00 2001 From: mmusich Date: Mon, 11 Oct 2021 16:03:12 +0200 Subject: [PATCH 05/12] clone instead of copy --- .../ALCARECOPromptCalibProdSiPixelLorentzAngle_Output_cff.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_Output_cff.py b/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_Output_cff.py index af2ec0ab2aa4b..d1894e8c425b2 100644 --- a/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_Output_cff.py +++ b/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_Output_cff.py @@ -9,7 +9,5 @@ 'keep *_MEtoEDMConvertSiPixelLorentzAngle_*_*', ) ) - -import copy -OutALCARECOPromptCalibProdSiPixelLorentzAngle=copy.deepcopy(OutALCARECOPromptCalibProdSiPixelLorentzAngle_noDrop) +OutALCARECOPromptCalibProdSiPixelLorentzAngle=OutALCARECOPromptCalibProdSiPixelLorentzAngle_noDrop.clone() OutALCARECOPromptCalibProdSiPixelLorentzAngle.outputCommands.insert(0, "drop *") From 5848729073f36a7b0e7a1bc3682361c88d1771a7 Mon Sep 17 00:00:00 2001 From: mmusich Date: Mon, 11 Oct 2021 16:03:32 +0200 Subject: [PATCH 06/12] update parseFwkJobReport to include SiPixelLA --- Calibration/TkAlCaRecoProducers/test/parseFwkJobReport.py | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/Calibration/TkAlCaRecoProducers/test/parseFwkJobReport.py b/Calibration/TkAlCaRecoProducers/test/parseFwkJobReport.py index 6d89d47f99823..6b1b78905d609 100644 --- a/Calibration/TkAlCaRecoProducers/test/parseFwkJobReport.py +++ b/Calibration/TkAlCaRecoProducers/test/parseFwkJobReport.py @@ -4,12 +4,11 @@ ## declare all constants here TARGET_LIST_OF_TAGS=['BeamSpotObject_ByLumi', 'BeamSpotObject_ByRun', 'BeamSpotObjectHP_ByLumi', 'BeamSpotObjectHP_ByRun', - 'SiPixelQualityFromDbRcd_other', 'SiPixelQualityFromDbRcd_prompt', 'SiPixelQualityFromDbRcd_stuckTBM', - 'SiStripApvGain_pcl', 'SiStripApvGainAAG_pcl', - 'SiStripBadStrip_pcl', 'SiPixelAli_pcl'] + 'SiPixelLA_pcl', 'SiPixelQualityFromDbRcd_other', 'SiPixelQualityFromDbRcd_prompt', 'SiPixelQualityFromDbRcd_stuckTBM', + 'SiStripApvGain_pcl', 'SiStripApvGainAAG_pcl', 'SiStripBadStrip_pcl', 'SiPixelAli_pcl'] TARGET_DQM_FILES=1 TARGET_DQM_FILENAME='./DQM_V0001_R000325022__Express__PCLTest__ALCAPROMPT.root' -TARGET_DB_FILES=11 +TARGET_DB_FILES=12 TARGET_DB_FILENAME='sqlite_file:promptCalibConditions.db' TOTAL_TARGET_FILES=TARGET_DQM_FILES+TARGET_DB_FILES From 6459eca2618c11ae8527c290ae1a7701cdc25428 Mon Sep 17 00:00:00 2001 From: mmusich Date: Mon, 11 Oct 2021 16:04:51 +0200 Subject: [PATCH 07/12] modify writing logic post PR #35556 --- .../src/SiPixelLorentzAnglePCLHarvester.cc | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc index 97c372bf6c293..4a8c0b6c71077 100644 --- a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc +++ b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc @@ -29,6 +29,7 @@ class SiPixelLorentzAnglePCLHarvester : public DQMEDHarvester { public: SiPixelLorentzAnglePCLHarvester(const edm::ParameterSet&); + ~SiPixelLorentzAnglePCLHarvester() override = default; void beginRun(const edm::Run&, const edm::EventSetup&) override; static void fillDescriptions(edm::ConfigurationDescriptions&); @@ -249,7 +250,7 @@ void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMS << "\t" << "newDetId" << std::endl; - SiPixelLorentzAngle* LorentzAngle = new SiPixelLorentzAngle(); + std::unique_ptr LorentzAngle = std::make_unique(); f1 = std::make_unique("f1", "[0] + [1]*x + [2]*x*x + [3]*x*x*x + [4]*x*x*x*x + [5]*x*x*x*x*x", 5., 280.); f1->SetParName(0, "offset"); @@ -444,7 +445,7 @@ void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMS edm::Service mydbservice; if (mydbservice.isAvailable()) { try { - mydbservice->writeOne(LorentzAngle, mydbservice->currentTime(), recordName_); + mydbservice->writeOneIOV(LorentzAngle.get(), mydbservice->currentTime(), recordName_); } catch (const cond::Exception& er) { edm::LogError("SiPixelLorentzAngleDB") << er.what() << std::endl; } catch (const std::exception& er) { @@ -453,7 +454,6 @@ void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMS } else { edm::LogError("SiPixelLorentzAngleDB") << "Service is unavailable" << std::endl; } - delete LorentzAngle; } void SiPixelLorentzAnglePCLHarvester::findMean(MonitorElement* h_drift_depth_adc_slice_, int i, int i_ring) { From 708afb2bf5219b41993fea2abb907d295e1508f7 Mon Sep 17 00:00:00 2001 From: mmusich Date: Mon, 11 Oct 2021 19:56:54 +0200 Subject: [PATCH 08/12] more clean-up --- .../SiPixelLorentzAnglePCLHarvester_cfi.py | 13 +-- .../SiPixelLorentzAnglePCLWorker_cfi.py | 17 ++-- .../src/SiPixelLorentzAnglePCLHarvester.cc | 91 ++++++------------ .../src/SiPixelLorentzAnglePCLWorker.cc | 92 ++++++++----------- 4 files changed, 80 insertions(+), 133 deletions(-) diff --git a/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py index 92766f0917920..2e38b6cda3282 100644 --- a/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py +++ b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py @@ -1,18 +1,11 @@ import FWCore.ParameterSet.Config as cms from DQMServices.Core.DQMEDHarvester import DQMEDHarvester +from CalibTracker.SiPixelLorentzAngle.SiPixelLorentzAnglePCLWorker_cfi import SiPixelLorentzAnglePCLWorker as worker SiPixelLorentzAnglePCLHarvester = DQMEDHarvester( "SiPixelLorentzAnglePCLHarvester", - newmodulelist = cms.vstring( - "BPix_BmI_SEC7_LYR2_LDR12F_MOD1", - "BPix_BmI_SEC8_LYR2_LDR14F_MOD1", - "BPix_BmO_SEC3_LYR2_LDR5F_MOD1", - "BPix_BmO_SEC3_LYR2_LDR5F_MOD2", - "BPix_BmO_SEC3_LYR2_LDR5F_MOD3", - "BPix_BpO_SEC1_LYR2_LDR1F_MOD1", - "BPix_BpO_SEC1_LYR2_LDR1F_MOD2", - "BPix_BpO_SEC1_LYR2_LDR1F_MOD3"), - dqmDir = cms.string('AlCaReco/SiPixelLorentzAngle'), + newmodulelist = cms.vstring(worker.newmodulelist.value()), # taken from worker configuration, need to stay in synch + dqmDir = cms.string(worker.folder.value()), # taken from worker configuration, need to stay in synch record = cms.string("SiPixelLorentzAngleRcd"), fitProbCut = cms.double(0.5) ) diff --git a/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py index 0ae9692e6b247..8e6bc0611eaac 100644 --- a/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py +++ b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py @@ -6,15 +6,14 @@ folder = cms.string('AlCaReco/SiPixelLorentzAngle'), notInPCL = cms.bool(False), fileName = cms.string('testrun.root'), - newmodulelist = cms.vstring( - "BPix_BmI_SEC7_LYR2_LDR12F_MOD1", - "BPix_BmI_SEC8_LYR2_LDR14F_MOD1", - "BPix_BmO_SEC3_LYR2_LDR5F_MOD1", - "BPix_BmO_SEC3_LYR2_LDR5F_MOD2", - "BPix_BmO_SEC3_LYR2_LDR5F_MOD3", - "BPix_BpO_SEC1_LYR2_LDR1F_MOD1", - "BPix_BpO_SEC1_LYR2_LDR1F_MOD2", - "BPix_BpO_SEC1_LYR2_LDR1F_MOD3"), + newmodulelist = cms.vstring("BPix_BmI_SEC7_LYR2_LDR12F_MOD1", + "BPix_BmI_SEC8_LYR2_LDR14F_MOD1", + "BPix_BmO_SEC3_LYR2_LDR5F_MOD1", + "BPix_BmO_SEC3_LYR2_LDR5F_MOD2", + "BPix_BmO_SEC3_LYR2_LDR5F_MOD3", + "BPix_BpO_SEC1_LYR2_LDR1F_MOD1", + "BPix_BpO_SEC1_LYR2_LDR1F_MOD2", + "BPix_BpO_SEC1_LYR2_LDR1F_MOD3"), src = cms.InputTag("TrackRefitter"), binsDepth = cms.int32(50), binsDrift = cms.int32(200), diff --git a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc index 4a8c0b6c71077..175cd60e398df 100644 --- a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc +++ b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc @@ -1,31 +1,26 @@ -#include "DQMServices/Core/interface/DQMEDHarvester.h" +#include +#include +#include +#include "CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h" +#include "CondCore/DBOutputService/interface/PoolDBOutputService.h" +#include "CondFormats/DataRecord/interface/SiPixelLorentzAngleRcd.h" +#include "CondFormats/SiPixelObjects/interface/SiPixelLorentzAngle.h" +#include "DQMServices/Core/interface/DQMEDHarvester.h" +#include "DataFormats/TrackerCommon/interface/PixelBarrelName.h" +#include "DataFormats/TrackerCommon/interface/PixelEndcapName.h" +#include "FWCore/Framework/interface/EventSetup.h" #include "FWCore/Framework/interface/Frameworkfwd.h" #include "FWCore/Framework/interface/MakerMacros.h" -#include "FWCore/Framework/interface/ESHandle.h" -#include "FWCore/Framework/interface/EventSetup.h" #include "FWCore/ParameterSet/interface/ParameterSet.h" #include "FWCore/ServiceRegistry/interface/Service.h" -#include "DataFormats/TrackerCommon/interface/PixelBarrelName.h" -#include "DataFormats/TrackerCommon/interface/PixelEndcapName.h" +#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" #include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyMap.h" #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" - -#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" -#include "CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h" -#include "CondCore/DBOutputService/interface/PoolDBOutputService.h" -#include "CondFormats/SiPixelObjects/interface/SiPixelLorentzAngle.h" -#include "CondFormats/DataRecord/interface/SiPixelLorentzAngleRcd.h" - #include "MagneticField/Engine/interface/MagneticField.h" #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h" -#include -#include -#include - //------------------------------------------------------------------------------ - class SiPixelLorentzAnglePCLHarvester : public DQMEDHarvester { public: SiPixelLorentzAnglePCLHarvester(const edm::ParameterSet&); @@ -56,7 +51,6 @@ class SiPixelLorentzAnglePCLHarvester : public DQMEDHarvester { }; //------------------------------------------------------------------------------ - SiPixelLorentzAnglePCLHarvester::SiPixelLorentzAnglePCLHarvester(const edm::ParameterSet& iConfig) : geomEsToken_(esConsumes()), topoEsToken_(esConsumes()), @@ -73,7 +67,6 @@ SiPixelLorentzAnglePCLHarvester::SiPixelLorentzAnglePCLHarvester(const edm::Para } //------------------------------------------------------------------------------ - void SiPixelLorentzAnglePCLHarvester::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) { // geometry const TrackerGeometry* geom = &iSetup.getData(geomEsToken_); @@ -129,7 +122,6 @@ void SiPixelLorentzAnglePCLHarvester::beginRun(const edm::Run& iRun, const edm:: } //------------------------------------------------------------------------------ - void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter) { // go in the right directory iGetter.cd(); @@ -212,43 +204,18 @@ void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMS MonitorElement* h_drift_depth_adc_slice_ = iBooker.book1D("h_drift_depth_adc_slice", "slice of adc histogram", hist_drift_, min_drift_, max_drift_); - edm::LogPrint("LorentzAngle") << "module" - << "\t" - << "layer" - << "\t" - << "offset" - << "\t" - << "e0" - << "\t" - << "slope" - << "\t" - << "e1" - << "\t" - "rel.err" - << "\t" - "pull" - << "\t" - << "p2" - << "\t" - << "e2" - << "\t" - << "p3" - << "\t" - << "e3" - << "\t" - << "p4" - << "\t" - << "e4" - << "\t" - << "p5" - << "\t" - << "e5" - << "\t" - << "chi2" - << "\t" - << "prob" - << "\t" + // clang-format off + edm::LogPrint("LorentzAngle") << "module" << "\t" << "layer" << "\t" + << "offset" << "\t" << "e0" << "\t" + << "slope" << "\t" << "e1" << "\t" + << "rel.err" << "\t" << "pull" << "\t" + << "p2" << "\t" << "e2" << "\t" + << "p3" << "\t" << "e3" << "\t" + << "p4" << "\t" << "e4" << "\t" + << "p5" << "\t" << "e5" << "\t" + << "chi2" << "\t" << "prob" << "\t" << "newDetId" << std::endl; + // clang-format on std::unique_ptr LorentzAngle = std::make_unique(); @@ -456,6 +423,7 @@ void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMS } } +//------------------------------------------------------------------------------ void SiPixelLorentzAnglePCLHarvester::findMean(MonitorElement* h_drift_depth_adc_slice_, int i, int i_ring) { double nentries = 0; h_drift_depth_adc_slice_->Reset(); @@ -495,12 +463,11 @@ void SiPixelLorentzAnglePCLHarvester::findMean(MonitorElement* h_drift_depth_adc //------------------------------------------------------------------------------ void SiPixelLorentzAnglePCLHarvester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { edm::ParameterSetDescription desc; - desc.add>( - "newmodulelist", - {"BPix_BmO_SEC3_LYR2_LDR5F_MOD2", "BPix_BpO_SEC1_LYR2_LDR1F_MOD1"}); // the layer 2 new module names - desc.add("dqmDir", "AlCaReco/SiPixelLorentzAngle"); // the directory of PCL Worker output - desc.add("fitProbCut", 0.5); - desc.add("record", "SiPixelLorentzAngleRcd"); + desc.setComment("Harvester module of the SiPixel Lorentz Angle PCL monitoring workflow"); + desc.add>("newmodulelist", {})->setComment("the list of DetIds for new sensors"); + desc.add("dqmDir", "AlCaReco/SiPixelLorentzAngle")->setComment("the directory of PCL Worker output"); + desc.add("fitProbCut", 0.5)->setComment("cut on fit chi2 probabiblity to accept measurement"); + desc.add("record", "SiPixelLorentzAngleRcd")->setComment("target DB record"); descriptions.addWithDefaultLabel(desc); } diff --git a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc index c4c8da7ccb9c5..d350e94d784b4 100644 --- a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc +++ b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc @@ -17,44 +17,40 @@ #include // user include files -#include "DQMServices/Core/interface/DQMStore.h" -#include "FWCore/Framework/interface/Frameworkfwd.h" +#include "CalibTracker/Records/interface/SiPixelTemplateDBObjectESProducerRcd.h" +#include "CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h" +#include "CondFormats/SiPixelObjects/interface/SiPixelTemplateDBObject.h" +#include "CondFormats/SiPixelTransient/interface/SiPixelTemplate.h" +#include "CondFormats/SiPixelTransient/interface/SiPixelTemplateDefs.h" #include "DQMServices/Core/interface/DQMEDAnalyzer.h" - -#include "FWCore/Framework/interface/Event.h" -#include "FWCore/Framework/interface/MakerMacros.h" - -#include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h" -#include "FWCore/Framework/interface/ESHandle.h" +#include "DQMServices/Core/interface/DQMStore.h" #include "DataFormats/GeometryVector/interface/GlobalVector.h" #include "DataFormats/GeometryVector/interface/LocalVector.h" -#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" -#include "Geometry/CommonDetUnit/interface/GeomDetType.h" -#include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h" -#include "DataFormats/TrackReco/interface/TrackExtra.h" #include "DataFormats/SiStripDetId/interface/StripSubdetector.h" -#include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h" -#include "Geometry/CommonTopologies/interface/StripTopology.h" -#include "DataFormats/TrackerCommon/interface/TrackerTopology.h" -#include "Geometry/Records/interface/IdealGeometryRecord.h" +#include "DataFormats/TrackReco/interface/TrackExtra.h" #include "DataFormats/TrackReco/interface/TrackFwd.h" #include "DataFormats/TrackerCommon/interface/PixelBarrelName.h" #include "DataFormats/TrackerCommon/interface/PixelEndcapName.h" +#include "DataFormats/TrackerCommon/interface/TrackerTopology.h" +#include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h" +#include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h" +#include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h" +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/Framework/interface/Frameworkfwd.h" +#include "FWCore/Framework/interface/MakerMacros.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "Geometry/CommonDetUnit/interface/GeomDetType.h" +#include "Geometry/CommonTopologies/interface/StripTopology.h" +#include "Geometry/Records/interface/IdealGeometryRecord.h" +#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" +#include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyMap.h" +#include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" +#include "Geometry/TrackerGeometryBuilder/interface/phase1PixelTopology.h" +#include "RecoTracker/TransientTrackingRecHit/interface/TkTransientTrackingRecHitBuilder.h" +#include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h" #include "TrackingTools/Records/interface/TransientRecHitRecord.h" #include "TrackingTools/TrackFitters/interface/TrajectoryStateCombiner.h" -#include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h" -#include "RecoTracker/TransientTrackingRecHit/interface/TkTransientTrackingRecHitBuilder.h" #include "TrackingTools/TransientTrack/interface/TransientTrack.h" -#include "CalibTracker/Records/interface/SiPixelTemplateDBObjectESProducerRcd.h" -#include "CondFormats/SiPixelObjects/interface/SiPixelTemplateDBObject.h" -#include "CondFormats/SiPixelTransient/interface/SiPixelTemplateDefs.h" -#include "CondFormats/SiPixelTransient/interface/SiPixelTemplate.h" -#include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyMap.h" -#include "Geometry/TrackerGeometryBuilder/interface/phase1PixelTopology.h" -#include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" - -#include "FWCore/ParameterSet/interface/ParameterSet.h" -#include "CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h" #include #include #include @@ -204,14 +200,6 @@ class SiPixelLorentzAnglePCLWorker : public DQMEDAnalyzer { edm::EDGetTokenT t_trajTrack; }; -// -// constants, enums and typedefs -// - -// -// static data member definitions -// - // // constructors and destructor // @@ -811,22 +799,22 @@ const std::pair SiPixelLorentzAnglePCLWorker::surface_de // ------------ method fills 'descriptions' with the allowed parameters for the module ------------ void SiPixelLorentzAnglePCLWorker::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { edm::ParameterSetDescription desc; - desc.add("folder", "AlCaReco/SiPixelLorentzAngle"); // directory of PCL Worker output - desc.add("notInPCL", false); // create TTree (true) or not (false) - desc.add("fileName", "testrun.root"); // name of the TTree file if notInPCL = true - desc.add>( - "newmodulelist", - {"BPix_BmO_SEC3_LYR2_LDR5F_MOD2", "BPix_BpO_SEC1_LYR2_LDR1F_MOD1"}); // list of layer 2 new module names - desc.add("src", edm::InputTag("TrackRefitter")); - desc.add("ptMin", 3.); // minimum pt on tracks - desc.add("normChi2Max", 2.); // maximum reduced chi squared - desc.add("clustSizeYMin", 4); // minimum cluster size on Y axis for Layer 1-3 - desc.add("clustSizeYMinL4", 3); // minimum cluster size on Y axis for Layer 4 - desc.add("clustSizeXMax", 5); // maximum cluster size on X axis - desc.add("residualMax", 0.005); // maximum residual - desc.add("clustChargeMaxPerLength", 50000); // maximum cluster charge per unit length of pixel depth (z) - desc.add("binsDepth", 50); // # bins for electron production depth axis - desc.add("binsDrift", 200); // # bins for electron drift axis + desc.setComment("Worker module of the SiPixel Lorentz Angle PCL monitoring workflow"); + desc.add("folder", "AlCaReco/SiPixelLorentzAngle")->setComment("directory of PCL Worker output"); + desc.add("notInPCL", false)->setComment("create TTree (true) or not (false)"); + desc.add("fileName", "testrun.root")->setComment("name of the TTree file if notInPCL = true"); + desc.add>("newmodulelist", {})->setComment("the list of DetIds for new sensors"); + desc.add("src", edm::InputTag("TrackRefitter"))->setComment("input track collections"); + desc.add("ptMin", 3.)->setComment("minimum pt on tracks"); + desc.add("normChi2Max", 2.)->setComment("maximum reduced chi squared"); + desc.add("clustSizeYMin", 4)->setComment("minimum cluster size on Y axis for Layer 1-3"); + desc.add("clustSizeYMinL4", 3)->setComment("minimum cluster size on Y axis for Layer 4"); + desc.add("clustSizeXMax", 5)->setComment("maximum cluster size on X axis"); + desc.add("residualMax", 0.005)->setComment("maximum residual"); + desc.add("clustChargeMaxPerLength", 50000) + ->setComment("maximum cluster charge per unit length of pixel depth (z)"); + desc.add("binsDepth", 50)->setComment("# bins for electron production depth axis"); + desc.add("binsDrift", 200)->setComment("# bins for electron drift axis"); descriptions.addWithDefaultLabel(desc); } From 60c65a393b689cbf39f67de1bbe2aa8016ad2da2 Mon Sep 17 00:00:00 2001 From: wweiphy <67766900+wweiphy@users.noreply.github.com> Date: Wed, 13 Oct 2021 12:03:04 -0700 Subject: [PATCH 09/12] Update SiPixelLorentzAnglePCLHarvester.cc --- .../src/SiPixelLorentzAnglePCLHarvester.cc | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc index 175cd60e398df..a4ee5adcf7cfc 100644 --- a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc +++ b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc @@ -281,15 +281,15 @@ void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMS for (int i_module = 1; i_module <= hists.nModules_[i_layer - 1]; i_module++) { if (i_layer == 1) p1_simul[i_layer - 1][i_module - 1] = 0.436848; - if (i_layer == 2) + else if (i_layer == 2) p1_simul[i_layer - 1][i_module - 1] = 0.25802; - if (i_layer == 3 && i_module <= 4) + else if (i_layer == 3 && i_module <= 4) p1_simul[i_layer - 1][i_module - 1] = 0.29374; - if (i_layer == 3 && i_module >= 5) + else if (i_layer == 3 && i_module >= 5) p1_simul[i_layer - 1][i_module - 1] = 0.31084; - if (i_layer == 4 && i_module <= 4) + else if (i_layer == 4 && i_module <= 4) p1_simul[i_layer - 1][i_module - 1] = 0.29944; - if (i_layer == 4 && i_module >= 5) + else p1_simul[i_layer - 1][i_module - 1] = 0.31426; } } From b96cd38a577f4bf1156e57ce718272843a3ef2a3 Mon Sep 17 00:00:00 2001 From: wweiphy <67766900+wweiphy@users.noreply.github.com> Date: Wed, 13 Oct 2021 12:05:25 -0700 Subject: [PATCH 10/12] Update ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py --- .../ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py b/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py index 2c723527933a3..d7de69d618d6d 100644 --- a/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py +++ b/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_cff.py @@ -29,9 +29,9 @@ from RecoTracker.TrackProducer.TrackRefitter_cfi import * from TrackingTools.TransientTrack.TransientTrackBuilder_cfi import * -ALCARECOPixelLACalibrationTracksRefit = TrackRefitter.clone(src = cms.InputTag("ALCARECOPixelLACalibrationTracks"), - TrajectoryInEvent = cms.bool(True), - NavigationSchool = cms.string("") +ALCARECOPixelLACalibrationTracksRefit = TrackRefitter.clone(src = "ALCARECOPixelLACalibrationTracks", + TrajectoryInEvent = True, + NavigationSchool = "" ) # refit and BS can be dropped if done together with RECO. From 0b4260661ad7650d4ebd99b0c5c94ef67ac6cfc7 Mon Sep 17 00:00:00 2001 From: wweiphy <67766900+wweiphy@users.noreply.github.com> Date: Wed, 13 Oct 2021 12:06:39 -0700 Subject: [PATCH 11/12] Update AlcaSiPixelLorentzAngleHarvester_cff.py --- .../python/AlcaSiPixelLorentzAngleHarvester_cff.py | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cff.py b/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cff.py index d1920180bfb91..65e12849ff48a 100644 --- a/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cff.py +++ b/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cff.py @@ -3,9 +3,10 @@ from Calibration.TkAlCaRecoProducers.AlcaSiPixelLorentzAngleHarvester_cfi import * from DQMServices.Components.EDMtoMEConverter_cfi import * -EDMtoMEConvertSiPixelLorentzAngle = EDMtoMEConverter.clone() -EDMtoMEConvertSiPixelLorentzAngle.lumiInputTag = cms.InputTag("MEtoEDMConvertSiPixelLorentzAngle","MEtoEDMConverterLumi") -EDMtoMEConvertSiPixelLorentzAngle.runInputTag = cms.InputTag("MEtoEDMConvertSiPixelLorentzAngle","MEtoEDMConverterRun") +EDMtoMEConvertSiPixelLorentzAngle = EDMtoMEConverter.clone( + lumiInputTag = ("MEtoEDMConvertSiPixelLorentzAngle","MEtoEDMConverterLumi"), + runInputTag = ("MEtoEDMConvertSiPixelLorentzAngle","MEtoEDMConverterRun") +) DQMStore = cms.Service("DQMStore") From a9c114b004376c03aa626daf6e3adf3d988add5d Mon Sep 17 00:00:00 2001 From: wweiphy <67766900+wweiphy@users.noreply.github.com> Date: Wed, 13 Oct 2021 12:08:52 -0700 Subject: [PATCH 12/12] Remove unnecessary class members --- .../src/SiPixelLorentzAnglePCLWorker.cc | 48 +++++-------------- 1 file changed, 12 insertions(+), 36 deletions(-) diff --git a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc index d350e94d784b4..9e2c352e21333 100644 --- a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc +++ b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc @@ -160,20 +160,6 @@ class SiPixelLorentzAnglePCLWorker : public DQMEDAnalyzer { float qScaleF_; float rQmQtF_; - // histogram etc - int hist_x_; - int hist_y_; - double min_x_; - double max_x_; - double min_y_; - double max_y_; - double width_; - double min_depth_; - double max_depth_; - double min_drift_; - double max_drift_; - int bufsize; - // parameters from config file double ptmin_; double normChi2Max_; @@ -225,19 +211,7 @@ SiPixelLorentzAnglePCLWorker::SiPixelLorentzAnglePCLWorker(const edm::ParameterS t_trajTrack = consumes(iConfig.getParameter("src")); // now do what ever initialization is needed - hist_x_ = 50; - hist_y_ = 100; - min_x_ = -500.; - max_x_ = 500.; - min_y_ = -1500.; - max_y_ = 500.; - width_ = 0.0285; - min_depth_ = -100.; - max_depth_ = 400.; - min_drift_ = -1000.; //-200.;(iConfig.getParameter("residualMax")) - max_drift_ = 1000.; //400.; - - bufsize = 64000; + int bufsize = 64000; // create tree structure // Barrel pixel @@ -405,6 +379,7 @@ void SiPixelLorentzAnglePCLWorker::analyze(edm::Event const& iEvent, edm::EventS const PixelTopology* topol = &(theGeomDet->specificTopology()); float ypitch_ = topol->pitch().second; + float width_ = 0.0285; if (!topol) continue; @@ -514,23 +489,19 @@ void SiPixelLorentzAnglePCLWorker::analyze(edm::Event const& iEvent, edm::EventS } } - double residual = TMath::Sqrt((trackhitCorrX_ - rechitCorr_.x) * (trackhitCorrX_ - rechitCorr_.x) + - (trackhitCorrY_ - rechitCorr_.y) * (trackhitCorrY_ - rechitCorr_.y)); + double residualsq = (trackhitCorrX_ - rechitCorr_.x) * (trackhitCorrX_ - rechitCorr_.x) + + (trackhitCorrY_ - rechitCorr_.y) * (trackhitCorrY_ - rechitCorr_.y); double xlim1 = trackhitCorrX_ - width_ * cotalpha / 2.; double hypitch_ = ypitch_ / 2.; double ylim1 = trackhitCorrY_ - width_ * cotbeta / 2.; double ylim2 = trackhitCorrY_ + width_ * cotbeta / 2.; - int clustSizeY_cut; - if (layer_ < 4) { - clustSizeY_cut = clustSizeYMin_; - } else { - clustSizeY_cut = clustSizeYMinL4_; - } + int clustSizeY_cut = layer_ < 4 ? clustSizeYMin_ : clustSizeYMinL4_; if (!large_pix && (chi2_ / ndof_) < normChi2Max_ && cluster->sizeY() >= clustSizeY_cut && - residual < residualMax_ && cluster->charge() < clusterCharge_cut && cluster->sizeX() < clustSizeXMax_) { + residualsq < residualMax_ * residualMax_ && cluster->charge() < clusterCharge_cut && + cluster->sizeX() < clustSizeXMax_) { // iterate over pixels in hit for (int j = 0; j < pixinfo_.npix; j++) { // use trackhits and include bowing correction @@ -701,6 +672,11 @@ void SiPixelLorentzAnglePCLWorker::bookHistograms(DQMStore::IBooker& iBooker, iBooker.setCurrentFolder(folder_); iHists.h_tracks_ = iBooker.book1D("h_tracks", "h_tracks", 2, 0., 2.); + double min_depth_ = -100.; + double max_depth_ = 400.; + double min_drift_ = -1000.; + double max_drift_ = 1000.; + //book histograms char name[128]; for (int i_layer = 1; i_layer <= iHists.nlay; i_layer++) {