eic-opticks/sysrap/srec.h

0001 #pragma once
0002 /**
0003 srec.h : highly domain compressed photon step records used for debugging only
0004 ==============================================================================
0005
0006 Domain compression means that must carry around domain metadata in order to
0007 encode or decode the arrays.
0008
0009 Principal user of srec.h is qevent::add_rec
0010
0011 NB seqhis seqmat histories are defined at photon level, so it does not
0012 make sense to include them here at step-record level
0013
0014 For persisting srec arrays use::
0015
0016    NP* rec = NP<short>::Make(num_rec, max_rec, 4, 2)
0017
0018
0019
0020 **/
0021
0022 #include "scuda.h"
0023
0024 #if defined(__CUDACC__) || defined(__CUDABE__)
0025 #    define SREC_METHOD __device__ __forceinline__
0026 #else
0027 #    define SREC_METHOD inline
0028 #endif
0029
0030
0031 #if defined(__CUDACC__) || defined(__CUDABE__)
0032 #else
0033 #include <string>
0034 #include <iostream>
0035 #include <sstream>
0036 #endif
0037
0038
0039 struct srec
0040 {
0041     short4 post ; // position and time into 4*int16 = 64 bits
0042     char4  polw ; // tightly packed, polarization and wavelength into 4*int8 = 32 bits
0043     uchar4 flag ; // 4*int8 = 32 bits
0044
0045     SREC_METHOD void zero();
0046
0047     SREC_METHOD void set_position(const float3& pos, const float4& ce );
0048     SREC_METHOD void get_position(      float3& pos, const float4& ce ) const ;
0049
0050     SREC_METHOD void set_time( const float  t, const float2& td );
0051     SREC_METHOD void get_time(       float& t, const float2& td ) const ;
0052
0053     SREC_METHOD void set_polarization( const float3& pol );
0054     SREC_METHOD void get_polarization(       float3& pol ) const ;
0055
0056     SREC_METHOD void set_wavelength( const float  w, const float2& wd );
0057     SREC_METHOD void get_wavelength(       float& w, const float2& wd ) const ;
0058
0059
0060 #if defined(__CUDACC__) || defined(__CUDABE__)
0061 #else
0062     SREC_METHOD std::string desc() const ;
0063 #endif
0064
0065 };
0066
0067
0068 SREC_METHOD void srec::zero()
0069 {
0070     post.x = 0 ;
0071     post.y = 0 ;
0072     post.z = 0 ;
0073     post.w = 0 ;
0074
0075     polw.x = 0 ;
0076     polw.y = 0 ;
0077     polw.z = 0 ;
0078     polw.w = 0 ;
0079
0080     flag.x = 0 ;
0081     flag.y = 0 ;
0082     flag.z = 0 ;
0083     flag.w = 0 ;
0084 }
0085
0086
0087 #if defined(__CUDACC__) || defined(__CUDABE__)
0088 #else
0089
0090 SREC_METHOD std::string srec::desc() const
0091 {
0092     std::stringstream ss ;
0093     ss
0094        << " post " << post
0095        << " polw " << polw
0096        << " flag " << flag
0097        ;
0098
0099     std::string s = ss.str();
0100     return s ;
0101 }
0102 #endif
0103
0104
0105
0106 #if defined(__CUDACC__) || defined(__CUDABE__)
0107 #define FLOAT2INT_RN(x) (__float2int_rn(x))
0108 #else
0109 #define FLOAT2INT_RN(x) (lrint(x))
0110 #endif
0111
0112
0113 /**
0114 srec::set_position
0115 -------------------
0116
0117 NB positions outside the domain defined by the center-extent *ce* will just cycle.
0118 That is an inherent limitation of the domain compression.  But as compressed records
0119 are just used for debug, that is fine. Just have to make sure the center-extent covers
0120 the region of interest.
0121
0122 **/
0123
0124
0125
0126 SREC_METHOD void srec::set_position( const float3& pos, const float4& ce )
0127 {
0128     post.x = FLOAT2INT_RN( 32767.0f * (pos.x - ce.x)/ce.w ) ;
0129     post.y = FLOAT2INT_RN( 32767.0f * (pos.y - ce.y)/ce.w ) ;
0130     post.z = FLOAT2INT_RN( 32767.0f * (pos.z - ce.z)/ce.w ) ;
0131 }
0132
0133 /**
0134 srec::set_time
0135 ---------------
0136
0137 Time domain is treated as a center at zero with an extent
0138 even though that is wasting half the bits  (no -ve times)
0139 because it simplifies analysis to treat times the same as positions.
0140
0141 **/
0142
0143 SREC_METHOD void srec::set_time( const float  t, const float2& td )
0144 {
0145     post.w = FLOAT2INT_RN(  32767.0f*(t - td.x)/td.y );
0146 }
0147
0148 SREC_METHOD void srec::get_position( float3& pos, const float4& ce ) const
0149 {
0150     pos.x = float(post.x)*ce.w/32767.0f + ce.x ;
0151     pos.y = float(post.y)*ce.w/32767.0f + ce.y ;
0152     pos.z = float(post.z)*ce.w/32767.0f + ce.z ;
0153 }
0154
0155 SREC_METHOD void srec::get_time(       float& t, const float2& td ) const
0156 {
0157     t = float(post.w)*td.y/32767.0f + td.x ;
0158 }
0159
0160
0161 /**
0162 srec::set_polarization
0163 -----------------------
0164
0165 Components of polarization have implicit domain of -1.f to 1.f
0166
0167 **/
0168 SREC_METHOD void srec::set_polarization( const float3& pol )
0169 {
0170     polw.x = FLOAT2INT_RN(pol.x*127.f );
0171     polw.y = FLOAT2INT_RN(pol.y*127.f );
0172     polw.z = FLOAT2INT_RN(pol.z*127.f );
0173 }
0174
0175 SREC_METHOD void srec::get_polarization( float3& pol ) const
0176 {
0177     pol.x = (float(polw.x)/127.f) ;
0178     pol.y = (float(polw.y)/127.f) ;
0179     pol.z = (float(polw.z)/127.f) ;
0180 }
0181
0182 /**
0183 srec::set_wavelength
0184 ----------------------
0185
0186 Wavelength domain was previous specificed by a range,
0187 are now using a symmetrical center-extent form for consistency with other
0188 domains, see qevent.h
0189
0190 **/
0191
0192 SREC_METHOD void srec::set_wavelength( const float w, const float2& wd )
0193 {
0194     polw.w = FLOAT2INT_RN(( w - wd.x)*127.f/ wd.y ) ;
0195 }
0196 SREC_METHOD void srec::get_wavelength(       float& w, const float2& wd ) const
0197 {
0198     w = float(polw.w)*wd.y/127.f + wd.x ;
0199 }
0200
0201