diff --git a/CalibTracker/SiPixelLorentzAngle/BuildFile.xml b/CalibTracker/SiPixelLorentzAngle/BuildFile.xml index 28ad08af08f7b..d01ba32f1258e 100644 --- a/CalibTracker/SiPixelLorentzAngle/BuildFile.xml +++ b/CalibTracker/SiPixelLorentzAngle/BuildFile.xml @@ -8,5 +8,9 @@ + + + + diff --git a/CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h b/CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h new file mode 100644 index 0000000000000..9f9be5aca4604 --- /dev/null +++ b/CalibTracker/SiPixelLorentzAngle/interface/SiPixelLorentzAngleCalibrationStruct.h @@ -0,0 +1,35 @@ +#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; + + int nlay; + 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_; + MonitorMap h_drift_depth_noadc_; + MonitorMap h_drift_depth_; + MonitorMap h_mean_; + + dqm::reco::MonitorElement* h_tracks_; +}; + +#endif diff --git a/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py new file mode 100644 index 0000000000000..2e38b6cda3282 --- /dev/null +++ b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLHarvester_cfi.py @@ -0,0 +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(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 new file mode 100644 index 0000000000000..8e6bc0611eaac --- /dev/null +++ b/CalibTracker/SiPixelLorentzAngle/python/SiPixelLorentzAnglePCLWorker_cfi.py @@ -0,0 +1,27 @@ +import FWCore.ParameterSet.Config as cms + +from DQMServices.Core.DQMEDAnalyzer import DQMEDAnalyzer +SiPixelLorentzAnglePCLWorker = DQMEDAnalyzer( + "SiPixelLorentzAnglePCLWorker", + 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"), + 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), + clustChargeMaxPerLength = cms.double(50000) +) diff --git a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc new file mode 100644 index 0000000000000..a4ee5adcf7cfc --- /dev/null +++ b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLHarvester.cc @@ -0,0 +1,474 @@ +#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/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ServiceRegistry/interface/Service.h" +#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" +#include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyMap.h" +#include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" +#include "MagneticField/Engine/interface/MagneticField.h" +#include "MagneticField/Records/interface/IdealMagneticFieldRecord.h" + +//------------------------------------------------------------------------------ +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&); + +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_; + 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; +}; + +//------------------------------------------------------------------------------ +SiPixelLorentzAnglePCLHarvester::SiPixelLorentzAnglePCLHarvester(const edm::ParameterSet& iConfig) + : geomEsToken_(esConsumes()), + topoEsToken_(esConsumes()), + 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()) + 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_); + + 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); + } + } +} + +//------------------------------------------------------------------------------ +void SiPixelLorentzAnglePCLHarvester::dqmEndJob(DQMStore::IBooker& iBooker, DQMStore::IGetter& iGetter) { + // go in the right directory + iGetter.cd(); + iGetter.setCurrentFolder(dqmDir_); + + 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)); + + 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)); + + 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 (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_; + 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_); + + // 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(); + + 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; + else if (i_layer == 2) + p1_simul[i_layer - 1][i_module - 1] = 0.25802; + else if (i_layer == 3 && i_module <= 4) + p1_simul[i_layer - 1][i_module - 1] = 0.29374; + else if (i_layer == 3 && i_module >= 5) + p1_simul[i_layer - 1][i_module - 1] = 0.31084; + else if (i_layer == 4 && i_module <= 4) + p1_simul[i_layer - 1][i_module - 1] = 0.29944; + else + 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++) { + findMean(h_drift_depth_adc_slice_, i, i_index); + } // end loop over bins in depth + + 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(); + + 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 + edm::Service mydbservice; + if (mydbservice.isAvailable()) { + try { + mydbservice->writeOneIOV(LorentzAngle.get(), 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; + } +} + +//------------------------------------------------------------------------------ +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) { + 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.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); +} + +DEFINE_FWK_MODULE(SiPixelLorentzAnglePCLHarvester); diff --git a/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc new file mode 100644 index 0000000000000..9e2c352e21333 --- /dev/null +++ b/CalibTracker/SiPixelLorentzAngle/src/SiPixelLorentzAnglePCLWorker.cc @@ -0,0 +1,798 @@ +// -*- 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 "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 "DQMServices/Core/interface/DQMStore.h" +#include "DataFormats/GeometryVector/interface/GlobalVector.h" +#include "DataFormats/GeometryVector/interface/LocalVector.h" +#include "DataFormats/SiStripDetId/interface/StripSubdetector.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/TransientTrack/interface/TransientTrack.h" +#include +#include +#include +// +// 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]; +}; + +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 = default; + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +private: + void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override; + + void analyze(edm::Event const&, edm::EventSetup const&) override; + + void dqmBeginRun(edm::Run const&, edm::EventSetup 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_; + + // parameters from config file + double ptmin_; + double normChi2Max_; + int clustSizeYMin_; + int clustSizeYMinL4_; + int clustSizeXMax_; + double residualMax_; + double clustChargeMaxPerLength_; + int hist_depth_; + int hist_drift_; + + std::unique_ptr hFile_; + std::unique_ptr SiPixelLorentzAngleTreeBarrel_; + std::unique_ptr SiPixelLorentzAngleTreeForward_; + + // es consumes + edm::ESGetToken geomEsToken_; + edm::ESGetToken topoEsToken_; + edm::ESGetToken siPixelTemplateEsToken_; + edm::ESGetToken topoPerEventEsToken_; + edm::ESGetToken geomPerEventEsToken_; + + // event consumes + edm::EDGetTokenT t_trajTrack; +}; + +// +// constructors and destructor +// +SiPixelLorentzAnglePCLWorker::SiPixelLorentzAnglePCLWorker(const edm::ParameterSet& iConfig) + : folder_(iConfig.getParameter("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")), + clustChargeMaxPerLength_(iConfig.getParameter("clustChargeMaxPerLength")), + 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 + int 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); + } +} + +// +// member functions +// + +// ------------ method called for each event ------------ + +void SiPixelLorentzAnglePCLWorker::analyze(edm::Event const& iEvent, edm::EventSetup const& iSetup) { + // 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); + + 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(); + ++it) { + const reco::Track& track = *it->val; + const Trajectory& traj = *it->key; + + // get the trajectory measurements + std::vector tmColl = traj.measurements(); + pt_ = track.pt(); + eta_ = track.eta(); + phi_ = track.phi(); + chi2_ = traj.chiSquared(); + 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()); + + float ypitch_ = topol->pitch().second; + float width_ = 0.0285; + + if (!topol) + continue; + + 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(); + + isflipped_ = (tmp2 < tmp1) ? 1 : 0; + + const SiPixelRecHit* recHitPix = dynamic_cast((*recHit).hit()); + if (!recHitPix) + continue; + rechit_.x = recHitPix->localPosition().x(); + rechit_.y = recHitPix->localPosition().y(); + SiPixelRecHit::ClusterRef const& cluster = recHitPix->cluster(); + + pixinfo_ = fillPix(*cluster, topol); + + // fill entries in clust_ + + 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(); + 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 + 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 = 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.; + float locBz = locBx; + if (cotalpha < 0.) + locBz = -locBx; + + auto detId = detIdObj.rawId(); + 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); + 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; + LocalPoint lp_rechit = lp_pair.second; + + 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; + 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 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 = layer_ < 4 ? clustSizeYMin_ : clustSizeYMinL4_; + + if (!large_pix && (chi2_ / ndof_) < normChi2Max_ && cluster->sizeY() >= clustSizeY_cut && + 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 + 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) { + 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) { + // 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) { + 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++) { + 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_); + } + } + + 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 { + 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)); + + 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.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); +} + +// 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..d1894e8c425b2 --- /dev/null +++ b/Calibration/TkAlCaRecoProducers/python/ALCARECOPromptCalibProdSiPixelLorentzAngle_Output_cff.py @@ -0,0 +1,13 @@ +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_*_*', + ) +) +OutALCARECOPromptCalibProdSiPixelLorentzAngle=OutALCARECOPromptCalibProdSiPixelLorentzAngle_noDrop.clone() +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..d7de69d618d6d --- /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 = "ALCARECOPixelLACalibrationTracks", + TrajectoryInEvent = True, + NavigationSchool = "" + ) + +# 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/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cff.py b/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cff.py new file mode 100644 index 0000000000000..65e12849ff48a --- /dev/null +++ b/Calibration/TkAlCaRecoProducers/python/AlcaSiPixelLorentzAngleHarvester_cff.py @@ -0,0 +1,18 @@ +import FWCore.ParameterSet.Config as cms + +from Calibration.TkAlCaRecoProducers.AlcaSiPixelLorentzAngleHarvester_cfi import * +from DQMServices.Components.EDMtoMEConverter_cfi import * + +EDMtoMEConvertSiPixelLorentzAngle = EDMtoMEConverter.clone( + lumiInputTag = ("MEtoEDMConvertSiPixelLorentzAngle","MEtoEDMConverterLumi"), + runInputTag = ("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/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 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/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/AlCa/python/autoPCL.py b/Configuration/AlCa/python/autoPCL.py index 5d5a06f7ecd14..0876ddef0fc4c 100644 --- a/Configuration/AlCa/python/autoPCL.py +++ b/Configuration/AlCa/python/autoPCL.py @@ -3,11 +3,12 @@ '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', 'PromptCalibProdPPSDiamondSampic' : 'PPSDiamondSampicTimingCalibration' } 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_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 5586e1657abfd..af8689b01b987 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', @@ -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) 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',