Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:29:16

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