g4simulation/g4main/PHG4ConsistencyCheck.cc

0001 #include "PHG4ConsistencyCheck.h"
0002
0003 #include "PHG4Hit.h"
0004 #include "PHG4HitContainer.h"
0005 #include "PHG4Particle.h"
0006 #include "PHG4TruthInfoContainer.h"
0007
0008 #include <fun4all/Fun4AllReturnCodes.h>
0009
0010 #include <phool/getClass.h>
0011
0012 #include <algorithm>
0013 #include <iostream>  // for operator<<, basic_ostream::opera...
0014 #include <map>       // for _Rb_tree_const_iterator, map<>::...
0015 #include <set>       // for set
0016 #include <utility>   // for pair
0017
0018 class PHCompositeNode;
0019
0020 PHG4ConsistencyCheck::PHG4ConsistencyCheck(const std::string &name)
0021   : SubsysReco(name)
0022 {
0023 }
0024
0025 int PHG4ConsistencyCheck::process_event(PHCompositeNode *topNode)
0026 {
0027   PHG4TruthInfoContainer *truthcont = findNode::getClass<PHG4TruthInfoContainer>(topNode, "G4TruthInfo");
0028   if (!truthcont)
0029   {
0030     return 0;
0031   }
0032   PHG4TruthInfoContainer::ConstRange trange = truthcont->GetParticleRange();
0033   PHG4TruthInfoContainer::ConstIterator titer;
0034   int imax = 1000000;
0035   for (titer = trange.first; titer != trange.second; ++titer)
0036   {
0037     imax = std::min(titer->first, imax);
0038   }
0039   std::cout << "min index: " << imax << std::endl;
0040   std::pair<std::map<int, int>::const_iterator, std::map<int, int>::const_iterator> embtrk_b_e = truthcont->GetEmbeddedTrkIds();
0041   std::map<int, int>::const_iterator embiter;
0042   for (embiter = embtrk_b_e.first; embiter != embtrk_b_e.second; ++embiter)
0043   {
0044     std::cout << "embedded trkid: " << embiter->first << std::endl;
0045   }
0046   PHG4HitContainer *ghit = findNode::getClass<PHG4HitContainer>(topNode, "G4HIT_CEMC_E");
0047   if (ghit)
0048   {
0049     PHG4HitContainer::ConstIterator hit;
0050     PHG4HitContainer::ConstRange hit_begin_end = ghit->getHits();
0051     std::set<int> printpart;
0052     for (hit = hit_begin_end.first; hit != hit_begin_end.second; ++hit)
0053     {
0054       int trkid = hit->second->get_trkid();
0055       PHG4Particle *part = truthcont->GetParticle(trkid);
0056       if (!part)
0057       {
0058         hit->second->identify();
0059         std::cout << "could not locate geant particle " << trkid << " in G4HIT_CEMC_E" << std::endl;
0060         errorcnt++;
0061       }
0062       else
0063       {
0064         int primary_id = part->get_primary_id();
0065         if (truthcont->isEmbeded(primary_id) > 0)
0066         {
0067           if (!printpart.contains(primary_id))
0068           {
0069             std::cout << "primary id " << primary_id << " is embedded" << std::endl;
0070             printpart.insert(primary_id);
0071             PHG4Particle *parta = truthcont->GetParticle(primary_id);
0072             parta->identify();
0073           }
0074         }
0075       }
0076     }
0077   }
0078   ghit = findNode::getClass<PHG4HitContainer>(topNode, "G4HIT_SVTX");
0079   if (ghit)
0080   {
0081     PHG4HitContainer::ConstIterator hit;
0082     PHG4HitContainer::ConstRange hit_begin_end = ghit->getHits();
0083     for (hit = hit_begin_end.first; hit != hit_begin_end.second; ++hit)
0084     {
0085       int trkid = hit->second->get_trkid();
0086       PHG4Particle *part = truthcont->GetParticle(trkid);
0087       if (!part)
0088       {
0089         std::cout << "could not locate geant particle " << trkid << " in G4HIT_SVTX" << std::endl;
0090         errorcnt++;
0091       }
0092     }
0093   }
0094
0095   return Fun4AllReturnCodes::EVENT_OK;
0096 }