Material and Surface Properties¶
Opticks handling of material and surface properties is all about preparing GPU textures : which means they have to be standardized to fit into arrays with specific properties and common wavelength domain.
Materials must have a non-zero REEMISSIONPROB for it to be regarded as a scintillator and collected by GScintillatorLib.
The below guide to the relevant classes is intended to show how material/surface properties get translated into GPU textures and how those are accessed GPU side.
GPU texture preparation¶
ggeo/GScintillatorLib
scintillator lib is populated by GGeo::prepareScintillatorLib (invoked by GGeo::prepare) from raw materials with the 3 properties : SLOWCOMPONENT,FASTCOMPONENT,REEMISSIONPROB
For a material to be regarded as a scintillator it must have all these three properties.
Also GScintillatorLib prepares the reemission inverse CDF array.
- optixrap/OScintillatorLib
- Converts the GScintillatorLib array into the reemission GPU texture. Lookups on the reemission texture with uniform randoms generates wavelengths according to the desired distribution.
ggeo/GMaterialLib,GMaterial
materials get standardized into 8 properties and a common wavelength domain, as they are destined to get interleaved into the boundary GPU texture.
For a material to be a scintillator must have non-zero “reemission_prob:REEMISSIONPROB,”
The 8 material properties (3 are spare).
63 const char* GMaterialLib::keyspec =
64 "refractive_index:RINDEX,"
65 "absorption_length:ABSLENGTH,"
66 "scattering_length:RAYLEIGH,"
67 "reemission_prob:REEMISSIONPROB,"
68 "group_velocity:GROUPVEL,"
69 "extra_y:EXTRA_Y,"
70 "extra_z:EXTRA_Z,"
71 "extra_w:EXTRA_W,"
ggeo/GSurfaceLib,GSurface
surface properties are standardized into 8 properties (4 are spare)
51 const char* GSurfaceLib::detect = "detect" ;
52 const char* GSurfaceLib::absorb = "absorb" ;
53 const char* GSurfaceLib::reflect_specular = "reflect_specular" ;
54 const char* GSurfaceLib::reflect_diffuse = "reflect_diffuse" ;
55
56 const char* GSurfaceLib::extra_x = "extra_x" ;
57 const char* GSurfaceLib::extra_y = "extra_y" ;
58 const char* GSurfaceLib::extra_z = "extra_z" ;
59 const char* GSurfaceLib::extra_w = "extra_w" ;
60
87 const char* GSurfaceLib::keyspec =
88 "detect:EFFICIENCY,"
89 "absorb:DUMMY,"
90 "reflect_specular:REFLECTIVITY,"
91 "reflect_diffuse:REFLECTIVITY,"
92 ;
- ggeo/GBnd
- A boundary is the 4 integers (omat,osur,isur,imat) representing materials and surfaces on either side of a transition been materials.
omat:outer-material
osur:outer-surface
isur:inner-surface
imat:inner-material
- ggeo/GBndLib
collects all unique boundaries GBnd for the entire geometry, boundary indices for all bits of the geometry are stored with the geometry
interleaves the standardized material and surface property arrays from GMaterialLib and GSurfaceLib into the boundary dynamic array
optixrap/OBndLib
converts GBndLib into the GPU boundary texture and uploads its content to GPU
GPU texture usage¶
- optixrap/cu/boundary_lookup.h
- static __device__ __inline__ float4 boundary_lookup( float nm, unsigned int line, unsigned int k)
optixrap/cu/state.h
fill_state does the boundary lookups
All pieces of geometry have a boundary index assigned to it, which means that on intersection with a ray you get the boundary index which via offsets into the boundary texture enable you to access the material and surface properties omat/osur/isur/imat relevant to that intersect.
This approach copies material and surface properties multiple times, but the boundary texture is very small compared to what GPUs are designed to handle so its not a problem and it simplifies property lookup.
23 struct State
24 {
25 unsigned int flag ;
26 float4 material1 ; // refractive_index/absorption_length/scattering_length/reemission_prob
27 float4 m1group2 ; // group_velocity/spare1/spare2/spare3
28 float4 material2 ;
29 float4 surface ; // detect/absorb/reflect_specular/reflect_diffuse
30 float3 surface_normal ;
31 float cos_theta ;
32 float distance_to_boundary ;
33 uint4 optical ; // x/y/z/w index/type/finish/value
34 uint4 index ; // indices of m1/m2/surf/sensor
35 uint4 identity ; // node/mesh/boundary/sensor indices of last intersection
36 float ureflectcheat ;
37 };
38
..
48 __device__ void fill_state( State& s, int boundary, uint4 identity, float wavelength )
49 {
50 // boundary : 1 based code, signed by cos_theta of photon direction to outward geometric normal
51 // >0 outward going photon
52 // <0 inward going photon
53 //
54 // NB the line is above the details of the payload (ie how many float4 per matsur)
55 // it is just
56 // boundaryIndex*4 + 0/1/2/3 for OMAT/OSUR/ISUR/IMAT
57 //
58
59 int line = boundary > 0 ? (boundary - 1)*BOUNDARY_NUM_MATSUR : (-boundary - 1)*BOUNDARY_NUM_MATSUR ;
60
61 // pick relevant lines depening on boundary sign, ie photon direction relative to normal
62 //
63 int m1_line = boundary > 0 ? line + IMAT : line + OMAT ;
64 int m2_line = boundary > 0 ? line + OMAT : line + IMAT ;
65 int su_line = boundary > 0 ? line + ISUR : line + OSUR ;
66
67 // consider photons arriving at PMT cathode surface
68 // geometry normals are expected to be out of the PMT
69 //
70 // boundary sign will be -ve : so line+3 outer-surface is the relevant one
71
72 s.material1 = boundary_lookup( wavelength, m1_line, 0);
73 s.m1group2 = boundary_lookup( wavelength, m1_line, 1);
74
75 s.material2 = boundary_lookup( wavelength, m2_line, 0);
76 s.surface = boundary_lookup( wavelength, su_line, 0);
77
78 s.optical = optical_buffer[su_line] ; // index/type/finish/value
79
optixrap/cu/generate.cu
fill_state is the first thing down after an OptiX ray trace intersection
564 rtTrace(top_object, optix::make_Ray(p.position, p.direction, propagate_ray_type, propagate_epsilon, RT_DEFAULT_MAX), prd );
565
566 if(prd.boundary == 0)
567 {
568 s.flag = MISS ; // overwrite CERENKOV/SCINTILLATION for the no hitters
569 // zero out no-hitters to avoid leftovers
570 s.index.x = 0 ;
571 s.index.y = 0 ;
572 s.index.z = 0 ;
573 s.index.w = 0 ;
574 break ;
575 }
576 // initial and CONTINUE-ing records
577
578 // use boundary index at intersection point to do optical constant + material/surface property lookups
579 fill_state(s, prd.boundary, prd.identity, p.wavelength );
optixrap/cu/wavelength_lookup.h
reemission_lookup is used depending on reemission_prob and a random throw, note use of the state s.material1.w
optixrap/cu/propagate.h
122
123 if (absorption_distance <= scattering_distance)
124 {
125 if (absorption_distance <= s.distance_to_boundary)
126 {
127 p.time += absorption_distance/speed ;
128 p.position += absorption_distance*p.direction;
129
130 const float& reemission_prob = s.material1.w ;
131 float u_reemit = reemission_prob == 0.f ? 2.f : curand_uniform(&rng); // avoid consumption at absorption when not scintillator
132
133 if (u_reemit < reemission_prob)
134 {
135 // no materialIndex input to reemission_lookup as both scintillators share same CDF
136 // non-scintillators have zero reemission_prob
137 p.wavelength = reemission_lookup(curand_uniform(&rng));
138 p.direction = uniform_sphere(&rng);
139 p.polarization = normalize(cross(uniform_sphere(&rng), p.direction));
140 p.flags.i.x = 0 ; // no-boundary-yet for new direction
141
142 s.flag = BULK_REEMIT ;
143 return CONTINUE;
144 }
145 else
146 {
147 s.flag = BULK_ABSORB ;
148 return BREAK;
149 }
150 }
151 // otherwise sail to boundary
152 }
153 else