plugins/dd4hep/DDCMSDetElementCreator.cc

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, totalDepth = 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     totalDepth += l.second.second;
0137     ++totalCount;
0138   }
0139   printout(INFO, "DDCMSDetElementCreator", "++ Summary:     %-24s  %d.", "Total DetElements:", totalCount);
0140   printout(INFO, "DDCMSDetElementCreator", "+++++++++++++++ Readout structure generation  ++++++++++++++++++++++++");
0141   str << endl;
0142   for (const auto& f : fields) {
0143     string roName = f.first.name() + string("Hits");
0144     int num_bits = 8;
0145     id_str.str("");
0146     id_str << "system:" << num_bits;
0147     for (const auto& q : f.second) {
0148       int bits = 0;
0149       if (q.second < 1 << 0)
0150         bits = 1;
0151       else if (q.second < 1 << 1)
0152         bits = 1;
0153       else if (q.second < 1 << 2)
0154         bits = 2;
0155       else if (q.second < 1 << 3)
0156         bits = 3;
0157       else if (q.second < 1 << 4)
0158         bits = 4;
0159       else if (q.second < 1 << 5)
0160         bits = 5;
0161       else if (q.second < 1 << 6)
0162         bits = 6;
0163       else if (q.second < 1 << 7)
0164         bits = 7;
0165       else if (q.second < 1 << 8)
0166         bits = 8;
0167       else if (q.second < 1 << 9)
0168         bits = 9;
0169       else if (q.second < 1 << 10)
0170         bits = 10;
0171       else if (q.second < 1 << 11)
0172         bits = 11;
0173       else if (q.second < 1 << 12)
0174         bits = 12;
0175       else if (q.second < 1 << 13)
0176         bits = 13;
0177       else if (q.second < 1 << 14)
0178         bits = 14;
0179       else if (q.second < 1 << 15)
0180         bits = 15;
0181       bits += 1;
0182       id_str << ",Lv" << q.first << ":" << bits;
0183       num_bits += bits;
0184     }
0185     string idspec = id_str.str();
0186     str << "<readout name=\"" << roName << "\">" << endl
0187         << "\t<id>" << idspec << "</id>  <!-- Number of bits: " << num_bits << " -->" << endl
0188         << "</readout>" << endl;
0189
0190     /// Create ID Descriptors and readout configurations
0191     IDDescriptor dsc(roName, idspec);
0192     m_description.addIDSpecification(dsc);
0193     Readout ro(roName);
0194     ro.setIDDescriptor(dsc);
0195     m_description.addReadout(ro);
0196     dd4hep::SensitiveDetector sd = m_description.sensitiveDetector(f.first.name());
0197     sd.setHitsCollection(ro.name());
0198     sd.setReadout(ro);
0199     printout(INFO,
0200              "DDCMSDetElementCreator",
0201              "++ Setting up readout for subdetector:%-24s id:%04X",
0202              f.first.name(),
0203              f.first.id());
0204   }
0205   printout(INFO, "DDCMSDetElementCreator", "+++++++++++++++ ID Descriptor generation  ++++++++++++++++++++++++++++");
0206   printout(INFO, "", str.str().c_str());
0207   char volId[32];
0208   for (auto& p : m_allPlacements) {
0209     dd4hep::PlacedVolume place = p.first;
0210     dd4hep::Volume volume = place.volume();
0211     ::snprintf(volId, sizeof(volId), "Lv%d", p.second.first);
0212     printout(DEBUG,
0213              "DDCMSDetElementCreator",
0214              "++ Set volid (%-24s): %-6s = %3d  -> %s  (%p)",
0215              volume.isSensitive() ? volume.sensitiveDetector().name() : "Not Sensitive",
0216              volId,
0217              p.second.second,
0218              place.name(),
0219              place.ptr());
0220     place.addPhysVolID(volId, p.second.second);
0221   }
0222   printout(ALWAYS,
0223            "DDCMSDetElementCreator",
0224            "++ Instrumented %ld subdetectors with %d DetElements %d sensitives out of %d volumes and %ld sensitive "
0225            "placements.",
0226            fields.size(),
0227            total.elements,
0228            total.sensitives,
0229            total.volumes,
0230            m_allPlacements.size());
0231 }
0232
0233 dd4hep::DetElement DDCMSDetElementCreator::createElement(const char*, PlacedVolume volume, int id) {
0234   string name = detElementName(volume);
0235   dd4hep::DetElement det(name, id);
0236   det.setPlacement(volume);
0237   return det;
0238 }
0239 void DDCMSDetElementCreator::createTopLevelDetectors(PlacedVolume volume) {
0240   auto& data = m_stack.back();
0241   if (m_stack.size() == 2) {  // Main subssystem: tracker:Tracker
0242     data.element = m_tracker = addSubdetector(cms::detElementName(volume), volume, false);
0243     m_tracker->SetTitle("compound");
0244   } else if (m_stack.size() == 3) {  // Main subsystem detector: TIB, TEC, ....
0245     data.element = m_currentDetector = addSubdetector(cms::detElementName(volume), volume, true);
0246   }
0247 }
0248
0249 dd4hep::DetElement DDCMSDetElementCreator::addSubdetector(const std::string& nam,
0250                                                           dd4hep::PlacedVolume volume,
0251                                                           bool valid) {
0252   auto idet = m_subdetectors.find(nam);
0253   if (idet == m_subdetectors.end()) {
0254     dd4hep::DetElement det(nam, m_subdetectors.size() + 1);
0255     det.setPlacement(volume);
0256     if (valid) {
0257       det.placement().addPhysVolID("system", det.id());
0258     }
0259     idet = m_subdetectors.insert(make_pair(nam, det)).first;
0260     m_description.add(det);
0261   }
0262   return idet->second;
0263 }
0264
0265 int DDCMSDetElementCreator::operator()(dd4hep::PlacedVolume volume, int volumeLevel) {
0266   double fracSi = volume.volume().material().fraction(m_silicon);
0267   if (fracSi > 90e-2) {
0268     Data& data = m_stack.back();
0269     data.sensitive = true;
0270     data.hasSensitive = true;
0271     ++data.volumeCount;
0272     int idx = volume->GetMotherVolume()->GetIndex(volume.ptr()) + 1;
0273     auto& cnt = m_leafCount[make_pair(m_currentDetector, volumeLevel)];
0274     cnt.first = std::max(cnt.first, idx);
0275     ++cnt.second;
0276     m_allPlacements[volume] = make_pair(volumeLevel, idx);
0277     return 1;
0278   }
0279   return 0;
0280 }
0281
0282 int DDCMSDetElementCreator::process(dd4hep::PlacedVolume volume, int level, bool recursive) {
0283   m_stack.push_back(Data(volume));
0284   if (m_stack.size() <= 3) {
0285     createTopLevelDetectors(volume);
0286   }
0287   int ret = dd4hep::PlacedVolumeProcessor::process(volume, level, recursive);
0288
0289   /// Complete structures if the m_stack size is > 3!
0290   if (m_stack.size() > 3) {
0291     // Note: short-cuts to entries in the m_stack MUST be local and
0292     // initialized AFTER the call to "process"! The vector may be resized!
0293     auto& data = m_stack.back();
0294     auto& parent = m_stack[m_stack.size() - 2];
0295     auto& counts = m_counters[m_currentDetector];
0296     if (data.sensitive) {
0297       /// If this volume is sensitve, we must attach a sensitive detector handle
0298       if (!m_currentSensitive.isValid()) {
0299         dd4hep::SensitiveDetector sd = m_description.sensitiveDetector(m_currentDetector.name());
0300         if (!sd.isValid()) {
0301           sd = dd4hep::SensitiveDetector(m_currentDetector.name(), "tracker");
0302           m_currentDetector->flag |= DetElement::Object::HAVE_SENSITIVE_DETECTOR;
0303           m_description.add(sd);
0304         }
0305         m_currentSensitive = sd;
0306       }
0307       volume.volume().setSensitiveDetector(m_currentSensitive);
0308       ++counts.sensitives;
0309     }
0310     ++counts.volumes;
0311     bool added = false;
0312     if (data.volumeCount > 0) {
0313       parent.daughterCount += data.volumeCount;
0314       parent.daughterCount += data.daughterCount;
0315       data.hasSensitive = true;
0316     } else {
0317       parent.daughterCount += data.daughterCount;
0318       data.hasSensitive = (data.daughterCount > 0);
0319     }
0320
0321     if (data.hasSensitive) {
0322       // If we have sensitive elements at this level or below,
0323       // we must complete the DetElement hierarchy
0324       if (!data.element.isValid()) {
0325         data.element = createElement("Element", data.volume, m_currentDetector.id());
0326         ++counts.elements;
0327       }
0328       if (!parent.element.isValid()) {
0329         parent.element = createElement("Parent ", parent.volume, m_currentDetector.id());
0330         ++counts.elements;
0331       }
0332       printout(DEBUG,
0333                "DDCMSDetElementCreator",
0334                "++ Assign detector element: %s (%p, %ld children) to %s (%p) with %ld vols",
0335                data.element.name(),
0336                data.element.ptr(),
0337                data.element.children().size(),
0338                parent.element.name(),
0339                parent.element.ptr(),
0340                data.volumeCount);
0341
0342       // Trickle up the tree only for sensitive pathes. Forget the passive rest
0343       // This should automatically omit non-sensitive pathes
0344       parent.hasSensitive = true;
0345       parent.element.add(data.element);
0346       added = true;
0347       // It is simpler to collect the volumes and later assign the volids
0348       // rather than checking if the volid already exists.
0349       int volumeLevel = level;
0350       int idx = data.volume->GetMotherVolume()->GetIndex(data.volume.ptr()) + 1;
0351       m_allPlacements[data.volume] = make_pair(volumeLevel, idx);  // 1...n
0352       // Update counters
0353       auto& cnt_det = m_leafCount[make_pair(m_currentDetector, volumeLevel)];
0354       cnt_det.first = std::max(cnt_det.first, idx);
0355       cnt_det.second += 1;
0356     }
0357     if (!added && data.element.isValid()) {
0358       printout(WARNING,
0359                "MEMORY-LEAK",
0360                "Level:%3d Orpahaned DetElement:%s Daugthers:%d Parent:%s",
0361                int(m_stack.size()),
0362                data.element.name(),
0363                data.volumeCount,
0364                parent.volume.name());
0365     }
0366   }
0367   /// Now the cleanup kicks in....
0368   if (m_stack.size() == 3) {
0369     m_currentSensitive = SensitiveDetector();
0370     m_currentDetector = DetElement();
0371     ret = 0;
0372   }
0373   m_stack.pop_back();
0374   return ret;
0375 }
0376
0377 static void* createObject(dd4hep::Detector& description, int /* argc */, char** /* argv */) {
0378   dd4hep::PlacedVolumeProcessor* proc = new DDCMSDetElementCreator(description);
0379   return (void*)proc;
0380 }
0381
0382 // first argument is the type from the xml file
0383 DECLARE_DD4HEP_CONSTRUCTOR(DDCMS_DetElementCreator, createObject);