diff --git a/DDG4/include/DDG4/Geant4Primary.h b/DDG4/include/DDG4/Geant4Primary.h index 86c1d45203f1b2b2e1acb5370cb0a57b861089ac..6ecae7b66b76f7020945dafadbda02922305a16f 100644 --- a/DDG4/include/DDG4/Geant4Primary.h +++ b/DDG4/include/DDG4/Geant4Primary.h @@ -100,11 +100,11 @@ namespace DD4hep { Geant4PrimaryInteraction& operator=(const Geant4PrimaryInteraction& c) = delete; public: - typedef Geant4Particle Particle; - typedef Geant4Vertex Vertex; - typedef std::map<int,Particle*> ParticleMap; - typedef std::map<int,Vertex*> VertexMap; - typedef dd4hep_ptr<PrimaryExtension> ExtensionHandle; + typedef Geant4Particle Particle; + typedef Geant4Vertex Vertex; + typedef std::map<int,Particle*> ParticleMap; + typedef std::map<int,std::vector<Vertex*>> VertexMap; + typedef dd4hep_ptr<PrimaryExtension> ExtensionHandle; /// The map of primary vertices for the particles. VertexMap vertices; diff --git a/DDG4/src/Geant4InputAction.cpp b/DDG4/src/Geant4InputAction.cpp index b9cbc71431c9270e6efcbbf15bf3e7ea702426db..d45f122c373a96b1820121790efd4079a9c777b4 100644 --- a/DDG4/src/Geant4InputAction.cpp +++ b/DDG4/src/Geant4InputAction.cpp @@ -207,7 +207,7 @@ void Geant4InputAction::operator()(G4Event* event) { for(size_t i=0; i<vertices.size(); ++i ) { - inter->vertices.insert(make_pair(m_mask,vertices[i])); + inter->vertices[m_mask].push_back( vertices[i] ); } // build collection of MCParticles diff --git a/DDG4/src/Geant4InputHandling.cpp b/DDG4/src/Geant4InputHandling.cpp index 3d1dd488999ff9c6095834c9d3615608dcdf2a50..fe174a7aa5e175dcdc9d95680384ac96a4a2c12f 100644 --- a/DDG4/src/Geant4InputHandling.cpp +++ b/DDG4/src/Geant4InputHandling.cpp @@ -109,7 +109,9 @@ static void collectPrimaries(Geant4PrimaryMap* pm, dv->mask = mask; dv->in.insert(p->id); - interaction->vertices.insert(make_pair(vid,dv)); + + interaction->vertices[vid].push_back(dv) ; + for(; dau; dau = dau->GetNext()) collectPrimaries(pm, interaction, dv, dau); } @@ -127,7 +129,8 @@ DD4hep::Simulation::createPrimary(int mask, interaction->locked = true; interaction->mask = mask; v->mask = mask; - interaction->vertices.insert(make_pair(vid,v)); + interaction->vertices[vid].push_back(v); + for (G4PrimaryParticle *gp = gv->GetPrimary(); gp; gp = gp->GetNext() ) collectPrimaries(pm, interaction, v, gp); return interaction; @@ -177,12 +180,13 @@ static void appendInteraction(const Geant4Action* caller, } Geant4PrimaryInteraction::VertexMap::iterator ivfnd, iv, ivend; for( iv=input->vertices.begin(), ivend=input->vertices.end(); iv != ivend; ++iv ) { - ivfnd = output->vertices.find((*iv).second->mask); + ivfnd = output->vertices.find((*iv).first) ; //(*iv).second->mask); if ( ivfnd != output->vertices.end() ) { caller->abortRun("Duplicate primary interaction identifier!", "Cannot handle 2 interactions with identical identifiers!"); } - output->vertices.insert(make_pair((*iv).second->mask,(*iv).second->addRef())); + for(Geant4Vertex* vtx : (*iv).second ) + output->vertices[(*iv).first].push_back( vtx->addRef() ); } } @@ -204,16 +208,16 @@ static void rebaseVertices(Geant4PrimaryInteraction::VertexMap& vertices, int pa set<int>::iterator i; // Now move begin and end-vertex of all primary vertices accordingly for(iv=vertices.begin(), ivend=vertices.end(); iv != ivend; ++iv) { - Geant4Vertex* v = (*iv).second; - in = v->in; - out = v->out; - for(in=v->in, i=in.begin(), v->in.clear(); i != in.end(); ++i) - v->in.insert((*i)+part_offset); - for(out=v->out, i=out.begin(), v->out.clear(); i != out.end(); ++i) - v->out.insert((*i)+part_offset); + for( Geant4Vertex* v : (*iv).second ){ + in = v->in; + out = v->out; + for(in=v->in, i=in.begin(), v->in.clear(); i != in.end(); ++i) + v->in.insert((*i)+part_offset); + for(out=v->out, i=out.begin(), v->out.clear(); i != out.end(); ++i) + v->out.insert((*i)+part_offset); + } } } - /// Merge all interactions present in the context int DD4hep::Simulation::mergeInteractions(const Geant4Action* caller, const Geant4Context* context) @@ -262,18 +266,18 @@ int DD4hep::Simulation::boostInteraction(const Geant4Action* caller, } else if ( alpha != 0.0 ) { // Now move begin and end-vertex of all primary vertices accordingly - for(iv=inter->vertices.begin(); iv != inter->vertices.end(); ++iv) { - Geant4Vertex* v = (*iv).second; - double t = gamma * v->time + betagamma * v->x / CLHEP::c_light; - double x = gamma * v->x + betagamma * CLHEP::c_light * v->time; - double y = v->y; - double z = v->z; - v->x = x; - v->y = y; - v->z = z; - v->time = t; + for(iv=inter->vertices.begin(); iv != inter->vertices.end(); ++iv){ + for( Geant4Vertex* v : (*iv).second ) { + double t = gamma * v->time + betagamma * v->x / CLHEP::c_light; + double x = gamma * v->x + betagamma * CLHEP::c_light * v->time; + double y = v->y; + double z = v->z; + v->x = x; + v->y = y; + v->z = z; + v->time = t; + } } - // Now move begin and end-vertex of all primary and generator particles accordingly for(ip=inter->particles.begin(); ip != inter->particles.end(); ++ip) { Geant4ParticleHandle p = (*ip).second; @@ -281,18 +285,18 @@ int DD4hep::Simulation::boostInteraction(const Geant4Action* caller, double x = gamma * p->vsx + betagamma * CLHEP::c_light * p->time; double y = p->vsy; double z = p->vsz; - + double m = p->mass; double e2 = SQR(p->psx)+SQR(p->psy)+SQR(p->psz)+SQR(m); double px = betagamma * std::sqrt(e2) + gamma * p->psx; double py = p->psy; double pz = p->psz; - + p->vsx = x; p->vsy = y; p->vsz = z; p->time = t; - + p->psx = px; p->psy = py; p->psz = pz; @@ -316,11 +320,12 @@ int DD4hep::Simulation::smearInteraction(const Geant4Action* caller, // Now move begin and end-vertex of all primary vertices accordingly for(iv=inter->vertices.begin(); iv != inter->vertices.end(); ++iv) { - Geant4Vertex* v = (*iv).second; - v->x += dx; - v->y += dy; - v->z += dz; - v->time += dt; + for( Geant4Vertex* v : (*iv).second ){ + v->x += dx; + v->y += dy; + v->z += dz; + v->time += dt; + } } // Now move begin and end-vertex of all primary and generator particles accordingly for(ip=inter->particles.begin(); ip != inter->particles.end(); ++ip) { @@ -431,33 +436,35 @@ int DD4hep::Simulation::generatePrimaries(const Geant4Action* caller, else { Geant4PrimaryInteraction::VertexMap::iterator ivfnd, iv, ivend; for(Interaction::VertexMap::const_iterator iend=vm.end(),i=vm.begin(); i!=iend; ++i) { - int num_part = 0; - Geant4Vertex* v = (*i).second; - G4PrimaryVertex* v4 = new G4PrimaryVertex(v->x,v->y,v->z,v->time); - event->AddPrimaryVertex(v4); - caller->print("+++++ G4PrimaryVertex at (%+.2e,%+.2e,%+.2e) [mm] %+.2e [ns]", - v->x/CLHEP::mm,v->y/CLHEP::mm,v->z/CLHEP::mm,v->time/CLHEP::ns); - for(Geant4Vertex::Particles::const_iterator ip=v->out.begin(); ip!=v->out.end(); ++ip) { - Geant4ParticleHandle p = pm[*ip]; - if ( p->daughters.size() > 0 ) { - PropertyMask mask(p->reason); - mask.set(G4PARTICLE_HAS_SECONDARIES); - } - if ( p->parents.size() == 0 ) { - Primaries relevant = getRelevant(visited,prim,pm,p); - for(Primaries::const_iterator j=relevant.begin(); j!= relevant.end(); ++j) { - Geant4ParticleHandle r = (*j).first; - G4PrimaryParticle* p4 = (*j).second; - PropertyMask reason(r->reason); - char text[64]; - - reason.set(G4PARTICLE_PRIMARY); - v4->SetPrimary(p4); - ::snprintf(text,sizeof(text),"-> G4Primary[%3d]",num_part); - r.dumpWithMomentum(caller->outputLevel()-1,caller->name(),text); - ++num_part; - } - } + for( Geant4Vertex* v : (*i).second ){ + + int num_part = 0; + G4PrimaryVertex* v4 = new G4PrimaryVertex(v->x,v->y,v->z,v->time); + event->AddPrimaryVertex(v4); + caller->print("+++++ G4PrimaryVertex at (%+.2e,%+.2e,%+.2e) [mm] %+.2e [ns]", + v->x/CLHEP::mm,v->y/CLHEP::mm,v->z/CLHEP::mm,v->time/CLHEP::ns); + for(Geant4Vertex::Particles::const_iterator ip=v->out.begin(); ip!=v->out.end(); ++ip) { + Geant4ParticleHandle p = pm[*ip]; + if ( p->daughters.size() > 0 ) { + PropertyMask mask(p->reason); + mask.set(G4PARTICLE_HAS_SECONDARIES); + } + if ( p->parents.size() == 0 ) { + Primaries relevant = getRelevant(visited,prim,pm,p); + for(Primaries::const_iterator j=relevant.begin(); j!= relevant.end(); ++j) { + Geant4ParticleHandle r = (*j).first; + G4PrimaryParticle* p4 = (*j).second; + PropertyMask reason(r->reason); + char text[64]; + + reason.set(G4PARTICLE_PRIMARY); + v4->SetPrimary(p4); + ::snprintf(text,sizeof(text),"-> G4Primary[%3d]",num_part); + r.dumpWithMomentum(caller->outputLevel()-1,caller->name(),text); + ++num_part; + } + } + } } } for(map<int,G4PrimaryParticle*>::iterator i=prim.begin(); i!=prim.end(); ++i) { diff --git a/DDG4/src/Geant4ParticleGenerator.cpp b/DDG4/src/Geant4ParticleGenerator.cpp index 17742a8d9dd04f253b3d7b3b3b9c0b565f3d5a19..2f822eb8b7c8648e099ac570489b1da4f900055d 100644 --- a/DDG4/src/Geant4ParticleGenerator.cpp +++ b/DDG4/src/Geant4ParticleGenerator.cpp @@ -83,14 +83,15 @@ void Geant4ParticleGenerator::printInteraction(Geant4PrimaryInteraction* inter) return; } for(_V::const_iterator iv=inter->vertices.begin(); iv!=inter->vertices.end(); ++iv) { - Geant4Vertex* v = (*iv).second; - print("+-> Interaction [%d] %.3f GeV %s pos:(%.3f %.3f %.3f)[mm]", - count, m_energy/CLHEP::GeV, m_particleName.c_str(), - v->x/CLHEP::mm, v->y/CLHEP::mm, v->z/CLHEP::mm); - ++count; - for(set<int>::const_iterator i=v->out.begin(); i!=v->out.end(); ++i) { - Geant4ParticleHandle p = inter->particles[*i]; - p.dumpWithVertex(outputLevel(),name()," +->"); + for( Geant4Vertex* v : (*iv).second ){ + print("+-> Interaction [%d] %.3f GeV %s pos:(%.3f %.3f %.3f)[mm]", + count, m_energy/CLHEP::GeV, m_particleName.c_str(), + v->x/CLHEP::mm, v->y/CLHEP::mm, v->z/CLHEP::mm); + ++count; + for(set<int>::const_iterator i=v->out.begin(); i!=v->out.end(); ++i) { + Geant4ParticleHandle p = inter->particles[*i]; + p.dumpWithVertex(outputLevel(),name()," +->"); + } } } } @@ -120,7 +121,7 @@ void Geant4ParticleGenerator::operator()(G4Event*) { vtx->x = position.X(); vtx->y = position.Y(); vtx->z = position.Z(); - inter->vertices.insert(make_pair(inter->vertices.size(),vtx)); + inter->vertices[m_mask].push_back( vtx ); for(int i=0; i<m_multiplicity; ++i) { double momentum = m_energy; Particle* p = new Particle(); diff --git a/DDG4/src/Geant4Primary.cpp b/DDG4/src/Geant4Primary.cpp index 45f307dd1168654e518d324d06a7fc69060aaba1..72ab573b75423e83f93fd926e7a15db42901f7e8 100644 --- a/DDG4/src/Geant4Primary.cpp +++ b/DDG4/src/Geant4Primary.cpp @@ -54,7 +54,12 @@ const Geant4Particle* Geant4PrimaryMap::get(const G4PrimaryParticle* particle) c /// Default destructor Geant4PrimaryInteraction::~Geant4PrimaryInteraction() { - releaseObjects(vertices); + + Geant4PrimaryInteraction::VertexMap::iterator iv, ivend; + for( iv=vertices.begin(), ivend=vertices.end(); iv != ivend; ++iv ){ + for( Geant4Vertex* vtx : (*iv).second ) + ReleaseObject<Geant4Vertex*>()( vtx ); + } releaseObjects(particles); } @@ -70,17 +75,19 @@ void Geant4PrimaryInteraction::setNextPID(int new_value) { /// Apply mask to all contained vertices and particles bool Geant4PrimaryInteraction::applyMask() { - if ( vertices.size() <= 1 ) { - Geant4PrimaryInteraction::ParticleMap::iterator ip, ipend; - for( ip=particles.begin(), ipend=particles.end(); ip != ipend; ++ip ) - (*ip).second->mask = mask; - - Geant4PrimaryInteraction::VertexMap::iterator iv, ivend; - for( iv=vertices.begin(), ivend=vertices.end(); iv != ivend; ++iv ) - (*iv).second->mask = mask; - return true; + + Geant4PrimaryInteraction::ParticleMap::iterator ip, ipend; + for( ip=particles.begin(), ipend=particles.end(); ip != ipend; ++ip ) + (*ip).second->mask = mask; + + Geant4PrimaryInteraction::VertexMap::iterator iv, ivend; + for( iv=vertices.begin(), ivend=vertices.end(); iv != ivend; ++iv ){ + for( auto vtx : (*iv).second ) + vtx->mask = mask; } - return false; + + return true; + } /// Default destructor