eic-opticks/qudarap/QPMT_MOCK.h

0001 #pragma once
0002 /**
0003 QPMT_MOCK.h
0004 =============
0005
0006 CPU mockup impl following GPU code : QPMT.cu
0007
0008 * HMM : this code is very similar to QPMT.cu
0009
0010 Q: How to rearrange such code for single source ?
0011 A: Presumably use a header with such an impl from the QPMT.cu
0012
0013    * but as this is used just for testing here.. its not so compelling
0014      to be worthy of the effort
0015
0016 Notice the game for maximum reuse
0017
0018 1. minimize un-sharable .cc and .cu, instead put almost everything
0019    into small headers
0020
0021 **/
0022
0023
0024 #include "qpmt.h"
0025
0026 template <typename F>
0027 void QPMT_lpmtcat_MOCK(
0028     qpmt<F>* pmt,
0029     int etype,
0030     F* lookup,
0031     const F* domain,
0032     unsigned domain_width )
0033 {
0034
0035     if( etype == qpmt_RINDEX )
0036     {
0037         for(unsigned ix=0 ; ix < domain_width ; ix++)
0038         {
0039             F energy_eV = domain[ix] ;
0040
0041             const int& ni = qpmt_NUM_CAT ;
0042             const int& nj = qpmt_NUM_LAYR ;
0043             const int& nk = qpmt_NUM_PROP ;
0044
0045             for(int i=0 ; i < ni ; i++)
0046             for(int j=0 ; j < nj ; j++)
0047             for(int k=0 ; k < nk ; k++)
0048             {
0049                 int iprop = i*nj*nk+j*nk+k ;            // linearized higher dimensions
0050                 int index = iprop * domain_width + ix ; // output index into lookup
0051
0052                 F value = pmt->rindex_prop->interpolate(iprop, energy_eV );
0053                 lookup[index] = value ;
0054             }
0055         }
0056     }
0057     else if( etype == qpmt_QESHAPE )
0058     {
0059         for(unsigned ix=0 ; ix < domain_width ; ix++)
0060         {
0061             F energy_eV = domain[ix] ;
0062             const int& ni = qpmt_NUM_CAT ;
0063
0064             for(int i=0 ; i < ni ; i++)
0065             {
0066                 int index = i * domain_width + ix ; // output index into lookup
0067                 F value = pmt->qeshape_prop->interpolate(i, energy_eV );
0068                 lookup[index] = value ;
0069             }
0070         }
0071     }
0072     else if( etype == qpmt_CATSPEC )
0073     {
0074         for(unsigned ix=0 ; ix < domain_width ; ix++)
0075         {
0076             F energy_eV = domain[ix] ;
0077
0078             const int& ni = qpmt_NUM_CAT ;
0079             const int& nj = domain_width ;
0080             const int  nk = 16 ;
0081             const int&  j = ix ;
0082
0083             F ss[nk] ;
0084
0085             for(int i=0 ; i < ni ; i++)  // over pmtcat
0086             {
0087                 int index = i*nj*nk + j*nk  ;
0088                 pmt->get_lpmtcat_stackspec(ss, i, energy_eV );
0089                 for( int k=0 ; k < nk ; k++) lookup[index+k] = ss[k] ;
0090             }
0091         }
0092     }
0093 }
0094
0095
0096
0097 #ifdef WITH_CUSTOM4
0098 // templated payload size P as it needs to be a compile time constant
0099 template <typename F, int P>
0100 void _QPMT_mct_lpmtid_MOCK(
0101     qpmt<F>* pmt,
0102     int etype,
0103     F* lookup ,
0104     const F* domain,
0105     unsigned domain_width,
0106     const int* lpmtid,
0107     unsigned num_lpmtid )
0108 {
0109     for(unsigned ix = 0 ; ix < domain_width ; ix++)
0110     {
0111         F minus_cos_theta = domain[ix] ;
0112         F wavelength_nm = 440.f ;
0113         F dot_pol_cross_mom_nrm = 0.f ; // SPOL zero is pure P polarized
0114
0115         const int& ni = num_lpmtid ;
0116         const int& nj = domain_width ;   // minus_cos_theta values "AOI"
0117         const int&  j = ix ;
0118
0119         F payload[P] ;
0120
0121         for(int i=0 ; i < ni ; i++)  // over num_lpmtid
0122         {
0123             int index = i*nj*P + j*P  ;
0124             int pmtid = lpmtid[i] ;
0125
0126             if( etype == qpmt_SPEC )
0127             {
0128                 pmt->get_lpmtid_SPEC(payload, pmtid, wavelength_nm, minus_cos_theta, dot_pol_cross_mom_nrm );
0129             }
0130             else if( etype == qpmt_LL )
0131             {
0132                 pmt->get_lpmtid_LL(payload, pmtid, wavelength_nm, minus_cos_theta, dot_pol_cross_mom_nrm );
0133             }
0134             else if( etype == qpmt_COMP )
0135             {
0136                 pmt->get_lpmtid_COMP(payload, pmtid, wavelength_nm, minus_cos_theta, dot_pol_cross_mom_nrm );
0137             }
0138             else if( etype == qpmt_ART )
0139             {
0140                 pmt->get_lpmtid_ART(payload, pmtid, wavelength_nm, minus_cos_theta, dot_pol_cross_mom_nrm );
0141             }
0142             else if( etype == qpmt_ARTE )
0143             {
0144                 pmt->get_lpmtid_ARTE(payload, pmtid, wavelength_nm, minus_cos_theta, dot_pol_cross_mom_nrm );
0145             }
0146
0147             for( int k=0 ; k < P ; k++) lookup[index+k] = payload[k] ;
0148         }
0149     }
0150 }
0151
0152
0153
0154 template <typename F>
0155 void QPMT_mct_lpmtid_MOCK(
0156     qpmt<F>* pmt,
0157     int etype,
0158     F* lookup,
0159     const F* domain,
0160     unsigned domain_width,
0161     const int* lpmtid,
0162     unsigned num_lpmtid
0163 )
0164 {
0165     switch(etype)
0166     {
0167         case qpmt_SPEC:
0168            _QPMT_mct_lpmtid_MOCK<F,16>(
0169               pmt, etype, lookup, domain, domain_width, lpmtid, num_lpmtid ) ;  break ;
0170
0171         case qpmt_ART:
0172            _QPMT_mct_lpmtid_MOCK<F,16>(
0173               pmt, etype, lookup, domain, domain_width, lpmtid, num_lpmtid ) ;  break ;
0174
0175         case qpmt_COMP:
0176            _QPMT_mct_lpmtid_MOCK<F,32>(
0177               pmt, etype, lookup, domain, domain_width, lpmtid, num_lpmtid ) ;  break ;
0178
0179         case qpmt_LL:
0180            _QPMT_mct_lpmtid_MOCK<F,128>(
0181               pmt, etype, lookup, domain, domain_width, lpmtid, num_lpmtid ) ;  break ;
0182
0183         case qpmt_ARTE:
0184            _QPMT_mct_lpmtid_MOCK<F,4>(
0185               pmt, etype, lookup, domain, domain_width, lpmtid, num_lpmtid ) ;  break ;
0186     }
0187 }
0188 #endif
0189
0190
0191
0192