Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-01-17 23:31:44

0001 // system includes
0002 #include <map>
0003 #include <vector>
0004 #include <algorithm>
0005 
0006 // user includes
0007 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0008 #include "CommonTools/Utils/interface/TFileDirectory.h"
0009 #include "DataFormats/Common/interface/DetSetVector.h"
0010 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0011 #include "DataFormats/Common/interface/Handle.h"
0012 #include "DataFormats/DetId/interface/DetId.h"
0013 #include "DataFormats/Phase2TrackerCluster/interface/Phase2TrackerCluster1D.h"
0014 #include "DataFormats/Phase2TrackerDigi/interface/Phase2TrackerDigi.h"
0015 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0016 #include "DataFormats/SiPixelDigi/interface/PixelDigi.h"
0017 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0018 #include "FWCore/Framework/interface/ConsumesCollector.h"
0019 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0020 #include "FWCore/Framework/interface/ESHandle.h"
0021 #include "FWCore/Framework/interface/Event.h"
0022 #include "FWCore/Framework/interface/Frameworkfwd.h"
0023 #include "FWCore/Framework/interface/MakerMacros.h"
0024 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0025 #include "FWCore/ServiceRegistry/interface/Service.h"
0026 #include "FWCore/ServiceRegistry/interface/Service.h"
0027 #include "FWCore/Utilities/interface/InputTag.h"
0028 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0029 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0030 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0031 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0032 #include "SimDataFormats/Track/interface/SimTrackContainer.h"
0033 #include "SimDataFormats/TrackerDigiSimLink/interface/PixelDigiSimLink.h"
0034 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
0035 #include "SimDataFormats/Vertex/interface/SimVertexContainer.h"
0036 
0037 // ROOT includes
0038 #include <TH2F.h>
0039 #include <TH1F.h>
0040 #include <THStack.h>
0041 
0042 struct ClusterHistos {
0043   // use TH1D instead of TH1F to avoid stauration at 2^31
0044   // above this increments with +1 don't work for float, need double
0045 
0046   TH1D* numberClusters[3];
0047   TH1D* clusterSize[3];
0048 
0049   TH2F* globalPosXY[3];
0050   TH2F* localPosXY[3];
0051 
0052   TH1F* deltaX[3];
0053   TH1F* deltaY[3];
0054   TH1F* deltaX_P[3];
0055   TH1F* deltaY_P[3];
0056 
0057   TH1D* primarySimHits[3];
0058   TH1D* otherSimHits[3];
0059 };
0060 
0061 class Phase2TrackerClusterizerValidation : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0062 public:
0063   typedef std::map<unsigned int, std::vector<PSimHit> > SimHitsMap;
0064   typedef std::map<unsigned int, SimTrack> SimTracksMap;
0065 
0066   explicit Phase2TrackerClusterizerValidation(const edm::ParameterSet&);
0067   ~Phase2TrackerClusterizerValidation();
0068   void beginJob() override;
0069   void analyze(const edm::Event&, const edm::EventSetup&) override;
0070 
0071 private:
0072   std::map<unsigned int, ClusterHistos>::iterator createLayerHistograms(unsigned int);
0073   std::vector<unsigned int> getSimTrackId(const edm::Handle<edm::DetSetVector<PixelDigiSimLink> >&,
0074                                           const DetId&,
0075                                           unsigned int);
0076 
0077   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> esTokenGeom_;
0078   const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> esTokenTopo_;
0079   const edm::EDGetTokenT<Phase2TrackerCluster1DCollectionNew> tokenClusters_;
0080   const edm::EDGetTokenT<edm::DetSetVector<PixelDigiSimLink> > tokenLinks_;
0081   const edm::EDGetTokenT<edm::PSimHitContainer> tokenSimHitsB_;
0082   const edm::EDGetTokenT<edm::PSimHitContainer> tokenSimHitsE_;
0083   const edm::EDGetTokenT<edm::SimTrackContainer> tokenSimTracks_;
0084 
0085   const bool catECasRings_;
0086   const double simtrackminpt_;
0087 
0088   TH2F* trackerLayout_;
0089   TH2F* trackerLayoutXY_;
0090   TH2F* trackerLayoutXYBar_;
0091   TH2F* trackerLayoutXYEC_;
0092 
0093   std::map<unsigned int, ClusterHistos> histograms_;
0094 };
0095 
0096 Phase2TrackerClusterizerValidation::Phase2TrackerClusterizerValidation(const edm::ParameterSet& conf)
0097     : esTokenGeom_(esConsumes()),
0098       esTokenTopo_(esConsumes()),
0099       tokenClusters_(consumes<Phase2TrackerCluster1DCollectionNew>(conf.getParameter<edm::InputTag>("src"))),
0100       tokenLinks_(consumes<edm::DetSetVector<PixelDigiSimLink> >(conf.getParameter<edm::InputTag>("links"))),
0101       tokenSimHitsB_(consumes<edm::PSimHitContainer>(conf.getParameter<edm::InputTag>("simhitsbarrel"))),
0102       tokenSimHitsE_(consumes<edm::PSimHitContainer>(conf.getParameter<edm::InputTag>("simhitsendcap"))),
0103       tokenSimTracks_(consumes<edm::SimTrackContainer>(conf.getParameter<edm::InputTag>("simtracks"))),
0104       catECasRings_(conf.getParameter<bool>("ECasRings")),
0105       simtrackminpt_(conf.getParameter<double>("SimTrackMinPt")) {
0106   usesResource(TFileService::kSharedResource);
0107 }
0108 
0109 Phase2TrackerClusterizerValidation::~Phase2TrackerClusterizerValidation() = default;
0110 
0111 void Phase2TrackerClusterizerValidation::beginJob() {
0112   edm::Service<TFileService> fs;
0113   fs->file().cd("/");
0114   TFileDirectory td = fs->mkdir("Common");
0115   // Create common histograms
0116   trackerLayout_ = td.make<TH2F>("RVsZ", "R vs. z position", 6000, -300.0, 300.0, 1200, 0.0, 120.0);
0117   trackerLayoutXY_ = td.make<TH2F>("XVsY", "x vs. y position", 2400, -120.0, 120.0, 2400, -120.0, 120.0);
0118   trackerLayoutXYBar_ = td.make<TH2F>("XVsYBar", "x vs. y position", 2400, -120.0, 120.0, 2400, -120.0, 120.0);
0119   trackerLayoutXYEC_ = td.make<TH2F>("XVsYEC", "x vs. y position", 2400, -120.0, 120.0, 2400, -120.0, 120.0);
0120 }
0121 
0122 void Phase2TrackerClusterizerValidation::analyze(const edm::Event& event, const edm::EventSetup& eventSetup) {
0123   /*
0124      * Get the needed objects
0125      */
0126 
0127   // Get the clusters
0128   edm::Handle<Phase2TrackerCluster1DCollectionNew> clusters;
0129   event.getByToken(tokenClusters_, clusters);
0130 
0131   // Get the PixelDigiSimLinks
0132   edm::Handle<edm::DetSetVector<PixelDigiSimLink> > pixelSimLinks;
0133   event.getByToken(tokenLinks_, pixelSimLinks);
0134 
0135   // Get the SimHits
0136   edm::Handle<edm::PSimHitContainer> simHitsRaw[2];
0137   event.getByToken(tokenSimHitsB_, simHitsRaw[0]);
0138   event.getByToken(tokenSimHitsE_, simHitsRaw[1]);
0139 
0140   // Get the SimTracks
0141   edm::Handle<edm::SimTrackContainer> simTracksRaw;
0142   event.getByToken(tokenSimTracks_, simTracksRaw);
0143 
0144   // Get the geometry
0145   const TrackerGeometry* tkGeom = &eventSetup.getData(esTokenGeom_);
0146   const TrackerTopology* tTopo = &eventSetup.getData(esTokenTopo_);
0147 
0148   /*
0149      * Rearrange the simTracks
0150      */
0151 
0152   // Rearrange the simTracks for ease of use <simTrackID, simTrack>
0153   SimTracksMap simTracks;
0154   for (edm::SimTrackContainer::const_iterator simTrackIt(simTracksRaw->begin()); simTrackIt != simTracksRaw->end();
0155        ++simTrackIt) {
0156     if (simTrackIt->momentum().pt() > simtrackminpt_) {
0157       simTracks.emplace(simTrackIt->trackId(), *simTrackIt);
0158     }
0159   }
0160 
0161   /*
0162      * Validation   
0163      */
0164 
0165   // Number of clusters
0166   std::map<unsigned int, unsigned int> nClusters[3];
0167   std::map<unsigned int, unsigned int> nPrimarySimHits[3];
0168   std::map<unsigned int, unsigned int> nOtherSimHits[3];
0169 
0170   // Loop over modules
0171   for (Phase2TrackerCluster1DCollectionNew::const_iterator DSViter = clusters->begin(); DSViter != clusters->end();
0172        ++DSViter) {
0173     // Get the detector unit's id
0174     unsigned int rawid(DSViter->detId());
0175     DetId detId(rawid);
0176     unsigned int layer = (tTopo->side(detId) != 0) * 1000;  // don't split up endcap sides
0177     if (!layer) {
0178       layer += tTopo->layer(detId);
0179     } else {
0180       layer += (catECasRings_ ? tTopo->tidRing(detId) * 10 : tTopo->layer(detId));
0181     }
0182     TrackerGeometry::ModuleType mType = tkGeom->getDetectorType(detId);
0183     unsigned int det = 0;
0184     if (mType == TrackerGeometry::ModuleType::Ph2PSP) {
0185       det = 1;
0186     } else if (mType == TrackerGeometry::ModuleType::Ph2PSS || mType == TrackerGeometry::ModuleType::Ph2SS) {
0187       det = 2;
0188     } else {
0189       std::cout << "UNKNOWN DETECTOR TYPE!" << std::endl;
0190     }
0191 
0192     // Get the geomdet
0193     const GeomDetUnit* geomDetUnit(tkGeom->idToDetUnit(detId));
0194     if (!geomDetUnit)
0195       continue;
0196 
0197     // initialize the nhit counters if they don't exist for this layer
0198     auto nhitit(nClusters[det].find(layer));
0199     if (nhitit == nClusters[det].end()) {
0200       nClusters[det].emplace(layer, 0);
0201       nPrimarySimHits[det].emplace(layer, 0);
0202       nOtherSimHits[det].emplace(layer, 0);
0203     }
0204 
0205     // Create histograms for the layer if they do not yet exist
0206     std::map<unsigned int, ClusterHistos>::iterator histogramLayer(histograms_.find(layer));
0207     if (histogramLayer == histograms_.end())
0208       histogramLayer = createLayerHistograms(layer);
0209 
0210     // Loop over the clusters in the detector unit
0211     for (edmNew::DetSet<Phase2TrackerCluster1D>::const_iterator clustIt = DSViter->begin(); clustIt != DSViter->end();
0212          ++clustIt) {
0213       // determine the position
0214       MeasurementPoint mpClu(clustIt->center(), clustIt->column() + 0.5);
0215       Local3DPoint localPosClu = geomDetUnit->topology().localPosition(mpClu);
0216       Global3DPoint globalPosClu = geomDetUnit->surface().toGlobal(localPosClu);
0217 
0218       // Get all the simTracks that form the cluster
0219       std::vector<unsigned int> clusterSimTrackIds;
0220       for (unsigned int i(0); i < clustIt->size(); ++i) {
0221         unsigned int channel(Phase2TrackerDigi::pixelToChannel(clustIt->firstRow() + i, clustIt->column()));
0222         std::vector<unsigned int> simTrackIds(getSimTrackId(pixelSimLinks, detId, channel));
0223         for (auto it : simTrackIds) {
0224           bool add = true;
0225           for (unsigned int j = 0; j < clusterSimTrackIds.size(); ++j) {
0226             // only save simtrackids that are not present yet
0227             if (it == clusterSimTrackIds.at(j))
0228               add = false;
0229           }
0230           if (add)
0231             clusterSimTrackIds.push_back(it);
0232         }
0233       }
0234       std::sort(clusterSimTrackIds.begin(), clusterSimTrackIds.end());
0235 
0236       // find the closest simhit
0237       // this is needed because otherwise you get cases with simhits and clusters being swapped
0238       // when there are more than 1 cluster with common simtrackids
0239       const PSimHit* simhit = 0;  // bad naming to avoid changing code below. This is the closest simhit in x
0240       float minx = 10000;
0241       for (unsigned int simhitidx = 0; simhitidx < 2; ++simhitidx) {  // loop over both barrel and endcap hits
0242         for (auto simhitIt : *simHitsRaw[simhitidx]) {
0243           if (rawid == simhitIt.detUnitId()) {
0244             //std::cout << "=== " << rawid << " " << &simhitIt << " " << simhitIt.trackId() << " " << simhitIt.localPosition().x() << " " << simhitIt.localPosition().y() << std::endl;
0245             auto it = std::lower_bound(clusterSimTrackIds.begin(), clusterSimTrackIds.end(), simhitIt.trackId());
0246             if (it != clusterSimTrackIds.end() && *it == simhitIt.trackId()) {
0247               if (!simhit || fabs(simhitIt.localPosition().x() - localPosClu.x()) < minx) {
0248                 minx = fabs(simhitIt.localPosition().x() - localPosClu.x());
0249                 simhit = &simhitIt;
0250               }
0251             }
0252           }
0253         }
0254       }
0255       if (!simhit)
0256         continue;
0257 
0258       // only look at simhits from highpT tracks
0259       auto simTrackIt(simTracks.find(simhit->trackId()));
0260       if (simTrackIt == simTracks.end())
0261         continue;
0262 
0263       /*
0264              * Cluster related variables
0265              */
0266 
0267       // cluster size
0268       ++(nClusters[det].at(layer));
0269       ++(nOtherSimHits[det].at(layer));
0270 
0271       // cluster size
0272       histogramLayer->second.clusterSize[det]->Fill(clustIt->size());
0273 
0274       // Fill the position histograms
0275       trackerLayout_->Fill(globalPosClu.z(), globalPosClu.perp());
0276       trackerLayoutXY_->Fill(globalPosClu.x(), globalPosClu.y());
0277       if (layer < 1000) {
0278         trackerLayoutXYBar_->Fill(globalPosClu.x(), globalPosClu.y());
0279       } else {
0280         trackerLayoutXYEC_->Fill(globalPosClu.x(), globalPosClu.y());
0281       }
0282 
0283       histogramLayer->second.localPosXY[det]->Fill(localPosClu.x(), localPosClu.y());
0284       if (layer < 1000) {
0285         histogramLayer->second.globalPosXY[det]->Fill(globalPosClu.z(), globalPosClu.perp());
0286       } else {
0287         histogramLayer->second.globalPosXY[det]->Fill(globalPosClu.x(), globalPosClu.y());
0288       }
0289 
0290       // now get the position of the closest hit
0291       Local3DPoint localPosHit(simhit->localPosition());
0292 
0293       histogramLayer->second.deltaX[det]->Fill(localPosClu.x() - localPosHit.x());
0294       histogramLayer->second.deltaY[det]->Fill(localPosClu.y() - localPosHit.y());
0295 
0296       // Primary particles only
0297       unsigned int procT(simhit->processType());
0298       if (simTrackIt->second.vertIndex() == 0 and
0299           (procT == 2 || procT == 7 || procT == 9 || procT == 11 || procT == 13 || procT == 15)) {
0300         ++(nPrimarySimHits[det].at(layer));
0301         --(nOtherSimHits[det].at(layer));  // avoid double counting
0302         histogramLayer->second.deltaX_P[det]->Fill(localPosClu.x() - localPosHit.x());
0303         histogramLayer->second.deltaY_P[det]->Fill(localPosClu.y() - localPosHit.y());
0304       }
0305     }
0306   }
0307 
0308   // fill the counter histos per layer
0309   for (unsigned int det = 1; det < 3; ++det) {
0310     for (auto it : nClusters[det]) {
0311       auto histogramLayer(histograms_.find(it.first));
0312       if (histogramLayer == histograms_.end())
0313         std::cout << "*** SL *** No histogram for an existing counter! This should not happen!" << std::endl;
0314       histogramLayer->second.numberClusters[det]->Fill(it.second);
0315     }
0316     for (auto it : nPrimarySimHits[det]) {
0317       auto histogramLayer(histograms_.find(it.first));
0318       if (histogramLayer == histograms_.end())
0319         std::cout << "*** SL *** No histogram for an existing counter! This should not happen!" << std::endl;
0320       histogramLayer->second.primarySimHits[det]->Fill(it.second);
0321     }
0322     for (auto it : nOtherSimHits[det]) {
0323       auto histogramLayer(histograms_.find(it.first));
0324       if (histogramLayer == histograms_.end())
0325         std::cout << "*** SL *** No histogram for an existing counter! This should not happen!" << std::endl;
0326       histogramLayer->second.otherSimHits[det]->Fill(it.second);
0327     }
0328   }
0329 }
0330 
0331 // Create the histograms
0332 std::map<unsigned int, ClusterHistos>::iterator Phase2TrackerClusterizerValidation::createLayerHistograms(
0333     unsigned int ival) {
0334   std::ostringstream fname1, fname2;
0335 
0336   edm::Service<TFileService> fs;
0337   fs->file().cd("/");
0338 
0339   std::string tag;
0340   unsigned int id;
0341   if (ival < 1000) {
0342     id = ival;
0343     fname1 << "Barrel";
0344     fname2 << "Layer_" << id;
0345     tag = "_layer_";
0346   } else {
0347     int side = ival / 1000;
0348     id = ival - side * 1000;
0349     if (ival > 10) {
0350       id /= 10;
0351       //        fname1 << "EndCap_Side_" << side;
0352       fname1 << "EndCap";
0353       fname2 << "Ring_" << id;
0354       tag = "_ring_";
0355     } else {
0356       id = ival;
0357       //        fname1 << "EndCap_Side_" << side;
0358       fname1 << "EndCap";
0359       fname2 << "Disk_" << id;
0360       tag = "_disk_";
0361     }
0362   }
0363 
0364   TFileDirectory td1 = fs->mkdir(fname1.str().c_str());
0365   TFileDirectory td = td1.mkdir(fname2.str().c_str());
0366 
0367   ClusterHistos local_histos;
0368 
0369   std::ostringstream histoName;
0370 
0371   /*
0372      * Number of clusters
0373      */
0374 
0375   histoName.str("");
0376   histoName << "Number_Clusters_Pixel" << tag.c_str() << id;
0377   local_histos.numberClusters[1] = td.make<TH1D>(histoName.str().c_str(), histoName.str().c_str(), 50, 0., 0.);
0378   histoName.str("");
0379   histoName << "Number_Clusters_Strip" << tag.c_str() << id;
0380   local_histos.numberClusters[2] = td.make<TH1D>(histoName.str().c_str(), histoName.str().c_str(), 50, 0., 0.);
0381 
0382   /*
0383      * Cluster size
0384      */
0385 
0386   histoName.str("");
0387   histoName << "Cluster_Size_Pixel" << tag.c_str() << id;
0388   local_histos.clusterSize[1] = td.make<TH1D>(histoName.str().c_str(), histoName.str().c_str(), 21, -0.5, 20.5);
0389   histoName.str("");
0390   histoName << "Cluster_Size_Strip" << tag.c_str() << id;
0391   local_histos.clusterSize[2] = td.make<TH1D>(histoName.str().c_str(), histoName.str().c_str(), 21, -0.5, 20.5);
0392 
0393   /*
0394      * Local and Global positions
0395      */
0396 
0397   histoName.str("");
0398   histoName << "Local_Position_XY_Pixel" << tag.c_str() << id;
0399   local_histos.localPosXY[1] =
0400       td.make<TH2F>(histoName.str().c_str(), histoName.str().c_str(), 500, 0., 0., 500, 0., 0.);
0401 
0402   histoName.str("");
0403   histoName << "Local_Position_XY_Strip" << tag.c_str() << id;
0404   local_histos.localPosXY[2] =
0405       td.make<TH2F>(histoName.str().c_str(), histoName.str().c_str(), 500, 0., 0., 500, 0., 0.);
0406 
0407   histoName.str("");
0408   histoName << "Global_Position_XY_Pixel" << tag.c_str() << id;
0409   local_histos.globalPosXY[1] =
0410       td.make<TH2F>(histoName.str().c_str(), histoName.str().c_str(), 500, -120.0, 120.0, 2400, -120.0, 120.0);
0411 
0412   histoName.str("");
0413   histoName << "Global_Position_XY_Strip" << tag.c_str() << id;
0414   local_histos.globalPosXY[2] =
0415       td.make<TH2F>(histoName.str().c_str(), histoName.str().c_str(), 500, -120.0, 120.0, 2400, -120.0, 120.0);
0416 
0417   /*
0418      * Delta positions with SimHits
0419      */
0420 
0421   histoName.str("");
0422   histoName << "Delta_X_Pixel" << tag.c_str() << id;
0423   local_histos.deltaX[1] = td.make<TH1F>(histoName.str().c_str(), histoName.str().c_str(), 100, -0.02, 0.02);
0424 
0425   histoName.str("");
0426   histoName << "Delta_X_Strip" << tag.c_str() << id;
0427   local_histos.deltaX[2] = td.make<TH1F>(histoName.str().c_str(), histoName.str().c_str(), 100, -0.02, 0.02);
0428 
0429   histoName.str("");
0430   histoName << "Delta_Y_Pixel" << tag.c_str() << id;
0431   local_histos.deltaY[1] = td.make<TH1F>(histoName.str().c_str(), histoName.str().c_str(), 100, -0.2, 0.2);
0432 
0433   histoName.str("");
0434   histoName << "Delta_Y_Strip" << tag.c_str() << id;
0435   local_histos.deltaY[2] = td.make<TH1F>(histoName.str().c_str(), histoName.str().c_str(), 100, -3., 3.);
0436 
0437   /*
0438      * Delta position with simHits for primary tracks only
0439      */
0440 
0441   histoName.str("");
0442   histoName << "Delta_X_P" << tag.c_str() << id;
0443   local_histos.deltaX_P[1] = td.make<TH1F>(histoName.str().c_str(), histoName.str().c_str(), 100, -0.02, 0.02);
0444 
0445   histoName.str("");
0446   histoName << "Delta_X_P" << tag.c_str() << id;
0447   local_histos.deltaX_P[2] = td.make<TH1F>(histoName.str().c_str(), histoName.str().c_str(), 100, -0.02, 0.02);
0448 
0449   histoName.str("");
0450   histoName << "Delta_Y_P" << tag.c_str() << id;
0451   local_histos.deltaY_P[1] = td.make<TH1F>(histoName.str().c_str(), histoName.str().c_str(), 100, -0.2, 0.2);
0452 
0453   histoName.str("");
0454   histoName << "Delta_Y_P" << tag.c_str() << id;
0455   local_histos.deltaY_P[2] = td.make<TH1F>(histoName.str().c_str(), histoName.str().c_str(), 100, -3., 3.);
0456 
0457   /*
0458      * Information on the Digis per cluster
0459      */
0460 
0461   histoName.str("");
0462   histoName << "Primary_Digis_Pixel" << tag.c_str() << id;
0463   local_histos.primarySimHits[1] = td.make<TH1D>(histoName.str().c_str(), histoName.str().c_str(), 100, 0., 0.);
0464   histoName.str("");
0465   histoName << "Primary_Digis_Strip" << tag.c_str() << id;
0466   local_histos.primarySimHits[2] = td.make<TH1D>(histoName.str().c_str(), histoName.str().c_str(), 100, 0., 0.);
0467 
0468   histoName.str("");
0469   histoName << "Other_Digis_Pixel" << tag.c_str() << id;
0470   local_histos.otherSimHits[1] = td.make<TH1D>(histoName.str().c_str(), histoName.str().c_str(), 100, 0., 0.);
0471   histoName.str("");
0472   histoName << "Other_Digis_Strip" << tag.c_str() << id;
0473   local_histos.otherSimHits[2] = td.make<TH1D>(histoName.str().c_str(), histoName.str().c_str(), 100, 0., 0.);
0474 
0475   /*
0476      * End
0477      */
0478 
0479   std::pair<std::map<unsigned int, ClusterHistos>::iterator, bool> insertedIt(
0480       histograms_.insert(std::make_pair(ival, local_histos)));
0481   fs->file().cd("/");
0482 
0483   return insertedIt.first;
0484 }
0485 
0486 std::vector<unsigned int> Phase2TrackerClusterizerValidation::getSimTrackId(
0487     const edm::Handle<edm::DetSetVector<PixelDigiSimLink> >& pixelSimLinks, const DetId& detId, unsigned int channel) {
0488   std::vector<unsigned int> retvec;
0489   edm::DetSetVector<PixelDigiSimLink>::const_iterator DSViter(pixelSimLinks->find(detId));
0490   if (DSViter == pixelSimLinks->end())
0491     return retvec;
0492   for (edm::DetSet<PixelDigiSimLink>::const_iterator it = DSViter->data.begin(); it != DSViter->data.end(); ++it) {
0493     if (channel == it->channel()) {
0494       retvec.push_back(it->SimTrackId());
0495     }
0496   }
0497   return retvec;
0498 }
0499 
0500 DEFINE_FWK_MODULE(Phase2TrackerClusterizerValidation);