Back to home page

Project CMSSW displayed by LXR

 
 

    


Warning, /RecoLocalCalo/HGCalRecProducers/plugins/HGCalRecHitKernelImpl.cu is written in an unsupported language. File is not indexed.

0001 #include <cuda.h>
0002 #include <cuda_runtime.h>
0003 #include <inttypes.h>
0004 #include "DataFormats/HGCRecHit/interface/HGCRecHit.h"
0005 #include "DataFormats/ForwardDetId/interface/HGCalDetId.h"
0006 #include "HGCalRecHitKernelImpl.cuh"
0007 
0008 __device__ float get_weight_from_layer(const int32_t& layer,
0009                                        const double (&weights)[HGChefUncalibRecHitConstantData::hef_weights]) {
0010   return (float)weights[layer];
0011 }
0012 
0013 __device__ void make_rechit_silicon(unsigned tid,
0014                                     HGCRecHitSoA dst_soa,
0015                                     HGCUncalibRecHitSoA src_soa,
0016                                     const float& weight,
0017                                     const float& rcorr,
0018                                     const float& cce_correction,
0019                                     const float& sigmaNoiseGeV,
0020                                     const float& xmin,
0021                                     const float& xmax,
0022                                     const float& aterm,
0023                                     const float& cterm) {
0024   dst_soa.id_[tid] = src_soa.id_[tid];
0025   dst_soa.energy_[tid] = src_soa.amplitude_[tid] * weight * 0.001f * __fdividef(rcorr, cce_correction);
0026   dst_soa.time_[tid] = src_soa.jitter_[tid];
0027 
0028   HeterogeneousHGCSiliconDetId detid(src_soa.id_[tid]);
0029   dst_soa.flagBits_[tid] = 0 | (0x1 << HGCRecHit::kGood);
0030   float son = __fdividef(dst_soa.energy_[tid], sigmaNoiseGeV);
0031   float son_norm = fminf(32.f, son) / 32.f * ((1 << 8) - 1);
0032   long int son_round = lroundf(son_norm);
0033   //there is an extra 0.125 factor in HGCRecHit::signalOverSigmaNoise(), which should not affect CPU/GPU comparison
0034   dst_soa.son_[tid] = static_cast<uint8_t>(son_round);
0035 
0036   //get time resolution
0037   //https://github.com/cms-sw/cmssw/blob/master/RecoLocalCalo/HGCalRecProducers/src/ComputeClusterTime.cc#L50
0038   /*Maxmin trick to avoid conditions within the kernel (having xmin < xmax)
0039     3 possibilities: 1) xval -> xmin -> xmax
0040     2) xmin -> xval -> xmax
0041     3) xmin -> xmax -> xval
0042     The time error is calculated with the number in the middle.
0043   */
0044   float denominator = fminf(fmaxf(son, xmin), xmax);
0045   float div_ = __fdividef(aterm, denominator);
0046   dst_soa.timeError_[tid] = dst_soa.time_[tid] < 0 ? -1 : __fsqrt_rn(div_ * div_ + cterm * cterm);
0047   //if dst_soa.time_[tid] < 1 always, then the above conditional expression can be replaced by
0048   //dst_soa.timeError_[tid] = fminf( fmaxf( dst_soa.time_[tid]-1, -1 ), sqrt( div_*div_ + cterm*cterm ) )
0049   //which is *not* conditional, and thus potentially faster; compare to HGCalRecHitWorkerSimple.cc
0050 }
0051 
0052 __device__ void make_rechit_scintillator(
0053     unsigned tid, HGCRecHitSoA dst_soa, HGCUncalibRecHitSoA src_soa, const float& weight, const float& sigmaNoiseGeV) {
0054   dst_soa.id_[tid] = src_soa.id_[tid];
0055   dst_soa.energy_[tid] = src_soa.amplitude_[tid] * weight * 0.001f;
0056   dst_soa.time_[tid] = src_soa.jitter_[tid];
0057 
0058   HeterogeneousHGCScintillatorDetId detid(src_soa.id_[tid]);
0059   dst_soa.flagBits_[tid] = 0 | (0x1 << HGCRecHit::kGood);
0060   float son = __fdividef(dst_soa.energy_[tid], sigmaNoiseGeV);
0061   float son_norm = fminf(32.f, son) / 32.f * ((1 << 8) - 1);
0062   long int son_round = lroundf(son_norm);
0063   //there is an extra 0.125 factor in HGCRecHit::signalOverSigmaNoise(), which should not affect CPU/GPU comparison
0064   dst_soa.son_[tid] = static_cast<uint8_t>(son_round);
0065   dst_soa.timeError_[tid] = -1;
0066 }
0067 
0068 __device__ float get_thickness_correction(const int& type,
0069                                           const double (&rcorr)[HGChefUncalibRecHitConstantData::hef_rcorr]) {
0070   return __fdividef(1.f, (float)rcorr[type]);
0071 }
0072 
0073 __device__ float get_noise(const int& type, const double (&noise_fC)[HGChefUncalibRecHitConstantData::hef_noise_fC]) {
0074   return (float)noise_fC[type];
0075 }
0076 
0077 __device__ float get_cce_correction(const int& type, const double (&cce)[HGChefUncalibRecHitConstantData::hef_cce]) {
0078   return (float)cce[type];
0079 }
0080 
0081 __device__ float get_fCPerMIP(const int& type,
0082                               const double (&fCPerMIP)[HGChefUncalibRecHitConstantData::hef_fCPerMIP]) {
0083   return (float)fCPerMIP[type];
0084 }
0085 
0086 __global__ void ee_to_rechit(HGCRecHitSoA dst_soa,
0087                              HGCUncalibRecHitSoA src_soa,
0088                              const HGCeeUncalibRecHitConstantData cdata,
0089                              int length) {
0090   unsigned tid = blockDim.x * blockIdx.x + threadIdx.x;
0091 
0092   for (unsigned i = tid; i < length; i += blockDim.x * gridDim.x) {
0093     HeterogeneousHGCSiliconDetId detid(src_soa.id_[i]);
0094     int32_t layer = detid.layer();
0095     float weight = get_weight_from_layer(layer, cdata.weights_);
0096     float rcorr = get_thickness_correction(detid.type(), cdata.rcorr_);
0097     float noise = get_noise(detid.type(), cdata.noise_fC_);
0098     float cce_correction = get_cce_correction(detid.type(), cdata.cce_);
0099     float fCPerMIP = get_fCPerMIP(detid.type(), cdata.fCPerMIP_);
0100     float sigmaNoiseGeV = 1e-3 * weight * rcorr * __fdividef(noise, fCPerMIP);
0101     make_rechit_silicon(i,
0102                         dst_soa,
0103                         src_soa,
0104                         weight,
0105                         rcorr,
0106                         cce_correction,
0107                         sigmaNoiseGeV,
0108                         cdata.xmin_,
0109                         cdata.xmax_,
0110                         cdata.aterm_,
0111                         cdata.cterm_);
0112   }
0113   __syncthreads();
0114 }
0115 
0116 __global__ void hef_to_rechit(HGCRecHitSoA dst_soa,
0117                               HGCUncalibRecHitSoA src_soa,
0118                               const HGChefUncalibRecHitConstantData cdata,
0119                               int length) {
0120   unsigned tid = blockDim.x * blockIdx.x + threadIdx.x;
0121   for (unsigned i = tid; i < length; i += blockDim.x * gridDim.x) {
0122     /*Uncomment the lines set to 1. as soon as those factors are centrally defined for the HSi.
0123         CUDADataFormats/HGCal/interface/HGCUncalibRecHitsToRecHitsConstants.h maxsizes_constanats will perhaps have to be changed (change some 3's to 6's) 
0124       */
0125     HeterogeneousHGCSiliconDetId detid(src_soa.id_[i]);
0126     int32_t layer = detid.layer() + cdata.layerOffset_;
0127     float weight = get_weight_from_layer(layer, cdata.weights_);
0128     float rcorr = 1.f;  //get_thickness_correction(detid.type(), cdata.rcorr_);
0129     float noise = get_noise(detid.type(), cdata.noise_fC_);
0130     float cce_correction = 1.f;  //get_cce_correction(detid.type(), cdata.cce_);
0131     float fCPerMIP = get_fCPerMIP(detid.type(), cdata.fCPerMIP_);
0132     float sigmaNoiseGeV = 1e-3 * weight * rcorr * __fdividef(noise, fCPerMIP);
0133     make_rechit_silicon(i,
0134                         dst_soa,
0135                         src_soa,
0136                         weight,
0137                         rcorr,
0138                         cce_correction,
0139                         sigmaNoiseGeV,
0140                         cdata.xmin_,
0141                         cdata.xmax_,
0142                         cdata.aterm_,
0143                         cdata.cterm_);
0144   }
0145   __syncthreads();
0146 }
0147 
0148 __global__ void heb_to_rechit(HGCRecHitSoA dst_soa,
0149                               HGCUncalibRecHitSoA src_soa,
0150                               const HGChebUncalibRecHitConstantData cdata,
0151                               int length) {
0152   unsigned tid = blockDim.x * blockIdx.x + threadIdx.x;
0153 
0154   for (unsigned i = tid; i < length; i += blockDim.x * gridDim.x) {
0155     HeterogeneousHGCScintillatorDetId detid(src_soa.id_[i]);
0156     int32_t layer = detid.layer() + cdata.layerOffset_;
0157     float weight = get_weight_from_layer(layer, cdata.weights_);
0158     float noise = cdata.noise_MIP_;
0159     float sigmaNoiseGeV = 1e-3 * noise * weight;
0160     make_rechit_scintillator(i, dst_soa, src_soa, weight, sigmaNoiseGeV);
0161   }
0162   __syncthreads();
0163 }