238 auto result = std::make_unique<TracksterCollection>();
239 auto output_mask = std::make_unique<std::vector<float>>();
240 auto result_fromCP = std::make_unique<TracksterCollection>();
241 auto resultPU = std::make_unique<TracksterCollection>();
242 auto output_mask_fromCP = std::make_unique<std::vector<float>>();
243 auto cpToSc_SimTrackstersMap = std::make_unique<std::map<uint, std::vector<uint>>>();
244 auto result_ticlCandidates = std::make_unique<std::vector<TICLCandidate>>();
270 const auto& recoTracks = *recoTracks_h;
275 result->reserve(num_simclusters);
277 result_fromCP->resize(num_caloparticles);
278 std::map<uint, uint> SimClusterToCaloParticleMap;
280 for (
const auto& [
key, lcVec] : caloParticlesToRecoColl) {
281 auto const&
cp = *(
key);
283 for (
const auto& scRef :
cp.simClusters()) {
284 auto const& sc = *(scRef);
286 SimClusterToCaloParticleMap[scIndex] = cpIndex;
289 auto regr_energy =
cp.energy();
290 std::vector<uint> scSimTracksterIdx;
291 scSimTracksterIdx.reserve(
cp.simClusters().size());
294 if (
cp.g4Tracks()[0].crossedBoundary()) {
295 regr_energy =
cp.g4Tracks()[0].getMomentumAtBoundary().energy();
296 float time =
cp.g4Tracks()[0].getPositionAtBoundary().t();
311 for (
const auto& scRef :
cp.simClusters()) {
312 const auto&
it = simClustersToRecoColl.find(scRef);
313 if (
it == simClustersToRecoColl.end())
315 const auto& lcVec =
it->val;
316 auto const& sc = *(scRef);
323 sc.g4Tracks()[0].getMomentumAtBoundary().energy(),
326 sc.g4Tracks()[0].getPositionAtBoundary().t(),
336 if (
std::find(scSimTracksterIdx.begin(), scSimTracksterIdx.end(),
index) == scSimTracksterIdx.end()) {
337 scSimTracksterIdx.emplace_back(
index);
340 scSimTracksterIdx.shrink_to_fit();
342 float time = simVertices[
cp.g4Tracks()[0].vertIndex()].position().t();
359 if (result_fromCP->empty())
361 const auto index = loop_index - 1;
362 if (cpToSc_SimTrackstersMap->find(
index) == cpToSc_SimTrackstersMap->end()) {
363 (*cpToSc_SimTrackstersMap)[
index] = scSimTracksterIdx;
382 makePUTrackster(inputClusterMask, *output_mask, *resultPU, caloParticles_h.
id(), 0);
384 auto simTrackToRecoTrack = [&](
UniqueSimTrackId simTkId) -> std::pair<int, float> {
387 auto ipos = simTrackToTPMap.mapping.find(simTkId);
388 if (ipos != simTrackToTPMap.mapping.end()) {
389 auto jpos = TPtoRecoTrackMap.find((ipos->second));
390 if (jpos != TPtoRecoTrackMap.end()) {
391 auto& associatedRecoTracks = jpos->val;
392 if (!associatedRecoTracks.empty()) {
395 trackIdx = &(*associatedRecoTracks[0].first) - &recoTracks[0];
396 quality = associatedRecoTracks[0].second;
405 auto& simTrackstersFromCP = *result_fromCP;
406 for (
unsigned int i = 0;
i < simTrackstersFromCP.size(); ++
i) {
407 if (simTrackstersFromCP[
i].
vertices().empty())
411 auto bestAssociatedRecoTrack = simTrackToRecoTrack(simTkIds);
412 if (bestAssociatedRecoTrack.first != -1 and bestAssociatedRecoTrack.second >
qualityCutTrack_) {
413 auto trackIndex = bestAssociatedRecoTrack.first;
414 simTrackstersFromCP[
i].setTrackIdx(trackIndex);
418 auto& simTracksters = *
result;
420 std::unordered_map<unsigned int, std::vector<unsigned int>> TPtoSimTracksterMap;
421 for (
unsigned int i = 0;
i < simTracksters.size(); ++
i) {
422 const auto&
simTrack = (simTracksters[
i].seedID() == caloParticles_h.
id())
424 :
simclusters[simTracksters[
i].seedIndex()].g4Tracks()[0];
426 auto bestAssociatedRecoTrack = simTrackToRecoTrack(simTkIds);
427 if (bestAssociatedRecoTrack.first != -1 and bestAssociatedRecoTrack.second >
qualityCutTrack_) {
428 auto trackIndex = bestAssociatedRecoTrack.first;
429 simTracksters[
i].setTrackIdx(trackIndex);
436 std::unordered_map<unsigned int, const MtdSimTrackster*> SimTrackToMtdST;
437 for (
unsigned int i = 0;
i < MTDSimTracksters_h->size(); ++
i) {
438 const auto&
simTrack = (*MTDSimTracksters_h)[
i].g4Tracks()[0];
439 SimTrackToMtdST[
simTrack.trackId()] = &((*MTDSimTracksters_h)[
i]);
442 result_ticlCandidates->resize(result_fromCP->size());
443 std::vector<int> toKeep;
444 for (
size_t i = 0;
i < simTracksters_h->size(); ++
i) {
445 const auto& simTrackster = (*simTracksters_h)[
i];
446 int cp_index = (simTrackster.seedID() == caloParticles_h.
id())
447 ? simTrackster.seedIndex()
448 : SimClusterToCaloParticleMap[simTrackster.seedIndex()];
449 auto const& tCP = (*result_fromCP)[cp_index];
450 if (!tCP.vertices().empty()) {
451 auto trackIndex = tCP.trackIdx();
453 auto&
cand = (*result_ticlCandidates)[cp_index];
457 toKeep.push_back(cp_index);
461 auto isHad = [](
int pdgId) {
468 for (
size_t i = 0;
i < result_ticlCandidates->size(); ++
i) {
469 auto cp_index = (*result_fromCP)[
i].seedIndex();
472 auto&
cand = (*result_ticlCandidates)[
i];
475 float regressedEnergy = 0.f;
478 auto pos = SimTrackToMtdST.find(
simTrack.trackId());
479 if (
pos != SimTrackToMtdST.end()) {
480 auto MTDst =
pos->second;
482 cand.setMTDTime(MTDst->time(), 0);
485 cand.setTime(simVertices[
cp.g4Tracks()[0].vertIndex()].position().t() *
pow(10, 9), 0);
487 for (
const auto& trackster :
cand.tracksters()) {
489 regressedEnergy += trackster->regressed_energy();
495 if (
cand.trackPtr().isNonnull()) {
496 auto const&
track =
cand.trackPtr().get();
504 regressedEnergy *
track->momentum().unit().y(),
505 regressedEnergy *
track->momentum().unit().z(),
513 else if (
pdgId == 111)
522 const auto& simTracksterFromCP = (*result_fromCP)[
i];
523 float regressedEnergy = simTracksterFromCP.regressed_energy();
525 regressedEnergy * simTracksterFromCP.barycenter().unit().y(),
526 regressedEnergy * simTracksterFromCP.barycenter().unit().z(),
532 std::vector<int> all_nums(result_fromCP->size());
533 std::iota(all_nums.begin(), all_nums.end(), 0);
535 std::vector<int> toRemove;
536 std::set_difference(all_nums.begin(), all_nums.end(), toKeep.begin(), toKeep.end(), std::back_inserter(toRemove));
537 std::sort(toRemove.begin(), toRemove.end(), [](
int x,
int y) {
return x >
y; });
538 for (
auto const&
r : toRemove) {
539 result_fromCP->erase(result_fromCP->begin() +
r);
540 result_ticlCandidates->erase(result_ticlCandidates->begin() +
r);
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
T const & getData(const ESGetToken< T, R > &iToken) const noexcept(false)
bool get(ProductID const &oid, Handle< PROD > &result) const
const edm::EDGetTokenT< reco::SimToRecoCollection > associatormapStRsToken_
std::pair< uint32_t, EncodedEventId > UniqueSimTrackId
const edm::EDGetTokenT< MtdSimTracksterCollection > MTDSimTrackstersToken_
const edm::EDGetTokenT< std::vector< SimVertex > > simVerticesToken_
const edm::EDGetTokenT< std::vector< float > > filtered_layerclusters_mask_token_
void assignPCAtoTracksters(std::vector< Trackster > &tracksters, const std::vector< reco::CaloCluster > &layerClusters, const edm::ValueMap< std::pair< float, float >> &layerClustersTime, double z_limit_em, hgcal::RecHitTools const &rhTools, bool computeLocalTime=false, bool energyWeight=true, bool clean=false, int minLayer=10, int maxLayer=10)
const edm::ESGetToken< CaloGeometry, CaloGeometryRecord > geom_token_
bool getByToken(EDGetToken token, Handle< PROD > &result) const
const edm::EDGetTokenT< SimTrackToTPMap > associationSimTrackToTPToken_
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
U second(std::pair< T, U > const &p)
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Trackster::ParticleType tracksterParticleTypeFromPdgId(int pdgId, int charge)
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
void addTrackster(const int index, const std::vector< std::pair< edm::Ref< reco::CaloClusterCollection >, std::pair< float, float >>> &lcVec, const std::vector< float > &inputClusterMask, const float fractionCut_, const float energy, const int pdgId, const int charge, const float time, const edm::ProductID seed, const Trackster::IterationIndex iter, std::vector< float > &output_mask, std::vector< Trackster > &result, int &loop_index, const bool add=false)
Abs< T >::type abs(const T &t)
key
prepare the HTCondor submission files and eventually submit them
const edm::EDGetTokenT< std::vector< CaloParticle > > caloparticles_token_
const edm::EDGetTokenT< std::vector< reco::CaloCluster > > clusters_token_
const edm::EDGetTokenT< edm::ValueMap< std::pair< float, float > > > clustersTime_token_
const edm::EDGetTokenT< ticl::SimToRecoCollectionWithSimClusters > associatorMapSimClusterToReco_token_
const edm::EDGetTokenT< std::vector< SimCluster > > simclusters_token_
hgcal::RecHitTools rhtools_
const edm::EDGetTokenT< std::vector< reco::Track > > recoTracksToken_
const edm::EDGetTokenT< ticl::SimToRecoCollection > associatorMapCaloParticleToReco_token_
const float qualityCutTrack_
const edm::EDGetTokenT< std::vector< TrackingParticle > > trackingParticleToken_
const bool computeLocalTime_
Power< A, B >::type pow(const A &a, const B &b)
void makePUTrackster(const std::vector< float > &inputClusterMask, std::vector< float > &output_mask, std::vector< Trackster > &result, const edm::ProductID seed, int loop_index)