Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-06-13 03:23:44

0001 #include <algorithm>
0002 #include <cassert>
0003 #include <vector>
0004 
0005 #include "FWCore/Sources/interface/DaqProvenanceHelper.h"
0006 
0007 #include "DataFormats/Provenance/interface/BranchChildren.h"
0008 #include "DataFormats/Provenance/interface/BranchIDList.h"
0009 #include "DataFormats/Provenance/interface/ProcessHistory.h"
0010 #include "DataFormats/Provenance/interface/ProcessHistoryRegistry.h"
0011 #include "DataFormats/Provenance/interface/ProductRegistry.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0013 #include "FWCore/Utilities/interface/EDMException.h"
0014 #include "FWCore/Utilities/interface/GetPassID.h"
0015 #include "FWCore/Reflection/interface/TypeWithDict.h"
0016 #include "FWCore/Version/interface/GetReleaseVersion.h"
0017 
0018 namespace {
0019   edm::BranchDescription makeDescriptionForDaqProvHelper(edm::TypeID const& rawDataType,
0020                                                          std::string const& collectionName,
0021                                                          std::string const& friendlyName,
0022                                                          std::string const& sourceLabel) {
0023     edm::BranchDescription desc(edm::InEvent,
0024                                 "rawDataCollector",
0025                                 // "source",
0026                                 "LHC",
0027                                 // "HLT",
0028                                 collectionName,
0029                                 friendlyName,
0030                                 "",
0031                                 sourceLabel,
0032                                 edm::ParameterSetID(),
0033                                 edm::TypeWithDict(rawDataType.typeInfo()),
0034                                 false);
0035     desc.setIsProvenanceSetOnRead();
0036     return desc;
0037   }
0038 }  // namespace
0039 
0040 namespace edm {
0041   DaqProvenanceHelper::DaqProvenanceHelper(TypeID const& rawDataType,
0042                                            std::string const& collectionName,
0043                                            std::string const& friendlyName,
0044                                            std::string const& sourceLabel)
0045       : constBranchDescription_(
0046             makeDescriptionForDaqProvHelper(rawDataType, collectionName, friendlyName, sourceLabel)),
0047         dummyProvenance_(constBranchDescription_.branchID()),
0048         processParameterSet_(),
0049         oldProcessName_(),
0050         oldBranchID_(),
0051         newBranchID_(),
0052         oldProcessHistoryID_(nullptr),
0053         phidMap_() {
0054     // Now we create a process parameter set for the "LHC" process.
0055     // We don't currently use the untracked parameters, However, we make them available, just in case.
0056     std::string const& moduleLabel = constBranchDescription_.moduleLabel();
0057     std::string const& processName = constBranchDescription_.processName();
0058     std::string const& moduleName = constBranchDescription_.moduleName();
0059     typedef std::vector<std::string> vstring;
0060     vstring empty;
0061 
0062     vstring modlbl;
0063     modlbl.reserve(1);
0064     modlbl.push_back(moduleLabel);
0065     processParameterSet_.addParameter("@all_sources", modlbl);
0066 
0067     ParameterSet triggerPaths;
0068     triggerPaths.addParameter<vstring>("@trigger_paths", empty);
0069     processParameterSet_.addParameter<ParameterSet>("@trigger_paths", triggerPaths);
0070 
0071     ParameterSet pseudoInput;
0072     pseudoInput.addParameter<std::string>("@module_edm_type", "Source");
0073     pseudoInput.addParameter<std::string>("@module_label", moduleLabel);
0074     pseudoInput.addParameter<std::string>("@module_type", moduleName);
0075     processParameterSet_.addParameter<ParameterSet>(moduleLabel, pseudoInput);
0076 
0077     processParameterSet_.addParameter<vstring>("@all_esmodules", empty);
0078     processParameterSet_.addParameter<vstring>("@all_esprefers", empty);
0079     processParameterSet_.addParameter<vstring>("@all_essources", empty);
0080     processParameterSet_.addParameter<vstring>("@all_loopers", empty);
0081     processParameterSet_.addParameter<vstring>("@all_modules", empty);
0082     processParameterSet_.addParameter<vstring>("@end_paths", empty);
0083     processParameterSet_.addParameter<vstring>("@paths", empty);
0084     processParameterSet_.addParameter<std::string>("@process_name", processName);
0085     // Now we register the process parameter set.
0086     processParameterSet_.registerIt();
0087 
0088     //std::cerr << processParameterSet_.dump() << std::endl;
0089   }
0090 
0091   //default
0092   DaqProvenanceHelper::DaqProvenanceHelper(TypeID const& rawDataType)
0093       : DaqProvenanceHelper(rawDataType, "FEDRawDataCollection", "FEDRawDataCollection", "FedRawDataInputSource") {}
0094 
0095   ProcessHistoryID DaqProvenanceHelper::daqInit(ProductRegistry& productRegistry,
0096                                                 ProcessHistoryRegistry& processHistoryRegistry) const {
0097     // Now we need to set all the metadata
0098     // Add the product to the product registry
0099     productRegistry.copyProduct(constBranchDescription_);
0100 
0101     // Insert an entry for this process in the process history registry
0102     ProcessHistory ph;
0103     ph.emplace_back(constBranchDescription_.processName(), processParameterSet_.id(), getReleaseVersion(), getPassID());
0104     processHistoryRegistry.registerProcessHistory(ph);
0105 
0106     // Save the process history ID for use every event.
0107     return ph.setProcessHistoryID();
0108   }
0109 
0110   bool DaqProvenanceHelper::matchProcesses(ProcessConfiguration const& newPC, ProcessHistory const& ph) const {
0111     for (auto const& pc : ph) {
0112       if (pc.processName() == oldProcessName_) {
0113         return (pc.releaseVersion() == newPC.releaseVersion() && pc.passID() == newPC.passID());
0114       }
0115     }
0116     return false;
0117   }
0118 
0119   void DaqProvenanceHelper::fixMetaData(std::vector<ProcessConfiguration>& pcv, std::vector<ProcessHistory>& phv) {
0120     phv.push_back(ProcessHistory());  // For new processHistory, containing only processConfiguration_
0121     std::vector<ProcessConfiguration> newPCs;
0122     for (auto const& pc : pcv) {
0123       if (pc.processName() == oldProcessName_) {
0124         newPCs.emplace_back(
0125             constBranchDescription_.processName(), processParameterSet_.id(), pc.releaseVersion(), pc.passID());
0126       }
0127     }
0128     if (newPCs.empty()) {
0129       throw Exception(errors::LogicError) << "\nFatal error in RootFile constructor. Most likely this is because\n"
0130                                           << "the input file contains a FEDRawDataCollection with module label\n"
0131                                           << "\"source\". This is against CMS naming conventions.\n"
0132                                           << "See GitHub Issue 45137 for related details.\n";
0133     }
0134     pcv.reserve(pcv.size() + newPCs.size());
0135     pcv.insert(pcv.end(), newPCs.begin(), newPCs.end());
0136     // update existing process histories
0137     for (auto& ph : phv) {
0138       for (auto const& newPC : newPCs) {
0139         if (ph.empty() || matchProcesses(newPC, ph)) {
0140           ProcessHistoryID oldPHID = ph.id();
0141           ph.push_front(newPC);
0142           ProcessHistoryID newPHID = ph.id();
0143           phidMap_.insert(std::make_pair(oldPHID, newPHID));
0144           break;
0145         }
0146       }
0147     }
0148     // For new process histories, containing only the new process configurations
0149     phv.reserve(phv.size() + newPCs.size());
0150     for (auto const& newPC : newPCs) {
0151       phv.emplace_back();
0152       phv.back().push_front(newPC);
0153     }
0154   }
0155 
0156   void DaqProvenanceHelper::fixMetaData(std::vector<BranchID>& branchID) const {
0157     std::replace(branchID.begin(), branchID.end(), oldBranchID_, newBranchID_);
0158   }
0159 
0160   void DaqProvenanceHelper::fixMetaData(BranchIDLists const& branchIDLists) const {
0161     BranchID::value_type oldID = oldBranchID_.id();
0162     BranchID::value_type newID = newBranchID_.id();
0163     // The const_cast is ugly, but it beats the alternatives.
0164     BranchIDLists& lists = const_cast<BranchIDLists&>(branchIDLists);
0165     for (auto& list : lists) {
0166       std::replace(list.begin(), list.end(), oldID, newID);
0167     }
0168   }
0169 
0170   void DaqProvenanceHelper::fixMetaData(BranchChildren& branchChildren) const {
0171     typedef std::map<BranchID, std::set<BranchID> > BCMap;
0172     // The const_cast is ugly, but it beats the alternatives.
0173     BCMap& childLookup = const_cast<BCMap&>(branchChildren.childLookup());
0174     // First fix any old branchID's in the key.
0175     {
0176       BCMap::iterator i = childLookup.find(oldBranchID_);
0177       if (i != childLookup.end()) {
0178         childLookup.insert(std::make_pair(newBranchID_, i->second));
0179         childLookup.erase(i);
0180       }
0181     }
0182     // Now fix any old branchID's in the sets;
0183     for (auto& child : childLookup) {
0184       if (child.second.erase(oldBranchID_) != 0) {
0185         child.second.insert(newBranchID_);
0186       }
0187     }
0188   }
0189 
0190   // Replace process history ID.
0191   ProcessHistoryID const& DaqProvenanceHelper::mapProcessHistoryID(ProcessHistoryID const& phid) {
0192     ProcessHistoryIDMap::const_iterator it = phidMap_.find(phid);
0193     assert(it != phidMap_.end());
0194     oldProcessHistoryID_ = &it->first;
0195     return it->second;
0196   }
0197 
0198   // Replace parentage ID.
0199   ParentageID const& DaqProvenanceHelper::mapParentageID(ParentageID const& parentageID) const {
0200     ParentageIDMap::const_iterator it = parentageIDMap_.find(parentageID);
0201     if (it == parentageIDMap_.end()) {
0202       return parentageID;
0203     }
0204     return it->second;
0205   }
0206 
0207   // Replace branch ID if necessary.
0208   BranchID const& DaqProvenanceHelper::mapBranchID(BranchID const& branchID) const {
0209     return (branchID == oldBranchID_ ? newBranchID_ : branchID);
0210   }
0211 
0212   void DaqProvenanceHelper::setOldParentageIDToNew(ParentageID const& iOld, ParentageID const& iNew) {
0213     parentageIDMap_.insert(std::make_pair(iOld, iNew));
0214   }
0215 
0216 }  // namespace edm