HepMC3 event record library
testThreadssearch.cc
1 // -*- C++ -*-
2 //
3 // This file is part of HepMC
4 // Copyright (C) 2014-2019 The HepMC collaboration (see AUTHORS for details)
5 //
6 #include "HepMC3/Attribute.h"
7 #include "HepMC3/GenEvent.h"
8 #include "HepMC3/GenParticle.h"
9 #include "HepMC3/GenVertex.h"
10 #include "HepMC3/Relatives.h"
11 #include "HepMC3TestUtils.h"
12 #include <thread>
13 #include <random>
14 #include <iterator>
15 using namespace HepMC3;
16 const int NinputCopies=4;
17 const int NmaxThreads=3;
18 std::shared_ptr<GenEvent> generate(const int Z) {
19  std::shared_ptr<GenEvent> evt = std::make_shared<GenEvent>();
20  std::shared_ptr<GenRunInfo> run = std::make_shared<GenRunInfo>();
21  evt->set_run_info(run);
22  GenParticlePtr b2 = std::make_shared<GenParticle>( FourVector( 0.0, 0.0, 7000.0, 7000.0 ),2212, 3 );
23  GenParticlePtr b1 = std::make_shared<GenParticle>( FourVector( 0.750, -1.569, 32.191, 32.238), 1, 3 );
24  GenParticlePtr b3 = std::make_shared<GenParticle>( FourVector( 0.750, -1.569, 32.191, -32.238), 1, 3 );
25  GenVertexPtr v1 = std::make_shared<GenVertex>();
26  v1->add_particle_in (b1);
27  v1->add_particle_in(b2);
28  v1->add_particle_out(b3);
29  evt->add_vertex(v1);
30  for (int z = 0; z < Z; z++) {
31  std::vector<GenParticlePtr> particles = evt->particles();
32  for (auto p: particles) {
33  if (p->end_vertex()) continue;
34  GenParticlePtr p2 = std::make_shared<GenParticle>( FourVector( 0.0, 0.0, 7000.0+0.01*evt->particles().size(), 7000.0 ),2212, 3 );
35  GenParticlePtr p1 = std::make_shared<GenParticle>( FourVector( 0.750, -1.569, 32.191+0.01*evt->particles().size(), 32.238), 1, 3 );
36  GenVertexPtr v = std::make_shared<GenVertex>();
37  v->add_particle_in (p);
38  v->add_particle_out(p1);
39  v->add_particle_out(p2);
40  evt->add_vertex(v);
41  }
42  }
43  return evt;
44 }
45 
46 void attribute_function1(const std::vector<std::shared_ptr<GenEvent>>& evts, std::map<int,int>& res)
47 {
48  for (size_t i = 0; i < evts.size(); i++)
49  res[evts.at(i)->event_number()] =
50  (Relatives::DESCENDANTS(evts.at(i)->particles().at(2))).size();
51 }
52 
53 int main()
54 {
55  std::vector<std::shared_ptr<GenEvent> > evts;
56  for (int i=0; i < 20; i++) {
57  evts.push_back(generate(5+i/3));
58  evts.back()->set_event_number(i);
59  }
60  std::random_device rd;
61  std::mt19937 g(rd());
62  std::vector<std::shared_ptr<GenEvent>> thr_evts[NinputCopies];
63  for (int i=0; i<NinputCopies; i++) {
64  thr_evts[i] = evts;
65  std::shuffle(thr_evts[i].begin(), thr_evts[i].end(), g);
66  }
67  std::map<int,int> res[NinputCopies];
68  std::vector<std::thread> threads;
69 
70  for (int i = 0; i < NinputCopies; i++)
71  threads.push_back(std::thread(attribute_function1,std::cref(thr_evts[i]), std::ref(res[i])));
72  for (auto& th : threads) th.join();
73  threads.clear();
74 
75  for (int k = 1; k < NinputCopies; k++)
76  {
77  if (!std::equal(res[k].begin(), res[k].end(), res[0].begin())) {
78  return 1;
79  };
80  }
81  return 0;
82 }
Definition of class Attribute, class IntAttribute and class StringAttribute.
Definition of class GenEvent.
Definition of class GenParticle.
Definition of class GenVertex.
Defines helper classes to extract relatives of an input GenParticle or GenVertex.
Generic 4-vector.
Definition: FourVector.h:36
static HEPMC3search_Relatives_EXPORT_API thread_local const Descendants DESCENDANTS
Descendants.
Definition: Relatives.h:204
HepMC3 main namespace.
int main(int argc, char **argv)