eic-opticks/sysrap/SIntersect.h

0001 #pragma once
0002
0003 #include <vector>
0004 #include <string>
0005
0006 #include "plog/Severity.h"
0007 #include "G4ThreeVector.hh"
0008 #include "geomdefs.hh"
0009
0010 struct float4 ;
0011 struct SCenterExtentGenstep ;
0012 class G4VSolid ;
0013 class G4MultiUnion ;
0014
0015 struct SIntersect
0016 {
0017     static constexpr const  plog::Severity LEVEL = info ;
0018     static constexpr const bool VERBOSE = true ;
0019
0020     static void Scan(const G4VSolid* solid, const char* name, const char* basedir );
0021
0022     SIntersect( const G4VSolid* solid_ );
0023     void init();
0024     const char* desc() const ;
0025
0026     void scan_();
0027     void scan();
0028
0029     const G4VSolid* solid ;
0030     float4*         ce ;
0031     SCenterExtentGenstep* cegs ;
0032 };
0033
0034 #include <iostream>
0035 #include <iomanip>
0036
0037 #include "G4ThreeVector.hh"
0038 #include "G4VSolid.hh"
0039 #include "G4MultiUnion.hh"
0040
0041 #include "scuda.h"
0042 #include "squad.h"
0043
0044 #include "SSys.hh"
0045 #include "SDirect.hh"   // cout_redirect cerr_redirect
0046 #include "SPath.hh"
0047 #include "SCenterExtentGenstep.hh"
0048
0049
0050 #include "sgeomdefs.h"
0051 #include "ssolid.h"
0052
0053 #include "SLOG.hh"
0054
0055
0056 /**
0057 SIntersect::Scan
0058 -------------------
0059
0060 Used from tests/SIntersectSolidTest.cc which is used by xxs.sh
0061
0062 **/
0063
0064 inline void SIntersect::Scan(const G4VSolid* solid, const char* name, const char* basedir )  // static
0065 {
0066     assert( solid && "SIntersect::Scan requires solid");
0067
0068     SIntersect* si = new SIntersect(solid);
0069     si->scan();
0070
0071     const std::string& solidname = solid->GetName() ;
0072
0073     const char* outdir = SPath::Resolve(basedir, name, "SIntersect", DIRPATH );
0074
0075     LOG(LEVEL)
0076         << "si.desc " << si->desc()
0077         << " solidname " << solidname.c_str()
0078         << " name " << name
0079         << " outdir " << outdir
0080         ;
0081
0082     SCenterExtentGenstep* cegs = si->cegs ;
0083     cegs->set_meta<std::string>("name", name);
0084     cegs->set_meta<int>("iidx", 0 );
0085     cegs->save(outdir);
0086 }
0087
0088
0089 inline SIntersect::SIntersect( const G4VSolid* solid_  )
0090     :
0091     solid(solid_),
0092     ce(nullptr),   // float4*
0093     cegs(nullptr)  // SCenterExtentGenstep*
0094 {
0095     init();
0096 }
0097
0098 /**
0099 SIntersect::init
0100 ------------------
0101
0102 Uses G4VSolid::BoundingLimits to determine *ce* center and extent,
0103 and uses *ce* to create *cegs* SCenterExtentGenstep instance.
0104
0105 **/
0106
0107 inline void SIntersect::init()
0108 {
0109     ce = new float4 ;
0110     ssolid::GetCenterExtent( *ce, solid );
0111     cegs = new SCenterExtentGenstep(ce) ;
0112 }
0113
0114 inline const char* SIntersect::desc() const
0115 {
0116     return cegs->desc() ;
0117 }
0118
0119
0120 /**
0121 SIntersect::scan_
0122 -------------------
0123
0124 Using the *cegs.pp* vector of "photon" positions and directions
0125 calulate distances to the solid.  Collect intersections
0126 into *cegs.ii* vector.
0127
0128 TODO:
0129
0130 * adopt simtrace layout, even although some aspects like surface normal
0131   will  be missing from it. This means cegs->pp and cegs->ii will kinda merge
0132
0133 **/
0134
0135 inline void SIntersect::scan_()
0136 {
0137     const std::vector<quad4>& pp = cegs->pp ;
0138     std::vector<quad4>& ii = cegs->ii ;
0139
0140     bool dump = false ;
0141     for(unsigned i=0 ; i < pp.size() ; i++)
0142     {
0143         const quad4& p = pp[i];
0144
0145         G4ThreeVector pos(p.q0.f.x, p.q0.f.y, p.q0.f.z);
0146         G4ThreeVector dir(p.q1.f.x, p.q1.f.y, p.q1.f.z);
0147
0148         G4double t = ssolid::Distance( solid, pos, dir, dump );
0149
0150         if( t == kInfinity ) continue ;
0151         G4ThreeVector ipos = pos + dir*t ;
0152
0153         quad4 isect ;
0154         isect.zero();
0155
0156         isect.q0.f.x = float(ipos.x()) ;
0157         isect.q0.f.y = float(ipos.y()) ;
0158         isect.q0.f.z = float(ipos.z()) ;
0159         isect.q0.f.w = float(t) ;
0160         // TODO: normals, flags, ...
0161
0162         ii.push_back(isect);
0163     }
0164 }
0165
0166
0167 inline void SIntersect::scan()
0168 {
0169     std::stringstream coutbuf;
0170     std::stringstream cerrbuf;
0171     {
0172        cout_redirect out(coutbuf.rdbuf());
0173        cerr_redirect err(cerrbuf.rdbuf());
0174
0175        scan_();
0176     }
0177     std::string cout_ = coutbuf.str() ;
0178     std::string cerr_ = cerrbuf.str() ;
0179
0180
0181     LOG(LEVEL)
0182         << "scan"
0183         << " cout " << strlen(cout_.c_str())
0184         << " cerr " << strlen(cerr_.c_str())
0185         ;
0186
0187     if(VERBOSE)
0188     {
0189         bool with_cout = cout_.size() > 0 ;
0190         bool with_cerr = cerr_.size() > 0 ;
0191         LOG_IF(LEVEL, with_cout) << "cout from scan " << std::endl << cout_ ;
0192         LOG_IF(LEVEL, with_cerr) << "cerr from scan " << std::endl << cerr_ ;
0193     }
0194 }
0195