File indexing completed on 2025-10-16 08:21:02
0001 #include "TruthNeutralMesonBuilder.h"
0002
0003 #include "TruthNeutralMesonContainer.h"
0004 #include "TruthNeutralMesonv1.h"
0005
0006 #include <fun4all/Fun4AllReturnCodes.h>
0007
0008 #include <g4main/PHG4Particle.h>
0009 #include <g4main/PHG4TruthInfoContainer.h>
0010
0011 #include <phool/PHCompositeNode.h>
0012 #include <phool/PHIODataNode.h>
0013 #include <phool/PHNodeIterator.h>
0014 #include <phool/getClass.h>
0015 #include <phool/phool.h>
0016
0017 #include <phhepmc/PHHepMCGenEvent.h>
0018 #include <phhepmc/PHHepMCGenEventMap.h>
0019
0020 #pragma GCC diagnostic push
0021 #pragma GCC diagnostic ignored "-Wdeprecated-declarations"
0022 #include <HepMC/GenEvent.h>
0023 #include <HepMC/GenVertex.h>
0024 #pragma GCC diagnostic pop
0025 #include <HepMC/GenParticle.h>
0026
0027 #include <cmath>
0028 #include <iostream>
0029 #include <utility>
0030 #include <vector>
0031
0032 namespace
0033 {
0034 inline float safe_pt(float px, float py)
0035 {
0036 return std::sqrt(px * px + py * py);
0037 }
0038 inline float safe_p(float px, float py, float pz)
0039 {
0040 return std::sqrt(px * px + py * py + pz * pz);
0041 }
0042 inline float safe_phi(float px, float py)
0043 {
0044 return std::atan2(py, px);
0045 }
0046 inline float safe_eta(float p, float pz)
0047 {
0048 const float denom = p - std::fabs(pz);
0049 if (p <= 0 || denom <= 0)
0050 {
0051 return 0;
0052 }
0053 const float num = p + std::fabs(pz);
0054 return 0.5 * std::log(num / denom) * (pz < 0 ? -1 : 1);
0055 }
0056 }
0057
0058 TruthNeutralMesonBuilder::TruthNeutralMesonBuilder(const std::string &name)
0059 : SubsysReco(name)
0060 {
0061 }
0062
0063 int TruthNeutralMesonBuilder::Init(PHCompositeNode *topNode)
0064 {
0065 CreateNodes(topNode);
0066
0067 return Fun4AllReturnCodes::EVENT_OK;
0068 }
0069
0070 int TruthNeutralMesonBuilder::CreateNodes(PHCompositeNode *topNode)
0071 {
0072 PHNodeIterator iter(topNode);
0073
0074 PHCompositeNode *dstNode = dynamic_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "DST"));
0075 if (!dstNode)
0076 {
0077 std::cout << PHWHERE << " DST node missing." << std::endl;
0078 dstNode = new PHCompositeNode("DST");
0079 topNode->addNode(dstNode);
0080 std::cout << "Added PHComposite node: DST" << std::endl;
0081 }
0082
0083 _container = findNode::getClass<TruthNeutralMesonContainer>(dstNode, m_truthnmeson_node_name);
0084 if (!_container)
0085 {
0086 _container = new TruthNeutralMesonContainer();
0087 PHIODataNode<PHObject> *node = new PHIODataNode<PHObject>(_container, m_truthnmeson_node_name, "PHObject");
0088 dstNode->addNode(node);
0089 std::cout << "Added TruthNeutralMesonContainer node" << std::endl;
0090 }
0091 return Fun4AllReturnCodes::EVENT_OK;
0092 }
0093
0094 int TruthNeutralMesonBuilder::process_event(PHCompositeNode *topNode)
0095 {
0096 genevtmap = findNode::getClass<PHHepMCGenEventMap>(topNode, "PHHepMCGenEventMap");
0097 truthinfo = findNode::getClass<PHG4TruthInfoContainer>(topNode, "G4TruthInfo");
0098 if (!genevtmap && !truthinfo)
0099 {
0100 std::cout << PHWHERE << " Missing both PHHepMCGenEventMap or G4TruthInfo. Abort run" << std::endl;
0101 return Fun4AllReturnCodes::ABORTRUN;
0102 }
0103
0104 _container = findNode::getClass<TruthNeutralMesonContainer>(topNode, m_truthnmeson_node_name);
0105 if (!_container)
0106 {
0107 std::cout << PHWHERE << " Missing output node " << m_truthnmeson_node_name << std::endl;
0108 return Fun4AllReturnCodes::ABORTRUN;
0109 }
0110 m_seen_mother_keys.clear();
0111 hepmc_by_barcode.clear();
0112 m_seen_barcodes.clear();
0113
0114 if (genevtmap)
0115 {
0116 for (auto &it : *genevtmap)
0117 {
0118 PHHepMCGenEvent *g = it.second;
0119 if (!g || !g->getEvent())
0120 {
0121 continue;
0122 }
0123 const int embedding_id = g->get_embedding_id();
0124 HepMC::GenEvent *evt = g->getEvent();
0125 for (HepMC::GenEvent::particle_const_iterator hepmcgenit = evt->particles_begin(); hepmcgenit != evt->particles_end(); ++hepmcgenit)
0126 {
0127 HepMC::GenParticle *p = *hepmcgenit;
0128 hepmc_by_barcode[p->barcode()] = std::make_pair(p, embedding_id);
0129 }
0130 }
0131 }
0132
0133 if (truthinfo)
0134 {
0135 g4_by_barcode.clear();
0136 g4_by_id.clear();
0137 g4_children.clear();
0138
0139 PHG4TruthInfoContainer::ConstRange range = truthinfo->GetParticleRange();
0140 for (auto it = range.first; it != range.second; ++it)
0141 {
0142 PHG4Particle *p = it->second;
0143 if (!p)
0144 {
0145 continue;
0146 }
0147 g4_by_id[p->get_track_id()] = p;
0148 g4_by_barcode[p->get_barcode()] = p;
0149 g4_children[p->get_parent_id()].push_back(p);
0150 }
0151 }
0152
0153 auto fill_kin_from_hepmc = [](HepMC::GenParticle *p, float &e, float &pt, float &eta, float &phi, float &magp)
0154 {
0155 const auto &m = p->momentum();
0156 const float px = m.px();
0157 const float py = m.py();
0158 const float pz = m.pz();
0159 const float ee = m.e();
0160 const float pp = safe_p(px, py, pz);
0161 e = ee;
0162 pt = safe_pt(px, py);
0163 phi = safe_phi(px, py);
0164 eta = safe_eta(pp, pz);
0165 magp = pp;
0166 };
0167
0168 auto fill_kin_from_g4 = [](PHG4Particle *p, float &e, float &pt, float &eta, float &phi, float &magp)
0169 {
0170 const float px = p->get_px();
0171 const float py = p->get_py();
0172 const float pz = p->get_pz();
0173 const float ee = p->get_e();
0174 const float pp = safe_p(px, py, pz);
0175 e = ee;
0176 pt = safe_pt(px, py);
0177 phi = safe_phi(px, py);
0178 eta = safe_eta(pp, pz);
0179 magp = pp;
0180 };
0181
0182 auto is_hadron = [](int id)
0183 {
0184 int a = std::abs(id);
0185 return (a > 100 && a < 100000);
0186 };
0187
0188 auto classify_origin = [&](int pid, int barcode, bool &is_prompt, int &parent_pid, bool &pi0_has_eta)
0189 {
0190 is_prompt = true;
0191 parent_pid = 0;
0192 pi0_has_eta = false;
0193
0194 HepMC::GenParticle *cur = nullptr;
0195 auto itH = hepmc_by_barcode.find(barcode);
0196 if (itH != hepmc_by_barcode.end())
0197 {
0198 cur = itH->second.first;
0199 }
0200
0201 while (cur)
0202 {
0203 HepMC::GenVertex *pv = cur->production_vertex();
0204 if (!pv)
0205 {
0206 break;
0207 }
0208
0209 if (pv->particles_in_const_begin() == pv->particles_in_const_end())
0210 {
0211 break;
0212 }
0213
0214 HepMC::GenParticle *par = *(pv->particles_in_const_begin());
0215 if (!par)
0216 {
0217 break;
0218 }
0219
0220 int tmp_parentid = par->pdg_id();
0221 if (!is_hadron(tmp_parentid))
0222 {
0223 break;
0224 }
0225
0226 if (is_prompt)
0227 {
0228 parent_pid = par->pdg_id();
0229 }
0230 is_prompt = false;
0231 if (pid == 111 && std::abs(parent_pid) == 221)
0232 {
0233 pi0_has_eta = true;
0234 }
0235 cur = par;
0236 }
0237
0238 if (!is_prompt)
0239 {
0240 return;
0241 }
0242
0243 PHG4Particle *g4 = nullptr;
0244 auto itG = g4_by_barcode.find(barcode);
0245 if (itG != g4_by_barcode.end())
0246 {
0247 g4 = itG->second;
0248 }
0249
0250 while (g4 && g4->get_parent_id() > 0)
0251 {
0252 PHG4Particle *par = nullptr;
0253 auto itP = g4_by_id.find(g4->get_parent_id());
0254 if (itP == g4_by_id.end())
0255 {
0256 break;
0257 }
0258 par = itP->second;
0259 if (is_prompt)
0260 {
0261 parent_pid = par->get_pid();
0262 }
0263 is_prompt = false;
0264 if (pid == 111 && std::abs(parent_pid) == 221)
0265 {
0266 pi0_has_eta = true;
0267 }
0268 g4 = par;
0269 }
0270 };
0271
0272 auto collect_photons_hepmc = [&](HepMC::GenParticle *mom, std::vector<HepMC::GenParticle *> &gammas)
0273 {
0274 gammas.clear();
0275 HepMC::GenVertex *vtx = mom->end_vertex();
0276 if (!vtx)
0277 {
0278 return;
0279 }
0280 for (HepMC::GenVertex::particles_out_const_iterator it = vtx->particles_out_const_begin(); it != vtx->particles_out_const_end(); ++it)
0281 {
0282 HepMC::GenParticle *ch = *it;
0283 if (!ch)
0284 {
0285 continue;
0286 }
0287 if (std::abs(ch->pdg_id()) == 22)
0288 {
0289 gammas.push_back(ch);
0290 }
0291 }
0292 };
0293
0294 auto collect_photons_g4 = [&](int mother_barcode, std::vector<PHG4Particle *> &gammas)
0295 {
0296 gammas.clear();
0297 auto itMom = g4_by_barcode.find(mother_barcode);
0298 if (itMom == g4_by_barcode.end())
0299 {
0300 return;
0301 }
0302 PHG4Particle *g4mom = itMom->second;
0303 auto itKids = g4_children.find(g4mom->get_track_id());
0304 if (itKids == g4_children.end())
0305 {
0306 return;
0307 }
0308 for (auto *c : itKids->second)
0309 {
0310 if (!c)
0311 {
0312 continue;
0313 }
0314 if (std::abs(c->get_pid()) == 22)
0315 {
0316 gammas.push_back(c);
0317 }
0318 }
0319 };
0320
0321 if (genevtmap)
0322 {
0323 for (auto &it : *genevtmap)
0324 {
0325 PHHepMCGenEvent *g = it.second;
0326 if (!g || !g->getEvent())
0327 {
0328 continue;
0329 }
0330 const int embedding = g->get_embedding_id();
0331 HepMC::GenEvent *evt = g->getEvent();
0332
0333 for (HepMC::GenEvent::particle_const_iterator hepmcgenit = evt->particles_begin(); hepmcgenit != evt->particles_end(); ++hepmcgenit)
0334 {
0335 HepMC::GenParticle *p = *hepmcgenit;
0336 if (!p)
0337 {
0338 continue;
0339 }
0340 const int pid = p->pdg_id();
0341 if ((pid != 111 && pid != 221) || (pid == 111 && !m_save_pi0) || (pid == 221 && !m_save_eta))
0342 {
0343 continue;
0344 }
0345
0346 const int bc = p->barcode();
0347 const std::pair<int, int> key(embedding, bc);
0348 if (m_seen_mother_keys.contains(key))
0349 {
0350 continue;
0351 }
0352 m_seen_mother_keys.insert(key);
0353 m_seen_barcodes.insert(bc);
0354
0355 auto *mes = new TruthNeutralMesonv1();
0356 mes->set_pid(pid);
0357
0358 float e;
0359 float pt;
0360 float eta;
0361 float phi;
0362 float magp;
0363 fill_kin_from_hepmc(p, e, pt, eta, phi, magp);
0364 mes->set_e(e);
0365 mes->set_pt(pt);
0366 mes->set_eta(eta);
0367 mes->set_phi(phi);
0368 mes->set_p(magp);
0369
0370 bool is_prompt = false;
0371 bool pi0_has_eta = false;
0372 int parent_pid = 0;
0373 classify_origin(pid, bc, is_prompt, parent_pid, pi0_has_eta);
0374 mes->set_prompt(is_prompt);
0375 mes->set_parent_pid(parent_pid);
0376 if (pid == 111)
0377 {
0378 mes->set_mother_is_eta(pi0_has_eta);
0379 }
0380 else if (pid == 221)
0381 {
0382 bool eta_3pi0 = false;
0383 bool eta_pi0pipm = false;
0384 classify_eta_decay(bc, p, eta_3pi0, eta_pi0pipm);
0385 mes->set_eta_3pi0(eta_3pi0);
0386 mes->set_eta_pi0pipm(eta_pi0pipm);
0387 }
0388
0389 std::vector<HepMC::GenParticle *> gamsH;
0390 collect_photons_hepmc(p, gamsH);
0391
0392 if (gamsH.size() == 2)
0393 {
0394 for (auto *gph : gamsH)
0395 {
0396 float ge;
0397 float gpt;
0398 float geta;
0399 float gphi;
0400 float gp;
0401 fill_kin_from_hepmc(gph, ge, gpt, geta, gphi, gp);
0402
0403 bool isconverted = false;
0404 auto itg = g4_by_barcode.find(gph->barcode());
0405 if (itg != g4_by_barcode.end() && itg->second)
0406 {
0407 int gphtrackid = itg->second->get_track_id();
0408 isconverted = FindConversion(truthinfo, gphtrackid, ge);
0409 }
0410 mes->add_photon(ge, gpt, geta, gphi, gp, isconverted);
0411 }
0412 }
0413 else if (truthinfo)
0414 {
0415 std::vector<PHG4Particle *> gamsG4;
0416 collect_photons_g4(bc, gamsG4);
0417 if (Verbosity() > 0)
0418 {
0419 if (gamsG4.size() > 2)
0420 {
0421 std::cout << " In main loop / Geant part more than 2 photon check..." << std::endl;
0422 for (auto *gph : gamsG4)
0423 {
0424 std::cout << " mother pid " << pid << " daughter pid " << gph->get_pid() << ", parent id " << gph->get_parent_id() << std::endl;
0425 }
0426 }
0427 }
0428
0429 if (gamsG4.size() == 2)
0430 {
0431 for (auto *gph : gamsG4)
0432 {
0433 float ge;
0434 float gpt;
0435 float geta;
0436 float gphi;
0437 float gp;
0438 fill_kin_from_g4(gph, ge, gpt, geta, gphi, gp);
0439 bool isconverted = false;
0440 int gphtrackid = gph->get_track_id();
0441 isconverted = FindConversion(truthinfo, gphtrackid, ge);
0442 mes->add_photon(ge, gpt, geta, gphi, gp, isconverted);
0443 }
0444 }
0445 }
0446 _container->AddMeson(mes);
0447 }
0448 }
0449 }
0450
0451 if (truthinfo)
0452 {
0453 PHG4TruthInfoContainer::ConstRange range = truthinfo->GetParticleRange();
0454 for (auto it = range.first; it != range.second; ++it)
0455 {
0456 PHG4Particle *p = it->second;
0457 if (!p)
0458 {
0459 continue;
0460 }
0461 const int pid = p->get_pid();
0462 if ((pid != 111 && pid != 221) || (pid == 111 && !m_save_pi0) || (pid == 221 && !m_save_eta))
0463 {
0464 continue;
0465 }
0466
0467 const int bc = p->get_barcode();
0468 if (m_seen_barcodes.contains(bc))
0469 {
0470 continue;
0471 }
0472 m_seen_barcodes.insert(bc);
0473 const std::pair<int, int> key(0, bc);
0474 if (m_seen_mother_keys.contains(key))
0475 {
0476 continue;
0477 }
0478 m_seen_mother_keys.insert(key);
0479
0480 auto *mes = new TruthNeutralMesonv1();
0481 mes->set_pid(pid);
0482
0483 float e;
0484 float pt;
0485 float eta;
0486 float phi;
0487 float magp;
0488 fill_kin_from_g4(p, e, pt, eta, phi, magp);
0489 mes->set_e(e);
0490 mes->set_pt(pt);
0491 mes->set_eta(eta);
0492 mes->set_phi(phi);
0493 mes->set_p(magp);
0494
0495 bool is_prompt = false;
0496 bool pi0_has_eta = false;
0497 int parent_pid = 0;
0498 classify_origin(pid, bc, is_prompt, parent_pid, pi0_has_eta);
0499 if (parent_pid == pid || parent_pid == 22)
0500 {
0501 continue;
0502 }
0503
0504 auto g4parents = g4_by_id.find(p->get_parent_id());
0505 if (g4parents == g4_by_id.end())
0506 {
0507 continue;
0508 }
0509
0510 auto *gp = g4parents->second;
0511 if (!gp)
0512 {
0513 continue;
0514 }
0515 if (gp->get_track_id() < 0)
0516 {
0517 continue;
0518 }
0519 int parentbc = gp->get_barcode();
0520
0521 auto itH = hepmc_by_barcode.find(parentbc);
0522 if (itH == hepmc_by_barcode.end())
0523 {
0524 continue;
0525 }
0526
0527 mes->set_prompt(is_prompt);
0528 mes->set_parent_pid(parent_pid);
0529 if (pid == 111)
0530 {
0531 mes->set_mother_is_eta(pi0_has_eta);
0532 }
0533 else if (pid == 221)
0534 {
0535 bool eta_3pi0 = false;
0536 bool eta_pi0pipm = false;
0537 classify_eta_decay(bc, nullptr, eta_3pi0, eta_pi0pipm);
0538 mes->set_eta_3pi0(eta_3pi0);
0539 mes->set_eta_pi0pipm(eta_pi0pipm);
0540 }
0541
0542 std::vector<PHG4Particle *> gamsG4;
0543 auto itKids = g4_children.find(p->get_track_id());
0544 if (itKids == g4_children.end())
0545 {
0546 continue;
0547 }
0548
0549 for (auto *c : itKids->second)
0550 {
0551 if (!c)
0552 {
0553 continue;
0554 }
0555 if (std::abs(c->get_pid()) != 22)
0556 {
0557 continue;
0558 }
0559 gamsG4.push_back(c);
0560 }
0561
0562 if (Verbosity() > 0)
0563 {
0564 if (gamsG4.size() > 2)
0565 {
0566 std::cout << " In G4 loop more than 2 photon decay check..." << std::endl;
0567 for (auto *gph : gamsG4)
0568 {
0569 std::cout << " pid " << gph->get_pid() << ", parent id " << pid << std::endl;
0570 }
0571 }
0572 }
0573
0574 if (gamsG4.size() == 2)
0575 {
0576 for (auto *gph : gamsG4)
0577 {
0578 float ge;
0579 float gpt;
0580 float geta;
0581 float gphi;
0582 float gpmom;
0583 fill_kin_from_g4(gph, ge, gpt, geta, gphi, gpmom);
0584 bool isconverted = false;
0585 int gphtrackid = gph->get_track_id();
0586 isconverted = FindConversion(truthinfo, gphtrackid, ge);
0587 mes->add_photon(ge, gpt, geta, gphi, gpmom, isconverted);
0588 }
0589 }
0590 _container->AddMeson(mes);
0591 }
0592 }
0593
0594 return Fun4AllReturnCodes::EVENT_OK;
0595 }
0596
0597 void TruthNeutralMesonBuilder::classify_eta_decay(int bc, HepMC::GenParticle *hepMom, bool &eta_3pi0, bool &eta_pi0pipm)
0598 {
0599 eta_3pi0 = false;
0600 eta_pi0pipm = false;
0601
0602 if (hepMom)
0603 {
0604 HepMC::GenVertex *vtx = hepMom->end_vertex();
0605 if (vtx)
0606 {
0607 int nPi0 = 0;
0608 int nPip = 0;
0609 int nPim = 0;
0610 for (auto it = vtx->particles_out_const_begin(); it != vtx->particles_out_const_end(); ++it)
0611 {
0612 const HepMC::GenParticle *ch = *it;
0613 if (!ch)
0614 {
0615 continue;
0616 }
0617 const int chId = ch->pdg_id();
0618 if (chId == 111)
0619 {
0620 ++nPi0;
0621 }
0622 else if (chId == 211)
0623 {
0624 ++nPip;
0625 }
0626 else if (chId == -211)
0627 {
0628 ++nPim;
0629 }
0630 }
0631 if (nPi0 == 3)
0632 {
0633 eta_3pi0 = true;
0634 }
0635 if (nPi0 == 1 && nPip == 1 && nPim == 1)
0636 {
0637 eta_pi0pipm = true;
0638 }
0639 }
0640 }
0641
0642 if (!eta_3pi0 && !eta_pi0pipm && truthinfo)
0643 {
0644 auto itMom = g4_by_barcode.find(bc);
0645 if (itMom != g4_by_barcode.end())
0646 {
0647 PHG4Particle *g4mom = itMom->second;
0648 auto itKids = g4_children.find(g4mom->get_track_id());
0649 if (itKids != g4_children.end())
0650 {
0651 int nPi0 = 0;
0652 int nPip = 0;
0653 int nPim = 0;
0654 for (auto *c : itKids->second)
0655 {
0656 if (!c)
0657 {
0658 continue;
0659 }
0660 const int chId = c->get_pid();
0661 if (chId == 111)
0662 {
0663 ++nPi0;
0664 }
0665 else if (chId == 211)
0666 {
0667 ++nPip;
0668 }
0669 else if (chId == -211)
0670 {
0671 ++nPim;
0672 }
0673 }
0674 if (nPi0 == 3)
0675 {
0676 eta_3pi0 = true;
0677 }
0678 if (nPi0 == 1 && nPip == 1 && nPim == 1)
0679 {
0680 eta_pi0pipm = true;
0681 }
0682 }
0683 }
0684 }
0685 }
0686
0687 bool TruthNeutralMesonBuilder::FindConversion(PHG4TruthInfoContainer *_truthinfo, int trackid, float energy)
0688 {
0689 PHG4Shower *shower = _truthinfo->GetShower(trackid);
0690 if (!shower)
0691 {
0692 return false;
0693 }
0694 bool foundconversion = false;
0695 auto g4particle_ids = shower->g4particle_ids();
0696 for (auto g4particle_id : g4particle_ids)
0697 {
0698 PHG4Particle *g4particle = _truthinfo->GetParticle(g4particle_id);
0699 if (!g4particle)
0700 {
0701 continue;
0702 }
0703 int vertexid = g4particle->get_vtx_id();
0704 PHG4VtxPoint *vtxp = _truthinfo->GetVtx(vertexid);
0705 if (!vtxp)
0706 {
0707 continue;
0708 }
0709 float vertexr = sqrt(vtxp->get_x() * vtxp->get_x() + vtxp->get_y() * vtxp->get_y());
0710 if (vertexr < m_conversion_radius_limit)
0711 {
0712 float momentum = sqrt(g4particle->get_px() * g4particle->get_px() + g4particle->get_py() * g4particle->get_py() + g4particle->get_pz() * g4particle->get_pz());
0713 if (momentum > 0.3 * energy)
0714 {
0715 int g4particlepid = g4particle->get_pid();
0716 if (abs(g4particlepid) == 11)
0717 {
0718 foundconversion = true;
0719 }
0720 }
0721 }
0722 }
0723 return foundconversion;
0724 }
0725
0726 bool TruthNeutralMesonBuilder::RejectShowerMeson(PHG4TruthInfoContainer *_truthinfo, int parent_trackid, int this_trackid)
0727 {
0728 PHG4Shower *shower = _truthinfo->GetShower(parent_trackid);
0729 if (!shower)
0730 {
0731 return false;
0732 }
0733 bool reject = false;
0734 auto g4particle_ids = shower->g4particle_ids();
0735 for (auto g4particle_id : g4particle_ids)
0736 {
0737 if (g4particle_id != this_trackid)
0738 {
0739 continue;
0740 }
0741 PHG4Particle *g4particle = _truthinfo->GetParticle(g4particle_id);
0742 if (!g4particle)
0743 {
0744 continue;
0745 }
0746 int vertexid = g4particle->get_vtx_id();
0747 PHG4VtxPoint *vtxp = _truthinfo->GetVtx(vertexid);
0748 if (!vtxp)
0749 {
0750 continue;
0751 }
0752
0753 float vertexr = sqrt(vtxp->get_x() * vtxp->get_x() + vtxp->get_y() * vtxp->get_y());
0754 if (vertexr > m_shower_reject_radius)
0755 {
0756 int g4particlepid = g4particle->get_pid();
0757 if (abs(g4particlepid) == 111 || abs(g4particlepid) == 221)
0758 {
0759 reject = true;
0760 }
0761 }
0762 }
0763 return reject;
0764 }
0765
0766 int TruthNeutralMesonBuilder::End(PHCompositeNode * )
0767 {
0768 return Fun4AllReturnCodes::EVENT_OK;
0769 }
