File indexing completed on 2025-06-03 00:12:18
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/TrackerTFP/interface/KalmanFilterFormats.h"
0019
0020 #include <TProfile.h>
0021 #include <TH1F.h>
0022 #include <TEfficiency.h>
0023
0024 #include <vector>
0025 #include <deque>
0026 #include <set>
0027 #include <cmath>
0028 #include <numeric>
0029 #include <sstream>
0030 #include <utility>
0031
0032 namespace trackerTFP {
0033
0034
0035
0036
0037
0038
0039 class AnalyzerKF : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::SharedResources> {
0040 public:
0041 AnalyzerKF(const edm::ParameterSet& iConfig);
0042 void beginJob() override {}
0043 void beginRun(const edm::Run& iEvent, const edm::EventSetup& iSetup) override;
0044 void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0045 void endRun(const edm::Run& iEvent, const edm::EventSetup& iSetup) override {}
0046 void endJob() override;
0047
0048 private:
0049
0050 void associate(const std::vector<TrackKF>& tracks,
0051 const std::vector<std::vector<StubKF*>>& stubs,
0052 int region,
0053 const tt::StubAssociation* ass,
0054 std::set<TPPtr>& tps,
0055 int& sum,
0056 const std::vector<TH1F*>& his,
0057 const std::vector<TProfile*>& prof,
0058 bool perfect = true) const;
0059
0060 edm::EDGetTokenT<tt::StreamsStub> edGetTokenStubs_;
0061
0062 edm::EDGetTokenT<tt::StreamsTrack> edGetTokenTracks_;
0063
0064 edm::EDGetTokenT<int> edGetTokenNumStatesAccepted_;
0065
0066 edm::EDGetTokenT<int> edGetTokenNumStatesTruncated_;
0067
0068 edm::EDGetTokenT<std::vector<std::pair<double, double>>> edGetTokenChi2s_;
0069
0070 edm::EDGetTokenT<tt::StubAssociation> edGetTokenSelection_;
0071
0072 edm::EDGetTokenT<tt::StubAssociation> edGetTokenReconstructable_;
0073
0074 edm::ESGetToken<tt::Setup, tt::SetupRcd> esGetTokenSetup_;
0075
0076 edm::ESGetToken<DataFormats, DataFormatsRcd> esGetTokenDataFormats_;
0077
0078 const tt::Setup* setup_ = nullptr;
0079
0080 const DataFormats* dataFormats_ = nullptr;
0081
0082 bool useMCTruth_;
0083
0084 int nEvents_ = 0;
0085
0086
0087
0088 TProfile* prof_;
0089 TProfile* profChannel_;
0090 TH1F* hisChannel_;
0091 std::vector<TH1F*> hisRes_;
0092 std::vector<TProfile*> profRes_;
0093 TH1F* hisEffEta_;
0094 TH1F* hisEffEtaTotal_;
0095 TEfficiency* effEta_;
0096 TH1F* hisEffZT_;
0097 TH1F* hisEffZTTotal_;
0098 TEfficiency* effZT_;
0099 TH1F* hisEffInv2R_;
0100 TH1F* hisEffInv2RTotal_;
0101 TEfficiency* effInv2R_;
0102 TH1F* hisEffPT_;
0103 TH1F* hisEffPTTotal_;
0104 TEfficiency* effPT_;
0105 TH1F* hisChi20s_;
0106 TH1F* hisChi21s_;
0107 TH1F* hisChi2s_;
0108 TH1F* hisTracks_;
0109 TH1F* hisLayers_;
0110 TH1F* hisNumLayers_;
0111 TProfile* profNumLayers_;
0112
0113
0114 std::stringstream log_;
0115 };
0116
0117 AnalyzerKF::AnalyzerKF(const edm::ParameterSet& iConfig)
0118 : useMCTruth_(iConfig.getParameter<bool>("UseMCTruth")), nEvents_(0), hisRes_(4), profRes_(4) {
0119 usesResource("TFileService");
0120
0121 const std::string& label = iConfig.getParameter<std::string>("OutputLabelKF");
0122 const std::string& branchStubs = iConfig.getParameter<std::string>("BranchStubs");
0123 const std::string& branchTracks = iConfig.getParameter<std::string>("BranchTracks");
0124 const std::string& branchTruncated = iConfig.getParameter<std::string>("BranchTruncated");
0125 edGetTokenStubs_ = consumes<tt::StreamsStub>(edm::InputTag(label, branchStubs));
0126 edGetTokenTracks_ = consumes<tt::StreamsTrack>(edm::InputTag(label, branchTracks));
0127 edGetTokenNumStatesAccepted_ = consumes<int>(edm::InputTag(label, branchTracks));
0128 edGetTokenNumStatesTruncated_ = consumes<int>(edm::InputTag(label, branchTruncated));
0129 edGetTokenChi2s_ = consumes<std::vector<std::pair<double, double>>>(edm::InputTag(label, branchTracks));
0130 if (useMCTruth_) {
0131 const auto& inputTagSelecttion = iConfig.getParameter<edm::InputTag>("InputTagSelection");
0132 const auto& inputTagReconstructable = iConfig.getParameter<edm::InputTag>("InputTagReconstructable");
0133 edGetTokenSelection_ = consumes<tt::StubAssociation>(inputTagSelecttion);
0134 edGetTokenReconstructable_ = consumes<tt::StubAssociation>(inputTagReconstructable);
0135 }
0136
0137 esGetTokenSetup_ = esConsumes<edm::Transition::BeginRun>();
0138 esGetTokenDataFormats_ = esConsumes<edm::Transition::BeginRun>();
0139
0140 log_.setf(std::ios::fixed, std::ios::floatfield);
0141 log_.precision(4);
0142 }
0143
0144 void AnalyzerKF::beginRun(const edm::Run& iEvent, const edm::EventSetup& iSetup) {
0145
0146 setup_ = &iSetup.getData(esGetTokenSetup_);
0147 dataFormats_ = &iSetup.getData(esGetTokenDataFormats_);
0148
0149 edm::Service<TFileService> fs;
0150 TFileDirectory dir;
0151 dir = fs->mkdir("KF");
0152 prof_ = dir.make<TProfile>("Counts", ";", 12, 0.5, 12.5);
0153 prof_->GetXaxis()->SetBinLabel(1, "Stubs");
0154 prof_->GetXaxis()->SetBinLabel(2, "Tracks");
0155 prof_->GetXaxis()->SetBinLabel(4, "Matched Tracks");
0156 prof_->GetXaxis()->SetBinLabel(5, "All Tracks");
0157 prof_->GetXaxis()->SetBinLabel(6, "Found TPs");
0158 prof_->GetXaxis()->SetBinLabel(7, "Found selected TPs");
0159 prof_->GetXaxis()->SetBinLabel(9, "All TPs");
0160 prof_->GetXaxis()->SetBinLabel(10, "states");
0161 prof_->GetXaxis()->SetBinLabel(12, "max tp");
0162
0163 constexpr int maxOcc = 180;
0164 const int numChannels = dataFormats_->numChannel(Process::kf);
0165 hisChannel_ = dir.make<TH1F>("His Channel Occupancy", ";", maxOcc, -.5, maxOcc - .5);
0166 profChannel_ = dir.make<TProfile>("Prof Channel Occupancy", ";", numChannels, -.5, numChannels - .5);
0167
0168 static const std::vector<std::string> names = {"phi0", "inv2R", "z0", "cot"};
0169 static const std::vector<double> ranges = {.01, .004, 20., .4};
0170 for (int i = 0; i < 4; i++) {
0171 const double range = ranges[i];
0172 hisRes_[i] = dir.make<TH1F>(("HisRes" + names[i]).c_str(), ";", 100, -range, range);
0173 profRes_[i] = dir.make<TProfile>(("ProfRes" + names[i]).c_str(), ";", 32, 0, 2.4);
0174 }
0175
0176 hisEffEtaTotal_ = dir.make<TH1F>("HisTPEtaTotal", ";", 48, -2.4, 2.4);
0177 hisEffEta_ = dir.make<TH1F>("HisTPEta", ";", 48, -2.4, 2.4);
0178 effEta_ = dir.make<TEfficiency>("EffEta", ";", 48, -2.4, 2.4);
0179 const int zTBins = setup_->gpNumBinsZT();
0180 hisEffZTTotal_ = dir.make<TH1F>("HisTPZTTotal", ";", zTBins, -zTBins / 2, zTBins / 2);
0181 hisEffZT_ = dir.make<TH1F>("HisTPZT", ";", zTBins, -zTBins / 2, zTBins / 2);
0182 effZT_ = dir.make<TEfficiency>("EffZT", ";", zTBins, -zTBins / 2, zTBins / 2);
0183 const double rangeInv2R = dataFormats_->format(Variable::inv2R, Process::dr).range();
0184 hisEffInv2R_ = dir.make<TH1F>("HisTPInv2R", ";", 32, -rangeInv2R / 2., rangeInv2R / 2.);
0185 hisEffInv2RTotal_ = dir.make<TH1F>("HisTPInv2RTotal", ";", 32, -rangeInv2R / 2., rangeInv2R / 2.);
0186 effInv2R_ = dir.make<TEfficiency>("EffInv2R", ";", 32, -rangeInv2R / 2., rangeInv2R / 2.);
0187 hisEffPT_ = dir.make<TH1F>("HisTPPT", ";", 100, 0, 100);
0188 hisEffPTTotal_ = dir.make<TH1F>("HisTPPTTotal", ";", 100, 0, 100);
0189 effPT_ = dir.make<TEfficiency>("EffPT", ";", 100, 0, 100);
0190
0191 hisChi20s_ = dir.make<TH1F>("HisChi20", ";", 128, 0., 10);
0192 hisChi21s_ = dir.make<TH1F>("HisChi21", ";", 128, 0., 10);
0193 hisChi2s_ = dir.make<TH1F>("HisChi2", ";", 128, 0., 10);
0194
0195 hisTracks_ = dir.make<TH1F>("HisTracks", ";", 40, 0., 400);
0196
0197 hisLayers_ = dir.make<TH1F>("HisLayers", ";", 8, 0, 8);
0198 hisNumLayers_ = dir.make<TH1F>("HisNumLayers", ";", 9, 0, 9);
0199 profNumLayers_ = dir.make<TProfile>("Prof NumLayers", ";", 32, 0, 2.4);
0200 }
0201
0202 void AnalyzerKF::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0203 static const int numChannel = dataFormats_->numChannel(Process::kf);
0204 static const int numLayers = setup_->numLayers();
0205 auto fill = [this](const TPPtr& tpPtr, TH1F* hisEta, TH1F* hisZT, TH1F* hisInv2R, TH1F* hisPT) {
0206 const double tpPhi0 = tpPtr->phi();
0207 const double tpCot = sinh(tpPtr->eta());
0208 const math::XYZPointD& v = tpPtr->vertex();
0209 const double tpZ0 = v.z() - tpCot * (v.x() * cos(tpPhi0) + v.y() * sin(tpPhi0));
0210 const double tpZT = tpZ0 + tpCot * setup_->chosenRofZ();
0211 hisEta->Fill(tpPtr->eta());
0212 hisZT->Fill(dataFormats_->format(Variable::zT, Process::gp).integer(tpZT));
0213 hisInv2R->Fill(tpPtr->charge() / tpPtr->pt() * setup_->invPtToDphi());
0214 hisPT->Fill(tpPtr->pt());
0215 };
0216
0217 edm::Handle<tt::StreamsStub> handleStubs;
0218 iEvent.getByToken<tt::StreamsStub>(edGetTokenStubs_, handleStubs);
0219 const tt::StreamsStub& allStubs = *handleStubs;
0220 edm::Handle<tt::StreamsTrack> handleTracks;
0221 iEvent.getByToken<tt::StreamsTrack>(edGetTokenTracks_, handleTracks);
0222 const tt::StreamsTrack& allTracks = *handleTracks;
0223 edm::Handle<int> handleNumStatesAccepted;
0224 iEvent.getByToken<int>(edGetTokenNumStatesAccepted_, handleNumStatesAccepted);
0225 edm::Handle<int> handleNumStatesTruncated;
0226 iEvent.getByToken<int>(edGetTokenNumStatesTruncated_, handleNumStatesTruncated);
0227 edm::Handle<std::vector<std::pair<double, double>>> handleChi2s;
0228 iEvent.getByToken<std::vector<std::pair<double, double>>>(edGetTokenChi2s_, handleChi2s);
0229
0230 const tt::StubAssociation* selection = nullptr;
0231 const tt::StubAssociation* reconstructable = nullptr;
0232 if (useMCTruth_) {
0233 edm::Handle<tt::StubAssociation> handleSelection;
0234 iEvent.getByToken<tt::StubAssociation>(edGetTokenSelection_, handleSelection);
0235 selection = handleSelection.product();
0236 prof_->Fill(9, selection->numTPs());
0237 edm::Handle<tt::StubAssociation> handleReconstructable;
0238 iEvent.getByToken<tt::StubAssociation>(edGetTokenReconstructable_, handleReconstructable);
0239 reconstructable = handleReconstructable.product();
0240 for (const auto& p : selection->getTrackingParticleToTTStubsMap())
0241 fill(p.first, hisEffEtaTotal_, hisEffZTTotal_, hisEffInv2RTotal_, hisEffPTTotal_);
0242 }
0243
0244 for (const std::pair<double, double>& chi2s : *handleChi2s) {
0245 hisChi20s_->Fill(chi2s.first * 2.);
0246 hisChi21s_->Fill(chi2s.second * 2.);
0247 hisChi2s_->Fill(chi2s.first + chi2s.second);
0248 }
0249
0250 std::set<TPPtr> tpPtrs;
0251 std::set<TPPtr> tpPtrsSelection;
0252 std::set<TPPtr> tpPtrsMax;
0253 int numMatched(0);
0254 int numTracks(0);
0255 for (int region = 0; region < setup_->numRegions(); region++) {
0256 int nRegionStubs(0);
0257 int nRegionTracks(0);
0258 for (int channel = 0; channel < numChannel; channel++) {
0259 const int index = region * numChannel + channel;
0260 const int offset = index * numLayers;
0261 const tt::StreamTrack& channelTracks = allTracks[index];
0262 hisChannel_->Fill(channelTracks.size());
0263 profChannel_->Fill(channel, channelTracks.size());
0264 std::vector<TrackKF> tracks;
0265 std::vector<StubKF> stubs;
0266 std::vector<std::vector<StubKF*>> tracksStubs(channelTracks.size(), std::vector<StubKF*>(numLayers, nullptr));
0267 tracks.reserve(channelTracks.size());
0268 stubs.reserve(channelTracks.size() * numLayers);
0269 for (int frame = 0; frame < static_cast<int>(channelTracks.size()); frame++) {
0270 tracks.emplace_back(channelTracks[frame], dataFormats_);
0271 const double cot = tracks.back().zT() / setup_->chosenRofZ();
0272 int nLs(0);
0273 for (int layer = 0; layer < numLayers; layer++) {
0274 const tt::FrameStub& fs = allStubs[offset + layer][frame];
0275 if (fs.first.isNull())
0276 continue;
0277 stubs.emplace_back(fs, dataFormats_);
0278 tracksStubs[frame][layer] = &stubs.back();
0279 hisLayers_->Fill(layer);
0280 nLs++;
0281 }
0282 hisNumLayers_->Fill(nLs);
0283 profNumLayers_->Fill(abs(sinh(cot)), nLs);
0284 }
0285 nRegionStubs += stubs.size();
0286 nRegionTracks += tracks.size();
0287 if (!useMCTruth_)
0288 continue;
0289 int tmp(0);
0290 associate(tracks, tracksStubs, region, selection, tpPtrsSelection, tmp, hisRes_, profRes_);
0291 associate(tracks,
0292 tracksStubs,
0293 region,
0294 reconstructable,
0295 tpPtrs,
0296 numMatched,
0297 std::vector<TH1F*>(),
0298 std::vector<TProfile*>(),
0299 false);
0300 associate(tracks,
0301 tracksStubs,
0302 region,
0303 selection,
0304 tpPtrsMax,
0305 tmp,
0306 std::vector<TH1F*>(),
0307 std::vector<TProfile*>(),
0308 false);
0309 }
0310 numTracks += nRegionTracks;
0311 prof_->Fill(1, nRegionStubs);
0312 prof_->Fill(2, nRegionTracks);
0313 }
0314 for (const TPPtr& tpPtr : tpPtrsSelection)
0315 fill(tpPtr, hisEffEta_, hisEffZT_, hisEffInv2R_, hisEffPT_);
0316 prof_->Fill(4, numMatched);
0317 prof_->Fill(5, numTracks);
0318 prof_->Fill(6, tpPtrs.size());
0319 prof_->Fill(7, tpPtrsSelection.size());
0320 prof_->Fill(10, *handleNumStatesAccepted);
0321 prof_->Fill(11, *handleNumStatesTruncated);
0322 prof_->Fill(12, tpPtrsMax.size());
0323 hisTracks_->Fill(numTracks);
0324 nEvents_++;
0325 }
0326
0327 void AnalyzerKF::endJob() {
0328 if (nEvents_ == 0)
0329 return;
0330
0331 effEta_->SetPassedHistogram(*hisEffEta_, "f");
0332 effEta_->SetTotalHistogram(*hisEffEtaTotal_, "f");
0333 effZT_->SetPassedHistogram(*hisEffZT_, "f");
0334 effZT_->SetTotalHistogram(*hisEffZTTotal_, "f");
0335 effInv2R_->SetPassedHistogram(*hisEffInv2R_, "f");
0336 effInv2R_->SetTotalHistogram(*hisEffInv2RTotal_, "f");
0337 effPT_->SetPassedHistogram(*hisEffPT_, "f");
0338 effPT_->SetTotalHistogram(*hisEffPTTotal_, "f");
0339
0340 const double totalTPs = prof_->GetBinContent(9);
0341 const double numStubs = prof_->GetBinContent(1);
0342 const double numTracks = prof_->GetBinContent(2);
0343 const double totalTracks = prof_->GetBinContent(5);
0344 const double numTracksMatched = prof_->GetBinContent(4);
0345 const double numTPsAll = prof_->GetBinContent(6);
0346 const double numTPsEff = prof_->GetBinContent(7);
0347 const double numTPsEffMax = prof_->GetBinContent(12);
0348 const double errStubs = prof_->GetBinError(1);
0349 const double errTracks = prof_->GetBinError(2);
0350 const double fracFake = (totalTracks - numTracksMatched) / totalTracks;
0351 const double fracDup = (numTracksMatched - numTPsAll) / totalTracks;
0352 const double eff = numTPsEff / totalTPs;
0353 const double errEff = sqrt(eff * (1. - eff) / totalTPs / nEvents_);
0354 const double effMax = numTPsEffMax / totalTPs;
0355 const double errEffMax = sqrt(effMax * (1. - effMax) / totalTPs / nEvents_);
0356 const int numStates = prof_->GetBinContent(10);
0357 const int numStatesLost = prof_->GetBinContent(11);
0358 const double fracSatest = numStates / (double)(numStates + numStatesLost);
0359 const std::vector<double> nums = {numStubs, numTracks};
0360 const std::vector<double> errs = {errStubs, errTracks};
0361 const int wNums = std::ceil(std::log10(*std::max_element(nums.begin(), nums.end()))) + 5;
0362 const int wErrs = std::ceil(std::log10(*std::max_element(errs.begin(), errs.end()))) + 5;
0363 log_ << " KF SUMMARY " << std::endl;
0364 log_ << "number of stubs per TFP = " << std::setw(wNums) << numStubs << " +- " << std::setw(wErrs) << errStubs
0365 << std::endl;
0366 log_ << "number of tracks per TFP = " << std::setw(wNums) << numTracks << " +- " << std::setw(wErrs)
0367 << errTracks << std::endl;
0368 log_ << " tracking efficiency = " << std::setw(wNums) << eff << " +- " << std::setw(wErrs) << errEff
0369 << std::endl;
0370 log_ << " max tracking efficiency = " << std::setw(wNums) << effMax << " +- " << std::setw(wErrs) << errEffMax
0371 << std::endl;
0372 log_ << " fake rate = " << std::setw(wNums) << fracFake << std::endl;
0373 log_ << " duplicate rate = " << std::setw(wNums) << fracDup << std::endl;
0374 log_ << " state assessment fraction = " << std::setw(wNums) << fracSatest << std::endl;
0375 log_ << " number of states per TFP = " << std::setw(wNums) << (numStates + numStatesLost) / setup_->numRegions()
0376 << std::endl;
0377 log_ << "=============================================================";
0378 edm::LogPrint(moduleDescription().moduleName()) << log_.str();
0379 }
0380
0381
0382 void AnalyzerKF::associate(const std::vector<TrackKF>& tracks,
0383 const std::vector<std::vector<StubKF*>>& tracksStubs,
0384 int region,
0385 const tt::StubAssociation* ass,
0386 std::set<TPPtr>& tps,
0387 int& sum,
0388 const std::vector<TH1F*>& his,
0389 const std::vector<TProfile*>& prof,
0390 bool perfect) const {
0391 for (int frame = 0; frame < static_cast<int>(tracks.size()); frame++) {
0392 const TrackKF& track = tracks[frame];
0393 const std::vector<StubKF*>& stubs = tracksStubs[frame];
0394 std::vector<TTStubRef> ttStubRefs;
0395 ttStubRefs.reserve(stubs.size());
0396 TTBV hitPattern(0, setup_->numLayers());
0397 int layer(-1);
0398 for (StubKF* stub : stubs) {
0399 layer++;
0400 if (!stub)
0401 continue;
0402 hitPattern.set(layer);
0403 ttStubRefs.push_back(stub->frame().first);
0404 }
0405 const std::vector<TPPtr>& tpPtrs = perfect ? ass->associateFinal(ttStubRefs) : ass->associate(ttStubRefs);
0406 if (tpPtrs.empty())
0407 continue;
0408 sum++;
0409 std::copy(tpPtrs.begin(), tpPtrs.end(), std::inserter(tps, tps.begin()));
0410 if (his.empty())
0411 continue;
0412 const double zT = dataFormats_->format(Variable::zT, Process::gp).digi(track.zT());
0413 const double cot = zT / setup_->chosenRofZ() + track.cot();
0414 const double z0 = track.zT() - setup_->chosenRofZ() * cot;
0415 const double inv2R = track.inv2R();
0416 const double phi0 = tt::deltaPhi(track.phiT() - setup_->chosenRofPhi() * inv2R +
0417 region * dataFormats_->format(Variable::phiT, Process::kf).range());
0418 for (const TPPtr& tpPtr : tpPtrs) {
0419 const double tpPhi0 = tpPtr->phi();
0420 const double tpCot = std::sinh(tpPtr->eta());
0421 const double tpInv2R = -setup_->invPtToDphi() * tpPtr->charge() / tpPtr->pt();
0422 const math::XYZPointD& v = tpPtr->vertex();
0423 const double tpZ0 = v.z() - tpCot * (v.x() * cos(tpPhi0) + v.y() * sin(tpPhi0));
0424 const double dCot = tpCot - cot;
0425 const double dZ0 = tpZ0 - z0;
0426 const double dInv2R = tpInv2R - inv2R;
0427 const double dPhi0 = tt::deltaPhi(tpPhi0 - phi0);
0428 const std::vector<double> ds = {dPhi0, dInv2R / setup_->invPtToDphi(), dZ0, dCot};
0429 for (int i = 0; i < static_cast<int>(ds.size()); i++) {
0430 his[i]->Fill(ds[i]);
0431 prof[i]->Fill(abs(tpPtr->eta()), abs(ds[i]));
0432 }
0433 }
0434 }
0435 }
0436
0437 }
0438
0439 DEFINE_FWK_MODULE(trackerTFP::AnalyzerKF);