Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-12-17 09:21:49

0001 #include "JetRecoEval.h"
0002 
0003 #include "CaloEvalStack.h"
0004 #include "CaloRawClusterEval.h"
0005 #include "CaloRawTowerEval.h"  // for CaloRawTowerEval
0006 #include "JetTruthEval.h"
0007 #include "SvtxEvalStack.h"
0008 #include "SvtxTrackEval.h"
0009 
0010 #include <calobase/RawCluster.h>
0011 #include <calobase/RawClusterContainer.h>
0012 #include <calobase/RawTower.h>
0013 #include <calobase/RawTowerContainer.h>
0014 
0015 #include <jetbase/Jet.h>
0016 #include <jetbase/JetContainer.h>
0017 
0018 #include <g4main/PHG4Particle.h>
0019 
0020 #include <trackbase_historic/SvtxTrack.h>
0021 #include <trackbase_historic/SvtxTrackMap.h>
0022 
0023 #include <phool/getClass.h>
0024 #include <phool/phool.h>
0025 
0026 #include <cassert>
0027 #include <cfloat>
0028 #include <cmath>
0029 #include <cstdlib>
0030 #include <iostream>
0031 #include <map>
0032 #include <set>
0033 #include <string>
0034 
0035 JetRecoEval::JetRecoEval(PHCompositeNode* topNode,
0036                          const std::string& recojetname,
0037                          const std::string& truthjetname)
0038   : _jettrutheval(topNode, truthjetname)
0039   , _recojetname(recojetname)
0040   , _truthjetname(truthjetname)
0041 {
0042   get_node_pointers(topNode);
0043 }
0044 
0045 JetRecoEval::~JetRecoEval()
0046 {
0047   if (_verbosity > 0)
0048   {
0049     if ((_errors > 0) || (_verbosity > 1))
0050     {
0051       std::cout << "JetRecoEval::~JetRecoEval() - Error Count: " << _errors << std::endl;
0052     }
0053   }
0054 }
0055 
0056 void JetRecoEval::next_event(PHCompositeNode* topNode)
0057 {
0058   _cache_all_truth_showers.clear();
0059   _cache_all_truth_particles.clear();
0060   _cache_all_truth_jets.clear();
0061   _cache_max_truth_jet_by_energy.clear();
0062   _cache_all_jets_from.clear();
0063   _cache_best_jet_from.clear();
0064   _cache_get_energy_contribution.clear();
0065   _cache_get_energy_contribution_src.clear();
0066   _cache_all_truth_hits.clear();
0067 
0068   _jettrutheval.next_event(topNode);
0069 
0070   get_node_pointers(topNode);
0071 }
0072 
0073 void JetRecoEval::set_track_nodename(const std::string& name)
0074 {
0075   m_TrackNodeName = name;
0076   _jettrutheval.set_track_nodename(name);
0077 }
0078 
0079 std::set<PHG4Shower*> JetRecoEval::all_truth_showers(Jet* recojet)
0080 {
0081   if (_strict)
0082   {
0083     assert(recojet);
0084   }
0085   else if (!recojet)
0086   {
0087     ++_errors;
0088     return std::set<PHG4Shower*>();
0089   }
0090 
0091   if (_do_cache)
0092   {
0093     std::map<Jet*, std::set<PHG4Shower*> >::iterator iter =
0094         _cache_all_truth_showers.find(recojet);
0095     if (iter != _cache_all_truth_showers.end())
0096     {
0097       return iter->second;
0098     }
0099   }
0100 
0101   std::set<PHG4Shower*> truth_showers;
0102 
0103   // loop over all the jet constituents, backtrack each reco object to the
0104   // truth hits and combine with other consituents
0105   for (auto jter : recojet->get_comp_vec())
0106   {
0107     Jet::SRC source = jter.first;
0108     unsigned int index = jter.second;
0109 
0110     std::set<PHG4Shower*> new_showers;
0111 
0112     if (source == Jet::TRACK)
0113     {
0114       if (!_trackmap)
0115       {
0116         std::cout << PHWHERE << "ERROR: can't find SvtxTrackMap" << std::endl;
0117         exit(-1);
0118       }
0119 
0120       // SvtxTrack* track = _trackmap->get(index);
0121 
0122       // if (_strict) {assert(track);}
0123       // else if (!track) {++_errors; continue;}
0124 
0125       // new_showers = get_svtx_eval_stack()->get_track_eval()->all_truth_showers(track);
0126     }
0127 
0128     else if (source == Jet::CEMC_TOWER)
0129     {
0130       if (!_cemctowers)
0131       {
0132         std::cout << PHWHERE << "ERROR: can't find TOWER_CEMC" << std::endl;
0133         exit(-1);
0134       }
0135 
0136       RawTower* tower = _cemctowers->getTower(index);
0137 
0138       if (_strict)
0139       {
0140         assert(tower);
0141       }
0142       else if (!tower)
0143       {
0144         ++_errors;
0145         continue;
0146       }
0147 
0148       new_showers = get_cemc_eval_stack()->get_rawtower_eval()->all_truth_primary_showers(tower);
0149     }
0150     else if (source == Jet::CEMC_CLUSTER)
0151     {
0152       if (!_cemcclusters)
0153       {
0154         std::cout << PHWHERE << "ERROR: can't find CLUSTER_CEMC" << std::endl;
0155         exit(-1);
0156       }
0157 
0158       RawCluster* cluster = _cemcclusters->getCluster(index);
0159 
0160       if (_strict)
0161       {
0162         assert(cluster);
0163       }
0164       else if (!cluster)
0165       {
0166         ++_errors;
0167         continue;
0168       }
0169 
0170       new_showers = get_cemc_eval_stack()->get_rawcluster_eval()->all_truth_primary_showers(cluster);
0171     }
0172     else if (source == Jet::EEMC_TOWER)
0173     {
0174       if (!_eemctowers)
0175       {
0176         std::cout << PHWHERE << "ERROR: can't find TOWER_EEMC" << std::endl;
0177         exit(-1);
0178       }
0179 
0180       RawTower* tower = _eemctowers->getTower(index);
0181 
0182       if (_strict)
0183       {
0184         assert(tower);
0185       }
0186       else if (!tower)
0187       {
0188         ++_errors;
0189         continue;
0190       }
0191 
0192       new_showers = get_eemc_eval_stack()->get_rawtower_eval()->all_truth_primary_showers(tower);
0193     }
0194     else if (source == Jet::EEMC_CLUSTER)
0195     {
0196       if (!_eemcclusters)
0197       {
0198         std::cout << PHWHERE << "ERROR: can't find CLUSTER_EEMC" << std::endl;
0199         exit(-1);
0200       }
0201 
0202       RawCluster* cluster = _eemcclusters->getCluster(index);
0203 
0204       if (_strict)
0205       {
0206         assert(cluster);
0207       }
0208       else if (!cluster)
0209       {
0210         ++_errors;
0211         continue;
0212       }
0213 
0214       new_showers = get_eemc_eval_stack()->get_rawcluster_eval()->all_truth_primary_showers(cluster);
0215     }
0216     else if (source == Jet::HCALIN_TOWER)
0217     {
0218       if (!_hcalintowers)
0219       {
0220         std::cout << PHWHERE << "ERROR: can't find TOWER_HCALIN" << std::endl;
0221         exit(-1);
0222       }
0223 
0224       RawTower* tower = _hcalintowers->getTower(index);
0225 
0226       if (_strict)
0227       {
0228         assert(tower);
0229       }
0230       else if (!tower)
0231       {
0232         ++_errors;
0233         continue;
0234       }
0235 
0236       new_showers = get_hcalin_eval_stack()->get_rawtower_eval()->all_truth_primary_showers(tower);
0237     }
0238     else if (source == Jet::HCALIN_CLUSTER)
0239     {
0240       if (!_hcalinclusters)
0241       {
0242         std::cout << PHWHERE << "ERROR: can't find CLUSTER_HCALIN" << std::endl;
0243         exit(-1);
0244       }
0245 
0246       RawCluster* cluster = _hcalinclusters->getCluster(index);
0247 
0248       if (_strict)
0249       {
0250         assert(cluster);
0251       }
0252       else if (!cluster)
0253       {
0254         ++_errors;
0255         continue;
0256       }
0257 
0258       new_showers = get_hcalin_eval_stack()->get_rawcluster_eval()->all_truth_primary_showers(cluster);
0259     }
0260     else if (source == Jet::HCALOUT_TOWER)
0261     {
0262       if (!_hcalouttowers)
0263       {
0264         std::cout << PHWHERE << "ERROR: can't find TOWER_HCALOUT" << std::endl;
0265         exit(-1);
0266       }
0267 
0268       RawTower* tower = _hcalouttowers->getTower(index);
0269 
0270       if (_strict)
0271       {
0272         assert(tower);
0273       }
0274       else if (!tower)
0275       {
0276         ++_errors;
0277         continue;
0278       }
0279 
0280       new_showers = get_hcalout_eval_stack()->get_rawtower_eval()->all_truth_primary_showers(tower);
0281     }
0282     else if (source == Jet::HCALOUT_CLUSTER)
0283     {
0284       if (!_hcaloutclusters)
0285       {
0286         std::cout << PHWHERE << "ERROR: can't find CLUSTER_HCALOUT" << std::endl;
0287         exit(-1);
0288       }
0289 
0290       RawCluster* cluster = _hcaloutclusters->getCluster(index);
0291 
0292       if (_strict)
0293       {
0294         assert(cluster);
0295       }
0296       else if (!cluster)
0297       {
0298         ++_errors;
0299         continue;
0300       }
0301 
0302       new_showers = get_hcalout_eval_stack()->get_rawcluster_eval()->all_truth_primary_showers(cluster);
0303     }
0304     else if (source == Jet::FEMC_TOWER)
0305     {
0306       if (!_femctowers)
0307       {
0308         std::cout << PHWHERE << "ERROR: can't find TOWER_FEMC" << std::endl;
0309         exit(-1);
0310       }
0311 
0312       RawTower* tower = _femctowers->getTower(index);
0313 
0314       if (_strict)
0315       {
0316         assert(tower);
0317       }
0318       else if (!tower)
0319       {
0320         ++_errors;
0321         continue;
0322       }
0323 
0324       new_showers = get_femc_eval_stack()->get_rawtower_eval()->all_truth_primary_showers(tower);
0325     }
0326     else if (source == Jet::FEMC_CLUSTER)
0327     {
0328       if (!_femcclusters)
0329       {
0330         std::cout << PHWHERE << "ERROR: can't find CLUSTER_FEMC" << std::endl;
0331         exit(-1);
0332       }
0333 
0334       RawCluster* cluster = _femcclusters->getCluster(index);
0335 
0336       if (_strict)
0337       {
0338         assert(cluster);
0339       }
0340       else if (!cluster)
0341       {
0342         ++_errors;
0343         continue;
0344       }
0345 
0346       new_showers = get_femc_eval_stack()->get_rawcluster_eval()->all_truth_primary_showers(cluster);
0347     }
0348     else if (source == Jet::FHCAL_TOWER)
0349     {
0350       if (!_fhcaltowers)
0351       {
0352         std::cout << PHWHERE << "ERROR: can't find TOWER_FHCAL" << std::endl;
0353         exit(-1);
0354       }
0355 
0356       RawTower* tower = _fhcaltowers->getTower(index);
0357 
0358       if (_strict)
0359       {
0360         assert(tower);
0361       }
0362       else if (!tower)
0363       {
0364         ++_errors;
0365         continue;
0366       }
0367 
0368       new_showers = get_fhcal_eval_stack()->get_rawtower_eval()->all_truth_primary_showers(tower);
0369     }
0370     else if (source == Jet::FHCAL_CLUSTER)
0371     {
0372       if (!_fhcalclusters)
0373       {
0374         std::cout << PHWHERE << "ERROR: can't find CLUSTER_FHCAL" << std::endl;
0375         exit(-1);
0376       }
0377 
0378       RawCluster* cluster = _fhcalclusters->getCluster(index);
0379 
0380       if (_strict)
0381       {
0382         assert(cluster);
0383       }
0384       else if (!cluster)
0385       {
0386         ++_errors;
0387         continue;
0388       }
0389 
0390       new_showers = get_fhcal_eval_stack()->get_rawcluster_eval()->all_truth_primary_showers(cluster);
0391     }
0392 
0393     for (auto* new_shower : new_showers)
0394     {
0395       truth_showers.insert(new_shower);
0396     }
0397   }
0398 
0399   if (_do_cache)
0400   {
0401     _cache_all_truth_showers.insert(std::make_pair(recojet, truth_showers));
0402   }
0403 
0404   return truth_showers;
0405 }
0406 
0407 std::set<PHG4Particle*> JetRecoEval::all_truth_particles(Jet* recojet)
0408 {
0409   if (_strict)
0410   {
0411     assert(recojet);
0412   }
0413   else if (!recojet)
0414   {
0415     ++_errors;
0416     return std::set<PHG4Particle*>();
0417   }
0418 
0419   if (_do_cache)
0420   {
0421     std::map<Jet*, std::set<PHG4Particle*> >::iterator iter =
0422         _cache_all_truth_particles.find(recojet);
0423     if (iter != _cache_all_truth_particles.end())
0424     {
0425       return iter->second;
0426     }
0427   }
0428 
0429   std::set<PHG4Particle*> truth_particles;
0430 
0431   // loop over all the jet constituents, backtrack each reco object to the
0432   // truth hits and combine with other consituents
0433   for (auto jter : recojet->get_comp_vec())
0434   {
0435     Jet::SRC source = jter.first;
0436     unsigned int index = jter.second;
0437 
0438     std::set<PHG4Particle*> new_particles;
0439 
0440     if (source == Jet::TRACK)
0441     {
0442       if (!_trackmap)
0443       {
0444         std::cout << PHWHERE << "ERROR: can't find TrackMap" << std::endl;
0445         exit(-1);
0446       }
0447 
0448       SvtxTrack* track = _trackmap->get(index);
0449 
0450       if (_strict)
0451       {
0452         assert(track);
0453       }
0454       else if (!track)
0455       {
0456         ++_errors;
0457         continue;
0458       }
0459 
0460       new_particles = get_svtx_eval_stack()->get_track_eval()->all_truth_particles(track);
0461     }
0462     else if (source == Jet::CEMC_TOWER)
0463     {
0464       if (!_cemctowers)
0465       {
0466         std::cout << PHWHERE << "ERROR: can't find TOWER_CEMC" << std::endl;
0467         exit(-1);
0468       }
0469 
0470       RawTower* tower = _cemctowers->getTower(index);
0471 
0472       if (_strict)
0473       {
0474         assert(tower);
0475       }
0476       else if (!tower)
0477       {
0478         ++_errors;
0479         continue;
0480       }
0481 
0482       new_particles = get_cemc_eval_stack()->get_rawtower_eval()->all_truth_primary_particles(tower);
0483     }
0484     else if (source == Jet::CEMC_CLUSTER)
0485     {
0486       if (!_cemcclusters)
0487       {
0488         std::cout << PHWHERE << "ERROR: can't find CLUSTER_CEMC" << std::endl;
0489         exit(-1);
0490       }
0491 
0492       RawCluster* cluster = _cemcclusters->getCluster(index);
0493 
0494       if (_strict)
0495       {
0496         assert(cluster);
0497       }
0498       else if (!cluster)
0499       {
0500         ++_errors;
0501         continue;
0502       }
0503 
0504       new_particles = get_cemc_eval_stack()->get_rawcluster_eval()->all_truth_primary_particles(cluster);
0505     }
0506     else if (source == Jet::EEMC_TOWER)
0507     {
0508       if (!_eemctowers)
0509       {
0510         std::cout << PHWHERE << "ERROR: can't find TOWER_EEMC" << std::endl;
0511         exit(-1);
0512       }
0513 
0514       RawTower* tower = _eemctowers->getTower(index);
0515 
0516       if (_strict)
0517       {
0518         assert(tower);
0519       }
0520       else if (!tower)
0521       {
0522         ++_errors;
0523         continue;
0524       }
0525 
0526       new_particles = get_eemc_eval_stack()->get_rawtower_eval()->all_truth_primary_particles(tower);
0527     }
0528     else if (source == Jet::EEMC_CLUSTER)
0529     {
0530       if (!_eemcclusters)
0531       {
0532         std::cout << PHWHERE << "ERROR: can't find CLUSTER_EEMC" << std::endl;
0533         exit(-1);
0534       }
0535 
0536       RawCluster* cluster = _eemcclusters->getCluster(index);
0537 
0538       if (_strict)
0539       {
0540         assert(cluster);
0541       }
0542       else if (!cluster)
0543       {
0544         ++_errors;
0545         continue;
0546       }
0547 
0548       new_particles = get_eemc_eval_stack()->get_rawcluster_eval()->all_truth_primary_particles(cluster);
0549     }
0550     else if (source == Jet::HCALIN_TOWER)
0551     {
0552       if (!_hcalintowers)
0553       {
0554         std::cout << PHWHERE << "ERROR: can't find TOWER_HCALIN" << std::endl;
0555         exit(-1);
0556       }
0557 
0558       RawTower* tower = _hcalintowers->getTower(index);
0559 
0560       if (_strict)
0561       {
0562         assert(tower);
0563       }
0564       else if (!tower)
0565       {
0566         ++_errors;
0567         continue;
0568       }
0569 
0570       new_particles = get_hcalin_eval_stack()->get_rawtower_eval()->all_truth_primary_particles(tower);
0571     }
0572     else if (source == Jet::HCALIN_CLUSTER)
0573     {
0574       if (!_hcalinclusters)
0575       {
0576         std::cout << PHWHERE << "ERROR: can't find CLUSTER_HCALIN" << std::endl;
0577         exit(-1);
0578       }
0579 
0580       RawCluster* cluster = _hcalinclusters->getCluster(index);
0581 
0582       if (_strict)
0583       {
0584         assert(cluster);
0585       }
0586       else if (!cluster)
0587       {
0588         ++_errors;
0589         continue;
0590       }
0591 
0592       new_particles = get_hcalin_eval_stack()->get_rawcluster_eval()->all_truth_primary_particles(cluster);
0593     }
0594     else if (source == Jet::HCALOUT_TOWER)
0595     {
0596       if (!_hcalouttowers)
0597       {
0598         std::cout << PHWHERE << "ERROR: can't find TOWER_HCALOUT" << std::endl;
0599         exit(-1);
0600       }
0601 
0602       RawTower* tower = _hcalouttowers->getTower(index);
0603 
0604       if (_strict)
0605       {
0606         assert(tower);
0607       }
0608       else if (!tower)
0609       {
0610         ++_errors;
0611         continue;
0612       }
0613 
0614       new_particles = get_hcalout_eval_stack()->get_rawtower_eval()->all_truth_primary_particles(tower);
0615     }
0616     else if (source == Jet::HCALOUT_CLUSTER)
0617     {
0618       if (!_hcaloutclusters)
0619       {
0620         std::cout << PHWHERE << "ERROR: can't find CLUSTER_HCALOUT" << std::endl;
0621         exit(-1);
0622       }
0623 
0624       RawCluster* cluster = _hcaloutclusters->getCluster(index);
0625 
0626       if (_strict)
0627       {
0628         assert(cluster);
0629       }
0630       else if (!cluster)
0631       {
0632         ++_errors;
0633         continue;
0634       }
0635 
0636       new_particles = get_hcalout_eval_stack()->get_rawcluster_eval()->all_truth_primary_particles(cluster);
0637     }
0638     else if (source == Jet::FEMC_TOWER)
0639     {
0640       if (!_femctowers)
0641       {
0642         std::cout << PHWHERE << "ERROR: can't find TOWER_FEMC" << std::endl;
0643         exit(-1);
0644       }
0645 
0646       RawTower* tower = _femctowers->getTower(index);
0647 
0648       if (_strict)
0649       {
0650         assert(tower);
0651       }
0652       else if (!tower)
0653       {
0654         ++_errors;
0655         continue;
0656       }
0657 
0658       new_particles = get_femc_eval_stack()->get_rawtower_eval()->all_truth_primary_particles(tower);
0659     }
0660     else if (source == Jet::FEMC_CLUSTER)
0661     {
0662       if (!_femcclusters)
0663       {
0664         std::cout << PHWHERE << "ERROR: can't find CLUSTER_FEMC" << std::endl;
0665         exit(-1);
0666       }
0667 
0668       RawCluster* cluster = _femcclusters->getCluster(index);
0669 
0670       if (_strict)
0671       {
0672         assert(cluster);
0673       }
0674       else if (!cluster)
0675       {
0676         ++_errors;
0677         continue;
0678       }
0679 
0680       new_particles = get_femc_eval_stack()->get_rawcluster_eval()->all_truth_primary_particles(cluster);
0681     }
0682     else if (source == Jet::FHCAL_TOWER)
0683     {
0684       if (!_fhcaltowers)
0685       {
0686         std::cout << PHWHERE << "ERROR: can't find TOWER_FHCAL" << std::endl;
0687         exit(-1);
0688       }
0689 
0690       RawTower* tower = _fhcaltowers->getTower(index);
0691 
0692       if (_strict)
0693       {
0694         assert(tower);
0695       }
0696       else if (!tower)
0697       {
0698         ++_errors;
0699         continue;
0700       }
0701 
0702       new_particles = get_fhcal_eval_stack()->get_rawtower_eval()->all_truth_primary_particles(tower);
0703     }
0704     else if (source == Jet::FHCAL_CLUSTER)
0705     {
0706       if (!_fhcalclusters)
0707       {
0708         std::cout << PHWHERE << "ERROR: can't find CLUSTER_FHCAL" << std::endl;
0709         exit(-1);
0710       }
0711 
0712       RawCluster* cluster = _fhcalclusters->getCluster(index);
0713 
0714       if (_strict)
0715       {
0716         assert(cluster);
0717       }
0718       else if (!cluster)
0719       {
0720         ++_errors;
0721         continue;
0722       }
0723 
0724       new_particles = get_fhcal_eval_stack()->get_rawcluster_eval()->all_truth_primary_particles(cluster);
0725     }
0726 
0727     for (auto* new_particle : new_particles)
0728     {
0729       truth_particles.insert(new_particle);
0730     }
0731   }
0732 
0733   if (_do_cache)
0734   {
0735     _cache_all_truth_particles.insert(std::make_pair(recojet, truth_particles));
0736   }
0737 
0738   return truth_particles;
0739 }
0740 
0741 std::set<Jet*> JetRecoEval::all_truth_jets(Jet* recojet)
0742 {
0743   if (_strict)
0744   {
0745     assert(recojet);
0746   }
0747   else if (!recojet)
0748   {
0749     ++_errors;
0750     return std::set<Jet*>();
0751   }
0752 
0753   if (_do_cache)
0754   {
0755     std::map<Jet*, std::set<Jet*> >::iterator iter =
0756         _cache_all_truth_jets.find(recojet);
0757     if (iter != _cache_all_truth_jets.end())
0758     {
0759       return iter->second;
0760     }
0761   }
0762 
0763   std::set<Jet*> truth_jets;
0764 
0765   // get all truth particles (this can include muons and other truth excludes)...
0766   std::set<PHG4Particle*> particles = all_truth_particles(recojet);
0767 
0768   // backtrack from the truth particles to the truth jets...
0769   for (auto* particle : particles)
0770   {
0771     if (_strict)
0772     {
0773       assert(particle);
0774     }
0775     else if (!particle)
0776     {
0777       ++_errors;
0778       continue;
0779     }
0780 
0781     Jet* truth_jet = _jettrutheval.get_truth_jet(particle);
0782     if (!truth_jet)
0783     {
0784       continue;
0785     }
0786 
0787     truth_jets.insert(truth_jet);
0788   }
0789 
0790   if (_do_cache)
0791   {
0792     _cache_all_truth_jets.insert(std::make_pair(recojet, truth_jets));
0793   }
0794 
0795   return truth_jets;
0796 }
0797 
0798 Jet* JetRecoEval::max_truth_jet_by_energy(Jet* recojet)
0799 {
0800   if (_strict)
0801   {
0802     assert(recojet);
0803   }
0804   else if (!recojet)
0805   {
0806     ++_errors;
0807     return nullptr;
0808   }
0809 
0810   if (_do_cache)
0811   {
0812     std::map<Jet*, Jet*>::iterator iter =
0813         _cache_max_truth_jet_by_energy.find(recojet);
0814     if (iter != _cache_max_truth_jet_by_energy.end())
0815     {
0816       return iter->second;
0817     }
0818   }
0819 
0820   Jet* truthjet = nullptr;
0821   float max_energy = std::numeric_limits<float>::min();
0822 
0823   std::set<Jet*> truthjets = all_truth_jets(recojet);
0824   for (auto* candidate : truthjets)
0825   {
0826     if (_strict)
0827     {
0828       assert(candidate);
0829     }
0830     else if (!candidate)
0831     {
0832       ++_errors;
0833       continue;
0834     }
0835 
0836     float energy = get_energy_contribution(recojet, candidate);
0837     if (energy > max_energy)
0838     {
0839       truthjet = candidate;
0840       max_energy = energy;
0841     }
0842   }
0843 
0844   if (_do_cache)
0845   {
0846     _cache_max_truth_jet_by_energy.insert(std::make_pair(recojet, truthjet));
0847   }
0848 
0849   return truthjet;
0850 }
0851 
0852 std::set<Jet*> JetRecoEval::all_jets_from(Jet* truthjet)
0853 {
0854   if (_strict)
0855   {
0856     assert(truthjet);
0857   }
0858   else if (!truthjet)
0859   {
0860     ++_errors;
0861     return std::set<Jet*>();
0862   }
0863 
0864   if (_do_cache)
0865   {
0866     std::map<Jet*, std::set<Jet*> >::iterator iter =
0867         _cache_all_jets_from.find(truthjet);
0868     if (iter != _cache_all_jets_from.end())
0869     {
0870       return iter->second;
0871     }
0872   }
0873 
0874   std::set<Jet*> recojets;
0875 
0876   // loop over all reco jets
0877   for (auto* recojet : *_recojets)
0878   {
0879     /* Jet* recojet = _recojet.second; */
0880 
0881     // if this jet back tracks to the truth jet
0882     std::set<Jet*> truthcandidates = all_truth_jets(recojet);
0883     for (auto* truthcandidate : truthcandidates)
0884     {
0885       if (_strict)
0886       {
0887         assert(truthcandidate);
0888       }
0889       else if (!truthcandidate)
0890       {
0891         ++_errors;
0892         continue;
0893       }
0894 
0895       if (truthcandidate->get_id() == truthjet->get_id())
0896       {
0897         recojets.insert(recojet);
0898       }
0899     }
0900   }
0901 
0902   if (_do_cache)
0903   {
0904     _cache_all_jets_from.insert(std::make_pair(truthjet, recojets));
0905   }
0906 
0907   return recojets;
0908 }
0909 
0910 Jet* JetRecoEval::best_jet_from(Jet* truthjet)
0911 {
0912   if (_strict)
0913   {
0914     assert(truthjet);
0915   }
0916   else if (!truthjet)
0917   {
0918     ++_errors;
0919     return nullptr;
0920   }
0921 
0922   if (_do_cache)
0923   {
0924     std::map<Jet*, Jet*>::iterator iter =
0925         _cache_best_jet_from.find(truthjet);
0926     if (iter != _cache_best_jet_from.end())
0927     {
0928       return iter->second;
0929     }
0930   }
0931 
0932   Jet* bestrecojet = nullptr;
0933   float max_energy = std::numeric_limits<float>::min();
0934 
0935   std::set<Jet*> recojets = all_jets_from(truthjet);
0936   for (auto* recojet : recojets)
0937   {
0938     if (_strict)
0939     {
0940       assert(recojet);
0941     }
0942     else if (!recojet)
0943     {
0944       ++_errors;
0945       continue;
0946     }
0947 
0948     float energy = get_energy_contribution(recojet, truthjet);
0949     if (energy > max_energy)
0950     {
0951       bestrecojet = recojet;
0952       max_energy = energy;
0953     }
0954   }
0955 
0956   if (_do_cache)
0957   {
0958     _cache_best_jet_from.insert(std::make_pair(truthjet, bestrecojet));
0959   }
0960   return bestrecojet;
0961 }
0962 
0963 Jet* JetRecoEval::unique_reco_jet_from_truth(Jet* truthjet)
0964 {
0965   if (_strict)
0966   {
0967     assert(truthjet);
0968   }
0969   else if (!truthjet)
0970   {
0971     ++_errors;
0972     return nullptr;
0973   }
0974 
0975   Jet* recojet = best_jet_from(truthjet);
0976 
0977   if (recojet)
0978   {
0979     Jet* back_matching = max_truth_jet_by_energy(recojet);
0980 
0981     if (back_matching->get_id() == truthjet->get_id())
0982     {
0983       return recojet;  // uniquely matched
0984     }
0985 
0986     return nullptr;
0987   }
0988 
0989   return nullptr;
0990 }
0991 
0992 Jet* JetRecoEval::unique_truth_jet_from_reco(Jet* recojet)
0993 {
0994   if (_strict)
0995   {
0996     assert(recojet);
0997   }
0998   else if (!recojet)
0999   {
1000     ++_errors;
1001     return nullptr;
1002   }
1003 
1004   Jet* truthjet = max_truth_jet_by_energy(recojet);
1005 
1006   if (truthjet)
1007   {
1008     Jet* back_matching = best_jet_from(truthjet);
1009 
1010     if (back_matching->get_id() == recojet->get_id())
1011     {
1012       return truthjet;  // uniquely matched
1013     }
1014 
1015     return nullptr;
1016   }
1017 
1018   return nullptr;
1019 }
1020 
1021 // overlap calculations
1022 float JetRecoEval::get_energy_contribution(Jet* recojet, Jet* truthjet)
1023 {
1024   if (_strict)
1025   {
1026     assert(recojet);
1027     assert(truthjet);
1028   }
1029   else if (!recojet || !truthjet)
1030   {
1031     ++_errors;
1032     return std::numeric_limits<float>::quiet_NaN();
1033   }
1034 
1035   if (_do_cache)
1036   {
1037     std::map<std::pair<Jet*, Jet*>, float>::iterator iter =
1038         _cache_get_energy_contribution.find(std::make_pair(recojet, truthjet));
1039     if (iter != _cache_get_energy_contribution.end())
1040     {
1041       return iter->second;
1042     }
1043   }
1044 
1045   float energy_contribution = 0.0;
1046 
1047   std::set<PHG4Particle*> truthjetcomp = get_truth_eval()->all_truth_particles(truthjet);
1048   // loop over all truthjet constituents
1049   for (auto* truthparticle : truthjetcomp)
1050   {
1051     if (_strict)
1052     {
1053       assert(truthparticle);
1054     }
1055     else if (!truthparticle)
1056     {
1057       ++_errors;
1058       continue;
1059     }
1060 
1061     for (auto jter : recojet->get_comp_vec())
1062     {
1063       Jet::SRC source = jter.first;
1064       unsigned int index = jter.second;
1065 
1066       float energy = 0.0;
1067 
1068       if (source == Jet::TRACK)
1069       {
1070         SvtxTrack* track = _trackmap->get(index);
1071 
1072         if (_strict)
1073         {
1074           assert(track);
1075         }
1076         else if (!track)
1077         {
1078           ++_errors;
1079           continue;
1080         }
1081 
1082         PHG4Particle* maxtruthparticle = get_svtx_eval_stack()->get_track_eval()->max_truth_particle_by_nclusters(track);
1083 
1084         if (maxtruthparticle == nullptr)
1085         {
1086           // in extreme rare cases, noise hits can make a track with no maxtruthparticle matched
1087           energy = 0;
1088         }
1089         else if (maxtruthparticle->get_track_id() == truthparticle->get_track_id())
1090         {
1091           energy = track->get_p();
1092         }
1093       }
1094       else if (source == Jet::CEMC_TOWER)
1095       {
1096         RawTower* tower = _cemctowers->getTower(index);
1097 
1098         if (_strict)
1099         {
1100           assert(tower);
1101         }
1102         else if (!tower)
1103         {
1104           ++_errors;
1105           continue;
1106         }
1107 
1108         energy = get_cemc_eval_stack()->get_rawtower_eval()->get_energy_contribution(tower, truthparticle);
1109       }
1110       else if (source == Jet::CEMC_CLUSTER)
1111       {
1112         RawCluster* cluster = _cemcclusters->getCluster(index);
1113 
1114         if (_strict)
1115         {
1116           assert(cluster);
1117         }
1118         else if (!cluster)
1119         {
1120           ++_errors;
1121           continue;
1122         }
1123 
1124         energy = get_cemc_eval_stack()->get_rawcluster_eval()->get_energy_contribution(cluster, truthparticle);
1125       }
1126       else if (source == Jet::EEMC_TOWER)
1127       {
1128         RawTower* tower = _eemctowers->getTower(index);
1129 
1130         if (_strict)
1131         {
1132           assert(tower);
1133         }
1134         else if (!tower)
1135         {
1136           ++_errors;
1137           continue;
1138         }
1139 
1140         energy = get_eemc_eval_stack()->get_rawtower_eval()->get_energy_contribution(tower, truthparticle);
1141       }
1142       else if (source == Jet::EEMC_CLUSTER)
1143       {
1144         RawCluster* cluster = _eemcclusters->getCluster(index);
1145 
1146         if (_strict)
1147         {
1148           assert(cluster);
1149         }
1150         else if (!cluster)
1151         {
1152           ++_errors;
1153           continue;
1154         }
1155 
1156         energy = get_eemc_eval_stack()->get_rawcluster_eval()->get_energy_contribution(cluster, truthparticle);
1157       }
1158       else if (source == Jet::HCALIN_TOWER)
1159       {
1160         RawTower* tower = _hcalintowers->getTower(index);
1161 
1162         if (_strict)
1163         {
1164           assert(tower);
1165         }
1166         else if (!tower)
1167         {
1168           ++_errors;
1169           continue;
1170         }
1171 
1172         energy = get_hcalin_eval_stack()->get_rawtower_eval()->get_energy_contribution(tower, truthparticle);
1173       }
1174       else if (source == Jet::HCALIN_CLUSTER)
1175       {
1176         RawCluster* cluster = _hcalinclusters->getCluster(index);
1177 
1178         if (_strict)
1179         {
1180           assert(cluster);
1181         }
1182         else if (!cluster)
1183         {
1184           ++_errors;
1185           continue;
1186         }
1187 
1188         energy = get_hcalin_eval_stack()->get_rawcluster_eval()->get_energy_contribution(cluster, truthparticle);
1189       }
1190       else if (source == Jet::HCALOUT_TOWER)
1191       {
1192         RawTower* tower = _hcalouttowers->getTower(index);
1193 
1194         if (_strict)
1195         {
1196           assert(tower);
1197         }
1198         else if (!tower)
1199         {
1200           ++_errors;
1201           continue;
1202         }
1203 
1204         energy = get_hcalout_eval_stack()->get_rawtower_eval()->get_energy_contribution(tower, truthparticle);
1205       }
1206       else if (source == Jet::HCALOUT_CLUSTER)
1207       {
1208         RawCluster* cluster = _hcaloutclusters->getCluster(index);
1209 
1210         if (_strict)
1211         {
1212           assert(cluster);
1213         }
1214         else if (!cluster)
1215         {
1216           ++_errors;
1217           continue;
1218         }
1219 
1220         energy = get_hcalout_eval_stack()->get_rawcluster_eval()->get_energy_contribution(cluster, truthparticle);
1221       }
1222       else if (source == Jet::FEMC_TOWER)
1223       {
1224         RawTower* tower = _femctowers->getTower(index);
1225 
1226         if (_strict)
1227         {
1228           assert(tower);
1229         }
1230         else if (!tower)
1231         {
1232           ++_errors;
1233           continue;
1234         }
1235 
1236         energy = get_femc_eval_stack()->get_rawtower_eval()->get_energy_contribution(tower, truthparticle);
1237       }
1238       else if (source == Jet::FEMC_CLUSTER)
1239       {
1240         RawCluster* cluster = _femcclusters->getCluster(index);
1241 
1242         if (_strict)
1243         {
1244           assert(cluster);
1245         }
1246         else if (!cluster)
1247         {
1248           ++_errors;
1249           continue;
1250         }
1251 
1252         energy = get_femc_eval_stack()->get_rawcluster_eval()->get_energy_contribution(cluster, truthparticle);
1253       }
1254       else if (source == Jet::FHCAL_TOWER)
1255       {
1256         RawTower* tower = _fhcaltowers->getTower(index);
1257 
1258         if (_strict)
1259         {
1260           assert(tower);
1261         }
1262         else if (!tower)
1263         {
1264           ++_errors;
1265           continue;
1266         }
1267 
1268         energy = get_fhcal_eval_stack()->get_rawtower_eval()->get_energy_contribution(tower, truthparticle);
1269       }
1270       else if (source == Jet::FHCAL_CLUSTER)
1271       {
1272         RawCluster* cluster = _fhcalclusters->getCluster(index);
1273 
1274         if (_strict)
1275         {
1276           assert(cluster);
1277         }
1278         else if (!cluster)
1279         {
1280           ++_errors;
1281           continue;
1282         }
1283 
1284         energy = get_fhcal_eval_stack()->get_rawcluster_eval()->get_energy_contribution(cluster, truthparticle);
1285       }
1286 
1287       energy_contribution += energy;
1288     }
1289   }
1290 
1291   if (_do_cache)
1292   {
1293     _cache_get_energy_contribution.insert(std::make_pair(std::make_pair(recojet, truthjet), energy_contribution));
1294   }
1295 
1296   return energy_contribution;
1297 }
1298 
1299 // overlap calculations
1300 float JetRecoEval::get_energy_contribution(Jet* recojet, Jet::SRC src)
1301 {
1302   if (_strict)
1303   {
1304     assert(recojet);
1305   }
1306   else if (!recojet)
1307   {
1308     ++_errors;
1309     return std::numeric_limits<float>::quiet_NaN();
1310   }
1311 
1312   if (_do_cache)
1313   {
1314     std::map<std::pair<Jet*, Jet::SRC>, float>::iterator iter =
1315         _cache_get_energy_contribution_src.find(std::make_pair(recojet, src));
1316     if (iter != _cache_get_energy_contribution_src.end())
1317     {
1318       return iter->second;
1319     }
1320   }
1321 
1322   float energy = 0.0;
1323 
1324   // loop over all recojet constituents
1325   for (Jet::ITER_comp_vec jter = recojet->comp_begin(src);
1326        jter != recojet->comp_end(src);
1327        ++jter)
1328   {
1329     Jet::SRC source = jter->first;
1330     assert(source == src);  // jet container consistency check
1331     unsigned int index = jter->second;
1332 
1333     if (source == Jet::TRACK)
1334     {
1335       SvtxTrack* track = _trackmap->get(index);
1336       energy += track->get_p();
1337     }
1338     else if (source == Jet::CEMC_TOWER)
1339     {
1340       RawTower* tower = _cemctowers->getTower(index);
1341 
1342       if (_strict)
1343       {
1344         assert(tower);
1345       }
1346       else if (!tower)
1347       {
1348         ++_errors;
1349         continue;
1350       }
1351 
1352       energy += tower->get_energy();
1353     }
1354     else if (source == Jet::CEMC_CLUSTER)
1355     {
1356       RawCluster* cluster = _cemcclusters->getCluster(index);
1357 
1358       if (_strict)
1359       {
1360         assert(cluster);
1361       }
1362       else if (!cluster)
1363       {
1364         ++_errors;
1365         continue;
1366       }
1367 
1368       energy += cluster->get_energy();
1369     }
1370     else if (source == Jet::EEMC_TOWER)
1371     {
1372       RawTower* tower = _eemctowers->getTower(index);
1373 
1374       if (_strict)
1375       {
1376         assert(tower);
1377       }
1378       else if (!tower)
1379       {
1380         ++_errors;
1381         continue;
1382       }
1383 
1384       energy += tower->get_energy();
1385     }
1386     else if (source == Jet::EEMC_CLUSTER)
1387     {
1388       RawCluster* cluster = _eemcclusters->getCluster(index);
1389 
1390       if (_strict)
1391       {
1392         assert(cluster);
1393       }
1394       else if (!cluster)
1395       {
1396         ++_errors;
1397         continue;
1398       }
1399 
1400       energy += cluster->get_energy();
1401     }
1402     else if (source == Jet::HCALIN_TOWER)
1403     {
1404       RawTower* tower = _hcalintowers->getTower(index);
1405 
1406       if (_strict)
1407       {
1408         assert(tower);
1409       }
1410       else if (!tower)
1411       {
1412         ++_errors;
1413         continue;
1414       }
1415 
1416       energy += tower->get_energy();
1417     }
1418     else if (source == Jet::HCALIN_CLUSTER)
1419     {
1420       RawCluster* cluster = _hcalinclusters->getCluster(index);
1421 
1422       if (_strict)
1423       {
1424         assert(cluster);
1425       }
1426       else if (!cluster)
1427       {
1428         ++_errors;
1429         continue;
1430       }
1431 
1432       energy += cluster->get_energy();
1433     }
1434     else if (source == Jet::HCALOUT_TOWER)
1435     {
1436       RawTower* tower = _hcalouttowers->getTower(index);
1437 
1438       if (_strict)
1439       {
1440         assert(tower);
1441       }
1442       else if (!tower)
1443       {
1444         ++_errors;
1445         continue;
1446       }
1447 
1448       energy += tower->get_energy();
1449     }
1450     else if (source == Jet::HCALOUT_CLUSTER)
1451     {
1452       RawCluster* cluster = _hcaloutclusters->getCluster(index);
1453 
1454       if (_strict)
1455       {
1456         assert(cluster);
1457       }
1458       else if (!cluster)
1459       {
1460         ++_errors;
1461         continue;
1462       }
1463 
1464       energy += cluster->get_energy();
1465     }
1466     else if (source == Jet::FEMC_TOWER)
1467     {
1468       RawTower* tower = _femctowers->getTower(index);
1469 
1470       if (_strict)
1471       {
1472         assert(tower);
1473       }
1474       else if (!tower)
1475       {
1476         ++_errors;
1477         continue;
1478       }
1479 
1480       energy += tower->get_energy();
1481     }
1482     else if (source == Jet::FEMC_CLUSTER)
1483     {
1484       RawCluster* cluster = _femcclusters->getCluster(index);
1485 
1486       if (_strict)
1487       {
1488         assert(cluster);
1489       }
1490       else if (!cluster)
1491       {
1492         ++_errors;
1493         continue;
1494       }
1495 
1496       energy += cluster->get_energy();
1497     }
1498     else if (source == Jet::FHCAL_TOWER)
1499     {
1500       RawTower* tower = _fhcaltowers->getTower(index);
1501 
1502       if (_strict)
1503       {
1504         assert(tower);
1505       }
1506       else if (!tower)
1507       {
1508         ++_errors;
1509         continue;
1510       }
1511 
1512       energy += tower->get_energy();
1513     }
1514     else if (source == Jet::FHCAL_CLUSTER)
1515     {
1516       RawCluster* cluster = _fhcalclusters->getCluster(index);
1517 
1518       if (_strict)
1519       {
1520         assert(cluster);
1521       }
1522       else if (!cluster)
1523       {
1524         ++_errors;
1525         continue;
1526       }
1527 
1528       energy += cluster->get_energy();
1529     }  // else if (source == Jet::FHCAL_CLUSTER)
1530 
1531   }  // for (Jet::ConstIter jter = recojet->lower_bound_comp(src);
1532 
1533   if (_do_cache)
1534   {
1535     _cache_get_energy_contribution_src.insert(std::make_pair(std::make_pair(recojet, src), energy));
1536   }
1537 
1538   return energy;
1539 }
1540 
1541 std::set<PHG4Hit*> JetRecoEval::all_truth_hits(Jet* recojet)
1542 {
1543   if (_strict)
1544   {
1545     assert(recojet);
1546   }
1547   else if (!recojet)
1548   {
1549     ++_errors;
1550     return std::set<PHG4Hit*>();
1551   }
1552 
1553   if (_do_cache)
1554   {
1555     std::map<Jet*, std::set<PHG4Hit*> >::iterator iter =
1556         _cache_all_truth_hits.find(recojet);
1557     if (iter != _cache_all_truth_hits.end())
1558     {
1559       return iter->second;
1560     }
1561   }
1562 
1563   std::set<PHG4Hit*> truth_hits;
1564 
1565   // loop over all the jet constituents, backtrack each reco object to the
1566   // truth hits and combine with other consituents
1567 
1568   for (auto jter : recojet->get_comp_vec())
1569   {
1570     Jet::SRC source = jter.first;
1571     unsigned int index = jter.second;
1572 
1573     std::set<PHG4Hit*> new_hits;
1574 
1575     if (source == Jet::TRACK)
1576     {
1577       if (!_trackmap)
1578       {
1579         std::cout << PHWHERE << "ERROR: can't find SvtxTrackMap" << std::endl;
1580         exit(-1);
1581       }
1582 
1583       SvtxTrack* track = _trackmap->get(index);
1584 
1585       if (_strict)
1586       {
1587         assert(track);
1588       }
1589       else if (!track)
1590       {
1591         ++_errors;
1592         continue;
1593       }
1594 
1595       new_hits = get_svtx_eval_stack()->get_track_eval()->all_truth_hits(track);
1596     }
1597     else if (source == Jet::CEMC_TOWER)
1598     {
1599       if (!_cemctowers)
1600       {
1601         std::cout << PHWHERE << "ERROR: can't find TOWER_CEMC" << std::endl;
1602         exit(-1);
1603       }
1604 
1605       RawTower* tower = _cemctowers->getTower(index);
1606 
1607       if (_strict)
1608       {
1609         assert(tower);
1610       }
1611       else if (!tower)
1612       {
1613         ++_errors;
1614         continue;
1615       }
1616 
1617       new_hits = get_cemc_eval_stack()->get_rawtower_eval()->all_truth_hits(tower);
1618     }
1619     else if (source == Jet::CEMC_CLUSTER)
1620     {
1621       if (!_cemcclusters)
1622       {
1623         std::cout << PHWHERE << "ERROR: can't find CLUSTER_CEMC" << std::endl;
1624         exit(-1);
1625       }
1626 
1627       RawCluster* cluster = _cemcclusters->getCluster(index);
1628 
1629       if (_strict)
1630       {
1631         assert(cluster);
1632       }
1633       else if (!cluster)
1634       {
1635         ++_errors;
1636         continue;
1637       }
1638 
1639       new_hits = get_cemc_eval_stack()->get_rawcluster_eval()->all_truth_hits(cluster);
1640     }
1641     else if (source == Jet::EEMC_TOWER)
1642     {
1643       if (!_eemctowers)
1644       {
1645         std::cout << PHWHERE << "ERROR: can't find TOWER_EEMC" << std::endl;
1646         exit(-1);
1647       }
1648 
1649       RawTower* tower = _eemctowers->getTower(index);
1650 
1651       if (_strict)
1652       {
1653         assert(tower);
1654       }
1655       else if (!tower)
1656       {
1657         ++_errors;
1658         continue;
1659       }
1660 
1661       new_hits = get_eemc_eval_stack()->get_rawtower_eval()->all_truth_hits(tower);
1662     }
1663     else if (source == Jet::EEMC_CLUSTER)
1664     {
1665       if (!_eemcclusters)
1666       {
1667         std::cout << PHWHERE << "ERROR: can't find CLUSTER_EEMC" << std::endl;
1668         exit(-1);
1669       }
1670 
1671       RawCluster* cluster = _eemcclusters->getCluster(index);
1672 
1673       if (_strict)
1674       {
1675         assert(cluster);
1676       }
1677       else if (!cluster)
1678       {
1679         ++_errors;
1680         continue;
1681       }
1682 
1683       new_hits = get_eemc_eval_stack()->get_rawcluster_eval()->all_truth_hits(cluster);
1684     }
1685     else if (source == Jet::HCALIN_TOWER)
1686     {
1687       if (!_hcalintowers)
1688       {
1689         std::cout << PHWHERE << "ERROR: can't find TOWER_HCALIN" << std::endl;
1690         exit(-1);
1691       }
1692 
1693       RawTower* tower = _hcalintowers->getTower(index);
1694 
1695       if (_strict)
1696       {
1697         assert(tower);
1698       }
1699       else if (!tower)
1700       {
1701         ++_errors;
1702         continue;
1703       }
1704 
1705       new_hits = get_hcalin_eval_stack()->get_rawtower_eval()->all_truth_hits(tower);
1706     }
1707     else if (source == Jet::HCALIN_CLUSTER)
1708     {
1709       if (!_hcalinclusters)
1710       {
1711         std::cout << PHWHERE << "ERROR: can't find CLUSTER_HCALIN" << std::endl;
1712         exit(-1);
1713       }
1714 
1715       RawCluster* cluster = _hcalinclusters->getCluster(index);
1716 
1717       if (_strict)
1718       {
1719         assert(cluster);
1720       }
1721       else if (!cluster)
1722       {
1723         ++_errors;
1724         continue;
1725       }
1726 
1727       new_hits = get_hcalin_eval_stack()->get_rawcluster_eval()->all_truth_hits(cluster);
1728     }
1729     else if (source == Jet::HCALOUT_TOWER)
1730     {
1731       if (!_hcalouttowers)
1732       {
1733         std::cout << PHWHERE << "ERROR: can't find TOWER_HCALOUT" << std::endl;
1734         exit(-1);
1735       }
1736 
1737       RawTower* tower = _hcalouttowers->getTower(index);
1738 
1739       if (_strict)
1740       {
1741         assert(tower);
1742       }
1743       else if (!tower)
1744       {
1745         ++_errors;
1746         continue;
1747       }
1748 
1749       new_hits = get_hcalout_eval_stack()->get_rawtower_eval()->all_truth_hits(tower);
1750     }
1751     else if (source == Jet::HCALOUT_CLUSTER)
1752     {
1753       if (!_hcaloutclusters)
1754       {
1755         std::cout << PHWHERE << "ERROR: can't find CLUSTER_HCALOUT" << std::endl;
1756         exit(-1);
1757       }
1758 
1759       RawCluster* cluster = _hcaloutclusters->getCluster(index);
1760 
1761       if (_strict)
1762       {
1763         assert(cluster);
1764       }
1765       else if (!cluster)
1766       {
1767         ++_errors;
1768         continue;
1769       }
1770 
1771       new_hits = get_hcalout_eval_stack()->get_rawcluster_eval()->all_truth_hits(cluster);
1772     }
1773     else if (source == Jet::FEMC_TOWER)
1774     {
1775       if (!_femctowers)
1776       {
1777         std::cout << PHWHERE << "ERROR: can't find TOWER_FEMC" << std::endl;
1778         exit(-1);
1779       }
1780 
1781       RawTower* tower = _femctowers->getTower(index);
1782 
1783       if (_strict)
1784       {
1785         assert(tower);
1786       }
1787       else if (!tower)
1788       {
1789         ++_errors;
1790         continue;
1791       }
1792 
1793       new_hits = get_femc_eval_stack()->get_rawtower_eval()->all_truth_hits(tower);
1794     }
1795     else if (source == Jet::FEMC_CLUSTER)
1796     {
1797       if (!_femcclusters)
1798       {
1799         std::cout << PHWHERE << "ERROR: can't find CLUSTER_FEMC" << std::endl;
1800         exit(-1);
1801       }
1802 
1803       RawCluster* cluster = _femcclusters->getCluster(index);
1804 
1805       if (_strict)
1806       {
1807         assert(cluster);
1808       }
1809       else if (!cluster)
1810       {
1811         ++_errors;
1812         continue;
1813       }
1814 
1815       new_hits = get_femc_eval_stack()->get_rawcluster_eval()->all_truth_hits(cluster);
1816     }
1817     else if (source == Jet::FHCAL_TOWER)
1818     {
1819       if (!_fhcaltowers)
1820       {
1821         std::cout << PHWHERE << "ERROR: can't find TOWER_FHCAL" << std::endl;
1822         exit(-1);
1823       }
1824 
1825       RawTower* tower = _fhcaltowers->getTower(index);
1826 
1827       if (_strict)
1828       {
1829         assert(tower);
1830       }
1831       else if (!tower)
1832       {
1833         ++_errors;
1834         continue;
1835       }
1836 
1837       new_hits = get_fhcal_eval_stack()->get_rawtower_eval()->all_truth_hits(tower);
1838     }
1839     else if (source == Jet::FHCAL_CLUSTER)
1840     {
1841       if (!_fhcalclusters)
1842       {
1843         std::cout << PHWHERE << "ERROR: can't find CLUSTER_FHCAL" << std::endl;
1844         exit(-1);
1845       }
1846 
1847       RawCluster* cluster = _fhcalclusters->getCluster(index);
1848 
1849       if (_strict)
1850       {
1851         assert(cluster);
1852       }
1853       else if (!cluster)
1854       {
1855         ++_errors;
1856         continue;
1857       }
1858 
1859       new_hits = get_fhcal_eval_stack()->get_rawcluster_eval()->all_truth_hits(cluster);
1860     }
1861 
1862     for (auto* new_hit : new_hits)
1863     {
1864       truth_hits.insert(new_hit);
1865     }
1866   }
1867 
1868   if (_do_cache)
1869   {
1870     _cache_all_truth_hits.insert(std::make_pair(recojet, truth_hits));
1871   }
1872 
1873   return truth_hits;
1874 }
1875 
1876 void JetRecoEval::get_node_pointers(PHCompositeNode* topNode)
1877 {
1878   // need things off of the DST...
1879   _recojets = findNode::getClass<JetContainer>(topNode, _recojetname);
1880   if (!_recojets)
1881   {
1882     std::cout << PHWHERE << " ERROR: Can't find " << _recojetname << std::endl;
1883     exit(-1);
1884   }
1885 
1886   _truthjets = findNode::getClass<JetContainer>(topNode, _truthjetname);
1887   if (!_truthjets)
1888   {
1889     std::cout << PHWHERE << " ERROR: Can't find " << _truthjetname << std::endl;
1890     exit(-1);
1891   }
1892 
1893   _trackmap = findNode::getClass<SvtxTrackMap>(topNode, m_TrackNodeName);
1894   if (!_trackmap)
1895   {
1896     _trackmap = findNode::getClass<SvtxTrackMap>(topNode, "TrackMap");
1897   }
1898   _cemctowers = findNode::getClass<RawTowerContainer>(topNode, "TOWER_CALIB_CEMC");
1899   _hcalintowers = findNode::getClass<RawTowerContainer>(topNode, "TOWER_CALIB_HCALIN");
1900   _hcalouttowers = findNode::getClass<RawTowerContainer>(topNode, "TOWER_CALIB_HCALOUT");
1901   _femctowers = findNode::getClass<RawTowerContainer>(topNode, "TOWER_CALIB_FEMC");
1902   _fhcaltowers = findNode::getClass<RawTowerContainer>(topNode, "TOWER_CALIB_FHCAL");
1903   _eemctowers = findNode::getClass<RawTowerContainer>(topNode, "TOWER_CALIB_EEMC");
1904   _cemcclusters = findNode::getClass<RawClusterContainer>(topNode, "CLUSTER_CEMC");
1905   _hcalinclusters = findNode::getClass<RawClusterContainer>(topNode, "CLUSTER_HCALIN");
1906   _hcaloutclusters = findNode::getClass<RawClusterContainer>(topNode, "CLUSTER_HCALOUT");
1907   _femcclusters = findNode::getClass<RawClusterContainer>(topNode, "CLUSTER_FEMC");
1908   _fhcalclusters = findNode::getClass<RawClusterContainer>(topNode, "CLUSTER_FHCAL");
1909   _eemcclusters = findNode::getClass<RawClusterContainer>(topNode, "CLUSTER_EEMC");
1910 
1911   return;
1912 }