Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:05:28

0001 #include "DD4hep/VolumeProcessor.h"
0002 #include "DD4hep/detail/DetectorInterna.h"
0003 #include "DD4hep/DetFactoryHelper.h"
0004 #include "DD4hep/DetectorHelper.h"
0005 #include "DD4hep/Printout.h"
0006 #include "DetectorDescription/DDCMS/interface/DDAlgoArguments.h"
0007 
0008 #include <sstream>
0009 
0010 using namespace std;
0011 using namespace cms;
0012 using namespace dd4hep;
0013 
0014 namespace cms {
0015 
0016   // Heuristically assign DetElement structures
0017   // to the sensitive volume pathes
0018   //
0019   class DDCMSDetElementCreator : public dd4hep::PlacedVolumeProcessor {
0020   public:
0021     DDCMSDetElementCreator(dd4hep::Detector&);
0022     ~DDCMSDetElementCreator() override;
0023 
0024     /// Callback to output PlacedVolume information of an single Placement
0025     int operator()(dd4hep::PlacedVolume volume, int level) override;
0026     /// Callback to output PlacedVolume information of an entire Placement
0027     int process(dd4hep::PlacedVolume volume, int level, bool recursive) override;
0028 
0029   private:
0030     dd4hep::DetElement addSubdetector(const std::string& nam, dd4hep::PlacedVolume volume, bool valid);
0031     dd4hep::DetElement createElement(const char* debugTag, dd4hep::PlacedVolume volume, int id);
0032     void createTopLevelDetectors(dd4hep::PlacedVolume volume);
0033 
0034     struct Data {
0035       Data() = default;
0036       Data(dd4hep::PlacedVolume v) : volume(v) {}
0037       Data(const Data& d) = default;
0038       Data& operator=(const Data& d) = default;
0039 
0040       dd4hep::PlacedVolume volume{nullptr};
0041       dd4hep::DetElement element{};
0042       bool sensitive = false;
0043       bool hasSensitive = false;
0044       int volumeCount = 0;
0045       int daughterCount = 0;
0046       int sensitiveCount = 0;
0047     };
0048 
0049     struct Count {
0050       Count() = default;
0051       Count(const Count&) = default;
0052       Count& operator=(const Count&) = default;
0053 
0054       int elements = 0;
0055       int volumes = 0;
0056       int sensitives = 0;
0057     };
0058 
0059     using Detectors = std::map<std::string, dd4hep::DetElement>;
0060     using Counters = std::map<dd4hep::DetElement, Count>;
0061     using LeafCount = std::map<std::pair<dd4hep::DetElement, int>, std::pair<int, int> >;
0062     using VolumeStack = std::vector<Data>;
0063 
0064     std::map<dd4hep::PlacedVolume, std::pair<int, int> > m_allPlacements;
0065 
0066     Counters m_counters;
0067     LeafCount m_leafCount;
0068     VolumeStack m_stack;
0069     Detectors m_subdetectors;
0070     dd4hep::DetElement m_tracker, m_currentDetector;
0071     dd4hep::SensitiveDetector m_currentSensitive;
0072     dd4hep::Detector& m_description;
0073     dd4hep::Atom m_silicon;
0074   };
0075 
0076   std::string detElementName(dd4hep::PlacedVolume volume);
0077 }  // namespace cms
0078 
0079 std::string cms::detElementName(dd4hep::PlacedVolume volume) {
0080   if (volume.isValid()) {
0081     std::string name = volume.name();
0082     std::string nnam = name.substr(name.find(NAMESPACE_SEP) + 1);
0083     return nnam;
0084   }
0085   except("DD4CMS", "++ Cannot deduce name from invalid PlacedVolume handle!");
0086   return std::string();
0087 }
0088 
0089 DDCMSDetElementCreator::DDCMSDetElementCreator(dd4hep::Detector& desc) : m_description(desc) {
0090   dd4hep::DetectorHelper helper(m_description);
0091   m_silicon = helper.element("SI");
0092   if (!m_silicon.isValid()) {
0093     except("DDCMSDetElementCreator", "++ Failed to extract SILICON from the element table.");
0094   }
0095   m_stack.reserve(32);
0096 }
0097 
0098 DDCMSDetElementCreator::~DDCMSDetElementCreator() {
0099   Count total;
0100   stringstream str, id_str;
0101 
0102   printout(INFO, "DDCMSDetElementCreator", "+++++++++++++++ Summary of sensitve elements  ++++++++++++++++++++++++");
0103   for (const auto& c : m_counters) {
0104     printout(INFO,
0105              "DDCMSDetElementCreator",
0106              "++ Summary: SD: %-24s %7d DetElements %7d sensitives out of %7d volumes",
0107              (c.first.name() + string(":")).c_str(),
0108              c.second.elements,
0109              c.second.sensitives,
0110              c.second.volumes);
0111     total.elements += c.second.elements;
0112     total.sensitives += c.second.sensitives;
0113     total.volumes += c.second.volumes;
0114   }
0115   printout(INFO,
0116            "DDCMSDetElementCreator",
0117            "++ Summary:     %-24s %7d DetElements %7d sensitives out of %7d volumes",
0118            "Grand Total:",
0119            total.elements,
0120            total.sensitives,
0121            total.volumes);
0122   printout(INFO, "DDCMSDetElementCreator", "+++++++++++++++ Summary of geometry depth analysis  ++++++++++++++++++");
0123   int totalCount = 0;
0124   map<dd4hep::DetElement, vector<pair<int, int> > > fields;
0125   for (const auto& l : m_leafCount) {
0126     dd4hep::DetElement de = l.first.first;
0127     printout(INFO,
0128              "DDCMSDetElementCreator",
0129              "++ Summary: SD: %-24s system:%04X Lvl:%3d Sensitives: %6d [Max: %6d].",
0130              (de.name() + string(":")).c_str(),
0131              de.id(),
0132              l.first.second,
0133              l.second.second,
0134              l.second.first);
0135     fields[de].push_back(make_pair(l.first.second, l.second.first));
0136     ++totalCount;
0137   }
0138   printout(INFO, "DDCMSDetElementCreator", "++ Summary:     %-24s  %d.", "Total DetElements:", totalCount);
0139   printout(INFO, "DDCMSDetElementCreator", "+++++++++++++++ Readout structure generation  ++++++++++++++++++++++++");
0140   str << endl;
0141   for (const auto& f : fields) {
0142     string roName = f.first.name() + string("Hits");
0143     int num_bits = 8;
0144     id_str.str("");
0145     id_str << "system:" << num_bits;
0146     for (const auto& q : f.second) {
0147       int bits = 0;
0148       if (q.second < 1 << 0)
0149         bits = 1;
0150       else if (q.second < 1 << 1)
0151         bits = 1;
0152       else if (q.second < 1 << 2)
0153         bits = 2;
0154       else if (q.second < 1 << 3)
0155         bits = 3;
0156       else if (q.second < 1 << 4)
0157         bits = 4;
0158       else if (q.second < 1 << 5)
0159         bits = 5;
0160       else if (q.second < 1 << 6)
0161         bits = 6;
0162       else if (q.second < 1 << 7)
0163         bits = 7;
0164       else if (q.second < 1 << 8)
0165         bits = 8;
0166       else if (q.second < 1 << 9)
0167         bits = 9;
0168       else if (q.second < 1 << 10)
0169         bits = 10;
0170       else if (q.second < 1 << 11)
0171         bits = 11;
0172       else if (q.second < 1 << 12)
0173         bits = 12;
0174       else if (q.second < 1 << 13)
0175         bits = 13;
0176       else if (q.second < 1 << 14)
0177         bits = 14;
0178       else if (q.second < 1 << 15)
0179         bits = 15;
0180       bits += 1;
0181       id_str << ",Lv" << q.first << ":" << bits;
0182       num_bits += bits;
0183     }
0184     string idspec = id_str.str();
0185     str << "<readout name=\"" << roName << "\">" << endl
0186         << "\t<id>" << idspec << "</id>  <!-- Number of bits: " << num_bits << " -->" << endl
0187         << "</readout>" << endl;
0188 
0189     /// Create ID Descriptors and readout configurations
0190     IDDescriptor dsc(roName, idspec);
0191     m_description.addIDSpecification(dsc);
0192     Readout ro(roName);
0193     ro.setIDDescriptor(dsc);
0194     m_description.addReadout(ro);
0195     dd4hep::SensitiveDetector sd = m_description.sensitiveDetector(f.first.name());
0196     sd.setHitsCollection(ro.name());
0197     sd.setReadout(ro);
0198     printout(INFO,
0199              "DDCMSDetElementCreator",
0200              "++ Setting up readout for subdetector:%-24s id:%04X",
0201              f.first.name(),
0202              f.first.id());
0203   }
0204   printout(INFO, "DDCMSDetElementCreator", "+++++++++++++++ ID Descriptor generation  ++++++++++++++++++++++++++++");
0205   printout(INFO, "", str.str().c_str());
0206   char volId[32];
0207   for (auto& p : m_allPlacements) {
0208     dd4hep::PlacedVolume place = p.first;
0209     dd4hep::Volume volume = place.volume();
0210     ::snprintf(volId, sizeof(volId), "Lv%d", p.second.first);
0211     printout(DEBUG,
0212              "DDCMSDetElementCreator",
0213              "++ Set volid (%-24s): %-6s = %3d  -> %s  (%p)",
0214              volume.isSensitive() ? volume.sensitiveDetector().name() : "Not Sensitive",
0215              volId,
0216              p.second.second,
0217              place.name(),
0218              place.ptr());
0219     place.addPhysVolID(volId, p.second.second);
0220   }
0221   printout(ALWAYS,
0222            "DDCMSDetElementCreator",
0223            "++ Instrumented %ld subdetectors with %d DetElements %d sensitives out of %d volumes and %ld sensitive "
0224            "placements.",
0225            fields.size(),
0226            total.elements,
0227            total.sensitives,
0228            total.volumes,
0229            m_allPlacements.size());
0230 }
0231 
0232 dd4hep::DetElement DDCMSDetElementCreator::createElement(const char*, PlacedVolume volume, int id) {
0233   string name = detElementName(volume);
0234   dd4hep::DetElement det(name, id);
0235   det.setPlacement(volume);
0236   return det;
0237 }
0238 void DDCMSDetElementCreator::createTopLevelDetectors(PlacedVolume volume) {
0239   auto& data = m_stack.back();
0240   if (m_stack.size() == 2) {  // Main subssystem: tracker:Tracker
0241     data.element = m_tracker = addSubdetector(cms::detElementName(volume), volume, false);
0242     m_tracker->SetTitle("compound");
0243   } else if (m_stack.size() == 3) {  // Main subsystem detector: TIB, TEC, ....
0244     data.element = m_currentDetector = addSubdetector(cms::detElementName(volume), volume, true);
0245   }
0246 }
0247 
0248 dd4hep::DetElement DDCMSDetElementCreator::addSubdetector(const std::string& nam,
0249                                                           dd4hep::PlacedVolume volume,
0250                                                           bool valid) {
0251   auto idet = m_subdetectors.find(nam);
0252   if (idet == m_subdetectors.end()) {
0253     dd4hep::DetElement det(nam, m_subdetectors.size() + 1);
0254     det.setPlacement(volume);
0255     if (valid) {
0256       det.placement().addPhysVolID("system", det.id());
0257     }
0258     idet = m_subdetectors.insert(make_pair(nam, det)).first;
0259     m_description.add(det);
0260   }
0261   return idet->second;
0262 }
0263 
0264 int DDCMSDetElementCreator::operator()(dd4hep::PlacedVolume volume, int volumeLevel) {
0265   double fracSi = volume.volume().material().fraction(m_silicon);
0266   if (fracSi > 90e-2) {
0267     Data& data = m_stack.back();
0268     data.sensitive = true;
0269     data.hasSensitive = true;
0270     ++data.volumeCount;
0271     int idx = volume->GetMotherVolume()->GetIndex(volume.ptr()) + 1;
0272     auto& cnt = m_leafCount[make_pair(m_currentDetector, volumeLevel)];
0273     cnt.first = std::max(cnt.first, idx);
0274     ++cnt.second;
0275     m_allPlacements[volume] = make_pair(volumeLevel, idx);
0276     return 1;
0277   }
0278   return 0;
0279 }
0280 
0281 int DDCMSDetElementCreator::process(dd4hep::PlacedVolume volume, int level, bool recursive) {
0282   m_stack.push_back(Data(volume));
0283   if (m_stack.size() <= 3) {
0284     createTopLevelDetectors(volume);
0285   }
0286   int ret = dd4hep::PlacedVolumeProcessor::process(volume, level, recursive);
0287 
0288   /// Complete structures if the m_stack size is > 3!
0289   if (m_stack.size() > 3) {
0290     // Note: short-cuts to entries in the m_stack MUST be local and
0291     // initialized AFTER the call to "process"! The vector may be resized!
0292     auto& data = m_stack.back();
0293     auto& parent = m_stack[m_stack.size() - 2];
0294     auto& counts = m_counters[m_currentDetector];
0295     if (data.sensitive) {
0296       /// If this volume is sensitve, we must attach a sensitive detector handle
0297       if (!m_currentSensitive.isValid()) {
0298         dd4hep::SensitiveDetector sd = m_description.sensitiveDetector(m_currentDetector.name());
0299         if (!sd.isValid()) {
0300           sd = dd4hep::SensitiveDetector(m_currentDetector.name(), "tracker");
0301           m_currentDetector->flag |= DetElement::Object::HAVE_SENSITIVE_DETECTOR;
0302           m_description.add(sd);
0303         }
0304         m_currentSensitive = sd;
0305       }
0306       volume.volume().setSensitiveDetector(m_currentSensitive);
0307       ++counts.sensitives;
0308     }
0309     ++counts.volumes;
0310     bool added = false;
0311     if (data.volumeCount > 0) {
0312       parent.daughterCount += data.volumeCount;
0313       parent.daughterCount += data.daughterCount;
0314       data.hasSensitive = true;
0315     } else {
0316       parent.daughterCount += data.daughterCount;
0317       data.hasSensitive = (data.daughterCount > 0);
0318     }
0319 
0320     if (data.hasSensitive) {
0321       // If we have sensitive elements at this level or below,
0322       // we must complete the DetElement hierarchy
0323       if (!data.element.isValid()) {
0324         data.element = createElement("Element", data.volume, m_currentDetector.id());
0325         ++counts.elements;
0326       }
0327       if (!parent.element.isValid()) {
0328         parent.element = createElement("Parent ", parent.volume, m_currentDetector.id());
0329         ++counts.elements;
0330       }
0331       printout(DEBUG,
0332                "DDCMSDetElementCreator",
0333                "++ Assign detector element: %s (%p, %ld children) to %s (%p) with %ld vols",
0334                data.element.name(),
0335                data.element.ptr(),
0336                data.element.children().size(),
0337                parent.element.name(),
0338                parent.element.ptr(),
0339                data.volumeCount);
0340 
0341       // Trickle up the tree only for sensitive pathes. Forget the passive rest
0342       // This should automatically omit non-sensitive pathes
0343       parent.hasSensitive = true;
0344       parent.element.add(data.element);
0345       added = true;
0346       // It is simpler to collect the volumes and later assign the volids
0347       // rather than checking if the volid already exists.
0348       int volumeLevel = level;
0349       int idx = data.volume->GetMotherVolume()->GetIndex(data.volume.ptr()) + 1;
0350       m_allPlacements[data.volume] = make_pair(volumeLevel, idx);  // 1...n
0351       // Update counters
0352       auto& cnt_det = m_leafCount[make_pair(m_currentDetector, volumeLevel)];
0353       cnt_det.first = std::max(cnt_det.first, idx);
0354       cnt_det.second += 1;
0355     }
0356     if (!added && data.element.isValid()) {
0357       printout(WARNING,
0358                "MEMORY-LEAK",
0359                "Level:%3d Orpahaned DetElement:%s Daugthers:%d Parent:%s",
0360                int(m_stack.size()),
0361                data.element.name(),
0362                data.volumeCount,
0363                parent.volume.name());
0364     }
0365   }
0366   /// Now the cleanup kicks in....
0367   if (m_stack.size() == 3) {
0368     m_currentSensitive = SensitiveDetector();
0369     m_currentDetector = DetElement();
0370     ret = 0;
0371   }
0372   m_stack.pop_back();
0373   return ret;
0374 }
0375 
0376 static void* createObject(dd4hep::Detector& description, int /* argc */, char** /* argv */) {
0377   dd4hep::PlacedVolumeProcessor* proc = new DDCMSDetElementCreator(description);
0378   return (void*)proc;
0379 }
0380 
0381 // first argument is the type from the xml file
0382 DECLARE_DD4HEP_CONSTRUCTOR(DDCMS_DetElementCreator, createObject);