From 44e848c77aa26d76b313a2ea55dd4a2bd2e3c3ba Mon Sep 17 00:00:00 2001
From: Markus FRANK <Markus.Frank@cern.ch>
Date: Tue, 3 Nov 2020 09:08:19 +0100
Subject: [PATCH] Try to fix reflections for volumes and assemblies
---
DDCore/include/DD4hep/GeoHandler.h | 13 +-
DDCore/src/DetectorImp.cpp | 9 +-
DDCore/src/GeoHandler.cpp | 39 +++--
DDCore/src/Volumes.cpp | 2 +
DDCore/src/plugins/Compact2Objects.cpp | 10 +-
DDEve/DDEve/DDG4IO.C | 20 +--
DDEve/include/DDEve/DDEveEventData.h | 6 +-
DDEve/include/DDEve/Display.h | 10 +-
DDEve/include/DDEve/DisplayConfiguration.h | 2 +-
DDEve/include/DDEve/HitActors.h | 13 +-
DDEve/include/DDEve/ParticleActors.h | 33 +++-
DDEve/src/Display.cpp | 73 ++++----
DDEve/src/DisplayConfigurationParser.cpp | 21 +--
DDEve/src/HitActors.cpp | 16 +-
DDEve/src/ParticleActors.cpp | 65 +++++--
DDG4/include/DDG4/Geant4Converter.h | 1 +
DDG4/include/DDG4/Geant4ParticleGun.h | 5 +-
DDG4/src/Geant4Converter.cpp | 159 +++++++++++++++---
DDG4/src/Geant4IsotropeGenerator.cpp | 2 +-
DDG4/src/Geant4ParticleGenerator.cpp | 13 +-
DDG4/src/Geant4ParticleGun.cpp | 24 +--
.../compact/NestedBoxReflection.xml | 2 +
.../ClientTests/eve/NestedBoxReflection.xml | 5 +-
.../ClientTests/scripts/Check_reflection.py | 16 +-
.../src/NestedBoxReflection_geo.cpp | 44 +++--
25 files changed, 418 insertions(+), 185 deletions(-)
diff --git a/DDCore/include/DD4hep/GeoHandler.h b/DDCore/include/DD4hep/GeoHandler.h
index 4cba4d5ef..3f3bb964b 100644
--- a/DDCore/include/DD4hep/GeoHandler.h
+++ b/DDCore/include/DD4hep/GeoHandler.h
@@ -67,15 +67,15 @@ namespace dd4hep {
*/
class GeometryInfo {
public:
- std::set<TGeoShape*> solids;
- std::set<Volume> volumeSet;
- std::vector<Volume> volumes;
- std::vector<std::pair<std::string, TGeoMatrix*> > trafos;
+ std::set<TGeoShape*> solids;
+ std::set<Volume> volumeSet;
+ std::vector<Volume> volumes;
std::set<VisAttr> vis;
std::set<Ref_t> fields;
std::set<Material> materials;
std::set<TGeoMedium*> media;
std::set<TGeoElement*> elements;
+ std::vector<std::pair<std::string, TGeoMatrix*> > trafos;
};
};
@@ -96,9 +96,12 @@ namespace dd4hep {
protected:
bool m_propagateRegions;
std::map<int, std::set<const TGeoNode*> >* m_data;
+ std::map<int, std::set<std::pair<const TGeoNode*,const TGeoNode*> > >* m_places;
/// Internal helper to collect geometry information from traversal
- GeoHandler& i_collect(const TGeoNode* node, int level, Region rg, LimitSet ls);
+ GeoHandler& i_collect(const TGeoNode* parent,
+ const TGeoNode* node,
+ int level, Region rg, LimitSet ls);
private:
/// Never call Copy constructor
diff --git a/DDCore/src/DetectorImp.cpp b/DDCore/src/DetectorImp.cpp
index 0d49b27aa..b02cf9efc 100644
--- a/DDCore/src/DetectorImp.cpp
+++ b/DDCore/src/DetectorImp.cpp
@@ -670,12 +670,17 @@ void DetectorImp::endDocument(bool close_geometry) {
add(trackingVis);
#endif
m_worldVol.solid()->ComputeBBox();
- /// Since we allow now for anonymous shapes,
- /// we will rename them to use the name of the volume they are assigned to
+ // Since we allow now for anonymous shapes,
+ // we will rename them to use the name of the volume they are assigned to
mgr->CloseGeometry();
}
+ // Patching shape names of anaonymous shapes
ShapePatcher patcher(m_volManager, m_world);
patcher.patchShapes();
+ // Propagating reflections
+ ReflectionBuilder rb(description);
+ rb.execute();
+
mapDetectorTypes();
m_state = READY;
}
diff --git a/DDCore/src/GeoHandler.cpp b/DDCore/src/GeoHandler.cpp
index 6b4aae5c5..590c5db8d 100644
--- a/DDCore/src/GeoHandler.cpp
+++ b/DDCore/src/GeoHandler.cpp
@@ -30,12 +30,16 @@ using namespace dd4hep;
using namespace std;
namespace {
- void collectSolid(GeoHandler::GeometryInfo& geo, const string& name, const string& node, TGeoShape* shape,
- TGeoMatrix* matrix) {
- if (0 == ::strncmp(shape->GetName(), "TGeo", 4)) {
+ void collectSolid(GeoHandler::GeometryInfo& geo,
+ const string& name,
+ const string& node,
+ TGeoShape* shape,
+ TGeoMatrix* matrix)
+ {
+ if ( 0 == ::strncmp(shape->GetName(), "TGeo", 4) ) {
shape->SetName(name.c_str());
}
- if (shape->IsA() == TGeoCompositeShape::Class()) {
+ if ( shape->IsA() == TGeoCompositeShape::Class() ) {
const TGeoCompositeShape* s = (const TGeoCompositeShape*) shape;
const TGeoBoolNode* boolean = s->GetBoolNode();
collectSolid(geo, name + "_left", name + "_left", boolean->GetLeftShape(), boolean->GetLeftMatrix());
@@ -49,23 +53,28 @@ namespace {
/// Default constructor
GeoHandler::GeoHandler() : m_propagateRegions(false) {
m_data = new map<int,set<const TGeoNode*> >();
+ m_places = new map<int, set<pair<const TGeoNode*,const TGeoNode*> > >();
}
/// Initializing constructor
GeoHandler::GeoHandler(map<int,set<const TGeoNode*> >* ptr)
: m_propagateRegions(false), m_data(ptr) {
+ m_places = new map<int, set<pair<const TGeoNode*,const TGeoNode*> > >();
}
/// Default destructor
GeoHandler::~GeoHandler() {
if (m_data)
delete m_data;
- m_data = 0;
+ m_data = nullptr;
+ if ( m_places )
+ delete m_places;
+ m_places = nullptr;
}
map<int,set<const TGeoNode*> >* GeoHandler::release() {
map<int,set<const TGeoNode*> >* d = m_data;
- m_data = 0;
+ m_data = nullptr;
return d;
}
@@ -77,13 +86,19 @@ bool GeoHandler::setPropagateRegions(bool value) {
}
GeoHandler& GeoHandler::collect(DetElement element) {
+ DetElement par = element.parent();
+ TGeoNode* par_node = par.isValid() ? par.placement().ptr() : nullptr;
m_data->clear();
- return i_collect(element.placement().ptr(), 0, Region(), LimitSet());
+ m_places->clear();
+ return i_collect(par_node, element.placement().ptr(), 0, Region(), LimitSet());
}
GeoHandler& GeoHandler::collect(DetElement element, GeometryInfo& info) {
+ DetElement par = element.parent();
+ TGeoNode* par_node = par.isValid() ? par.placement().ptr() : nullptr;
m_data->clear();
- i_collect(element.placement().ptr(), 0, Region(), LimitSet());
+ m_places->clear();
+ i_collect(par_node, element.placement().ptr(), 0, Region(), LimitSet());
for (auto i = m_data->rbegin(); i != m_data->rend(); ++i) {
const auto& mapped = (*i).second;
for (const TGeoNode* n : mapped ) {
@@ -117,7 +132,10 @@ GeoHandler& GeoHandler::collect(DetElement element, GeometryInfo& info) {
return *this;
}
-GeoHandler& GeoHandler::i_collect(const TGeoNode* current, int level, Region rg, LimitSet ls) {
+GeoHandler& GeoHandler::i_collect(const TGeoNode* parent,
+ const TGeoNode* current,
+ int level, Region rg, LimitSet ls)
+{
TGeoVolume* volume = current->GetVolume();
TObjArray* nodes = volume->GetNodes();
int num_children = nodes ? nodes->GetEntriesFast() : 0;
@@ -136,11 +154,12 @@ GeoHandler& GeoHandler::i_collect(const TGeoNode* current, int level, Region rg,
}
}
(*m_data)[level].emplace(current);
+ (*m_places)[level].emplace(make_pair(parent,current));
//printf("GeoHandler: collect level:%d %s\n",level,current->GetName());
if (num_children > 0) {
for (int i = 0; i < num_children; ++i) {
TGeoNode* node = (TGeoNode*) nodes->At(i);
- i_collect(node, level + 1, region, limits);
+ i_collect(current, node, level + 1, region, limits);
}
}
return *this;
diff --git a/DDCore/src/Volumes.cpp b/DDCore/src/Volumes.cpp
index 69f246aba..dd23f31ea 100644
--- a/DDCore/src/Volumes.cpp
+++ b/DDCore/src/Volumes.cpp
@@ -208,12 +208,14 @@ namespace {
v->CloneNodesAndConnect(vol);
// The volume is now properly cloned, but with the same shape.
// Reflect the shape (if any) and connect it.
+#if 0
if (v->GetShape()) {
TGeoScale* scale = new TGeoScale( 1., 1.,-1.);
TGeoShape* reflected_shape =
TGeoScaledShape::MakeScaledShape((nam+"_shape_refl").c_str(), v->GetShape(), scale);
vol->SetShape(reflected_shape);
}
+#endif
// Reflect the daughters.
Int_t nd = vol->GetNdaughters();
if ( !nd ) return vol;
diff --git a/DDCore/src/plugins/Compact2Objects.cpp b/DDCore/src/plugins/Compact2Objects.cpp
index f5044b1c5..cceffd1c2 100644
--- a/DDCore/src/plugins/Compact2Objects.cpp
+++ b/DDCore/src/plugins/Compact2Objects.cpp
@@ -1474,8 +1474,8 @@ template <> void Converter<Compact>::operator()(xml_h element) const {
open_geometry = steer.attr<bool>(_U(open));
if ( steer.hasAttr(_U(close)) )
close_document = steer.attr<bool>(_U(close));
- if ( steer.hasAttr(_U(close_geometry)) )
- close_geometry = steer.attr<bool>(_U(close_geometry));
+ if ( steer.hasAttr(_U(reflect)) )
+ build_reflections = steer.attr<bool>(_U(reflect));
for (xml_coll_t clr(steer, _U(clear)); clr; ++clr) {
string nam = clr.hasAttr(_U(name)) ? clr.attr<string>(_U(name)) : string();
if ( nam.substr(0,6) == "elemen" ) {
@@ -1563,12 +1563,6 @@ template <> void Converter<Compact>::operator()(xml_h element) const {
}
#ifdef _WIN32
- void buildReflections();
-
-/// Build reflections the ROOT way. To be called once the geometry is closed
-void DetectorImp::buildReflections() {
-}
-
template Converter<Plugin>;
template Converter<Constant>;
template Converter<Material>;
diff --git a/DDEve/DDEve/DDG4IO.C b/DDEve/DDEve/DDG4IO.C
index 0d283807c..a9501f284 100644
--- a/DDEve/DDEve/DDG4IO.C
+++ b/DDEve/DDEve/DDG4IO.C
@@ -83,16 +83,16 @@ namespace {
void* _convertParticleFunc(void* source, dd4hep::DDEveParticle* p) {
if (source ) {
dd4hep::sim::Geant4Particle* s = (dd4hep::sim::Geant4Particle*)source;
- p->id = s->id;
- p->vsx = s->vsx;
- p->vsy = s->vsy;
- p->vsz = s->vsz;
- p->vex = s->vex;
- p->vey = s->vey;
- p->vez = s->vez;
- p->psx = s->psx;
- p->psy = s->psy;
- p->psz = s->psz;
+ p->id = s->id;
+ p->vsx = s->vsx;
+ p->vsy = s->vsy;
+ p->vsz = s->vsz;
+ p->vex = s->vex;
+ p->vey = s->vey;
+ p->vez = s->vez;
+ p->psx = s->psx;
+ p->psy = s->psy;
+ p->psz = s->psz;
p->pdgID = s->pdgID;
p->parent = s->parents.empty() ? -1 : *(s->parents.begin());
p->energy = std::sqrt(s->vsx*s->vsx + s->vsy*s->vsy + s->vsz*s->vsz + s->mass*s->mass);
diff --git a/DDEve/include/DDEve/DDEveEventData.h b/DDEve/include/DDEve/DDEveEventData.h
index cd9ca2107..7fb213d1e 100644
--- a/DDEve/include/DDEve/DDEveEventData.h
+++ b/DDEve/include/DDEve/DDEveEventData.h
@@ -29,11 +29,11 @@ namespace dd4hep {
class DDEveHit {
public:
/// Track/Particle, which produced this hit
- int particle;
+ int particle {0};
/// Hit position
- float x,y,z;
+ float x = 0e0, y = 0e0, z = 0e0;
/// Energy deposit
- float deposit;
+ float deposit {0};
/// Default constructor
DDEveHit();
/// Initializing constructor
diff --git a/DDEve/include/DDEve/Display.h b/DDEve/include/DDEve/Display.h
index c39513bc8..dc3eb505e 100644
--- a/DDEve/include/DDEve/Display.h
+++ b/DDEve/include/DDEve/Display.h
@@ -55,11 +55,11 @@ namespace dd4hep {
*/
class Display : public EventConsumer {
public:
- typedef DisplayConfiguration::ViewConfig ViewConfig;
- typedef DisplayConfiguration::Config DataConfig;
- typedef std::set<View*> Views;
- typedef std::set<DisplayConfiguration*> Configurations;
- typedef std::set<PopupMenu*> Menus;
+ typedef DisplayConfiguration::ViewConfig ViewConfig;
+ typedef DisplayConfiguration::Config DataConfig;
+ typedef std::set<View*> Views;
+ typedef std::set<DisplayConfiguration*> Configurations;
+ typedef std::set<PopupMenu*> Menus;
typedef std::map<std::string, TEveElementList*> Topics;
typedef std::map<std::string, ViewConfig> ViewConfigurations;
typedef std::map<std::string, DataConfig> DataConfigurations;
diff --git a/DDEve/include/DDEve/DisplayConfiguration.h b/DDEve/include/DDEve/DisplayConfiguration.h
index 4b38a86c2..81af0b57c 100644
--- a/DDEve/include/DDEve/DisplayConfiguration.h
+++ b/DDEve/include/DDEve/DisplayConfiguration.h
@@ -102,7 +102,7 @@ namespace dd4hep {
} data;
std::string name;
std::string hits;
- std::string use;
+ std::string use {"true"};
int type;
/// Default constructor
Config();
diff --git a/DDEve/include/DDEve/HitActors.h b/DDEve/include/DDEve/HitActors.h
index 2307f3cb6..599abba23 100644
--- a/DDEve/include/DDEve/HitActors.h
+++ b/DDEve/include/DDEve/HitActors.h
@@ -47,9 +47,10 @@ namespace dd4hep {
* \ingroup DD4HEP_EVE
*/
struct PointsetCreator : public DDEveHitActor {
- TEvePointSet* pointset;
- float deposit;
- int count;
+ TEvePointSet* pointset {nullptr};
+ float threshold {0};
+ float deposit {0};
+ int count {0};
/// Standard initializing constructor
PointsetCreator(const std::string& collection, size_t length);
/// Standard initializing constructor
@@ -69,9 +70,9 @@ namespace dd4hep {
* \ingroup DD4HEP_EVE
*/
struct BoxsetCreator : public DDEveHitActor {
- TEveBoxSet* boxset;
- float emax, towerH, deposit;
- int count;
+ TEveBoxSet* boxset {0};
+ float emax = 1e12, towerH = 1e12, deposit = 0e0;
+ int count {0};
/// Standard initializing constructor
BoxsetCreator(const std::string& collection, size_t length);
/// Standard initializing constructor
diff --git a/DDEve/include/DDEve/ParticleActors.h b/DDEve/include/DDEve/ParticleActors.h
index 612e51c29..7655f0861 100644
--- a/DDEve/include/DDEve/ParticleActors.h
+++ b/DDEve/include/DDEve/ParticleActors.h
@@ -22,6 +22,7 @@
// Forward declarations
class TEveTrackPropagator;
+class TEvePointSet;
class TEveCompound;
class TEveElement;
class TEveLine;
@@ -37,11 +38,13 @@ namespace dd4hep {
*/
struct MCParticleCreator : public DDEveParticleActor {
typedef std::map<std::string,TEveCompound*> Compounds;
- TEveTrackPropagator* propagator;
- TEveCompound* particles;
+ TEveTrackPropagator* propagator {0};
+ TEveCompound* particles {0};
Compounds types;
- int count;
- int lineWidth;
+ double threshold {10e0}; // 10 MeV
+ int count {0};
+ int lineWidth {4};
+
/// Standard initializing constructor
MCParticleCreator(TEveTrackPropagator* p, TEveCompound* ps, const DisplayConfiguration::Config* cfg);
/// Access sub-compound by name
@@ -54,6 +57,28 @@ namespace dd4hep {
void close();
};
+ /// Fill a 3D pointset with particle start vertices
+ /*
+ * \author M.Frank
+ * \version 1.0
+ * \ingroup DD4HEP_EVE
+ */
+ struct StartVertexCreator : public DDEveParticleActor {
+ TEvePointSet* pointset;
+ float deposit;
+ int count;
+ /// Standard initializing constructor
+ StartVertexCreator(const std::string& collection, size_t length);
+ /// Standard initializing constructor
+ StartVertexCreator(const std::string& collection, size_t length, const DisplayConfiguration::Config& cfg);
+ /// Standard destructor
+ virtual ~StartVertexCreator();
+ /// Return eve element
+ TEveElement* element() const;
+ /// Action callback of this functor:
+ virtual void operator()(const DDEveParticle& particle);
+ };
+
} /* End namespace dd4hep */
diff --git a/DDEve/src/Display.cpp b/DDEve/src/Display.cpp
index 7383a4919..1f98aa84c 100644
--- a/DDEve/src/Display.cpp
+++ b/DDEve/src/Display.cpp
@@ -368,34 +368,37 @@ void Display::OnNewEvent(EventHandler& handler ) {
const Collections& colls = (*ityp).second;
for(Collections::const_iterator j=colls.begin(); j!=colls.end(); ++j) {
size_t len = (*j).second;
- const char* nam = (*j).first;
if ( len > 0 ) {
+ const char* nam = (*j).first;
+ DataConfigurations::const_iterator icfg = m_collectionsConfigs.find(nam);
+ DataConfigurations::const_iterator cfgend = m_collectionsConfigs.end();
EventHandler::CollectionType typ = handler.collectionType(nam);
if ( typ == EventHandler::CALO_HIT_COLLECTION ||
typ == EventHandler::TRACKER_HIT_COLLECTION ) {
- const DataConfigurations::const_iterator i=m_collectionsConfigs.find(nam);
- if ( i != m_collectionsConfigs.end() ) {
- const DataConfig& cfg = (*i).second;
- if ( cfg.hits == "PointSet" ) {
- PointsetCreator cr(nam,len,cfg);
- handler.collectionLoop((*j).first, cr);
- ImportEvent(cr.element());
- }
- else if ( cfg.hits == "BoxSet" ) {
- BoxsetCreator cr(nam,len,cfg);
- handler.collectionLoop((*j).first, cr);
- ImportEvent(cr.element());
- }
- else if ( cfg.hits == "TowerSet" ) {
- TowersetCreator cr(nam,len,cfg);
- handler.collectionLoop((*j).first, cr);
- ImportEvent(cr.element());
- }
- else { // Default is point set
- PointsetCreator cr(nam,len);
- handler.collectionLoop((*j).first, cr);
- ImportEvent(cr.element());
- }
+ if ( icfg != cfgend ) {
+ const DataConfig& cfg = (*icfg).second;
+ if ( ::toupper(cfg.use[0]) == 'T' || ::toupper(cfg.use[0]) == 'Y' ) {
+ if ( cfg.hits == "PointSet" ) {
+ PointsetCreator cr(nam,len,cfg);
+ handler.collectionLoop((*j).first, cr);
+ ImportEvent(cr.element());
+ }
+ else if ( cfg.hits == "BoxSet" ) {
+ BoxsetCreator cr(nam,len,cfg);
+ handler.collectionLoop((*j).first, cr);
+ ImportEvent(cr.element());
+ }
+ else if ( cfg.hits == "TowerSet" ) {
+ TowersetCreator cr(nam,len,cfg);
+ handler.collectionLoop((*j).first, cr);
+ ImportEvent(cr.element());
+ }
+ else { // Default is point set
+ PointsetCreator cr(nam,len);
+ handler.collectionLoop((*j).first, cr);
+ ImportEvent(cr.element());
+ }
+ }
}
else {
PointsetCreator cr(nam,len);
@@ -409,13 +412,21 @@ void Display::OnNewEvent(EventHandler& handler ) {
// last track is gone ie. when we re-initialize the event scene
// $$$ Do not know exactly what the field parameters mean
- const DataConfigurations::const_iterator i=m_collectionsConfigs.find(nam);
- const DataConfig* cfg = (i==m_collectionsConfigs.end()) ? 0 : &((*i).second);
- MCParticleCreator cr(new TEveTrackPropagator("","",new TEveMagFieldDuo(350, -3.5, 2.0)),
- new TEveCompound("MC_Particles","MC_Particles"),cfg);
- handler.collectionLoop((*j).first, cr);
- cr.close();
- particles = cr.particles;
+ if ( (icfg=m_collectionsConfigs.find("StartVertexPoints")) != cfgend ) {
+ StartVertexCreator cr("StartVertexPoints", len, (*icfg).second);
+ handler.collectionLoop((*j).first, cr);
+ printout(INFO,"Display","+++ StartVertexPoints: Filled %d start vertex points.....",cr.count);
+ ImportEvent(cr.element());
+ }
+ if ( (icfg=m_collectionsConfigs.find("MCParticles")) != cfgend ) {
+ MCParticleCreator cr(new TEveTrackPropagator("","",new TEveMagFieldDuo(350, -3.5, 2.0)),
+ new TEveCompound("MC_Particles","MC_Particles"),
+ icfg == cfgend ? 0 : &((*icfg).second));
+ handler.collectionLoop((*j).first, cr);
+ printout(INFO,"Display","+++ StartVertexPoints: Filled %d patricle tracks.....",cr.count);
+ cr.close();
+ particles = cr.particles;
+ }
}
}
}
diff --git a/DDEve/src/DisplayConfigurationParser.cpp b/DDEve/src/DisplayConfigurationParser.cpp
index 07bdde59f..c2e39253e 100644
--- a/DDEve/src/DisplayConfigurationParser.cpp
+++ b/DDEve/src/DisplayConfigurationParser.cpp
@@ -128,10 +128,10 @@ template <> void Converter<collection_configs>::operator()(xml_h e) const {
DisplayConfiguration::Config c;
c.name = e.attr<string>(_U(name));
c.type = DisplayConfiguration::COLLECTION;
- c.data.hits.color = e.hasAttr(_U(color)) ? e.attr<int>(_U(color)) : 0xBBBBBB;
+ c.data.hits.color = e.hasAttr(_U(color)) ? e.attr<int>(_U(color)) : 0xBBBBBB;
c.data.hits.alpha = e.hasAttr(_U(alpha)) ? e.attr<float>(_U(alpha)) : -1.0;
- c.data.hits.emax = e.hasAttr(u_emax) ? e.attr<float>(u_emax) : 25.0;
- c.data.hits.towerH = e.hasAttr(u_towerH) ? e.attr<float>(u_towerH) : 25.0;
+ c.data.hits.emax = e.hasAttr(u_emax) ? e.attr<float>(u_emax) : 25.0;
+ c.data.hits.towerH = e.hasAttr(u_towerH) ? e.attr<float>(u_towerH) : 25.0;
c.data.hits.threshold = e.hasAttr(_U(threshold)) ? e.attr<float>(_U(threshold)) : 0.0;
if ( e.hasAttr(u_hits) ) c.hits = e.attr<string>(u_hits);
if ( e.hasAttr(u_use) ) c.use = e.attr<string>(u_use);
@@ -153,10 +153,10 @@ template <> void Converter<view>::operator()(xml_h e) const {
ViewConfigurations* configs = (ViewConfigurations*)param;
DisplayConfiguration::ViewConfig c;
extract(c,e,DisplayConfiguration::VIEW);
+ c.name = e.attr<string>(_U(name));
c.type = e.attr<string>(_U(type));
c.show_structure = e.hasAttr(_U(structure)) ? e.attr<bool>(_U(structure)) : true;
c.show_sensitive = e.hasAttr(_U(sensitive)) ? e.attr<bool>(_U(sensitive)) : true;
- c.name = e.attr<string>(_U(name));
printout(INFO,"DisplayConfiguration","+++ View: %s sensitive:%d structure:%d.",
c.name.c_str(), c.show_sensitive, c.show_structure);
xml_coll_t(e,_Unicode(panel)).for_each(Converter<panel>(description,&c.subdetectors));
@@ -180,10 +180,10 @@ template <> void Converter<view>::operator()(xml_h e) const {
template <> void Converter<calodata>::operator()(xml_h e) const {
Configurations* configs = (Configurations*)param;
DisplayConfiguration::Config c;
- c.name = e.attr<string>(_U(name));
- c.type = DisplayConfiguration::CALODATA;
+ c.name = e.attr<string>(_U(name));
+ c.type = DisplayConfiguration::CALODATA;
if ( e.hasAttr(u_use) ) {
- c.use = e.attr<string>(u_use);
+ c.use = e.attr<string>(u_use);
c.hits = e.attr<string>(u_hits);
}
else {
@@ -220,8 +220,9 @@ template <> void Converter<collection>::operator()(xml_h e) const {
Configurations* configs = (Configurations*)param;
DisplayConfiguration::Config c;
c.name = e.attr<string>(_U(name));
- c.type = DisplayConfiguration::COLLECTION;
c.hits = e.attr<string>(u_hits);
+ c.type = DisplayConfiguration::COLLECTION;
+ c.use = e.hasAttr(u_use) ? e.attr<string>(u_use) : string("true");
c.data.hits.size = e.attr<float>(_U(size));
c.data.hits.type = e.attr<float>(_U(type));
c.data.hits.color = e.hasAttr(_U(color)) ? e.attr<int>(_U(color)) : kRed;
@@ -269,9 +270,9 @@ template <> void Converter<include>::operator()(xml_h e) const {
*/
template <> void Converter<display>::operator()(xml_h e) const {
Display* d = (Display*)param;
- if ( e.hasAttr(_Unicode(visLevel)) ) d->setVisLevel(e.attr<int>(_Unicode(visLevel)));
+ if ( e.hasAttr(_Unicode(visLevel)) ) d->setVisLevel(e.attr<int>(_Unicode(visLevel)));
if ( e.hasAttr(_Unicode(eventHandler)) ) d->setEventHandlerName(e.attr<std::string>(_Unicode(eventHandler)));
- if ( e.hasAttr(_Unicode(loadLevel)) ) d->setLoadLevel(e.attr<int>(_Unicode(loadLevel)));
+ if ( e.hasAttr(_Unicode(loadLevel)) ) d->setLoadLevel(e.attr<int>(_Unicode(loadLevel)));
}
/** Convert display configuration elements of tag type ddeve
diff --git a/DDEve/src/HitActors.cpp b/DDEve/src/HitActors.cpp
index 646cd1b83..dde055102 100644
--- a/DDEve/src/HitActors.cpp
+++ b/DDEve/src/HitActors.cpp
@@ -27,9 +27,11 @@ using namespace std;
using namespace dd4hep;
#ifdef HAVE_GEANT4_UNITS
-#define MM_2_CM 1.0
+#define MM_2_CM 1.0
+#define MEV_TO_GEV 1000.0
#else
-#define MM_2_CM 0.1
+#define MM_2_CM 0.1
+#define MEV_TO_GEV 1.0
#endif
/// Action callback of this functor:
@@ -40,7 +42,6 @@ void EtaPhiHistogramActor::operator()(const DDEveHit& hit) {
/// Standard initializing constructor
PointsetCreator::PointsetCreator(const std::string& collection, size_t length)
- : pointset(0), deposit(0), count(0)
{
pointset = new TEvePointSet(collection.c_str(),length);
pointset->SetMarkerSize(0.2);
@@ -48,13 +49,13 @@ PointsetCreator::PointsetCreator(const std::string& collection, size_t length)
/// Standard initializing constructor
PointsetCreator::PointsetCreator(const std::string& collection, size_t length, const DisplayConfiguration::Config& cfg)
- : pointset(0), deposit(0), count(0)
{
pointset = new TEvePointSet(collection.c_str(),length);
pointset->SetMarkerSize(cfg.data.hits.size);
pointset->SetMarkerStyle(cfg.data.hits.type);
//pointset->SetMarkerAlpha(cfg.data.hits.alpha);
pointset->SetMainColor(cfg.data.hits.color);
+ threshold = cfg.data.hits.threshold * MEV_TO_GEV;
}
/// Return eve element
TEveElement* PointsetCreator::element() const {
@@ -74,12 +75,14 @@ PointsetCreator::~PointsetCreator() {
/// Action callback of this functor:
void PointsetCreator::operator()(const DDEveHit& hit) {
- pointset->SetPoint(count++, hit.x*MM_2_CM, hit.y*MM_2_CM, hit.z*MM_2_CM);
+ if ( hit.deposit > threshold ) {
+ deposit += hit.deposit;
+ pointset->SetPoint(count++, hit.x*MM_2_CM, hit.y*MM_2_CM, hit.z*MM_2_CM);
+ }
}
/// Standard initializing constructor
BoxsetCreator::BoxsetCreator(const std::string& collection, size_t /*length */, const DisplayConfiguration::Config& cfg)
- : boxset(0), emax(1e12), towerH(1e12), deposit(0.0), count(0)
{
emax = cfg.data.hits.emax;
towerH = cfg.data.hits.towerH;
@@ -94,7 +97,6 @@ BoxsetCreator::BoxsetCreator(const std::string& collection, size_t /*length */,
/// Standard initializing constructor
BoxsetCreator::BoxsetCreator(const std::string& collection, size_t /*length */)
- : boxset(0), emax(1e12), towerH(1e12), deposit(0.0), count(0)
{
boxset = new TEveBoxSet(collection.c_str());
boxset->SetMainTransparency(0);
diff --git a/DDEve/src/ParticleActors.cpp b/DDEve/src/ParticleActors.cpp
index 56aaea7d7..6bc62890a 100644
--- a/DDEve/src/ParticleActors.cpp
+++ b/DDEve/src/ParticleActors.cpp
@@ -13,11 +13,14 @@
// Framework include files
#include "DDEve/ParticleActors.h"
-#include "DD4hep/Objects.h"
#include "DD4hep/DD4hepUnits.h"
+#include "DD4hep/Printout.h"
+#include "DD4hep/Objects.h"
-#include "TEveCompound.h"
#include "TEveTrack.h"
+#include "TEveBoxSet.h"
+#include "TEvePointSet.h"
+#include "TEveCompound.h"
#include "TEveTrackPropagator.h"
#include "TParticle.h"
@@ -28,11 +31,13 @@ using namespace std;
using namespace dd4hep;
#ifdef HAVE_GEANT4_UNITS
-#define CM_2_MM 1.0
-#define MM_2_CM 1.0
+#define CM_2_MM 1.0
+#define MM_2_CM 1.0
+#define MEV_TO_GEV 1000.0
#else
-#define CM_2_MM 10.0
-#define MM_2_CM 0.1
+#define CM_2_MM 10.0
+#define MM_2_CM 0.1
+#define MEV_TO_GEV 1.0
#endif
namespace {
@@ -43,7 +48,7 @@ namespace {
/// Standard initializing constructor
MCParticleCreator::MCParticleCreator(TEveTrackPropagator* p, TEveCompound* ps, const DisplayConfiguration::Config* cfg)
- : propagator(p), particles(ps), count(0), lineWidth(4)
+ : propagator(p), particles(ps)
{
propagator->SetName("Track propagator for charged particles");
propagator->SetMaxR(1000);
@@ -55,8 +60,10 @@ MCParticleCreator::MCParticleCreator(TEveTrackPropagator* p, TEveCompound* ps, c
propagator->RefPMAtt().SetMarkerSize(1.0);
if ( cfg ) {
lineWidth = cfg->data.hits.width;
+ threshold = cfg->data.hits.threshold * MEV_TO_GEV;
propagator->RefPMAtt().SetMarkerSize(cfg->data.hits.size);
propagator->RefPMAtt().SetMarkerStyle(cfg->data.hits.type);
+ printout(ALWAYS,"MCParticleCreator","+++ Minimal particle energy: %8.3g [GeV]",threshold);
}
}
@@ -119,16 +126,16 @@ void MCParticleCreator::operator()(const DDEveParticle& p) {
TEveVector dir = end-start;
// Tracks longer than 100 micron and energy > 100 MeV
- if ( dir.R()*CM_2_MM > 100e-3 && p.energy > 10e0 ) {
+ if ( p.energy > 10e0 && p.energy > threshold && dir.R()*CM_2_MM > 100e-3 ) {
TDatabasePDG* db = TDatabasePDG::Instance();
TParticlePDG* def = db->GetParticle(p.pdgID);
TParticle part(p.pdgID,
0,0,0,0,0,
p.psx*MEV_2_GEV, p.psy*MEV_2_GEV, p.psz*MEV_2_GEV, p.energy*MEV_2_GEV,
p.vsx*MM_2_CM, p.vsy*MM_2_CM, p.vsz*MM_2_CM, p.time);
-
- TEveTrack* t = new TEveTrack(&part,p.id,propagator);
+ TEveTrack* t = new TEveTrack(&part,p.id,propagator);
+ ++count;
// Add start-vertex as path mark
//t->AddPathMark(TEvePathMark(TEvePathMark::kDecay,start));
// Add end-vertex as path mark (kDecay)
@@ -149,7 +156,7 @@ void MCParticleCreator::operator()(const DDEveParticle& p) {
p.vex*MM_2_CM, p.vey*MM_2_CM, p.vez*MM_2_CM,
dir.R(),
p.psx*MEV_2_GEV, p.psy*MEV_2_GEV, p.psz*MEV_2_GEV, p.energy*MEV_2_GEV));
-
+
// Add element to collection
int pdg = abs(p.pdgID);
if ( pdg == 11 )
@@ -170,9 +177,41 @@ void MCParticleCreator::operator()(const DDEveParticle& p) {
addCompound("Protons", t);
else
addCompound("Other", t);
- //cout << "Add particle " << p.id << " to compound." << endl;
}
else {
- cout << "SKIP particle " << p.id << "." << endl;
+ printout(ALWAYS,"MCParticleCreator","+++ SKIP particle %4d. Energy: %8.3g [MeV]",p.id,p.energy);
}
}
+
+/// Standard initializing constructor
+StartVertexCreator::StartVertexCreator(const std::string& collection, size_t length)
+ : pointset(0), deposit(0), count(0)
+{
+ pointset = new TEvePointSet(collection.c_str(),length);
+ pointset->SetMarkerSize(0.2);
+}
+
+/// Standard initializing constructor
+StartVertexCreator::StartVertexCreator(const std::string& collection, size_t length, const DisplayConfiguration::Config& cfg)
+ : pointset(0), deposit(0), count(0)
+{
+ pointset = new TEvePointSet(collection.c_str(),length);
+ pointset->SetMarkerSize(cfg.data.hits.size);
+ pointset->SetMarkerStyle(cfg.data.hits.type);
+ //pointset->SetMarkerAlpha(cfg.data.hits.alpha);
+ pointset->SetMainColor(cfg.data.hits.color);
+}
+/// Return eve element
+TEveElement* StartVertexCreator::element() const {
+ return pointset;
+}
+
+/// Standard destructor
+StartVertexCreator::~StartVertexCreator() {
+}
+
+/// Action callback of this functor:
+void StartVertexCreator::operator()(const DDEveParticle& p) {
+ pointset->SetPoint(count++, p.vsx*MM_2_CM, p.vsy*MM_2_CM, p.vsz*MM_2_CM);
+}
+
diff --git a/DDG4/include/DDG4/Geant4Converter.h b/DDG4/include/DDG4/Geant4Converter.h
index d5e097178..67f3d0989 100644
--- a/DDG4/include/DDG4/Geant4Converter.h
+++ b/DDG4/include/DDG4/Geant4Converter.h
@@ -101,6 +101,7 @@ namespace dd4hep {
/// Convert the geometry type volume placement into the corresponding Geant4 object(s).
virtual void* handlePlacement(const std::string& name, const TGeoNode* node) const;
+ virtual void* handlePlacement2(const std::pair<const TGeoNode*,const TGeoNode*>& node) const;
virtual void* handleAssembly(const std::string& name, const TGeoNode* node) const;
/// Convert the geometry type field into the corresponding Geant4 object(s).
diff --git a/DDG4/include/DDG4/Geant4ParticleGun.h b/DDG4/include/DDG4/Geant4ParticleGun.h
index 529d679e7..ff633b08a 100644
--- a/DDG4/include/DDG4/Geant4ParticleGun.h
+++ b/DDG4/include/DDG4/Geant4ParticleGun.h
@@ -64,10 +64,7 @@ namespace dd4hep {
/// Shot number in sequence
int m_shotNo;
/// Particle modification. Caller presets defaults to: ( direction = m_direction, momentum = m_energy)
- virtual void getParticleDirection(int, ROOT::Math::XYZVector& direction, double& momentum) const {
- direction = m_direction;
- momentum = m_energy;
- }
+ virtual void getParticleDirection(int, ROOT::Math::XYZVector& direction, double& momentum) const;
public:
/// Standard constructor
Geant4ParticleGun(Geant4Context* context, const std::string& name);
diff --git a/DDG4/src/Geant4Converter.cpp b/DDG4/src/Geant4Converter.cpp
index 5fe2d9b87..dcb7a1a4d 100644
--- a/DDG4/src/Geant4Converter.cpp
+++ b/DDG4/src/Geant4Converter.cpp
@@ -847,31 +847,18 @@ void* Geant4Converter::handlePlacement(const string& name, const TGeoNode* node)
printout(lvl, "Geant4Converter", "++ Placement %s not converted [Veto'ed for simulation]",node->GetName());
return 0;
}
-#if 0
- if ( g4 ) {
- const TGeoNode* old_node = (*g4it).first.ptr();
- const TGeoNode* new_node = node;
- TGeoVolume* old_vol = old_node->GetVolume();
- TGeoVolume* new_vol = new_node->GetVolume();
- TGeoMatrix* old_tr = old_node->GetMatrix();
- TGeoMatrix* new_tr = new_node->GetMatrix();
- printout(debugReflections ? ALWAYS : lvl, "Geant4Converter",
- "+++ Placement: **** Reuse: OLD: %s vol: %s %s NEW: %s vol: %s %s",
- old_node->GetName(), old_vol->GetName(), is_left_handed(old_tr) ? " REFLECTED" : "NOT REFLECTED",
- new_node->GetName(), new_vol->GetName(), is_left_handed(new_tr) ? " REFLECTED" : "NOT REFLECTED");
- }
-#endif
g4 = nullptr;
if (!g4) {
TGeoVolume* mot_vol = node->GetMotherVolume();
TGeoMatrix* tr = node->GetMatrix();
if (!tr) {
- except("Geant4Converter", "++ Attempt to handle placement without transformation:%p %s of type %s vol:%p", node,
- node->GetName(), node->IsA()->GetName(), vol);
+ except("Geant4Converter",
+ "++ Attempt to handle placement without transformation:%p %s of type %s vol:%p",
+ node, node->GetName(), node->IsA()->GetName(), vol);
}
else if (0 == vol) {
- except("Geant4Converter", "++ Unknown G4 volume:%p %s of type %s ptr:%p", node, node->GetName(),
- node->IsA()->GetName(), vol);
+ except("Geant4Converter", "++ Unknown G4 volume:%p %s of type %s ptr:%p",
+ node, node->GetName(), node->IsA()->GetName(), vol);
}
else {
int copy = node->GetNumber();
@@ -894,14 +881,14 @@ void* Geant4Converter::handlePlacement(const string& name, const TGeoNode* node)
transform.dx(), transform.dy(), transform.dz());
return 0;
}
- else if ( node_is_assembly ) {
+ G4Scale3D scale;
+ G4Rotate3D rot;
+ G4Translate3D trans;
+ if ( node_is_assembly ) {
//
// Node is an assembly:
// Imprint the assembly. The mother MUST already be transformed.
//
- G4Scale3D scale;
- G4Rotate3D rot;
- G4Translate3D trans;
transform.getDecomposition(scale, rot, trans);
printout(lvl, "Geant4Converter", "+++ Assembly: makeImprint: dau:%-12s %s in mother %-12s "
"Tr:x=%8.1f y=%8.1f z=%8.1f Scale:x=%4.2f y=%4.2f z=%4.2f",
@@ -918,9 +905,6 @@ void* Geant4Converter::handlePlacement(const string& name, const TGeoNode* node)
else if ( node != info.manager->GetTopNode() && volIt == info.g4Volumes.end() ) {
throw logic_error("Geant4Converter: Invalid mother volume found!");
}
- G4Scale3D scale;
- G4Rotate3D rot;
- G4Translate3D trans;
G4LogicalVolume* g4vol = info.g4Volumes[vol];
G4LogicalVolume* g4mot = info.g4Volumes[mot_vol];
G4PhysicalVolumesPair pvPlaced =
@@ -962,6 +946,120 @@ void* Geant4Converter::handlePlacement(const string& name, const TGeoNode* node)
return g4;
}
+void* Geant4Converter::handlePlacement2(const std::pair<const TGeoNode*,const TGeoNode*>& n) const {
+ const TGeoNode* par_node = n.first;
+ const TGeoNode* node = n.second;
+ Geant4GeometryInfo& info = data();
+ PrintLevel lvl = debugPlacements ? ALWAYS : outputLevel;
+ Geant4GeometryMaps::PlacementMap::const_iterator g4it = info.g4Placements.find(node);
+ G4VPhysicalVolume* g4 = (g4it == info.g4Placements.end()) ? 0 : (*g4it).second;
+ TGeoVolume* vol = node->GetVolume();
+ Volume _v(vol);
+
+ if ( _v.testFlagBit(Volume::VETO_SIMU) ) {
+ printout(lvl, "Geant4Converter", "++ Placement %s not converted [Veto'ed for simulation]",node->GetName());
+ return 0;
+ }
+ g4 = nullptr;
+ if (!g4) {
+ TGeoVolume* mot_vol = par_node ? par_node->GetVolume() : node->GetMotherVolume();
+ TGeoMatrix* tr = node->GetMatrix();
+ if (!tr) {
+ except("Geant4Converter",
+ "++ Attempt to handle placement without transformation:%p %s of type %s vol:%p",
+ node, node->GetName(), node->IsA()->GetName(), vol);
+ }
+ else if (0 == vol) {
+ except("Geant4Converter", "++ Unknown G4 volume:%p %s of type %s ptr:%p",
+ node, node->GetName(), node->IsA()->GetName(), vol);
+ }
+ else {
+ int copy = node->GetNumber();
+ bool node_is_reflected = is_left_handed(tr);
+ bool node_is_assembly = vol->IsA() == TGeoVolumeAssembly::Class();
+ bool mother_is_assembly = mot_vol ? mot_vol->IsA() == TGeoVolumeAssembly::Class() : false;
+ MyTransform3D transform(tr->GetTranslation(),tr->IsRotation() ? tr->GetRotationMatrix() : s_identity_rot);
+ Geant4GeometryMaps::VolumeMap::const_iterator volIt = info.g4Volumes.find(mot_vol);
+
+ if ( mother_is_assembly ) {
+ //
+ // Mother is an assembly:
+ // Nothing to do here, because:
+ // -- placed volumes were already added before in "handleAssembly"
+ // -- imprint cannot be made, because this requires a logical volume as a mother
+ //
+ printout(lvl, "Geant4Converter", "+++ Assembly: **** : dau:%s "
+ "to mother %s Tr:x=%8.3f y=%8.3f z=%8.3f",
+ vol->GetName(), mot_vol->GetName(),
+ transform.dx(), transform.dy(), transform.dz());
+ return 0;
+ }
+ G4Scale3D scale;
+ G4Rotate3D rot;
+ G4Translate3D trans;
+ if ( node_is_assembly ) {
+ //
+ // Node is an assembly:
+ // Imprint the assembly. The mother MUST already be transformed.
+ //
+ transform.getDecomposition(scale, rot, trans);
+ printout(lvl, "Geant4Converter", "+++ Assembly: makeImprint: dau:%-12s %s in mother %-12s "
+ "Tr:x=%8.1f y=%8.1f z=%8.1f Scale:x=%4.2f y=%4.2f z=%4.2f",
+ node->GetName(), node_is_reflected ? "(REFLECTED)" : "",
+ mot_vol ? mot_vol->GetName() : "<unknown>",
+ transform.dx(), transform.dy(), transform.dz(),
+ scale.xx(), scale.yy(), scale.zz());
+ Geant4AssemblyVolume* ass = (Geant4AssemblyVolume*)info.g4AssemblyVolumes[node];
+ Geant4AssemblyVolume::Chain chain;
+ chain.emplace_back(node);
+ ass->imprint(info,node,chain,ass,(*volIt).second, transform, copy, checkOverlaps);
+ return 0;
+ }
+ else if ( node != info.manager->GetTopNode() && volIt == info.g4Volumes.end() ) {
+ throw logic_error("Geant4Converter: Invalid mother volume found!");
+ }
+ G4LogicalVolume* g4vol = info.g4Volumes[vol];
+ G4LogicalVolume* g4mot = info.g4Volumes[mot_vol];
+ G4PhysicalVolumesPair pvPlaced =
+ G4ReflectionFactory::Instance()->Place(transform, // no rotation
+ node->GetName(), // its name
+ g4vol, // its logical volume
+ g4mot, // its mother (logical) volume
+ false, // no boolean operations
+ copy, // its copy number
+ checkOverlaps);
+ transform.getDecomposition(scale, rot, trans);
+ printout(debugReflections ? ALWAYS : lvl, "Geant4Converter",
+ "+++ Place %svolume %-12s in mother %-12s "
+ "Tr:x=%8.1f y=%8.1f z=%8.1f Scale:x=%4.2f y=%4.2f z=%4.2f",
+ node_is_reflected ? "REFLECTED " : "", _v.name(),
+ mot_vol ? mot_vol->GetName() : "<unknown>",
+ transform.dx(), transform.dy(), transform.dz(),
+ scale.xx(), scale.yy(), scale.zz());
+ // First 2 cases can be combined.
+ // Leave them separated for debugging G4ReflectionFactory for now...
+ if ( node_is_reflected && !pvPlaced.second )
+ return info.g4Placements[node] = pvPlaced.first;
+ else if ( !node_is_reflected && !pvPlaced.second )
+ return info.g4Placements[node] = pvPlaced.first;
+ //G4LogicalVolume* g4refMoth = G4ReflectionFactory::Instance()->GetReflectedLV(g4mot);
+ // Now deal with valid pvPlaced.second ...
+ if ( node_is_reflected )
+ return info.g4Placements[node] = pvPlaced.first;
+ else if ( !node_is_reflected )
+ return info.g4Placements[node] = pvPlaced.first;
+ g4 = pvPlaced.second ? pvPlaced.second : pvPlaced.first;
+ }
+ info.g4Placements[node] = g4;
+ printout(ERROR, "Geant4Converter", "++ DEAD code. Should not end up here!");
+ }
+ else {
+ printout(ERROR, "Geant4Converter", "++ Attempt to DOUBLE-place physical volume: %s No:%d",
+ node->GetName(), node->GetNumber());
+ }
+ return g4;
+}
+
/// Convert the geometry type region into the corresponding Geant4 object(s).
void* Geant4Converter::handleRegion(Region region, const set<const TGeoVolume*>& /* volumes */) const {
G4Region* g4 = data().g4Regions[region];
@@ -1458,6 +1556,14 @@ namespace {
handle(o, (*i).second, pmf);
}
}
+ template <typename O, typename C, typename F> void handleRMap_(const O* o, const C& c, F pmf) {
+ for (typename C::const_iterator i = c.begin(); i != c.end(); ++i) {
+ const auto& cc = (*i).second;
+ for (const auto& j : cc) {
+ (o->*pmf)(j);
+ }
+ }
+ }
}
/// Create geometry conversion
@@ -1477,7 +1583,7 @@ Geant4Converter& Geant4Converter::create(DetElement top) {
//debugMaterials = true;
//debugElements = true;
debugReflections = true;
-
+ debugPlacements = true;
#if ROOT_VERSION_CODE >= ROOT_VERSION(6,17,0)
handleArray(this, geo.manager->GetListOfGDMLMatrices(), &Geant4Converter::handleMaterialProperties);
handleArray(this, geo.manager->GetListOfOpticalSurfaces(), &Geant4Converter::handleOpticalSurface);
@@ -1497,6 +1603,7 @@ Geant4Converter& Geant4Converter::create(DetElement top) {
handleRMap(this, *m_data, &Geant4Converter::handleAssembly);
// Now place all this stuff appropriately
handleRMap(this, *m_data, &Geant4Converter::handlePlacement);
+ //handleRMap_(this, *m_places, &Geant4Converter::handlePlacement2);
#if ROOT_VERSION_CODE >= ROOT_VERSION(6,17,0)
/// Handle concrete surfaces
handleArray(this, geo.manager->GetListOfSkinSurfaces(), &Geant4Converter::handleSkinSurface);
diff --git a/DDG4/src/Geant4IsotropeGenerator.cpp b/DDG4/src/Geant4IsotropeGenerator.cpp
index 1a8aaddd8..54588aafd 100644
--- a/DDG4/src/Geant4IsotropeGenerator.cpp
+++ b/DDG4/src/Geant4IsotropeGenerator.cpp
@@ -133,5 +133,5 @@ void Geant4IsotropeGenerator::getParticleDirection(int num, ROOT::Math::XYZVecto
break;
}
except("Unknown distribution densitiy: %s. Cannot generate primaries.",
- m_distribution.c_str());
+ m_distribution.c_str());
}
diff --git a/DDG4/src/Geant4ParticleGenerator.cpp b/DDG4/src/Geant4ParticleGenerator.cpp
index 898c2ac7d..d16cf736c 100644
--- a/DDG4/src/Geant4ParticleGenerator.cpp
+++ b/DDG4/src/Geant4ParticleGenerator.cpp
@@ -42,8 +42,8 @@ Geant4ParticleGenerator::Geant4ParticleGenerator(Geant4Context* ctxt, const stri
declareProperty("Energy", m_energy = 50 * CLHEP::MeV);
declareProperty("Multiplicity", m_multiplicity = 1);
declareProperty("Mask", m_mask = 0);
- declareProperty("Position", m_position);
- declareProperty("Direction", m_direction);
+ declareProperty("Position", m_position = ROOT::Math::XYZVector(0.,0.,0.));
+ declareProperty("Direction", m_direction = ROOT::Math::XYZVector(1.,1.,1.));
}
/// Default destructor
@@ -117,7 +117,7 @@ void Geant4ParticleGenerator::operator()(G4Event*) {
Geant4Vertex* vtx = new Geant4Vertex();
int multiplicity = m_multiplicity;
- ROOT::Math::XYZVector unit_direction, direction, position = m_position;
+ ROOT::Math::XYZVector unit_direction, position = m_position;
getVertexPosition(position);
getParticleMultiplicity(multiplicity);
vtx->mask = m_mask;
@@ -127,8 +127,8 @@ void Geant4ParticleGenerator::operator()(G4Event*) {
inter->vertices[m_mask].emplace_back( vtx );
for(int i=0; i<m_multiplicity; ++i) {
double momentum = m_energy;
+ ROOT::Math::XYZVector direction = m_direction;
Particle* p = new Particle();
- direction = m_direction;
getParticleDirection(i, direction, momentum);
unit_direction = direction.unit();
p->id = inter->nextPID();
@@ -155,8 +155,9 @@ void Geant4ParticleGenerator::operator()(G4Event*) {
// p->vez = vtx->z;
inter->particles.emplace(p->id,p);
vtx->out.insert(p->id);
- printout(INFO,name(),"Particle [%d] %s %.3f GeV direction:(%6.3f %6.3f %6.3f)",
- p->id, m_particleName.c_str(), momentum/CLHEP::GeV,
+ printout(INFO,name(),"Particle [%d] %-12s Mom:%.3f GeV vertex:(%6.3f %6.3f %6.3f)[mm] direction:(%6.3f %6.3f %6.3f)",
+ p->id, m_particleName.c_str(), momentum/CLHEP::GeV,
+ vtx->x/CLHEP::mm, vtx->y/CLHEP::mm, vtx->z/CLHEP::mm,
direction.X(), direction.Y(), direction.Z());
}
diff --git a/DDG4/src/Geant4ParticleGun.cpp b/DDG4/src/Geant4ParticleGun.cpp
index 1bfd9b11d..b4bdea577 100644
--- a/DDG4/src/Geant4ParticleGun.cpp
+++ b/DDG4/src/Geant4ParticleGun.cpp
@@ -32,10 +32,13 @@ Geant4ParticleGun::Geant4ParticleGun(Geant4Context* ctxt, const string& nam)
{
InstanceCount::increment(this);
m_needsControl = true;
- declareProperty("isotrop", m_isotrop = false);
declareProperty("Standalone", m_standalone = true);
+ declareProperty("mask", m_mask);
+ declareProperty("isotrop", m_isotrop = false);
+ declareProperty("standalone", m_standalone);
// Backwards compatibility: Un-capitalize
declareProperty("position", m_position);
+ declareProperty("distribution", m_distribution);
declareProperty("direction", m_direction);
declareProperty("energy", m_energy);
declareProperty("particle", m_particleName);
@@ -48,17 +51,18 @@ Geant4ParticleGun::~Geant4ParticleGun() {
InstanceCount::decrement(this);
}
+/// Particle modification. Caller presets defaults to: ( direction = m_direction, momentum = m_energy)
+void Geant4ParticleGun::getParticleDirection(int num, ROOT::Math::XYZVector& direction, double& momentum) const {
+ ( m_isotrop )
+ ? this->Geant4IsotropeGenerator::getParticleDirection(num, direction, momentum)
+ : this->Geant4ParticleGenerator::getParticleDirection(num, direction, momentum);
+}
+
/// Callback to generate primary particles
void Geant4ParticleGun::operator()(G4Event* event) {
- if ( m_isotrop ) {
- double mom = 0.0;
- this->Geant4IsotropeGenerator::getParticleDirection(0,m_direction,mom);
- }
- else {
- double r = m_direction.R(), eps = numeric_limits<float>::epsilon();
- if ( r > eps ) {
- m_direction.SetXYZ(m_direction.X()/r, m_direction.Y()/r, m_direction.Z()/r);
- }
+ double r = m_direction.R(), eps = numeric_limits<float>::epsilon();
+ if ( r > eps ) {
+ m_direction.SetXYZ(m_direction.X()/r, m_direction.Y()/r, m_direction.Z()/r);
}
if ( m_standalone ) {
diff --git a/examples/ClientTests/compact/NestedBoxReflection.xml b/examples/ClientTests/compact/NestedBoxReflection.xml
index 745c5bcc0..af560b628 100644
--- a/examples/ClientTests/compact/NestedBoxReflection.xml
+++ b/examples/ClientTests/compact/NestedBoxReflection.xml
@@ -22,6 +22,7 @@
<vis name="VisibleBlue" alpha="0.6" r="0.0" g="0.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="VisibleYellow" alpha="0.6" r="1.0" g="1.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="VisibleGreen" alpha="0.6" r="0.0" g="1.0" b="0.0" showDaughters="true" visible="true"/>
+ <vis name="VisibleCyan" alpha="0.6" r="0.0" g="1.0" b="1.0" showDaughters="true" visible="true"/>
</display>
<limits>
@@ -30,6 +31,7 @@
</limitset>
</limits>
+ <geometry reflect="true"/>
<detectors>
<detector id="1" name="NestedBox" type="NestedBoxReflection" readout="NestedBoxHits" vis="VisibleGreen" limits="cal_limits">
<comment>A box with 3 boxes inside spanning a coordinate system</comment>
diff --git a/examples/ClientTests/eve/NestedBoxReflection.xml b/examples/ClientTests/eve/NestedBoxReflection.xml
index 86b1fa628..264df55ec 100644
--- a/examples/ClientTests/eve/NestedBoxReflection.xml
+++ b/examples/ClientTests/eve/NestedBoxReflection.xml
@@ -1,7 +1,8 @@
<ddeve>
<display visLevel="7" loadLevel="3"/>
- <collection name="NestedBoxHits" hits="PointSet" color="kMagenta" size="0.3" type="20"/>
- <collection name="MCParticles" hits="Particles" size="0.2" width="1" type="kCircle"/>
+ <collection name="NestedBoxHits" hits="PointSet" color="kMagenta" size="0.3" type="20" threshold="1000*eV"/>
+ <collectiona name="StartVertexPoints" hits="ParticleOrigin" color="kCyan" size="2.3" type="20"/>
+ <collectiona name="MCParticles" hits="Particles" width="0.01" size="0.01" type="kCircle" threshold="1.5*GeV" use="false"/>
<!--
<view name="3D Trackers R-Phi (Global)" type="RhoPhiProjection">
<detelement name="NestedBox" load_geo="-1" alpha="0.5"/>
diff --git a/examples/ClientTests/scripts/Check_reflection.py b/examples/ClientTests/scripts/Check_reflection.py
index 20f5184be..ab7681495 100644
--- a/examples/ClientTests/scripts/Check_reflection.py
+++ b/examples/ClientTests/scripts/Check_reflection.py
@@ -7,12 +7,13 @@
"""
from __future__ import absolute_import, unicode_literals
import logging
+import math
import time
import sys
import os
from g4units import rad, GeV, MeV, mm, m
from ddsix.moves import range
-from math import pi
+
logging.basicConfig(format='%(levelname)s: %(message)s', level=logging.INFO)
logger = logging.getLogger(__name__)
@@ -92,13 +93,14 @@ def run():
# Setup particle gun
geant4.setupGun(name="Gun",
particle='e-',
- energy=1000*GeV,
+ energy=1000 * GeV,
+ isotrop=True,
multiplicity=1,
- position=(0*m, 0*m, 0*m),
- PhiMin=0.0*rad,
- PhiMax=math.pi*2.0*rad,
- ThetaMin=0.0*rad,
- ThetaMax=math.pi*rad)
+ position=(0 * m, 0 * m, 0 * m),
+ PhiMin=0.0 * rad,
+ PhiMax=math.pi * 2.0 * rad,
+ ThetaMin=0.0 * rad,
+ ThetaMax=math.pi * rad)
logger.info("# ....and handle the simulation particles.")
part = DDG4.GeneratorAction(kernel, str('Geant4ParticleHandler/ParticleHandler'))
diff --git a/examples/ClientTests/src/NestedBoxReflection_geo.cpp b/examples/ClientTests/src/NestedBoxReflection_geo.cpp
index 1897890ab..c8cf2264d 100644
--- a/examples/ClientTests/src/NestedBoxReflection_geo.cpp
+++ b/examples/ClientTests/src/NestedBoxReflection_geo.cpp
@@ -82,28 +82,31 @@ namespace {
double by = box.y();
double bz = box.z();
Material mat = description.material("Si");
- Box small_box(bx*0.2, by*0.2, bz*0.2);
- const char* cols[4] = {"VisibleRed","VisibleBlue","VisibleGreen","VisibleYellow"};
+ bool sens = level == 2 || level == 0;
+ Box long_box(bx*0.2, by*0.2, bz*0.4);
+ Box flat_box(bx*1.0, by*1.0, bz*0.1);
+ Box mini_box(bx*0.2, by*0.2, bz*0.2);
+ const char* cols[5] = {"VisibleRed","VisibleBlue","VisibleGreen","VisibleYellow","VisibleCyan"};
const char* c;
PlacedVolume pv;
Volume v;
- c = cols[(0+level)%4];
- v = Volume(_toString(1,"box%d"), small_box, mat);
+ c = cols[(0+level)%5];
+ v = Volume(_toString(1,"box%d"), mini_box, mat);
v.setRegion(vol.region());
- v.setSensitiveDetector(sensitive);
v.setLimitSet(vol.limitSet());
+ if ( !sens ) v.setSensitiveDetector(sensitive);
v.setVisAttributes(description, c);
pv = vol.placeVolume(v, Position(0,0,0));
pv.addPhysVolID(_toString(level,"lvl%d"), 1);
printout(INFO,"NestedBoxReflection","++ Volume: %s Color: %s", v.name(), c);
place_boxes(level-1, v);
- c = cols[(1+level)%4];
- v = Volume(_toString(2,"box%d"), small_box, mat);
+ c = cols[(1+level)%5];
+ v = Volume(_toString(2,"box%d"), long_box, mat);
v.setRegion(vol.region());
- v.setSensitiveDetector(sensitive);
v.setLimitSet(vol.limitSet());
+ if ( !sens ) v.setSensitiveDetector(sensitive);
v.setVisAttributes(description, c);
pv = vol.placeVolume(v, Position(0.8*bx,0,0));
pv.addPhysVolID(_toString(level,"lvl%d"), 2);
@@ -115,11 +118,11 @@ namespace {
pv = vol.placeVolume(v, Position(0.4*bx,0,0));
printout(INFO,"NestedBoxReflection","++ Volume: %s Color: %s", v.name(), c);
- c = cols[(2+level)%4];
- v = Volume(_toString(3,"box%d"), small_box, mat);
+ c = cols[(2+level)%5];
+ v = Volume(_toString(3,"box%d"), mini_box, mat);
v.setRegion(vol.region());
- v.setSensitiveDetector(sensitive);
v.setLimitSet(vol.limitSet());
+ if ( !sens ) v.setSensitiveDetector(sensitive);
v.setVisAttributes(description, c);
pv = vol.placeVolume(v, Position(0,0.8*by,0));
pv.addPhysVolID(_toString(level,"lvl%d"), 3);
@@ -131,11 +134,11 @@ namespace {
pv = vol.placeVolume(v, Position(0,0.4*by,0));
printout(INFO,"NestedBoxReflection","++ Volume: %s Color: %s", v.name(), c);
- c = cols[(3+level)%4];
- v = Volume(_toString(4,"box%d"), small_box, mat);
+ c = cols[(3+level)%5];
+ v = Volume(_toString(4,"box%d"), mini_box, mat);
v.setRegion(vol.region());
- v.setSensitiveDetector(sensitive);
v.setLimitSet(vol.limitSet());
+ if ( !sens ) v.setSensitiveDetector(sensitive);
v.setVisAttributes(description, c);
pv = vol.placeVolume(v, Position(0,0,0.8*bz));
pv.addPhysVolID(_toString(level,"lvl%d"), 4);
@@ -146,6 +149,19 @@ namespace {
v.setVisAttributes(description, c);
pv = vol.placeVolume(v, Position(0,0,0.4*bz));
printout(INFO,"NestedBoxReflection","++ Volume: %s Color: %s", v.name(), c);
+
+ c = cols[(4+level)%5];
+ v = Volume(_toString(5,"box%d"), flat_box, mat);
+ v.setRegion(vol.region());
+ v.setLimitSet(vol.limitSet());
+ v.setVisAttributes(description, c);
+ if ( !sens ) v.setSensitiveDetector(sensitive);
+ TGeoHMatrix* mtx = detail::matrix::_transform(Position(0,0,0.9*bz));
+ mtx->ReflectZ(kTRUE, kFALSE);
+ pv = vol.placeVolume(v, mtx);
+ pv.addPhysVolID(_toString(level,"lvl%d"), 5);
+ printout(INFO,"NestedBoxReflection","++ Volume: %s Color: %s", v.name(), c);
+ place_boxes(level-1, v);
}
}
--
GitLab