File indexing completed on 2023-03-17 11:21:39
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010 #include <memory>
0011
0012 #include "FWCore/Framework/interface/Frameworkfwd.h"
0013 #include "FWCore/Framework/interface/stream/EDProducer.h"
0014 #include "FWCore/Framework/interface/Event.h"
0015 #include "FWCore/Framework/interface/MakerMacros.h"
0016
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018 #include "FWCore/Utilities/interface/StreamID.h"
0019
0020 #include "DataFormats/Common/interface/DetSetVector.h"
0021
0022 #include "DataFormats/CTPPSDetId/interface/CTPPSDiamondDetId.h"
0023 #include "DataFormats/CTPPSReco/interface/CTPPSDiamondRecHit.h"
0024 #include "DataFormats/CTPPSReco/interface/CTPPSDiamondLocalTrack.h"
0025
0026 #include "RecoPPS/Local/interface/CTPPSDiamondTrackRecognition.h"
0027
0028 class CTPPSDiamondLocalTrackFitter : public edm::stream::EDProducer<> {
0029 public:
0030 explicit CTPPSDiamondLocalTrackFitter(const edm::ParameterSet&);
0031
0032 static void fillDescriptions(edm::ConfigurationDescriptions&);
0033
0034 private:
0035 void produce(edm::Event&, const edm::EventSetup&) override;
0036
0037 edm::EDGetTokenT<edm::DetSetVector<CTPPSDiamondRecHit> > recHitsToken_;
0038 const edm::ParameterSet trk_algo_params_;
0039 std::unordered_map<CTPPSDetId, std::unique_ptr<CTPPSDiamondTrackRecognition> > trk_algo_;
0040 };
0041
0042 CTPPSDiamondLocalTrackFitter::CTPPSDiamondLocalTrackFitter(const edm::ParameterSet& iConfig)
0043 : recHitsToken_(
0044 consumes<edm::DetSetVector<CTPPSDiamondRecHit> >(iConfig.getParameter<edm::InputTag>("recHitsTag"))),
0045 trk_algo_params_(iConfig.getParameter<edm::ParameterSet>("trackingAlgorithmParams")) {
0046 produces<edm::DetSetVector<CTPPSDiamondLocalTrack> >();
0047 }
0048
0049 void CTPPSDiamondLocalTrackFitter::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0050
0051 auto pOut = std::make_unique<edm::DetSetVector<CTPPSDiamondLocalTrack> >();
0052
0053 edm::Handle<edm::DetSetVector<CTPPSDiamondRecHit> > recHits;
0054 iEvent.getByToken(recHitsToken_, recHits);
0055
0056
0057 for (auto& algo_vs_id : trk_algo_)
0058 algo_vs_id.second->clear();
0059
0060
0061 for (const auto& vec : *recHits) {
0062 const CTPPSDiamondDetId raw_detid(vec.detId()), detid(raw_detid.arm(), raw_detid.station(), raw_detid.rp());
0063
0064 if (trk_algo_.count(detid) == 0)
0065 trk_algo_[detid] = std::make_unique<CTPPSDiamondTrackRecognition>(trk_algo_params_);
0066 for (const auto& hit : vec)
0067
0068 if (hit.ootIndex() != CTPPSDiamondRecHit::TIMESLICE_WITHOUT_LEADING)
0069 trk_algo_[detid]->addHit(hit);
0070 }
0071
0072
0073 for (auto& algo_vs_id : trk_algo_) {
0074 auto& tracks = pOut->find_or_insert(algo_vs_id.first);
0075 algo_vs_id.second->produceTracks(tracks);
0076 }
0077
0078 iEvent.put(std::move(pOut));
0079 }
0080
0081 void CTPPSDiamondLocalTrackFitter::fillDescriptions(edm::ConfigurationDescriptions& descr) {
0082 edm::ParameterSetDescription desc;
0083 desc.add<edm::InputTag>("recHitsTag", edm::InputTag("ctppsDiamondRecHits"))
0084 ->setComment("input rechits collection to retrieve");
0085
0086 edm::ParameterSetDescription trackingAlgoParams;
0087 trackingAlgoParams.add<double>("threshold", 1.5)
0088 ->setComment("minimal number of rechits to be observed before launching the track recognition algorithm");
0089 trackingAlgoParams.add<double>("thresholdFromMaximum", 0.5);
0090 trackingAlgoParams.add<double>("resolution", 0.01 )
0091 ->setComment("spatial resolution on the horizontal coordinate (in mm)");
0092 trackingAlgoParams.add<double>("sigma", 0.1);
0093 trackingAlgoParams.add<double>("startFromX", -0.5 )
0094 ->setComment("starting horizontal coordinate of rechits for the track recognition");
0095 trackingAlgoParams.add<double>("stopAtX", 19.5 )
0096 ->setComment("ending horizontal coordinate of rechits for the track recognition");
0097 trackingAlgoParams.add<double>("tolerance", 0.1 )
0098 ->setComment("tolerance used for checking if the track contains certain hit");
0099
0100 trackingAlgoParams.add<std::string>("pixelEfficiencyFunction", "(x>[0]-0.5*[1])*(x<[0]+0.5*[1])+0*[2]")
0101 ->setComment(
0102 "efficiency function for single pixel\n"
0103 "can be defined as:\n"
0104 " * Precise: (TMath::Erf((x-[0]+0.5*[1])/([2]/4)+2)+1)*TMath::Erfc((x-[0]-0.5*[1])/([2]/4)-2)/4\n"
0105 " * Fast: "
0106 "(x>[0]-0.5*[1])*(x<[0]+0.5*[1])+((x-[0]+0.5*[1]+[2])/"
0107 "[2])*(x>[0]-0.5*[1]-[2])*(x<[0]-0.5*[1])+(2-(x-[0]-0.5*[1]+[2])/[2])*(x>[0]+0.5*[1])*(x<[0]+0.5*[1]+[2])\n"
0108 " * Legacy: (1/(1+exp(-(x-[0]+0.5*[1])/[2])))*(1/(1+exp((x-[0]-0.5*[1])/[2])))\n"
0109 " * Default (sigma ignored): (x>[0]-0.5*[1])*(x<[0]+0.5*[1])+0*[2]\n"
0110 "with:\n"
0111 " [0]: centre of pad\n"
0112 " [1]: width of pad\n"
0113 " [2]: sigma: distance between efficiency ~100 -> 0 outside width");
0114
0115 trackingAlgoParams.add<double>("yPosition", 0.0)->setComment("vertical offset of the outcoming track centre");
0116 trackingAlgoParams.add<double>("yWidth", 0.0)->setComment("vertical track width");
0117 trackingAlgoParams.add<bool>("excludeSingleEdgeHits", true)
0118 ->setComment("exclude rechits with missing leading/trailing edge");
0119
0120 desc.add<edm::ParameterSetDescription>("trackingAlgorithmParams", trackingAlgoParams)
0121 ->setComment("list of parameters associated to the track recognition algorithm");
0122
0123 descr.add("ctppsDiamondLocalTracks", desc);
0124 }
0125
0126 DEFINE_FWK_MODULE(CTPPSDiamondLocalTrackFitter);