eic-opticks/qudarap/QPMT.hh

0001 #pragma once
0002 /**
0003 QPMT.hh : projecting PMT properties onto device using qpmt.h
0004 ==============================================================
0005
0006 * narrowing (or widening, or copying) inputs to template type done in ctor
0007
0008 * QPMT.hh does not directly depend on the JUNO specific JPMT.h, instead
0009   the QPMT ctor accepts the arrays that JPMT.h collects from property files
0010
0011   * hence JPMT.h (via PMTSim) only used in tests qudarap/tests/QPMTTest.cc not in qudarap
0012
0013 * QPMT.hh/qpmt.h layout has some similarities to QCerenkov.hh/qcerenkov.h
0014
0015 **/
0016
0017 #include "plog/Severity.h"
0018 #include "NPFold.h"
0019
0020 #include "sproc.h"
0021 #include "qpmt_enum.h"
0022 #include "qpmt.h"
0023 #include "s_pmt.h"
0024 #include "QProp.hh"
0025
0026 #if defined(MOCK_CURAND) || defined(MOCK_CUDA)
0027 #else
0028 #include "QU.hh"
0029 #endif
0030
0031 #include "QUDARAP_API_EXPORT.hh"
0032
0033 template<typename T>
0034 struct QUDARAP_API QPMT
0035 {
0036     static const plog::Severity LEVEL ;
0037     static const QPMT<T>*    INSTANCE ;
0038     static const QPMT<T>*    Get();
0039
0040     static std::string Desc();
0041
0042     const char* ExecutableName ;
0043
0044     const NP* src_rindex ;    // (NUM_PMTCAT, NUM_LAYER, NUM_PROP, NEN, 2:[energy,value] )
0045     const NP* src_thickness ; // (NUM_PMTCAT, NUM_LAYER, 1:value )
0046     const NP* src_qeshape ;   // (NUM_PMTCAT, NUM_SAMPLES~44, 2:[energy,value] )
0047     const NP* src_cetheta ;   // (NUM_PMTCAT, NUM_SAMPLES~9, 2:[theta,value] )
0048     const NP* src_cecosth ;   // (NUM_PMTCAT, NUM_SAMPLES~9, 2:[costh,value] )
0049     const NP* src_lcqs ;      // (NUM_LPMT, 2:[cat,qescale])
0050
0051     const NP* rindex3 ;       // (NUM_PMTCAT*NUM_LAYER*NUM_PROP,  NEN, 2:[energy,value] )
0052     const NP* rindex ;
0053     const QProp<T>* rindex_prop ;
0054
0055     const NP* qeshape ;
0056     const QProp<T>* qeshape_prop ;
0057
0058     const NP* cetheta ;
0059     const QProp<T>* cetheta_prop ;
0060
0061     const NP* cecosth ;
0062     const QProp<T>* cecosth_prop ;
0063
0064
0065     const NP* thickness ;
0066     const NP* lcqs ;
0067     const int* i_lcqs ;  // CPU side lpmtid -> lpmtcat 0/1/2
0068
0069
0070     const NP* src_s_qeshape ;
0071     const NP* src_s_qescale ;
0072
0073     const NP* s_qeshape ;
0074     const QProp<T>* s_qeshape_prop ;
0075
0076     const NP* s_qescale ;
0077
0078
0079     qpmt<T>* pmt ;
0080     qpmt<T>* d_pmt ;
0081
0082     // .h
0083     QPMT(const NPFold* pf);
0084
0085     // .cc
0086     void init();
0087     void init_prop();
0088     void init_thickness();
0089     void init_lcqs();
0090     void init_s_qescale();
0091
0092     // .h
0093     NPFold* serialize() const ;  // formerly get_fold
0094     std::string desc() const ;
0095
0096     // .h : CPU side lookups
0097
0098
0099     int  get_lpmtcat_from_lpmtidx( int lpmtidx ) const ;
0100     int  get_lpmtcat_from_lpmtid(  int lpmtid  ) const ;
0101     int  get_lpmtcat_from_lpmtid(  int* lpmtcat, const int* lpmtid , int num ) const ;
0102
0103     int  get_lpmtidx_from_lpmtid(  int* lpmtidx, const int* lpmtid , int num ) const ;
0104     int  get_spmtidx_from_spmtid(  int* spmtidx, const int* spmtid , int num ) const ;
0105
0106     static NP* MakeArray_pmtcat(int etype, unsigned num_domain );
0107     static NP* MakeArray_pmtid( int etype, unsigned num_domain, unsigned num_pmtid );
0108
0109     // .cc
0110     void pmtcat_check_domain_lookup_shape( int etype, const NP* domain, const NP* lookup) const ;
0111
0112     static const T* Upload(const NP* arr, const char* label);
0113     static T* Alloc(NP* out, const char* label);
0114
0115     NP*  pmtcat_scan(     int etype, const NP* domain) const ;
0116     NP*  mct_lpmtid_scan(  int etype, const NP* domain, const NP* lpmtid) const ;
0117
0118     NP*  spmtid_scan(int etype, const NP* spmtid ) const ;
0119
0120
0121 };
0122
0123
0124 /**
0125 QPMT::QPMT
0126 ------------
0127
0128 1. copy rindex_ into 3D in rindex3 then narrows rindex3 into rindex,
0129    NB this order preserves last prop column integer annotations
0130 2. creates rindex_prop from rindex
0131 3. narrows src_qeshape into qeshape
0132 4. creates qeshape_prop from qeshape
0133 5. creates cetheta_prop from cetheta
0134 5. narrows src_thickness into thickness
0135 6. narrows src_lcqs into lcqs
0136
0137 NB the jpmt argument is the NPFold provided by SPMT::CreateFromJPMTAndSerialize
0138 not the raw fold from _PMTSimParamData. So all the data preparation done
0139 in SPMT.h is accessible from here.
0140
0141 **/
0142
0143 template<typename T>
0144 inline QPMT<T>::QPMT(const NPFold* jpmt )
0145     :
0146     ExecutableName(sproc::ExecutableName()),
0147     src_rindex(   jpmt->get("rindex")),
0148     src_thickness(jpmt->get("thickness")),
0149     src_qeshape(  jpmt->get("qeshape")),
0150     src_cetheta(  jpmt->get("cetheta")),
0151     src_cecosth(  jpmt->get("cecosth")),
0152     src_lcqs(     jpmt->get("lcqs")),
0153     rindex3(  NP::MakeCopy3D(src_rindex)),   // make copy and change shape to 3D
0154     rindex(   NP::MakeWithType<T>(rindex3)), // adopt template type, potentially narrowing
0155     rindex_prop(new QProp<T>(rindex)),
0156     qeshape(   NP::MakeWithType<T>(src_qeshape)), // adopt template type, potentially narrowing
0157     qeshape_prop(new QProp<T>(qeshape)),
0158     cetheta(   NP::MakeWithType<T>(src_cetheta)), // adopt template type, potentially narrowing
0159     cetheta_prop(new QProp<T>(cetheta)),
0160     cecosth(   NP::MakeWithType<T>(src_cecosth)), // adopt template type, potentially narrowing
0161     cecosth_prop(new QProp<T>(cecosth)),
0162     thickness(NP::MakeWithType<T>(src_thickness)),
0163     lcqs(src_lcqs ? NP::MakeWithType<T>(src_lcqs) : nullptr),
0164     i_lcqs( lcqs ? (int*)lcqs->cvalues<T>() : nullptr ),    // CPU side lookup lpmtidx->lpmtcat 0/1/2
0165     src_s_qeshape(  jpmt->get("s_qeshape")),
0166     src_s_qescale(  jpmt->get("s_qescale")),
0167     s_qeshape(   NP::MakeWithType<T>(src_s_qeshape)), // adopt template type, potentially narrowing
0168     s_qeshape_prop(new QProp<T>(s_qeshape)),
0169     s_qescale(src_s_qescale ? NP::MakeWithType<T>(src_s_qescale) : nullptr),
0170     pmt(new qpmt<T>()),                    // host-side qpmt.h instance
0171     d_pmt(nullptr)                         // device-side pointer set at upload in init
0172 {
0173     init();
0174 }
0175
0176 // init in .cc
0177 template<typename T>
0178 inline NPFold* QPMT<T>::serialize() const  // formerly get_fold
0179 {
0180     NPFold* fold = new NPFold ;
0181
0182     fold->add("src_thickness", src_thickness );
0183     fold->add("src_rindex", src_rindex );
0184     fold->add("src_qeshape", src_qeshape );
0185     fold->add("src_cetheta", src_cetheta );
0186     fold->add("src_cecosth", src_cecosth );
0187     fold->add("src_lcqs", src_lcqs );
0188
0189     fold->add("thickness", thickness );
0190     fold->add("rindex", rindex );
0191     fold->add("qeshape", qeshape );
0192     fold->add("cetheta", cetheta );
0193     fold->add("cecosth", cecosth );
0194     fold->add("lcqs", lcqs );
0195
0196     fold->add("rindex_prop_a", rindex_prop->a );
0197     fold->add("qeshape_prop_a", qeshape_prop->a );
0198     fold->add("cetheta_prop_a", cetheta_prop->a );
0199     fold->add("cecosth_prop_a", cecosth_prop->a );
0200
0201     fold->add("s_qeshape", s_qeshape );
0202     fold->add("s_qeshape_prop_a", s_qeshape_prop->a );
0203     fold->add("s_qescale", s_qescale );
0204
0205     return fold ;
0206 }
0207
0208 template<typename T>
0209 inline std::string QPMT<T>::desc() const
0210 {
0211     int w = 30 ;
0212     std::stringstream ss ;
0213     ss
0214        << "QPMT::desc"
0215        << std::endl
0216        << std::setw(w) << "rindex "    << rindex->sstr() << std::endl
0217        << std::setw(w) << "qeshape " << qeshape->sstr() << std::endl
0218        << std::setw(w) << "cetheta " << cetheta->sstr() << std::endl
0219        << std::setw(w) << "cecosth " << cecosth->sstr() << std::endl
0220        << std::setw(w) << "thickness " << thickness->sstr() << std::endl
0221        << std::setw(w) << "lcqs " << lcqs->sstr() << std::endl
0222        << std::setw(w) << "s_qeshape " << s_qeshape->sstr() << std::endl
0223        << std::setw(w) << "s_qescale " << s_qescale->sstr() << std::endl
0224        << std::setw(w) << " pmt.rindex_prop " << pmt->rindex_prop  << std::endl
0225        << std::setw(w) << " pmt.qeshape_prop " << pmt->qeshape_prop  << std::endl
0226        << std::setw(w) << " pmt.cetheta_prop " << pmt->cetheta_prop  << std::endl
0227        << std::setw(w) << " pmt.cecosth_prop " << pmt->cecosth_prop  << std::endl
0228        << std::setw(w) << " pmt.thickness " << pmt->thickness  << std::endl
0229        << std::setw(w) << " pmt.lcqs " << pmt->lcqs  << std::endl
0230        << std::setw(w) << " d_pmt " << d_pmt   << std::endl
0231        ;
0232     std::string s = ss.str();
0233     return s ;
0234 }
0235
0236 /**
0237 QPMT::get_lpmtcat_from_lpmtidx
0238 --------------------------------
0239
0240 CPU side lookup of lpmtcat from lpmtidx using i_lcqs array.
0241
0242 **/
0243
0244 template<typename T>
0245 inline int QPMT<T>::get_lpmtcat_from_lpmtidx( int lpmtidx ) const
0246 {
0247     assert( lpmtidx >= 0 && lpmtidx < s_pmt::NUM_LPMTIDX );  // extended from NUM_CD_LPMT_AND_WP
0248     const int& lpmtcat = i_lcqs[lpmtidx*2+0] ;
0249     return lpmtcat ;
0250 }
0251
0252 template<typename T>
0253 inline int QPMT<T>::get_lpmtcat_from_lpmtid( int lpmtid ) const
0254 {
0255     int lpmtidx = s_pmt::lpmtidx_from_pmtid(lpmtid);
0256     return get_lpmtcat_from_lpmtidx(lpmtidx);
0257 }
0258
0259
0260 template<typename T>
0261 inline int QPMT<T>::get_lpmtcat_from_lpmtid( int* lpmtcat_, const int* lpmtid_, int num_lpmtid ) const
0262 {
0263     for(int i=0 ; i < num_lpmtid ; i++)
0264     {
0265         int lpmtid = lpmtid_[i] ;
0266         int lpmtcat = get_lpmtcat_from_lpmtid(lpmtid) ;
0267         lpmtcat_[i] = lpmtcat ;
0268     }
0269     return num_lpmtid ;
0270 }
0271
0272 template<typename T>
0273 inline int QPMT<T>::get_lpmtidx_from_lpmtid( int* lpmtidx_, const int* lpmtid_, int num_lpmtid ) const
0274 {
0275     for(int i=0 ; i < num_lpmtid ; i++)
0276     {
0277         int lpmtid = lpmtid_[i] ;
0278         int lpmtidx = s_pmt::lpmtidx_from_pmtid(lpmtid);
0279         lpmtidx_[i] = lpmtidx ;
0280     }
0281     return num_lpmtid ;
0282 }
0283
0284 template<typename T>
0285 inline int QPMT<T>::get_spmtidx_from_spmtid( int* spmtidx_, const int* spmtid_, int num_spmtid ) const
0286 {
0287     for(int i=0 ; i < num_spmtid ; i++)
0288     {
0289         int spmtid = spmtid_[i] ;
0290         int spmtidx = s_pmt::spmtidx_from_pmtid(spmtid);
0291         spmtidx_[i] = spmtidx ;
0292     }
0293     return num_spmtid ;
0294 }
0295
0296
0297
0298
0299
0300
0301 /**
0302 QPMT::MakeArray_pmtcat
0303 -------------------------
0304
0305 HMM: this is mainly for testing, perhaps put in QPMTTest ?
0306
0307 **/
0308
0309 template<typename T>
0310 inline NP* QPMT<T>::MakeArray_pmtcat(int etype, unsigned num_domain )   // static
0311 {
0312     const int& ni = s_pmt::NUM_CAT ;
0313     const int& nj = s_pmt::NUM_LAYR ;
0314     const int& nk = s_pmt::NUM_PROP ;
0315     NP* lookup = nullptr ;
0316     switch(etype)
0317     {
0318        case qpmt_RINDEX:  lookup = NP::Make<T>( ni, nj, nk, num_domain ) ; break ;
0319        case qpmt_CATSPEC: lookup = NP::Make<T>( ni, num_domain, 4, 4  )  ; break ;
0320        case qpmt_QESHAPE: lookup = NP::Make<T>( ni,         num_domain ) ; break ;
0321        case qpmt_CETHETA: lookup = NP::Make<T>( ni,         num_domain ) ; break ;
0322        case qpmt_CECOSTH: lookup = NP::Make<T>( ni,         num_domain ) ; break ;
0323        case qpmt_S_QESHAPE: lookup = NP::Make<T>( 1,        num_domain ) ; break ;
0324     }
0325     return lookup ;
0326 }
0327
0328
0329 template<typename T>
0330 inline NP* QPMT<T>::MakeArray_pmtid(int etype, unsigned num_domain, unsigned num_pmtid )   // static
0331 {
0332     const int ni = num_pmtid ;
0333     const int nj = num_domain ;
0334
0335     NP* lookup = nullptr ;
0336     switch(etype)
0337     {
0338        case qpmt_SPEC:    lookup = NP::Make<T>( ni, nj, 4, 4  )       ; break ;
0339        case qpmt_SPEC_ce: lookup = NP::Make<T>( ni, nj, 4, 4  )       ; break ;
0340        case qpmt_ART:     lookup = NP::Make<T>( ni, nj, 4, 4  )       ; break ;
0341        case qpmt_COMP:    lookup = NP::Make<T>( ni, nj, 1, 4, 4, 2 )  ; break ;
0342        case qpmt_LL:      lookup = NP::Make<T>( ni, nj, 4, 4, 4, 2 )  ; break ;
0343        case qpmt_ARTE:    lookup = NP::Make<T>( ni, nj, 4  )          ; break ;
0344        case qpmt_ATQC:    lookup = NP::Make<T>( ni, nj, 4  )          ; break ;
0345        case qpmt_S_QESCALE:  lookup = NP::Make<T>( num_pmtid )        ; break ;
0346     }
0347     return lookup ;
0348 }
0349
0350