File indexing completed on 2024-04-06 12:32:14
0001 #include "FWCore/Framework/interface/Event.h"
0002 #include "FWCore/Framework/interface/EventSetup.h"
0003 #include "FWCore/Framework/interface/ESTransientHandle.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "FWCore/ServiceRegistry/interface/Service.h"
0007 #include "FWCore/Utilities/interface/Exception.h"
0008
0009 #include "SimG4Core/Notification/interface/BeginOfEvent.h"
0010 #include "SimG4Core/Notification/interface/BeginOfTrack.h"
0011 #include "SimG4Core/Notification/interface/EndOfTrack.h"
0012 #include "SimG4Core/Watcher/interface/SimProducer.h"
0013 #include "SimG4Core/Notification/interface/Observer.h"
0014
0015 #include "DataFormats/Math/interface/GeantUnits.h"
0016 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0017 #include "SimDataFormats/ValidationFormats/interface/MaterialAccountingCalo.h"
0018
0019 #include "DetectorDescription/Core/interface/DDCompactView.h"
0020 #include "DetectorDescription/Core/interface/DDFilteredView.h"
0021 #include "DetectorDescription/Core/interface/DDFilter.h"
0022 #include "DetectorDescription/Core/interface/DDLogicalPart.h"
0023 #include "DetectorDescription/Core/interface/DDSplit.h"
0024 #include "DetectorDescription/Core/interface/DDValue.h"
0025 #include "DetectorDescription/Core/interface/DDsvalues.h"
0026 #include "DetectorDescription/DDCMS/interface/DDCompactView.h"
0027 #include "DetectorDescription/DDCMS/interface/DDFilteredView.h"
0028
0029 #include <CLHEP/Vector/LorentzVector.h>
0030
0031 #include "DD4hep/Filter.h"
0032
0033 #include "G4Step.hh"
0034 #include "G4Track.hh"
0035
0036 #include <iostream>
0037 #include <memory>
0038 #include <string>
0039 #include <vector>
0040
0041 using namespace geant_units::operators;
0042
0043 class MaterialBudgetHcalProducer : public SimProducer,
0044 public Observer<const BeginOfEvent*>,
0045 public Observer<const BeginOfTrack*>,
0046 public Observer<const G4Step*>,
0047 public Observer<const EndOfTrack*> {
0048 public:
0049 MaterialBudgetHcalProducer(const edm::ParameterSet&);
0050 ~MaterialBudgetHcalProducer() override = default;
0051
0052 MaterialBudgetHcalProducer(const MaterialBudgetHcalProducer&) = delete;
0053 const MaterialBudgetHcalProducer& operator=(const MaterialBudgetHcalProducer&) = delete;
0054
0055 void registerConsumes(edm::ConsumesCollector) override;
0056 void produce(edm::Event&, const edm::EventSetup&) override;
0057 void beginRun(edm::EventSetup const&) override;
0058
0059 private:
0060 void update(const BeginOfEvent*) override;
0061 void update(const BeginOfTrack*) override;
0062 void update(const G4Step*) override;
0063 void update(const EndOfTrack*) override;
0064
0065 bool stopAfter(const G4Step*);
0066 std::vector<std::string> getNames(DDFilteredView& fv);
0067 std::vector<std::string> getNames(cms::DDFilteredView& fv);
0068 std::vector<double> getDDDArray(const std::string& str, const DDsvalues_type& sv);
0069 bool isSensitive(const std::string&);
0070 bool isItHF(const G4VTouchable*);
0071 bool isItEC(const std::string&);
0072
0073 static const int maxSet_ = 25, maxSet2_ = 9;
0074 double rMax_, zMax_, etaMinP_, etaMaxP_;
0075 bool printSum_, fromdd4hep_;
0076 std::string name_;
0077
0078 edm::ESGetToken<DDCompactView, IdealGeometryRecord> cpvTokenDDD_;
0079 edm::ESGetToken<cms::DDCompactView, IdealGeometryRecord> cpvTokenDD4hep_;
0080
0081 std::vector<std::string> sensitives_, hfNames_, sensitiveEC_;
0082 std::vector<int> hfLevels_;
0083 MaterialAccountingCaloCollection matcoll_;
0084 double stepLens_[maxSet_], radLens_[maxSet_], intLens_[maxSet_];
0085 std::vector<std::string> matList_;
0086 std::vector<double> stepLength_, radLength_, intLength_;
0087 int id_, layer_, steps_;
0088 double radLen_, intLen_, stepLen_;
0089 double eta_, phi_;
0090 int nlayHB_, nlayHE_, nlayHO_, nlayHF_;
0091 };
0092
0093 MaterialBudgetHcalProducer::MaterialBudgetHcalProducer(const edm::ParameterSet& p) {
0094 edm::ParameterSet m_p = p.getParameter<edm::ParameterSet>("MaterialBudgetHcalProducer");
0095 rMax_ = m_p.getUntrackedParameter<double>("RMax", 4.5) * CLHEP::m;
0096 zMax_ = m_p.getUntrackedParameter<double>("ZMax", 13.0) * CLHEP::m;
0097 etaMinP_ = m_p.getUntrackedParameter<double>("EtaMinP", 5.2);
0098 etaMaxP_ = m_p.getUntrackedParameter<double>("EtaMaxP", 0.0);
0099 fromdd4hep_ = m_p.getUntrackedParameter<bool>("Fromdd4hep", false);
0100 printSum_ = m_p.getUntrackedParameter<bool>("PrintSummary", false);
0101 name_ = m_p.getUntrackedParameter<std::string>("Name", "Hcal");
0102 edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalProducer initialized with rMax " << rMax_ << " mm and zMax "
0103 << zMax_ << " mm printSummary is set to " << printSum_ << " and Fromdd4hep to "
0104 << fromdd4hep_;
0105
0106 produces<MaterialAccountingCaloCollection>(Form("%sMatBCalo", name_.c_str()));
0107 }
0108
0109 void MaterialBudgetHcalProducer::registerConsumes(edm::ConsumesCollector cc) {
0110 if (fromdd4hep_)
0111 cpvTokenDD4hep_ = cc.esConsumes<edm::Transition::BeginRun>();
0112 else
0113 cpvTokenDDD_ = cc.esConsumes<edm::Transition::BeginRun>();
0114 edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalProducer: Initialize the token for CompactView";
0115 }
0116
0117 void MaterialBudgetHcalProducer::produce(edm::Event& e, const edm::EventSetup&) {
0118 std::unique_ptr<MaterialAccountingCaloCollection> hgc(new MaterialAccountingCaloCollection);
0119 for (auto const& mbc : matcoll_) {
0120 hgc->emplace_back(mbc);
0121 }
0122 e.put(std::move(hgc), Form("%sMatBCalo", name_.c_str()));
0123 }
0124
0125 void MaterialBudgetHcalProducer::beginRun(edm::EventSetup const& es) {
0126
0127
0128 if (fromdd4hep_) {
0129 const cms::DDCompactView cpv = es.getData(cpvTokenDD4hep_);
0130 constexpr int32_t addLevel = 1;
0131 std::string attribute = "ReadOutName";
0132 std::string value = "HcalHits";
0133 const cms::DDFilter filter1(attribute, value);
0134 cms::DDFilteredView fv1(cpv, filter1);
0135 std::vector<std::string> names = getNames(fv1);
0136 for (auto& name : names) {
0137 std::string namx = (name.find('_') == std::string::npos) ? name : name.substr(0, name.find('_'));
0138 if (std::find(sensitives_.begin(), sensitives_.end(), namx) == sensitives_.end())
0139 sensitives_.emplace_back(namx);
0140 }
0141 edm::LogVerbatim("MaterialBudgetFull") << "MaterialBudgetHcalProducer: Names to be tested for " << attribute
0142 << " = " << value << " has " << sensitives_.size() << " elements";
0143 for (unsigned int i = 0; i < sensitives_.size(); i++)
0144 edm::LogVerbatim("MaterialBudgetFull")
0145 << "MaterialBudgetHcalProducer: sensitives[" << i << "] = " << sensitives_[i];
0146 attribute = "Volume";
0147 value = "HF";
0148 const cms::DDFilter filter2(attribute, value);
0149 cms::DDFilteredView fv2(cpv, filter2);
0150 std::vector<int> temp = fv2.get<std::vector<int> >("hf", "Levels");
0151 hfNames_ = getNames(fv2);
0152 edm::LogVerbatim("MaterialBudgetFull") << "MaterialBudgetHcalProducer: Names to be tested for " << attribute
0153 << " = " << value << " has " << hfNames_.size() << " elements";
0154 for (unsigned int i = 0; i < hfNames_.size(); i++) {
0155 hfLevels_.push_back(temp[i] + addLevel);
0156 edm::LogVerbatim("MaterialBudgetFull")
0157 << "MaterialBudgetHcalProducer: HF[" << i << "] = " << hfNames_[i] << " at level " << hfLevels_[i];
0158 }
0159
0160 const std::string ecalRO[2] = {"EcalHitsEB", "EcalHitsEE"};
0161 attribute = "ReadOutName";
0162 for (int k = 0; k < 2; k++) {
0163 value = ecalRO[k];
0164 const cms::DDFilter filter(attribute, value);
0165 cms::DDFilteredView fv(cpv, filter);
0166 std::vector<std::string> senstmp = getNames(fv);
0167 edm::LogVerbatim("MaterialBudgetFull") << "MaterialBudgetHcalProducer: Names to be tested for " << attribute
0168 << " = " << value << " has " << senstmp.size() << " elements";
0169 for (unsigned int i = 0; i < senstmp.size(); i++) {
0170 std::string name = senstmp[i].substr(0, 4);
0171 if (std::find(sensitiveEC_.begin(), sensitiveEC_.end(), name) == sensitiveEC_.end())
0172 sensitiveEC_.push_back(name);
0173 }
0174 }
0175 for (unsigned int i = 0; i < sensitiveEC_.size(); i++)
0176 edm::LogVerbatim("MaterialBudgetFull")
0177 << "MaterialBudgetHcalProducer:sensitiveEC[" << i << "] = " << sensitiveEC_[i];
0178 } else {
0179 const DDCompactView& cpv = es.getData(cpvTokenDDD_);
0180 constexpr int32_t addLevel = 0;
0181 std::string attribute = "ReadOutName";
0182 std::string value = "HcalHits";
0183 DDSpecificsMatchesValueFilter filter1{DDValue(attribute, value, 0)};
0184 DDFilteredView fv1(cpv, filter1);
0185 std::vector<std::string> names = getNames(fv1);
0186 for (auto& name : names) {
0187 std::string namx = (name.find('_') == std::string::npos) ? name : name.substr(0, name.find('_'));
0188 if (std::find(sensitives_.begin(), sensitives_.end(), namx) == sensitives_.end())
0189 sensitives_.emplace_back(namx);
0190 }
0191 edm::LogVerbatim("MaterialBudgetFull") << "MaterialBudgetHcalProducer: Names to be tested for " << attribute
0192 << " = " << value << " has " << sensitives_.size() << " elements";
0193 for (unsigned int i = 0; i < sensitives_.size(); i++)
0194 edm::LogVerbatim("MaterialBudgetFull")
0195 << "MaterialBudgetHcalProducer: sensitives[" << i << "] = " << sensitives_[i];
0196 attribute = "Volume";
0197 value = "HF";
0198 DDSpecificsMatchesValueFilter filter2{DDValue(attribute, value, 0)};
0199 DDFilteredView fv2(cpv, filter2);
0200 hfNames_ = getNames(fv2);
0201 fv2.firstChild();
0202 DDsvalues_type sv(fv2.mergedSpecifics());
0203 std::vector<double> temp = getDDDArray("Levels", sv);
0204 edm::LogVerbatim("MaterialBudgetFull") << "MaterialBudgetHcalProducer: Names to be tested for " << attribute
0205 << " = " << value << " has " << hfNames_.size() << " elements";
0206 for (unsigned int i = 0; i < hfNames_.size(); i++) {
0207 int level = static_cast<int>(temp[i]);
0208 hfLevels_.push_back(level + addLevel);
0209 edm::LogVerbatim("MaterialBudgetFull")
0210 << "MaterialBudgetHcalProducer: HF[" << i << "] = " << hfNames_[i] << " at level " << hfLevels_[i];
0211 }
0212
0213 const std::string ecalRO[2] = {"EcalHitsEB", "EcalHitsEE"};
0214 attribute = "ReadOutName";
0215 for (int k = 0; k < 2; k++) {
0216 value = ecalRO[k];
0217 DDSpecificsMatchesValueFilter filter3{DDValue(attribute, value, 0)};
0218 DDFilteredView fv3(cpv, filter3);
0219 std::vector<std::string> senstmp = getNames(fv3);
0220 edm::LogVerbatim("MaterialBudgetFull") << "MaterialBudgetHcalProducer: Names to be tested for " << attribute
0221 << " = " << value << " has " << senstmp.size() << " elements";
0222 for (unsigned int i = 0; i < senstmp.size(); i++) {
0223 std::string name = senstmp[i].substr(0, 4);
0224 if (std::find(sensitiveEC_.begin(), sensitiveEC_.end(), name) == sensitiveEC_.end())
0225 sensitiveEC_.push_back(name);
0226 }
0227 }
0228 for (unsigned int i = 0; i < sensitiveEC_.size(); i++)
0229 edm::LogVerbatim("MaterialBudgetFull")
0230 << "MaterialBudgetHcalProducer:sensitiveEC[" << i << "] = " << sensitiveEC_[i];
0231 }
0232 }
0233
0234 void MaterialBudgetHcalProducer::update(const BeginOfEvent*) { matcoll_.clear(); }
0235
0236 void MaterialBudgetHcalProducer::update(const BeginOfTrack* trk) {
0237 const G4Track* aTrack = (*trk)();
0238
0239 id_ = layer_ = steps_ = 0;
0240 radLen_ = intLen_ = stepLen_ = 0;
0241 nlayHB_ = nlayHE_ = nlayHF_ = nlayHO_ = 0;
0242 for (int i = 0; i < maxSet_; ++i)
0243 stepLens_[i] = radLens_[i] = intLens_[i] = 0;
0244
0245 const G4ThreeVector& dir = aTrack->GetMomentum();
0246 if (dir.theta() != 0) {
0247 eta_ = dir.eta();
0248 } else {
0249 eta_ = -99;
0250 }
0251 phi_ = dir.phi();
0252 edm::LogVerbatim("MaterialBudget") << "Start track with (eta, phi) = (" << eta_ << ", " << phi_ << " Energy "
0253 << aTrack->GetTotalEnergy() << " and ID "
0254 << (aTrack->GetDefinition()->GetPDGEncoding());
0255
0256 if (printSum_) {
0257 matList_.clear();
0258 stepLength_.clear();
0259 radLength_.clear();
0260 intLength_.clear();
0261 }
0262 }
0263
0264 void MaterialBudgetHcalProducer::update(const G4Step* aStep) {
0265
0266 G4Material* material = aStep->GetPreStepPoint()->GetMaterial();
0267 double step = aStep->GetStepLength();
0268 double radl = material->GetRadlen();
0269 double intl = material->GetNuclearInterLength();
0270 double density = convertUnitsTo(1._g_per_cm3, material->GetDensity());
0271
0272 int idOld = id_;
0273 const G4VTouchable* touch = aStep->GetPreStepPoint()->GetTouchable();
0274 std::string name = (static_cast<std::string>(dd4hep::dd::noNamespace(touch->GetVolume(0)->GetName())));
0275 std::string matName = (static_cast<std::string>(dd4hep::dd::noNamespace(material->GetName())));
0276 if (printSum_) {
0277 bool found = false;
0278 for (unsigned int ii = 0; ii < matList_.size(); ii++) {
0279 if (matList_[ii] == matName) {
0280 stepLength_[ii] += step;
0281 radLength_[ii] += (step / radl);
0282 intLength_[ii] += (step / intl);
0283 found = true;
0284 break;
0285 }
0286 }
0287 if (!found) {
0288 matList_.push_back(matName);
0289 stepLength_.push_back(step);
0290 radLength_.push_back(step / radl);
0291 intLength_.push_back(step / intl);
0292 }
0293 if ((std::abs(eta_) >= etaMinP_) && (std::abs(eta_) <= etaMaxP_))
0294 edm::LogVerbatim("MaterialBudget") << "Volume " << name << " id " << id_ << ":" << idOld << " Step " << step
0295 << " Material " << matName << " Old Length " << stepLen_ << " X0 "
0296 << step / radl << ":" << radLen_ << " Lambda " << step / intl << ":"
0297 << intLen_;
0298 } else {
0299 if ((std::abs(eta_) >= etaMinP_) && (std::abs(eta_) <= etaMaxP_))
0300 edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalProducer: Step at " << name << " id " << id_ << ":"
0301 << idOld << " Length " << step << " in " << matName << " of density "
0302 << density << " g/cc; Radiation Length " << radl << " mm; Interaction Length "
0303 << intl << " mm\n Position "
0304 << aStep->GetPreStepPoint()->GetPosition() << " Cylindrical R "
0305 << aStep->GetPreStepPoint()->GetPosition().perp() << " Length (so far) "
0306 << stepLen_ << " L/X0 " << step / radl << "/" << radLen_ << " L/Lambda "
0307 << step / intl << "/" << intLen_;
0308 }
0309
0310 int det = 0, lay = 0;
0311 double abseta = std::abs(eta_);
0312 edm::LogVerbatim("MaterialBudgetFull") << "Volume " << name << ":" << matName << " EC:Sensitive:HF " << isItEC(name)
0313 << ":" << isSensitive(name) << ":" << isItHF(touch) << " Eta " << abseta
0314 << " HC " << ((touch->GetReplicaNumber(1)) / 1000) << ":"
0315 << ((touch->GetReplicaNumber(0) / 10) % 100 + 3) << " X0 "
0316 << (radLen_ + (step / radl)) << " Lambda " << (intLen_ + (step / intl));
0317
0318 if (isItEC(name)) {
0319 det = 1;
0320 lay = 1;
0321 } else {
0322 if (isSensitive(name)) {
0323 if (isItHF(touch)) {
0324 det = 5;
0325 lay = 21;
0326 if (lay != layer_)
0327 ++nlayHF_;
0328 } else {
0329 det = (touch->GetReplicaNumber(1)) / 1000;
0330 lay = (touch->GetReplicaNumber(0) / 10) % 100 + 3;
0331 if (det == 4) {
0332 if (abseta < 1.479)
0333 lay = layer_ + 1;
0334 else
0335 lay--;
0336 if (lay < 3)
0337 lay = 3;
0338 if (lay == layer_)
0339 lay++;
0340 if (lay > 20)
0341 lay = 20;
0342 if (lay != layer_)
0343 ++nlayHE_;
0344 } else if (lay != layer_) {
0345 if (lay >= 20)
0346 ++nlayHO_;
0347 else
0348 ++nlayHB_;
0349 }
0350 }
0351 edm::LogVerbatim("MaterialBudgetFull") << "MaterialBudgetHcalProducer: Det " << det << " Layer " << lay << " Eta "
0352 << eta_ << " Phi " << convertRadToDeg(phi_);
0353 } else if (layer_ == 1) {
0354 det = -1;
0355 lay = 2;
0356 }
0357 }
0358 if (det != 0) {
0359 if (lay != layer_) {
0360 id_ = lay;
0361 layer_ = lay;
0362 }
0363 }
0364
0365 if (id_ > idOld) {
0366 if ((abseta >= etaMinP_) && (abseta <= etaMaxP_))
0367 edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalProducer: Step at " << name << " calls filHisto with "
0368 << (id_ - 1);
0369 stepLens_[id_ - 1] = stepLen_;
0370 radLens_[id_ - 1] = radLen_;
0371 intLens_[id_ - 1] = intLen_;
0372 }
0373
0374 stepLen_ += step;
0375 radLen_ += (step / radl);
0376 intLen_ += (step / intl);
0377 if (id_ == 21) {
0378 if (!isItHF(aStep->GetPostStepPoint()->GetTouchable())) {
0379 if ((abseta >= etaMinP_) && (abseta <= etaMaxP_))
0380 edm::LogVerbatim("MaterialBudget") << "MaterialBudgetHcalProducer: After HF in " << name << ":"
0381 << static_cast<std::string>(dd4hep::dd::noNamespace(
0382 aStep->GetPostStepPoint()->GetTouchable()->GetVolume(0)->GetName()))
0383 << " calls fillHisto with " << id_;
0384 stepLens_[idOld] = stepLen_;
0385 radLens_[idOld] = radLen_;
0386 intLens_[idOld] = intLen_;
0387 ++id_;
0388 layer_ = 0;
0389 }
0390 }
0391
0392
0393 if (stopAfter(aStep)) {
0394 G4Track* track = aStep->GetTrack();
0395 track->SetTrackStatus(fStopAndKill);
0396 }
0397 }
0398
0399 void MaterialBudgetHcalProducer::update(const EndOfTrack* trk) {
0400 MaterialAccountingCalo matCalo;
0401 matCalo.m_eta = eta_;
0402 matCalo.m_phi = phi_;
0403 for (int k = 0; k < maxSet_; ++k) {
0404 matCalo.m_stepLen.emplace_back(stepLens_[k]);
0405 matCalo.m_radLen.emplace_back(radLens_[k]);
0406 matCalo.m_intLen.emplace_back(intLens_[k]);
0407 }
0408 matCalo.m_layers.emplace_back(nlayHB_);
0409 matCalo.m_layers.emplace_back(nlayHE_);
0410 matCalo.m_layers.emplace_back(nlayHO_);
0411 matCalo.m_layers.emplace_back(nlayHF_);
0412 matcoll_.emplace_back(matCalo);
0413 }
0414
0415 bool MaterialBudgetHcalProducer::stopAfter(const G4Step* aStep) {
0416 G4ThreeVector hitPoint = aStep->GetPreStepPoint()->GetPosition();
0417 double rr = hitPoint.perp();
0418 double zz = std::abs(hitPoint.z());
0419
0420 if (rr > rMax_ || zz > zMax_) {
0421 edm::LogVerbatim("MaterialBudget") << " MaterialBudgetHcalProducer::StopAfter R = " << rr << " and Z = " << zz;
0422 return true;
0423 } else {
0424 return false;
0425 }
0426 }
0427
0428 std::vector<std::string> MaterialBudgetHcalProducer::getNames(DDFilteredView& fv) {
0429 std::vector<std::string> tmp;
0430 bool dodet = fv.firstChild();
0431 while (dodet) {
0432 const DDLogicalPart& log = fv.logicalPart();
0433 std::string namx = log.name().name();
0434 if (std::find(tmp.begin(), tmp.end(), namx) == tmp.end())
0435 tmp.push_back(namx);
0436 dodet = fv.next();
0437 }
0438 return tmp;
0439 }
0440
0441 std::vector<std::string> MaterialBudgetHcalProducer::getNames(cms::DDFilteredView& fv) {
0442 std::vector<std::string> tmp;
0443 const std::vector<std::string> notIn = {
0444 "CALO", "HCal", "MBBTL", "MBBTR", "MBBTC", "MBAT", "MBBT_R1M", "MBBT_R1P", "MBBT_R1MX", "MBBT_R1PX", "VCAL"};
0445 while (fv.firstChild()) {
0446 const std::string n{fv.name().data(), fv.name().size()};
0447 if (std::find(notIn.begin(), notIn.end(), n) == notIn.end()) {
0448 std::string::size_type pos = n.find(':');
0449 const std::string namx = (pos == std::string::npos) ? n : std::string(n, pos + 1, n.size() - 1);
0450 if (std::find(tmp.begin(), tmp.end(), namx) == tmp.end())
0451 tmp.push_back(namx);
0452 }
0453 }
0454 return tmp;
0455 }
0456
0457 std::vector<double> MaterialBudgetHcalProducer::getDDDArray(const std::string& str, const DDsvalues_type& sv) {
0458 edm::LogVerbatim("MaterialBudgetFull") << "MaterialBudgetHcalProducer:getDDDArray called for " << str;
0459 DDValue value(str);
0460 if (DDfetch(&sv, value)) {
0461 edm::LogVerbatim("MaterialBudgetFull") << value;
0462 const std::vector<double>& fvec = value.doubles();
0463 int nval = fvec.size();
0464 if (nval < 1) {
0465 throw cms::Exception("MaterialBudgetHcalProducer") << "nval = " << nval << " < 1 for array " << str << "\n";
0466 }
0467
0468 return fvec;
0469 } else {
0470 throw cms::Exception("MaterialBudgetHcalProducer") << "cannot get array " << str << "\n";
0471 }
0472 }
0473
0474 bool MaterialBudgetHcalProducer::isSensitive(const std::string& name) {
0475 std::vector<std::string>::const_iterator it = sensitives_.begin();
0476 std::vector<std::string>::const_iterator itEnd = sensitives_.end();
0477 std::string namx = (name.find('_') == std::string::npos) ? name : name.substr(0, name.find('_'));
0478 for (; it != itEnd; ++it)
0479 if (namx == *it)
0480 return true;
0481 return false;
0482 }
0483
0484 bool MaterialBudgetHcalProducer::isItHF(const G4VTouchable* touch) {
0485 int levels = ((touch->GetHistoryDepth()) + 1);
0486 for (unsigned int it = 0; it < hfNames_.size(); it++) {
0487 if (levels >= hfLevels_[it]) {
0488 std::string name =
0489 (static_cast<std::string>(dd4hep::dd::noNamespace(touch->GetVolume(levels - hfLevels_[it])->GetName())))
0490 .substr(0, 4);
0491 if (name == hfNames_[it]) {
0492 return true;
0493 }
0494 }
0495 }
0496 return false;
0497 }
0498
0499 bool MaterialBudgetHcalProducer::isItEC(const std::string& name) {
0500 std::vector<std::string>::const_iterator it = sensitiveEC_.begin();
0501 std::vector<std::string>::const_iterator itEnd = sensitiveEC_.end();
0502 for (; it != itEnd; ++it)
0503 if (name.substr(0, 4) == *it)
0504 return true;
0505 return false;
0506 }
0507
0508 #include "SimG4Core/Watcher/interface/SimWatcherFactory.h"
0509 #include "FWCore/PluginManager/interface/ModuleDef.h"
0510
0511 DEFINE_SIMWATCHER(MaterialBudgetHcalProducer);
