Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-10-14 01:44:14

0001 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0002 #include "FWCore/Framework/interface/Run.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Framework/interface/EventSetup.h"
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "FWCore/ServiceRegistry/interface/Service.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/Utilities/interface/EDGetToken.h"
0010 #include "FWCore/Utilities/interface/InputTag.h"
0011 #include "FWCore/Utilities/interface/Exception.h"
0012 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0013 #include "DataFormats/Common/interface/Handle.h"
0014 
0015 #include "SimTracker/TrackTriggerAssociation/interface/StubAssociation.h"
0016 #include "L1Trigger/TrackTrigger/interface/Setup.h"
0017 #include "L1Trigger/TrackerTFP/interface/DataFormats.h"
0018 #include "L1Trigger/TrackFindingTracklet/interface/ChannelAssignment.h"
0019 #include "L1Trigger/TrackFindingTracklet/interface/Settings.h"
0020 
0021 #include <TProfile.h>
0022 #include <TProfile2D.h>
0023 #include <TH1F.h>
0024 #include <TH2F.h>
0025 
0026 #include <vector>
0027 #include <deque>
0028 #include <set>
0029 #include <cmath>
0030 #include <numeric>
0031 #include <sstream>
0032 
0033 using namespace std;
0034 using namespace edm;
0035 using namespace trackerTFP;
0036 using namespace tt;
0037 
0038 namespace trklet {
0039 
0040   // stub resolution plots helper
0041   enum Resolution { Phi, Z, NumResolution };
0042   constexpr initializer_list<Resolution> AllResolution = {Phi, Z};
0043   constexpr auto NameResolution = {"Phi", "Z"};
0044   inline string name(Resolution r) { return string(*(NameResolution.begin() + r)); }
0045 
0046   /*! \class  trklet::AnalyzerTBout
0047    *  \brief  Class to analyze hardware like structured TTStub Collection generated by Tracklet
0048    *  \author Thomas Schuh
0049    *  \date   2021, Nov
0050    */
0051   class AnalyzerTBout : public one::EDAnalyzer<one::WatchRuns, one::SharedResources> {
0052   public:
0053     AnalyzerTBout(const ParameterSet& iConfig);
0054     void beginJob() override {}
0055     void beginRun(const Run& iEvent, const EventSetup& iSetup) override;
0056     void analyze(const Event& iEvent, const EventSetup& iSetup) override;
0057     void endRun(const Run& iEvent, const EventSetup& iSetup) override {}
0058     void endJob() override;
0059 
0060   private:
0061     //
0062     void formTracks(const StreamsTrack& streamsTrack,
0063                     const StreamsStub& streamsStubs,
0064                     vector<vector<TTStubRef>>& tracks,
0065                     int channel);
0066     //
0067     void associate(const vector<vector<TTStubRef>>& tracks, const StubAssociation* ass, set<TPPtr>& tps, int& sum) const;
0068     //
0069     void fill(const FrameTrack& frameTrack, const FrameStub& frameStub);
0070 
0071     // ED input token of dtc stubs
0072     EDGetTokenT<TTDTC> edGetTokenTTDTC_;
0073     // ED input token of stubs
0074     EDGetTokenT<StreamsStub> edGetTokenAcceptedStubs_;
0075     // ED input token of tracks
0076     EDGetTokenT<StreamsTrack> edGetTokenAcceptedTracks_;
0077     // ED input token of lost stubs
0078     EDGetTokenT<StreamsStub> edGetTokenLostStubs_;
0079     // ED input token of lost tracks
0080     EDGetTokenT<StreamsTrack> edGetTokenLostTracks_;
0081     // ED input token of TTStubRef to TPPtr association for tracking efficiency
0082     EDGetTokenT<StubAssociation> edGetTokenSelection_;
0083     // ED input token of TTStubRef to recontructable TPPtr association
0084     EDGetTokenT<StubAssociation> edGetTokenReconstructable_;
0085     // Setup token
0086     ESGetToken<Setup, SetupRcd> esGetTokenSetup_;
0087     // DataFormats token
0088     ESGetToken<DataFormats, DataFormatsRcd> esGetTokenDataFormats_;
0089     // ChannelAssignment token
0090     ESGetToken<ChannelAssignment, ChannelAssignmentRcd> esGetTokenChannelAssignment_;
0091     // stores, calculates and provides run-time constants
0092     const Setup* setup_ = nullptr;
0093     //
0094     const Settings settings_;
0095     // helper class to extract structured data from tt::Frames
0096     const DataFormats* dataFormats_ = nullptr;
0097     // helper class to assign tracklet track to channel
0098     const ChannelAssignment* channelAssignment_ = nullptr;
0099     // enables analyze of TPs
0100     bool useMCTruth_;
0101     //
0102     int nEvents_ = 0;
0103     //
0104     vector<deque<FrameStub>> regionStubs_;
0105     //
0106     int region_;
0107     //
0108     vector<int> nOverflows_;
0109 
0110     // Histograms
0111 
0112     TProfile* prof_;
0113     TProfile* profChannel_;
0114     TH1F* hisChannel_;
0115     vector<TH1F*> hisResolution_;
0116     vector<TProfile2D*> profResolution_;
0117     vector<TH1F*> hisResolutionMe_;
0118     vector<TH1F*> hisResolutionThey_;
0119     vector<TH2F*> his2Resolution_;
0120 
0121     // printout
0122     stringstream log_;
0123   };
0124 
0125   AnalyzerTBout::AnalyzerTBout(const ParameterSet& iConfig) : useMCTruth_(iConfig.getParameter<bool>("UseMCTruth")) {
0126     usesResource("TFileService");
0127     // book in- and output ED products
0128     const InputTag& inputTag = iConfig.getParameter<InputTag>("InputTagDTC");
0129     const string& label = iConfig.getParameter<string>("LabelTBout");
0130     const string& branchAcceptedStubs = iConfig.getParameter<string>("BranchAcceptedStubs");
0131     const string& branchAcceptedTracks = iConfig.getParameter<string>("BranchAcceptedTracks");
0132     const string& branchLostStubs = iConfig.getParameter<string>("BranchLostStubs");
0133     const string& branchLostTracks = iConfig.getParameter<string>("BranchLostTracks");
0134     edGetTokenTTDTC_ = consumes<TTDTC>(inputTag);
0135     edGetTokenAcceptedStubs_ = consumes<StreamsStub>(InputTag(label, branchAcceptedStubs));
0136     edGetTokenAcceptedTracks_ = consumes<StreamsTrack>(InputTag(label, branchAcceptedTracks));
0137     edGetTokenLostStubs_ = consumes<StreamsStub>(InputTag(label, branchLostStubs));
0138     edGetTokenLostTracks_ = consumes<StreamsTrack>(InputTag(label, branchLostTracks));
0139     if (useMCTruth_) {
0140       const auto& inputTagSelecttion = iConfig.getParameter<InputTag>("InputTagSelection");
0141       const auto& inputTagReconstructable = iConfig.getParameter<InputTag>("InputTagReconstructable");
0142       edGetTokenSelection_ = consumes<StubAssociation>(inputTagSelecttion);
0143       edGetTokenReconstructable_ = consumes<StubAssociation>(inputTagReconstructable);
0144     }
0145     // book ES products
0146     esGetTokenSetup_ = esConsumes<Setup, SetupRcd, Transition::BeginRun>();
0147     esGetTokenDataFormats_ = esConsumes<DataFormats, DataFormatsRcd, Transition::BeginRun>();
0148     esGetTokenChannelAssignment_ = esConsumes<ChannelAssignment, ChannelAssignmentRcd, Transition::BeginRun>();
0149     //
0150     nOverflows_ = vector<int>(2, 0);
0151     // log config
0152     log_.setf(ios::fixed, ios::floatfield);
0153     log_.precision(4);
0154   }
0155 
0156   void AnalyzerTBout::beginRun(const Run& iEvent, const EventSetup& iSetup) {
0157     // helper class to store configurations
0158     setup_ = &iSetup.getData(esGetTokenSetup_);
0159     // helper class to extract structured data from tt::Frames
0160     dataFormats_ = &iSetup.getData(esGetTokenDataFormats_);
0161     // helper class to assign tracklet track to channel
0162     channelAssignment_ = &iSetup.getData(esGetTokenChannelAssignment_);
0163     // book histograms
0164     Service<TFileService> fs;
0165     TFileDirectory dir;
0166     dir = fs->mkdir("TBout");
0167     prof_ = dir.make<TProfile>("Counts", ";", 9, 0.5, 9.5);
0168     prof_->GetXaxis()->SetBinLabel(1, "Stubs");
0169     prof_->GetXaxis()->SetBinLabel(2, "Tracks");
0170     prof_->GetXaxis()->SetBinLabel(3, "Lost Tracks");
0171     prof_->GetXaxis()->SetBinLabel(4, "Matched Tracks");
0172     prof_->GetXaxis()->SetBinLabel(5, "All Tracks");
0173     prof_->GetXaxis()->SetBinLabel(6, "Found TPs");
0174     prof_->GetXaxis()->SetBinLabel(7, "Found selected TPs");
0175     prof_->GetXaxis()->SetBinLabel(8, "Lost TPs");
0176     prof_->GetXaxis()->SetBinLabel(9, "All TPs");
0177     // channel occupancy
0178     constexpr int maxOcc = 180;
0179     const int numChannels = channelAssignment_->numChannelsStub() * setup_->numRegions();
0180     hisChannel_ = dir.make<TH1F>("His Channel Occupancy", ";", maxOcc, -.5, maxOcc - .5);
0181     profChannel_ = dir.make<TProfile>("Prof Channel Occupancy", ";", numChannels, -.5, numChannels - .5);
0182     // stub parameter resolutions
0183     constexpr int bins = 400;
0184     constexpr int binsHis = 100;
0185     constexpr double maxZ = 300.;
0186     constexpr double maxR = 120.;
0187     constexpr array<double, NumResolution> ranges{{.01, 5.}};
0188     hisResolution_.reserve(NumResolution);
0189     profResolution_.reserve(NumResolution);
0190     for (Resolution r : AllResolution) {
0191       hisResolution_.emplace_back(dir.make<TH1F>(("HisRes" + name(r)).c_str(), ";", binsHis, -ranges[r], ranges[r]));
0192       profResolution_.emplace_back(
0193           dir.make<TProfile2D>(("ProfRes" + name(r)).c_str(), ";;", bins, -maxZ, maxZ, bins, 0., maxR));
0194       hisResolutionMe_.emplace_back(
0195           dir.make<TH1F>(("HisResMe" + name(r)).c_str(), ";", binsHis, -ranges[r], ranges[r]));
0196       hisResolutionThey_.emplace_back(
0197           dir.make<TH1F>(("HisResThey" + name(r)).c_str(), ";", binsHis, -ranges[r], ranges[r]));
0198       his2Resolution_.emplace_back(
0199           dir.make<TH2F>(("His2" + name(r)).c_str(), ";;", bins, -ranges[r], ranges[r], bins, -ranges[r], ranges[r]));
0200     }
0201     regionStubs_ = vector<deque<FrameStub>>(setup_->numRegions());
0202   }
0203 
0204   void AnalyzerTBout::analyze(const Event& iEvent, const EventSetup& iSetup) {
0205     // read in TTDTC
0206     Handle<TTDTC> handleTTDTC;
0207     iEvent.getByToken<TTDTC>(edGetTokenTTDTC_, handleTTDTC);
0208     const TTDTC& ttDTC = *handleTTDTC;
0209     for (deque<FrameStub>& region : regionStubs_)
0210       region.clear();
0211     for (int r : ttDTC.tfpRegions()) {
0212       for (int c : ttDTC.tfpChannels()) {
0213         const StreamStub& s = ttDTC.stream(r, c);
0214         copy_if(
0215             s.begin(), s.end(), back_inserter(regionStubs_[r]), [](const FrameStub& f) { return f.first.isNonnull(); });
0216       }
0217     }
0218     // read in TBout products
0219     Handle<StreamsStub> handleAcceptedStubs;
0220     iEvent.getByToken<StreamsStub>(edGetTokenAcceptedStubs_, handleAcceptedStubs);
0221     const StreamsStub& acceptedStubs = *handleAcceptedStubs;
0222     Handle<StreamsTrack> handleAcceptedTracks;
0223     iEvent.getByToken<StreamsTrack>(edGetTokenAcceptedTracks_, handleAcceptedTracks);
0224     const StreamsTrack& acceptedTracks = *handleAcceptedTracks;
0225     Handle<StreamsStub> handleLostStubs;
0226     iEvent.getByToken<StreamsStub>(edGetTokenLostStubs_, handleLostStubs);
0227     const StreamsStub& lostStubs = *handleLostStubs;
0228     Handle<StreamsTrack> handleLostTracks;
0229     iEvent.getByToken<StreamsTrack>(edGetTokenLostTracks_, handleLostTracks);
0230     const StreamsTrack& lostTracks = *handleLostTracks;
0231     // read in MCTruth
0232     const StubAssociation* selection = nullptr;
0233     const StubAssociation* reconstructable = nullptr;
0234     if (useMCTruth_) {
0235       Handle<StubAssociation> handleSelection;
0236       iEvent.getByToken<StubAssociation>(edGetTokenSelection_, handleSelection);
0237       selection = handleSelection.product();
0238       prof_->Fill(9, selection->numTPs());
0239       Handle<StubAssociation> handleReconstructable;
0240       iEvent.getByToken<StubAssociation>(edGetTokenReconstructable_, handleReconstructable);
0241       reconstructable = handleReconstructable.product();
0242     }
0243     // analyze ht products and associate found tracks with reconstrucable TrackingParticles
0244     set<TPPtr> tpPtrs;
0245     set<TPPtr> tpPtrsSelection;
0246     set<TPPtr> tpPtrsLost;
0247     int allMatched(0);
0248     int allTracks(0);
0249     for (region_ = 0; region_ < setup_->numRegions(); region_++) {
0250       const int offset = region_ * channelAssignment_->numChannelsTrack();
0251       int nStubs(0);
0252       int nTracks(0);
0253       int nLost(0);
0254       for (int channel = 0; channel < channelAssignment_->numChannelsTrack(); channel++) {
0255         vector<vector<TTStubRef>> tracks;
0256         formTracks(acceptedTracks, acceptedStubs, tracks, offset + channel);
0257         vector<vector<TTStubRef>> lost;
0258         formTracks(lostTracks, lostStubs, lost, offset + channel);
0259         nTracks += tracks.size();
0260         nStubs +=
0261             accumulate(tracks.begin(), tracks.end(), 0, [](int& sum, const auto& v) { return sum += (int)v.size(); });
0262         nLost += lost.size();
0263         allTracks += tracks.size();
0264         if (!useMCTruth_)
0265           continue;
0266         int tmp(0);
0267         associate(tracks, selection, tpPtrsSelection, tmp);
0268         associate(lost, selection, tpPtrsLost, tmp);
0269         associate(tracks, reconstructable, tpPtrs, allMatched);
0270       }
0271       prof_->Fill(1, nStubs);
0272       prof_->Fill(2, nTracks);
0273       prof_->Fill(3, nLost);
0274     }
0275     vector<TPPtr> recovered;
0276     recovered.reserve(tpPtrsLost.size());
0277     set_intersection(tpPtrsLost.begin(), tpPtrsLost.end(), tpPtrs.begin(), tpPtrs.end(), back_inserter(recovered));
0278     for (const TPPtr& tpPtr : recovered)
0279       tpPtrsLost.erase(tpPtr);
0280     prof_->Fill(4, allMatched);
0281     prof_->Fill(5, allTracks);
0282     prof_->Fill(6, tpPtrs.size());
0283     prof_->Fill(7, tpPtrsSelection.size());
0284     prof_->Fill(8, tpPtrsLost.size());
0285     nEvents_++;
0286   }
0287 
0288   void AnalyzerTBout::endJob() {
0289     if (nEvents_ == 0)
0290       return;
0291     // printout TBout summary
0292     const double totalTPs = prof_->GetBinContent(9);
0293     const double numStubs = prof_->GetBinContent(1);
0294     const double numTracks = prof_->GetBinContent(2);
0295     const double numTracksLost = prof_->GetBinContent(3);
0296     const double totalTracks = prof_->GetBinContent(5);
0297     const double numTracksMatched = prof_->GetBinContent(4);
0298     const double numTPsAll = prof_->GetBinContent(6);
0299     const double numTPsEff = prof_->GetBinContent(7);
0300     const double numTPsLost = prof_->GetBinContent(8);
0301     const double errStubs = prof_->GetBinError(1);
0302     const double errTracks = prof_->GetBinError(2);
0303     const double errTracksLost = prof_->GetBinError(3);
0304     const double fracFake = (totalTracks - numTracksMatched) / totalTracks;
0305     const double fracDup = (numTracksMatched - numTPsAll) / totalTracks;
0306     const double eff = numTPsEff / totalTPs;
0307     const double errEff = sqrt(eff * (1. - eff) / totalTPs / nEvents_);
0308     const double effLoss = numTPsLost / totalTPs;
0309     const double errEffLoss = sqrt(effLoss * (1. - effLoss) / totalTPs / nEvents_);
0310     const vector<double> nums = {numStubs, numTracks, numTracksLost};
0311     const vector<double> errs = {errStubs, errTracks, errTracksLost};
0312     const int wNums = ceil(log10(*max_element(nums.begin(), nums.end()))) + 5;
0313     const int wErrs = ceil(log10(*max_element(errs.begin(), errs.end()))) + 5;
0314     log_ << "                        TBout SUMMARY                        " << endl;
0315     log_ << "number of stubs       per TFP = " << setw(wNums) << numStubs << " +- " << setw(wErrs) << errStubs << endl;
0316     log_ << "number of tracks      per TFP = " << setw(wNums) << numTracks << " +- " << setw(wErrs) << errTracks
0317          << endl;
0318     log_ << "number of lost tracks per TFP = " << setw(wNums) << numTracksLost << " +- " << setw(wErrs) << errTracksLost
0319          << endl;
0320     log_ << "     max  tracking efficiency = " << setw(wNums) << eff << " +- " << setw(wErrs) << errEff << endl;
0321     log_ << "     lost tracking efficiency = " << setw(wNums) << effLoss << " +- " << setw(wErrs) << errEffLoss << endl;
0322     log_ << "                    fake rate = " << setw(wNums) << fracFake << endl;
0323     log_ << "               duplicate rate = " << setw(wNums) << fracDup << endl;
0324     log_ << "number of overflowed phi residuals: " << nOverflows_[0] << " and z: " << nOverflows_[1] << endl;
0325     log_ << "=============================================================";
0326     LogPrint("L1Trigger/TrackerTFP") << log_.str();
0327   }
0328 
0329   //
0330   void AnalyzerTBout::formTracks(const StreamsTrack& streamsTrack,
0331                                  const StreamsStub& streamsStubs,
0332                                  vector<vector<TTStubRef>>& tracks,
0333                                  int channel) {
0334     const int seedType = channel % channelAssignment_->numChannelsTrack();
0335     const int offset = channelAssignment_->offsetStub(channel);
0336     const StreamTrack& streamTrack = streamsTrack[channel];
0337     const int numTracks = accumulate(streamTrack.begin(), streamTrack.end(), 0, [](int& sum, const FrameTrack& frame) {
0338       return sum += (frame.first.isNonnull() ? 1 : 0);
0339     });
0340     tracks.reserve(numTracks);
0341     for (int frame = 0; frame < (int)streamTrack.size(); frame++) {
0342       const FrameTrack& frameTrack = streamTrack[frame];
0343       if (frameTrack.first.isNull())
0344         continue;
0345       vector<TTStubRef> ttStubRefs;
0346       ttStubRefs.reserve(channelAssignment_->numProjectionLayers(seedType) +
0347                          channelAssignment_->seedingLayers(seedType).size());
0348       for (int layer = 0; layer < channelAssignment_->numProjectionLayers(seedType); layer++) {
0349         const FrameStub& stub = streamsStubs[offset + layer][frame];
0350         if (stub.first.isNonnull()) {
0351           this->fill(frameTrack, stub);
0352           ttStubRefs.push_back(stub.first);
0353         }
0354       }
0355       for (const TTStubRef& ttStubRef : frameTrack.first->getStubRefs()) {
0356         int layer(0);
0357         if (!channelAssignment_->layerId(seedType, ttStubRef, layer))
0358           ttStubRefs.push_back(ttStubRef);
0359       }
0360       tracks.push_back(ttStubRefs);
0361     }
0362   }
0363 
0364   //
0365   void AnalyzerTBout::associate(const vector<vector<TTStubRef>>& tracks,
0366                                 const StubAssociation* ass,
0367                                 set<TPPtr>& tps,
0368                                 int& sum) const {
0369     for (const vector<TTStubRef>& ttStubRefs : tracks) {
0370       const vector<TPPtr>& tpPtrs = ass->associate(ttStubRefs);
0371       if (tpPtrs.empty())
0372         continue;
0373       sum++;
0374       copy(tpPtrs.begin(), tpPtrs.end(), inserter(tps, tps.begin()));
0375     }
0376   }
0377 
0378   //
0379   void AnalyzerTBout::fill(const FrameTrack& frameTrack, const FrameStub& frameStub) {
0380     // get dtc stub frame
0381     const deque<FrameStub>& region = regionStubs_[region_];
0382     const auto it =
0383         find_if(region.begin(), region.end(), [&frameStub](const FrameStub& f) { return f.first == frameStub.first; });
0384     if (it == region.end())
0385       throw cms::Exception("LgociError.") << "Stub on track was not in DTC collection.";
0386     const GlobalPoint ttPos = setup_->stubPos(frameStub.first);
0387     const GlobalPoint pos = setup_->stubPos(true, *it, region_);
0388     static constexpr int widthPhi = 12;
0389     static constexpr int widthZ = 9;
0390     static constexpr int widthR = 7;
0391     const bool barrel = setup_->barrel(frameStub.first);
0392     const int layerIdTracklet = setup_->trackletLayerId(frameStub.first);
0393     static const double baseR = settings_.kz();
0394     const double basePhi = barrel ? settings_.kphi1() : settings_.kphi(layerIdTracklet);
0395     const double baseZ = settings_.kz(layerIdTracklet);
0396     static const double baseInvR = settings_.kphi1() / settings_.kr() * pow(2, settings_.rinv_shift());
0397     static const double basePhi0 = settings_.kphi1() * pow(2, settings_.phi0_shift());
0398     static const double baseZ0 = settings_.kz() * pow(2, settings_.z0_shift());
0399     static const double baseTanL = settings_.kz() / settings_.kr() * pow(2, settings_.t_shift());
0400     const int widthRZ = barrel ? widthZ : widthR;
0401     const double baseRZ = barrel ? baseZ : baseR;
0402     // calc residuals
0403     const double rInv = (frameTrack.first->rInv() / baseInvR + .5) * baseInvR;
0404     const double phi0 = (frameTrack.first->phi() / basePhi0 + .5) * basePhi0;
0405     const double z0 = (frameTrack.first->z0() / baseZ0 + .5) * baseZ0;
0406     const double tanL = (frameTrack.first->tanL() / baseTanL + .5) * baseTanL;
0407     const double phi = deltaPhi(pos.phi() - (phi0 - rInv * pos.perp() / 2.));
0408     const double r = pos.perp() - (pos.z() - z0) / tanL;
0409     const double z = pos.z() - (z0 + tanL * pos.perp());
0410     const double rz = barrel ? z : r;
0411     const int phii = floor(phi / basePhi);
0412     const int rzi = floor(rz / baseRZ);
0413     const double phid = (phii + .5) * basePhi;
0414     const double rzd = (rzi + .5) * baseRZ;
0415     // parse residuals
0416     TTBV hw(frameStub.second);
0417     const TTBV hwRZ(hw, widthRZ, 0, true);
0418     hw >>= widthRZ;
0419     const TTBV hwPhi(hw, widthPhi, 0, true);
0420     bool overflowPhi = abs(phii - hwPhi.val()) > pow(2, widthPhi) * 7. / 8.;
0421     bool overflowZ = abs(rzi - hwRZ.val()) > pow(2, widthRZ) * 7. / 8.;
0422     const double hwPhid = hwPhi.val(basePhi);
0423     const double hwRZd = hwRZ.val(baseRZ);
0424     const vector<double> resolutions = {phid - hwPhid, rzd - hwRZd};
0425     const vector<bool> overflows = {overflowPhi, overflowZ};
0426     for (Resolution r : AllResolution) {
0427       if (overflows[r]) {
0428         cout << rzi << " " << hwRZ.val() << " " << barrel << " " << setup_->psModule(frameStub.first) << " "
0429              << frameStub.second << endl;
0430         nOverflows_[r]++;
0431         continue;
0432       }
0433       hisResolution_[r]->Fill(resolutions[r]);
0434       profResolution_[r]->Fill(ttPos.z(), ttPos.perp(), abs(resolutions[r]));
0435     }
0436     hisResolutionMe_[0]->Fill(phid);
0437     hisResolutionMe_[1]->Fill(rzd);
0438     hisResolutionThey_[0]->Fill(hwPhid);
0439     hisResolutionThey_[1]->Fill(hwRZd);
0440     his2Resolution_[0]->Fill(phid, hwPhid);
0441     his2Resolution_[1]->Fill(rzd, hwRZd);
0442   }
0443 
0444 }  // namespace trklet
0445 
0446 DEFINE_FWK_MODULE(trklet::AnalyzerTBout);