File indexing completed on 2024-04-06 11:59:33
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021 #include "CalibTracker/SiPixelESProducers/interface/SiPixelDetInfoFileReader.h"
0022 #include "CalibTracker/SiPixelESProducers/interface/SiPixelFakeLorentzAngleESSource.h"
0023 #include "CalibTracker/StandaloneTrackerTopology/interface/StandaloneTrackerTopology.h"
0024 #include "DataFormats/TrackerCommon/interface/PixelBarrelName.h"
0025 #include "DataFormats/TrackerCommon/interface/PixelEndcapName.h"
0026 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0027 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0028 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0029 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0030 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0031
0032
0033
0034
0035 SiPixelFakeLorentzAngleESSource::SiPixelFakeLorentzAngleESSource(const edm::ParameterSet& conf_)
0036 : fp_(conf_.getParameter<edm::FileInPath>("file")),
0037 t_topo_fp_(conf_.getParameter<edm::FileInPath>("topologyInput")),
0038 myLabel_(conf_.getParameter<std::string>("appendToDataLabel")),
0039 BPixParameters_(conf_.getParameter<Parameters>("BPixParameters")),
0040 FPixParameters_(conf_.getParameter<Parameters>("FPixParameters")),
0041 ModuleParameters_(conf_.getParameter<Parameters>("ModuleParameters")),
0042 bPixLorentzAnglePerTesla_((float)conf_.getUntrackedParameter<double>("bPixLorentzAnglePerTesla", -9999.)),
0043 fPixLorentzAnglePerTesla_((float)conf_.getUntrackedParameter<double>("fPixLorentzAnglePerTesla", -9999.)) {
0044 edm::LogInfo("SiPixelFakeLorentzAngleESSource::SiPixelFakeLorentzAngleESSource");
0045
0046 setWhatProduced(this);
0047 findingRecord<SiPixelLorentzAngleRcd>();
0048 }
0049
0050 std::unique_ptr<SiPixelLorentzAngle> SiPixelFakeLorentzAngleESSource::produce(const SiPixelLorentzAngleRcd&) {
0051 using namespace edm::es;
0052 unsigned int nmodules = 0;
0053 SiPixelLorentzAngle* obj = new SiPixelLorentzAngle();
0054 SiPixelDetInfoFileReader reader(fp_.fullPath());
0055 const std::vector<uint32_t>& DetIds = reader.getAllDetIds();
0056
0057 TrackerTopology tTopo = StandaloneTrackerTopology::fromTrackerParametersXMLFile(t_topo_fp_.fullPath());
0058
0059
0060 for (const auto& detit : DetIds) {
0061 nmodules++;
0062 const DetId detid(detit);
0063 auto rawId = detid.rawId();
0064 int found = 0;
0065 int side = tTopo.side(detid);
0066
0067
0068 if (detid.subdetId() == static_cast<int>(PixelSubdetector::PixelBarrel)) {
0069 int layer = tTopo.pxbLayer(detid);
0070
0071 int ladder = tTopo.pxbLadder(detid);
0072
0073 int module = tTopo.pxbModule(detid);
0074 if (module < 5) {
0075 side = 1;
0076 } else {
0077 side = 2;
0078 }
0079
0080 LogDebug("SiPixelFakeLorentzAngleESSource") << " pixel barrel:"
0081 << " layer=" << layer << " ladder=" << ladder << " module=" << module
0082 << " rawId=" << rawId << " side " << side;
0083
0084
0085 if (bPixLorentzAnglePerTesla_ != -9999.) {
0086 edm::LogInfo("SiPixelFakeLorentzAngleESSource")
0087 << " LA = " << bPixLorentzAnglePerTesla_ << " common for all bpix" << std::endl;
0088 if (!obj->putLorentzAngle(detid.rawId(), bPixLorentzAnglePerTesla_))
0089 edm::LogError("SiPixelFakeLorentzAngleESSource")
0090 << "ERROR!: detid " << rawId << " already exists" << std::endl;
0091 } else {
0092
0093 for (const auto& it : ModuleParameters_) {
0094 if (it.getParameter<unsigned int>("rawid") == detid.rawId()) {
0095 float lorentzangle = (float)it.getParameter<double>("angle");
0096 if (!found) {
0097 obj->putLorentzAngle(detid.rawId(), lorentzangle);
0098 edm::LogInfo("SiPixelFakeLorentzAngleESSource")
0099 << " LA= " << lorentzangle << " individual value " << detid.rawId() << std::endl;
0100 found = 1;
0101 } else
0102 edm::LogError("SiPixelFakeLorentzAngleESSource") << "ERROR!: detid already exists" << std::endl;
0103 }
0104 }
0105
0106
0107 for (const auto& it : BPixParameters_) {
0108 if (it.exists("layer"))
0109 if (it.getParameter<int>("layer") != layer)
0110 continue;
0111 if (it.exists("ladder"))
0112 if (it.getParameter<int>("ladder") != ladder)
0113 continue;
0114 if (it.exists("module"))
0115 if (it.getParameter<int>("module") != module)
0116 continue;
0117 if (it.exists("side"))
0118 if (it.getParameter<int>("side") != side)
0119 continue;
0120 if (!found) {
0121 float lorentzangle = (float)it.getParameter<double>("angle");
0122 obj->putLorentzAngle(detid.rawId(), lorentzangle);
0123 edm::LogInfo("SiPixelFakeLorentzAngleESSource") << " LA= " << lorentzangle << std::endl;
0124 found = 2;
0125 } else if (found == 1) {
0126 edm::LogWarning("SiPixelFakeLorentzAngleESSource")
0127 << "detid already given in ModuleParameters, skipping ..." << std::endl;
0128 } else
0129 edm::LogError("SiPixelFakeLorentzAngleESSource")
0130 << " ERROR!: detid " << rawId << " already exists" << std::endl;
0131 }
0132 }
0133 } else if (detid.subdetId() == static_cast<int>(PixelSubdetector::PixelEndcap)) {
0134
0135
0136
0137
0138 PixelEndcapName pen(detid, &tTopo, true);
0139
0140
0141
0142
0143 int disk = pen.diskName();
0144 int blade = pen.bladeName();
0145 int panel = pen.pannelName();
0146 int ring = pen.ringName();
0147
0148 LogDebug("SiPixelFakeLorentzAngleESSource") << " pixel endcap:"
0149 << " side=" << side << " disk=" << disk << " blade =" << blade
0150 << " pannel=" << panel << " ring=" << ring << " rawId=" << rawId;
0151
0152
0153 if (fPixLorentzAnglePerTesla_ != -9999.) {
0154 edm::LogInfo("SiPixelFakeLorentzAngleESSource")
0155 << " LA = " << fPixLorentzAnglePerTesla_ << " common for all fpix" << std::endl;
0156 if (!obj->putLorentzAngle(detid.rawId(), fPixLorentzAnglePerTesla_))
0157 edm::LogError("SiPixelFakeLorentzAngleESSource") << " ERROR! detid already exists" << std::endl;
0158 } else {
0159
0160 for (const auto& it : ModuleParameters_) {
0161 if (it.getParameter<unsigned int>("rawid") == detid.rawId()) {
0162 float lorentzangle = (float)it.getParameter<double>("angle");
0163 if (!found) {
0164 obj->putLorentzAngle(detid.rawId(), lorentzangle);
0165 edm::LogInfo("SiPixelFakeLorentzAngleESSource")
0166 << " LA= " << lorentzangle << " individual value " << detid.rawId() << std::endl;
0167 found = 1;
0168 } else
0169 edm::LogError("SiPixelFakeLorentzAngleESSource")
0170 << "ERROR!: detid " << rawId << " already exists" << std::endl;
0171 }
0172 }
0173
0174
0175 for (const auto& it : FPixParameters_) {
0176 if (it.exists("side"))
0177 if (it.getParameter<int>("side") != side)
0178 continue;
0179 if (it.exists("disk"))
0180 if (it.getParameter<int>("disk") != disk)
0181 continue;
0182 if (it.exists("ring"))
0183 if (it.getParameter<int>("ring") != ring)
0184 continue;
0185 if (it.exists("blade"))
0186 if (it.getParameter<int>("blade") != blade)
0187 continue;
0188 if (it.exists("panel"))
0189 if (it.getParameter<int>("panel") != panel)
0190 continue;
0191 if (it.exists("HVgroup"))
0192 if (it.getParameter<int>("HVgroup") != HVgroup(panel, ring))
0193 continue;
0194 if (!found) {
0195 float lorentzangle = (float)it.getParameter<double>("angle");
0196 obj->putLorentzAngle(detid.rawId(), lorentzangle);
0197 edm::LogInfo("SiPixelFakeLorentzAngleESSource") << " LA= " << lorentzangle << std::endl;
0198 found = 2;
0199 } else if (found == 1) {
0200 edm::LogWarning("SiPixelFakeLorentzAngleESSource")
0201 << " detid " << rawId << " already given in ModuleParameters, skipping ..." << std::endl;
0202 } else
0203 edm::LogError("SiPixelFakeLorentzAngleESSource")
0204 << "ERROR!: detid" << rawId << "already exists" << std::endl;
0205 }
0206 }
0207 }
0208 }
0209
0210 edm::LogInfo("SiPixelFakeLorentzAngleESSource") << "Modules = " << nmodules << std::endl;
0211
0212 return std::unique_ptr<SiPixelLorentzAngle>(obj);
0213 }
0214
0215 void SiPixelFakeLorentzAngleESSource::setIntervalFor(const edm::eventsetup::EventSetupRecordKey&,
0216 const edm::IOVSyncValue& iosv,
0217 edm::ValidityInterval& oValidity) {
0218 edm::ValidityInterval infinity(iosv.beginOfTime(), iosv.endOfTime());
0219 oValidity = infinity;
0220 }
0221
0222 int SiPixelFakeLorentzAngleESSource::HVgroup(int panel, int module) {
0223 if (1 == panel && (1 == module || 2 == module)) {
0224 return 1;
0225 } else if (1 == panel && (3 == module || 4 == module)) {
0226 return 2;
0227 } else if (2 == panel && 1 == module) {
0228 return 1;
0229 } else if (2 == panel && (2 == module || 3 == module)) {
0230 return 2;
0231 } else {
0232 edm::LogError("SiPixelFakeLorentzAngleESSource")
0233 << " ERROR! in SiPixelFakeLorentzAngleESSource::HVgroup(...), panel = " << panel << ", module = " << module
0234 << std::endl;
0235 return 0;
0236 }
0237 }
0238
0239 void SiPixelFakeLorentzAngleESSource::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0240 edm::ParameterSetDescription desc;
0241 desc.setComment("ESSource to supply per-module SiPixelLorentzAngle payloads in the EventSetup");
0242
0243 desc.add<edm::FileInPath>(
0244 "file", edm::FileInPath("SLHCUpgradeSimulations/Geometry/data/PhaseI/PixelSkimmedGeometry_phase1.txt"))
0245 ->setComment("Tracker skimmed geometry");
0246 desc.add<edm::FileInPath>("topologyInput",
0247 edm::FileInPath("Geometry/TrackerCommonData/data/PhaseI/trackerParameters.xml"))
0248 ->setComment("Tracker Topology");
0249
0250 desc.add<std::string>("appendToDataLabel", "")->setComment("label to which write the data");
0251 desc.addUntracked<double>("bPixLorentzAnglePerTesla", -9999.)->setComment("LA value for all BPix");
0252 desc.addUntracked<double>("fPixLorentzAnglePerTesla", -9999.)->setComment("LA value for all FPix");
0253
0254 edm::ParameterSetDescription desc_BPixParameters;
0255 desc_BPixParameters.addOptional<int>("layer");
0256 desc_BPixParameters.addOptional<int>("ladder");
0257 desc_BPixParameters.addOptional<int>("module");
0258 desc_BPixParameters.addOptional<int>("side");
0259 desc_BPixParameters.add<double>("angle");
0260 std::vector<edm::ParameterSet> default_BPixParametersToAdd;
0261 desc.addVPSet("BPixParameters", desc_BPixParameters, default_BPixParametersToAdd)
0262 ->setComment("LA values for given BPix regions");
0263
0264 edm::ParameterSetDescription desc_FPixParameters;
0265 desc_FPixParameters.addOptional<int>("side");
0266 desc_FPixParameters.addOptional<int>("disk");
0267 desc_FPixParameters.addOptional<int>("ring");
0268 desc_FPixParameters.addOptional<int>("blade");
0269 desc_FPixParameters.addOptional<int>("panel");
0270 desc_FPixParameters.addOptional<int>("HVgroup");
0271 desc_FPixParameters.add<double>("angle");
0272 std::vector<edm::ParameterSet> default_FPixParametersToAdd;
0273 desc.addVPSet("FPixParameters", desc_FPixParameters, default_FPixParametersToAdd)
0274 ->setComment("LA values for given FPix regions");
0275
0276 edm::ParameterSetDescription desc_ModuleParameters;
0277 desc_ModuleParameters.add<unsigned int>("rawid");
0278 desc_ModuleParameters.add<double>("angle");
0279 std::vector<edm::ParameterSet> default_ModuleParametersToAdd;
0280 desc.addVPSet("ModuleParameters", desc_ModuleParameters, default_ModuleParametersToAdd)
0281 ->setComment("LA values for given modules");
0282
0283 descriptions.addWithDefaultLabel(desc);
0284 }