eic-opticks/u4/U4Engine.h

0001 #pragma once
0002 /**
0003 U4Engine.h
0004 ============
0005
0006 Hmm not convenient to persist state into files when have lots of them.
0007 put/get interface looks better::
0008
0009     108   virtual std::vector<unsigned long> put () const;
0010     109   virtual bool get (const std::vector<unsigned long> & v);
0011
0012 Note the inverted API naming, "get" is restoring and "put" is saving.
0013 For MixMaxRng the vector length is 17*2+4 = 38 and 32 bit unsigned long is OK
0014
0015 **/
0016
0017 #include <csignal>
0018 #include "CLHEP/Random/Randomize.h"
0019 #include "CLHEP/Random/RandomEngine.h"
0020 #include "U4_API_EXPORT.hh"
0021 #include "spath.h"
0022 #include "NP.hh"
0023
0024 struct U4_API U4Engine
0025 {
0026     static std::string Desc();
0027     static void ShowState();
0028     static std::string DescStateArray();
0029     static std::string DescState();
0030
0031     static std::string ConfPath(const char* fold, const char* name=nullptr);
0032     static void SaveState(   const char* fold, const char* name=nullptr) ;
0033     static void RestoreState(const char* fold, const char* name=nullptr) ;
0034
0035     static void SaveState(          NP* states, int idx );
0036     static void RestoreState( const NP* states, int idx );
0037
0038 };
0039
0040 /**
0041 U4Engine::Desc
0042 -----------------
0043
0044 Introspect the engine::
0045
0046     g4-cls Randomize
0047     g4-cls MixMaxRng   ## by inspection the default engine in 10.4.2
0048     g4-cls RandomEngine
0049
0050
0051 NB for breakpointing to work with gdb need the CLHEP it seems::
0052
0053    BP=CLHEP::MixMaxRng::flat
0054
0055 **/
0056
0057 inline std::string U4Engine::Desc() // static
0058 {
0059     CLHEP::HepRandomEngine* engine = CLHEP::HepRandom::getTheEngine() ;
0060     std::stringstream ss ;
0061     ss << "U4Engine::Desc"
0062        << " engine " << ( engine ? "Y" : "N" )
0063        << " engine.name " << ( engine ? engine->name() : "-" )
0064        ;
0065     std::string s = ss.str() ;
0066     return s ;
0067 }
0068 inline void U4Engine::ShowState() // static
0069 {
0070     CLHEP::HepRandomEngine* engine = CLHEP::HepRandom::getTheEngine() ;
0071     engine->showStatus();
0072 }
0073
0074 inline std::string U4Engine::DescStateArray() // static
0075 {
0076     CLHEP::HepRandomEngine* engine = CLHEP::HepRandom::getTheEngine() ;
0077     assert( engine );
0078     std::vector<unsigned long> state = engine->put();
0079
0080     std::stringstream ss ;
0081     ss << "U4Engine::DescStateArray" << std::endl << std::endl  ;
0082     ss << "state = " << NP::ArrayString<unsigned long>( state, 10 ) << std::endl ;
0083
0084     std::string s = ss.str();
0085     return s ;
0086 }
0087
0088 inline std::string U4Engine::DescState() // static
0089 {
0090     CLHEP::HepRandomEngine* engine = CLHEP::HepRandom::getTheEngine() ;
0091     assert( engine );
0092     std::stringstream ss ;
0093     ss << "U4Engine::DescState" << std::endl ;
0094     engine->put(ss);
0095     std::string s = ss.str();
0096     return s ;
0097 }
0098
0099
0100
0101 inline std::string U4Engine::ConfPath(const char* fold, const char* name) // static
0102 {
0103     CLHEP::HepRandomEngine* engine = CLHEP::HepRandom::getTheEngine() ;
0104     std::string defconf = engine->name() + "State.conf" ;
0105     if(name == nullptr) name = defconf.c_str() ;
0106     std::string path = spath::Join(fold, name) ;
0107     return path ;
0108 }
0109
0110 inline void U4Engine::SaveState(const char* fold, const char* name) // status
0111 {
0112     std::string conf = ConfPath(fold, name) ;
0113     std::cout << "U4Engine::SaveStatus " << conf << std::endl ;
0114     CLHEP::HepRandomEngine* engine = CLHEP::HepRandom::getTheEngine() ;
0115     engine->saveStatus(conf.c_str()) ;
0116 }
0117
0118 inline void U4Engine::RestoreState(const char* fold, const char* name) // status
0119 {
0120     std::string conf = ConfPath(fold, name) ;
0121     std::cout << "U4Engine::RestoreStatus " << conf << std::endl ;
0122     CLHEP::HepRandomEngine* engine = CLHEP::HepRandom::getTheEngine() ;
0123     engine->restoreStatus(conf.c_str()) ;
0124 }
0125
0126 /**
0127 U4Engine::SaveState
0128 ----------------------
0129
0130 * currently idx is treated as an index (not an id)
0131   with the requirement that the idx is less than the
0132   number of states items in the array
0133
0134 * currently saving 1000 states takes 320K::
0135
0136     In [2]: 8*38*1000
0137     Out[2]: 304000
0138
0139     640 -rw-r--r--   1 blyth  wheel  304128 Nov 23 22:46 g4state.npy
0140
0141     epsilon:ALL blyth$ du -h g4state.npy
0142     320K    g4state.npy
0143
0144 * NB could half the size as the g4state values would fit into 32 bit
0145   but are stored into 64 bit : need to type shift from
0146   "unsigned long" to "unsigned"
0147
0148 * so that gets expensive when want to to save/restore the g4state
0149   for one in a million photon, as would need 320M bytes just to access 304 bytes
0150
0151 * TODO: sparse (by id) approach
0152   (eg by adding idset/idget accessors to NP using id string vec metadata
0153   to give the mapping between id and storage index)
0154
0155 **/
0156
0157 inline void U4Engine::SaveState( NP* states, int idx ) // static
0158 {
0159     assert( states && states->shape.size() > 0 );
0160     if( idx >= states->shape[0] ) return ;
0161
0162     CLHEP::HepRandomEngine* engine = CLHEP::HepRandom::getTheEngine() ;
0163     std::vector<unsigned long> state = engine->put() ; // NB: back-to-front API name
0164     unsigned niv = states->num_itemvalues() ;
0165     unsigned size = states->item_bytes();
0166
0167     bool niv_expect = niv*sizeof(unsigned long) == size ;
0168     assert( niv_expect );
0169     if(!niv_expect) std::raise(SIGINT);
0170
0171     assert( state.size() == niv );
0172     memcpy( states->bytes() + size*idx,  state.data(), size  );
0173 }
0174
0175 inline void U4Engine::RestoreState( const NP* states, int idx ) // static
0176 {
0177     assert( states && states->shape.size() > 0 );
0178     if( idx >= states->shape[0] ) return ;
0179
0180     unsigned niv = states->num_itemvalues() ;
0181     std::vector<unsigned long> state ;
0182     states->slice(state, idx, -1 );  // -1: terminates the dimensions to select on
0183
0184     bool niv_expect = state.size() == niv ;
0185     assert( niv_expect );
0186     if(!niv_expect) std::raise(SIGINT);
0187
0188     CLHEP::HepRandomEngine* engine = CLHEP::HepRandom::getTheEngine() ;
0189
0190     bool rc = engine->get(state);  // HUH: back-to-front API name ?
0191     assert( rc );
0192     if(!rc) std::raise(SIGINT) ;
0193 }
0194