eic-opticks/u4/U4Touchable.h

0001 #pragma once
0002 /**
0003 U4Touchable.h
0004 ================
0005
0006 The touchable history is from inside to outside, so "deeper" corresponds to wider volumes.
0007
0008 A version of this is also available within the monolith at
0009
0010     $JUNOTOP/junosw/Simulation/DetSimV2/SimUtil/SimUtil/S4Touchable.h
0011
0012 1042
0013     G4VTouchable is an ordinary class
0014 1120
0015     G4VTouchable is using "alias" for G4TouchableHistory
0016
0017 **/
0018
0019 #include <string>
0020 #include <cstring>
0021 #include <csignal>
0022 #include <sstream>
0023 #include <iomanip>
0024
0025 #include "G4VPhysicalVolume.hh"
0026 #include "G4LogicalVolume.hh"
0027 #include "G4TouchableHistory.hh"
0028 #include "G4VSolid.hh"
0029
0030 struct U4Touchable
0031 {
0032     enum { MATCH_ALL, MATCH_START, MATCH_END } ;
0033     static bool Match(      const char* s, const char* q, int mode) ;
0034     static bool MatchAll(   const char* s, const char* q) ;
0035     static bool MatchStart( const char* s, const char* q) ;
0036     static bool MatchEnd(   const char* s, const char* q) ;
0037
0038     template<typename T>
0039     static const G4VPhysicalVolume* FindPV( const T* touch, const char* qname, int mode=MATCH_ALL );
0040
0041     template<typename T>
0042     static int ImmediateReplicaNumber(const T* touch );
0043
0044     template<typename T>
0045     static int AncestorReplicaNumber(const T* touch, int d=1 );
0046
0047
0048     template<typename T>
0049     static int ReplicaNumber(const T* touch, const char* replica_name_select ) ;
0050
0051     template<typename T>
0052     static int ReplicaDepth(const T* touch, const char* replica_name_select ) ;
0053
0054     template<typename T>
0055     static int TouchDepth(const T* touch );
0056
0057     static bool HasMoreThanOneDaughterWithName( const G4LogicalVolume* lv, const char* name);
0058
0059     template<typename T>
0060     static std::string Brief(const T* touch );
0061
0062     template<typename T>
0063     static std::string Desc(const T* touch, int so_wid=20, int pv_wid=20);
0064 };
0065
0066
0067 inline bool U4Touchable::Match( const char* s, const char* q, int mode )
0068 {
0069     bool ret = false ;
0070     switch(mode)
0071     {
0072         case MATCH_ALL:    ret = MatchAll(  s, q) ; break ;
0073         case MATCH_START:  ret = MatchStart(s, q) ; break ;
0074         case MATCH_END:    ret = MatchEnd(  s, q) ; break ;
0075     }
0076     return ret ;
0077 }
0078
0079 inline bool U4Touchable::MatchAll( const char* s, const char* q)
0080 {
0081     return s && q && strcmp(s, q) == 0 ;
0082 }
0083 inline bool U4Touchable::MatchStart( const char* s, const char* q)
0084 {
0085     return s && q && strlen(q) <= strlen(s) && strncmp(s, q, strlen(q)) == 0 ;
0086 }
0087 inline bool U4Touchable::MatchEnd( const char* s, const char* q)
0088 {
0089     int pos = strlen(s) - strlen(q) ;
0090     return pos > 0 && strncmp(s + pos, q, strlen(q)) == 0 ;
0091 }
0092
0093 /**
0094 U4Touchable::FindPV
0095 ---------------------
0096
0097 Find a PV by name in the touch stack, this is much quicker than the recursive U4Volume::FindPV
0098
0099 **/
0100
0101 template<typename T>
0102 inline const G4VPhysicalVolume* U4Touchable::FindPV( const T* touch, const char* qname, int mode )
0103 {
0104     int nd = touch->GetHistoryDepth();
0105     int count = 0 ;
0106     const G4VPhysicalVolume* qpv = nullptr ;
0107     for (int d=0 ; d < nd ; ++d )
0108     {
0109         const G4VPhysicalVolume* dpv = touch->GetVolume(d);
0110         const char* dpv_name = dpv->GetName().c_str() ;
0111         if(Match(dpv_name, qname, mode))
0112         {
0113             qpv = dpv ;
0114             count += 1 ;
0115         }
0116     }
0117     return qpv ;
0118 }
0119
0120 /**
0121 U4Touchable::ImmediateReplicaNumber
0122 -------------------------------------
0123
0124 This is used from U4Recorder::UserSteppingAction_Optical
0125 for step points classified as SURFACE_DETECT "SD".
0126
0127 Gets ReplicaNumber of parent or grandparent volume
0128 in the touch history. Using the (sometimes incorrect)
0129 heuristic that ReplicaNumber zero is not valid.
0130 However that would normally result in getting a zero
0131 in anycase.
0132
0133 Actually it would be better to arrange the ReplicaNumbers
0134 of singleton volumes that are sensitive or containers of sensitive
0135 volumes to adopt ReplicaNumber -1 in order to distinguish from
0136 valid zeros.
0137
0138 Calling this on the touchable from a G4Track that has just stepped
0139 onto a sensitive volume, eg within ProcessHits, would be
0140 expected to obtain the ReplicaNumber (aka pmtID) for several
0141 common ways to organize PMT geometry hierarchy.
0142
0143 **/
0144
0145 template<typename T>
0146 inline int U4Touchable::ImmediateReplicaNumber(const T* touch )
0147 {
0148     int copyNo = touch->GetReplicaNumber(1);
0149     if(copyNo <= 0) copyNo = touch->GetReplicaNumber(2);
0150     return copyNo ;
0151 }
0152
0153 /**
0154 U4Touchable::AncestorReplicaNumber
0155 -----------------------------------
0156
0157 Loops over ancestors starting from d (default d=1 corresponds to parent)
0158 looking for ReplicaNumber > 0 in the history to return.
0159
0160 * NB AGAIN NASTY PRAGMATIC ASSUMPTION THAT copyNo ZERO IS INVALID
0161 * IT WOULD BE AN ADVANTAGE FOR DEFAULT COPYNUMBER TO BE -1 (NOT 0)
0162
0163 **/
0164
0165 template<typename T>
0166 inline int U4Touchable::AncestorReplicaNumber(const T* touch, int d )
0167 {
0168     int depth = touch->GetHistoryDepth();
0169     int copyNo = -1 ;
0170     while( copyNo <= 0 && d < depth )
0171     {
0172         copyNo = touch->GetReplicaNumber(d);
0173         d++ ;
0174     }
0175     return copyNo ;
0176 }
0177
0178
0179 template<typename T>
0180 inline int U4Touchable::ReplicaNumber(const T* touch, const char* replica_name_select )  // static
0181 {
0182     int d = ReplicaDepth(touch, replica_name_select);
0183     bool found = d > -1 ;
0184     int repno = found ? touch->GetReplicaNumber(d) : d  ;
0185
0186 #ifdef U4TOUCHABLE_DEBUG
0187     if(found) std::cerr
0188         << "U4Touchable::ReplicaNumber"
0189         << " found " << found
0190         << " repno " << repno
0191         << std::endl
0192         ;
0193 #endif
0194
0195     return repno ;
0196 }
0197
0198 /**
0199 U4Touchable::ReplicaDepth
0200 ---------------------------
0201
0202 Starting from touch depth look outwards at (volume, mother_volume)
0203 pairs checking for a depth at which the mother_volume has more than one
0204 daughter with the same name as the volume.
0205
0206 This means the volume has
0207 at least one same named sibling making this the replica depth.
0208 When no such depth is found return -1.
0209
0210 For non-null replica_name_select the name search is restricted to names
0211 of logical volumes that contain the replica_name_select string.
0212 Depending on the the names of the replica logical volumes
0213 a suitable "replica_name_select" string (eg "PMT")
0214 may dramatically speedup the search as pointless searches
0215 over thousands of volumes are avoided. Of course this
0216 depends on suitable naming of replica logical volumes.
0217
0218 **/
0219
0220 template<typename T>
0221 inline int U4Touchable::ReplicaDepth(const T* touch, const char* replica_name_select )   // static
0222 {
0223     int nd = touch->GetHistoryDepth();
0224     int t = TouchDepth(touch);
0225     bool expected = t > -1 && t < nd ;
0226
0227     /*
0228     if(!expected) std::cerr
0229         << "U4Touchable::ReplicaDepth"
0230         << " UNEXPECTED "
0231         << " t " << t
0232         << " nd " << nd
0233         << std::endl
0234         ;
0235
0236     assert(expected);
0237     */
0238
0239     if(!expected) return -2 ;
0240
0241     for (int d=t ; d < nd-1; ++d )
0242     {
0243         const G4VPhysicalVolume* dpv = touch->GetVolume(d);
0244         const G4VPhysicalVolume* mpv = touch->GetVolume(d+1);
0245
0246         const G4LogicalVolume* dlv = dpv->GetLogicalVolume();
0247         const G4LogicalVolume* mlv = mpv->GetLogicalVolume();
0248
0249         //const G4VSensitiveDetector* dsd = dlv->GetSensitiveDetector();
0250         //const G4VSensitiveDetector* msd = mlv->GetSensitiveDetector();
0251         //bool sd_skip = dsd == nullptr && msd == nullptr ;
0252
0253         bool hierarchy = dpv->GetMotherLogical() == mlv ;
0254         assert(hierarchy);
0255         if(!hierarchy) std::raise(SIGINT);
0256
0257         const char* dlv_name = dlv->GetName().c_str() ;
0258         bool name_skip = replica_name_select && strstr(dlv_name, replica_name_select) == nullptr ;
0259         //bool skip = name_skip || sd_skip ;
0260
0261 #ifdef U4TOUCHABLE_DEBUG
0262         std::cerr
0263             << "U4Touchable::ReplicaDepth"
0264             << " d " << d
0265             << " nd " << nd
0266             << " dlv_name " << dlv_name
0267             << " replica_name_select " << ( replica_name_select ? replica_name_select : "-" )
0268             << " name_skip " << name_skip
0269             << " skip " << skip
0270             << std::endl
0271             ;
0272 #endif
0273
0274         // skip:true when replica_name_select is provided but the string is not found within the dlv name
0275         // HMM: thats a negative way of doing things, positive approach would be more restrictive so faster
0276         if(name_skip) continue ;
0277         bool found = HasMoreThanOneDaughterWithName(mlv, dlv_name) ;
0278
0279 #ifdef U4TOUCHABLE_DEBUG
0280         if(found)std::cerr
0281             << "U4Touchable::ReplicaDepth"
0282             << " d " << d
0283             << " dlv_name " << dlv_name
0284             << " found " << found
0285             << std::endl
0286             ;
0287 #endif
0288
0289         if(found) return d ;
0290     }
0291     return -1 ;
0292 }
0293
0294 /**
0295 U4Touchable::TouchDepth
0296 -------------------------
0297
0298 Depth of touch volume, -1 if not found.
0299
0300 **/
0301
0302 template<typename T>
0303 inline int U4Touchable::TouchDepth(const T* touch ) // static
0304 {
0305     const G4VPhysicalVolume* tpv = touch->GetVolume() ;
0306     int t = -1 ;
0307     for(int i=0 ; i < touch->GetHistoryDepth() ; i++)
0308     {
0309         const G4VPhysicalVolume* ipv = touch->GetVolume(i) ;
0310         if(ipv == tpv)
0311         {
0312             t = i ;
0313             break ;
0314         }
0315     }
0316     return t ;
0317 }
0318
0319
0320 /**
0321 U4Touchable::HasMoreThanOneDaughterWithName
0322 ---------------------------------------------
0323
0324 When called with name lAcrylic which has num_dau 46276
0325 this is real expensive.
0326
0327 **/
0328
0329 inline bool U4Touchable::HasMoreThanOneDaughterWithName( const G4LogicalVolume* lv, const char* name)  // static
0330 {
0331     int num_dau = lv->GetNoDaughters();
0332     if(num_dau <= 1) return false ;
0333
0334 #ifdef U4TOUCHABLE_DEBUG
0335     bool heavy = num_dau > 45000 ;
0336     if(heavy) std::cerr
0337         << "U4Touchable::HasMoreThanOneDaughterWithName"
0338         << " num_dau " << num_dau
0339         << " name " << name
0340         << " lv.name " << lv->GetName()
0341         << std::endl
0342         ;
0343 #endif
0344
0345
0346     int count = 0;
0347     for (int k=0; k < num_dau ; ++k)
0348     {
0349         const G4VPhysicalVolume* kpv = lv->GetDaughter(k) ;
0350         const G4LogicalVolume*   klv = kpv->GetLogicalVolume() ;
0351         const char* klv_name = klv->GetName().c_str() ;
0352
0353 #ifdef U4TOUCHABLE_DEBUG
0354         if(heavy) std::cerr
0355            << "U4Touchable::HasMoreThanOneDaughterWithName"
0356            << " k " << k
0357            << " kpv.name " << kpv->GetName()
0358            << " klv_name " << klv_name
0359            << " count " << count
0360            << std::endl
0361            ;
0362 #endif
0363
0364         if(strcmp(name, klv_name)==0) count += 1 ;
0365         if(count > 1) return true ;
0366     }
0367     return false ;
0368 }
0369
0370
0371
0372
0373 template<typename T>
0374 inline std::string U4Touchable::Brief(const T* touch )
0375 {
0376     std::stringstream ss ;
0377     ss << "U4Touchable::Brief"
0378        << " HistoryDepth " << std::setw(2) <<  touch->GetHistoryDepth()
0379        << " TouchDepth " << std::setw(2) << TouchDepth(touch)
0380        << " ReplicaDepth " << std::setw(2) << ReplicaDepth(touch, nullptr)
0381        << " ReplicaNumber " << std::setw(6) << ReplicaNumber(touch, nullptr)
0382        ;
0383     return ss.str();
0384 }
0385 template<typename T>
0386 inline std::string U4Touchable::Desc(const T* touch, int so_wid, int pv_wid )
0387 {
0388     int history_depth = touch->GetHistoryDepth();
0389     int touch_depth = TouchDepth(touch);
0390     int replica_depth = ReplicaDepth(touch, nullptr);
0391     int replica_number = ReplicaNumber(touch, nullptr);
0392     int immediate_replica_number = ImmediateReplicaNumber(touch);
0393     int ancestor_replica_number = AncestorReplicaNumber(touch);
0394
0395     std::stringstream ss ;
0396     ss << "U4Touchable::Desc"
0397        << " HistoryDepth " << std::setw(2) << history_depth
0398        << " TouchDepth " << std::setw(2) << touch_depth
0399        << " ReplicaDepth " << std::setw(2) << replica_depth
0400        << " ReplicaNumber " << std::setw(6) << replica_number
0401        << " ImmediateReplicaNumber " << std::setw(6) << immediate_replica_number
0402        << " AncestorReplicaNumber " << std::setw(6) << ancestor_replica_number
0403        << std::endl
0404        ;
0405
0406     for(int i=0 ; i< history_depth ; i++)
0407     {
0408         G4VPhysicalVolume* pv = touch->GetVolume(i);
0409         G4LogicalVolume* lv = pv->GetLogicalVolume();
0410         G4int nd = lv->GetNoDaughters();
0411         G4VSolid* so = touch->GetSolid(i);
0412         G4int cp = touch->GetReplicaNumber(i);
0413
0414         ss << " i " << std::setw(2) << i
0415            << " cp " << std::setw(6)  << cp
0416            << " nd " << std::setw(6)  << nd
0417            << " so " << std::setw(so_wid) << so->GetName()
0418            << " pv " << std::setw(pv_wid) << pv->GetName()
0419            << std::endl ;
0420     }
0421     std::string s = ss.str();
0422     return s ;
0423 }
0424
0425