Xmipp  v3.23.11-Nereus
Public Member Functions | Public Attributes | Static Public Attributes | List of all members
Gpu::CDF< T > Struct Template Reference

#include <cuda_cdf.h>

Collaboration diagram for Gpu::CDF< T >:
Collaboration graph
[legend]

Public Member Functions

 CDF (size_t volume_size, T multConst=1.0, T probStep=0.005)
 
 ~CDF ()
 
void calculateCDF (const T *d_filtered1, const T *d_filtered2)
 
void calculateCDF (const T *d_S)
 
void _calculateDifference (const T *__restrict__ d_filtered1, const T *__restrict__ d_filtered2)
 
void _calculateSquare (const T *__restrict__ d_S)
 
void _updateProbabilities ()
 

Public Attributes

T * d_V
 
T * d_x
 
T * d_probXLessThanx
 
size_t volume_size
 
probStep
 
multConst
 
Nsteps
 

Static Public Attributes

static constexpr size_t type_size = sizeof(T)
 

Detailed Description

template<typename T>
struct Gpu::CDF< T >

Gpu version of Cumulative density function. This function computes a table with the cumulative density function

Definition at line 33 of file cuda_cdf.h.

Constructor & Destructor Documentation

◆ CDF()

template<typename T >
CDF::CDF ( size_t  volume_size,
multConst = 1.0,
probStep = 0.005 
)

Definition at line 41 of file cuda_cdf.cpp.

45 , Nsteps(round(1.0/probStep)) {
46  gpuErrchk( cudaMalloc((void**)&d_V, volume_size * type_size) );
47  gpuErrchk( cudaMalloc((void**)&d_x, Nsteps * type_size) );
48  gpuErrchk( cudaMalloc((void**)&d_probXLessThanx, Nsteps * type_size) );
49 }
T probStep
Definition: cuda_cdf.h:41
#define gpuErrchk(code)
Definition: cuda_asserts.h:31
T * d_x
Definition: cuda_cdf.h:37
static constexpr size_t type_size
Definition: cuda_cdf.h:34
size_t volume_size
Definition: cuda_cdf.h:40
T * d_probXLessThanx
Definition: cuda_cdf.h:38
T multConst
Definition: cuda_cdf.h:42
int round(double x)
Definition: ap.cpp:7245
T * d_V
Definition: cuda_cdf.h:36
T Nsteps
Definition: cuda_cdf.h:43

◆ ~CDF()

template<typename T >
CDF::~CDF ( )

Definition at line 52 of file cuda_cdf.cpp.

52  {
53  gpuErrchk( cudaFree(d_V) );
54  gpuErrchk( cudaFree(d_x) );
55  gpuErrchk( cudaFree(d_probXLessThanx) );
56 }
#define gpuErrchk(code)
Definition: cuda_asserts.h:31
T * d_x
Definition: cuda_cdf.h:37
T * d_probXLessThanx
Definition: cuda_cdf.h:38
T * d_V
Definition: cuda_cdf.h:36

Member Function Documentation

◆ _calculateDifference()

template<typename T >
void CDF::_calculateDifference ( const T *__restrict__  d_filtered1,
const T *__restrict__  d_filtered2 
)

Definition at line 73 of file cuda_cdf.cpp.

73  {
74  T* __restrict__ d_V = this->d_V;
75  const T multConst = this->multConst;
76 
77  auto compute_diff = [=] __device__ (int index) {
78  T diff = d_filtered1[index] - d_filtered2[index];
79  d_V[index] = multConst * diff * diff;
80  };
81 
82  thrust::for_each_n(thrust::device, thrust::counting_iterator<int>(0), volume_size, compute_diff);
83 }
viol index
size_t volume_size
Definition: cuda_cdf.h:40
T multConst
Definition: cuda_cdf.h:42
T * d_V
Definition: cuda_cdf.h:36

◆ _calculateSquare()

template<typename T >
void CDF::_calculateSquare ( const T *__restrict__  d_S)

Definition at line 86 of file cuda_cdf.cpp.

86  {
87  T* __restrict__ d_V = this->d_V;
88 
89  auto k = [=] __host__ __device__ (int index) {
90  T val = d_S[index];
91  d_V[index] = val * val;
92  };
93 
94  thrust::for_each_n(thrust::device, thrust::counting_iterator<int>(0), volume_size, k);
95 }
ql0001_ & k(htemp+1),(cvec+1),(atemp+1),(bj+1),(bl+1),(bu+1),(x+1),(clamda+1), &iout, infoqp, &zero,(w+1), &lenw,(iw+1), &leniw, &glob_grd.epsmac
viol index
size_t volume_size
Definition: cuda_cdf.h:40
T * d_V
Definition: cuda_cdf.h:36

◆ _updateProbabilities()

template<typename T >
void CDF::_updateProbabilities ( )

Definition at line 103 of file cuda_cdf.cpp.

103  {
104  T* __restrict__ d_V = this->d_V;
105  T* __restrict__ kd_x = this->d_x;
106  T* __restrict__ kd_probXLessThanx = this->d_probXLessThanx;
107  const T kprobStep = this->probStep;
108  const size_t k_volume_size = this->volume_size;
109 
110  auto k = [d_V, kd_x, kd_probXLessThanx, kprobStep, k_volume_size] __device__ (int index) {
111  int i = 0;
112  for (T p = kprobStep / 2; p < 1; p += kprobStep, i++) {
113  int idx = static_cast<int>(round(p * k_volume_size));
114  kd_probXLessThanx[i] = p;
115  kd_x[i] = d_V[idx];
116  }
117  };
118 
119  thrust::for_each_n(thrust::device, thrust::counting_iterator<int>(0), 1, k);
120 }
T probStep
Definition: cuda_cdf.h:41
#define i
ql0001_ & k(htemp+1),(cvec+1),(atemp+1),(bj+1),(bl+1),(bu+1),(x+1),(clamda+1), &iout, infoqp, &zero,(w+1), &lenw,(iw+1), &leniw, &glob_grd.epsmac
T * d_x
Definition: cuda_cdf.h:37
viol index
size_t volume_size
Definition: cuda_cdf.h:40
T * d_probXLessThanx
Definition: cuda_cdf.h:38
int round(double x)
Definition: ap.cpp:7245
T * d_V
Definition: cuda_cdf.h:36

◆ calculateCDF() [1/2]

template<typename T >
void CDF::calculateCDF ( const T *  d_filtered1,
const T *  d_filtered2 
)

Definition at line 59 of file cuda_cdf.cpp.

59  {
60  _calculateDifference(d_filtered1, d_filtered2);
61  sort();
63 }
void _updateProbabilities()
Definition: cuda_cdf.cpp:103
void _calculateDifference(const T *__restrict__ d_filtered1, const T *__restrict__ d_filtered2)
Definition: cuda_cdf.cpp:73

◆ calculateCDF() [2/2]

template<typename T >
void CDF::calculateCDF ( const T *  d_S)

Definition at line 66 of file cuda_cdf.cpp.

66  {
67  _calculateSquare(d_S);
68  sort();
70 }
void _updateProbabilities()
Definition: cuda_cdf.cpp:103
void _calculateSquare(const T *__restrict__ d_S)
Definition: cuda_cdf.cpp:86

Member Data Documentation

◆ d_probXLessThanx

template<typename T>
T* Gpu::CDF< T >::d_probXLessThanx

Definition at line 38 of file cuda_cdf.h.

◆ d_V

template<typename T>
T* Gpu::CDF< T >::d_V

Definition at line 36 of file cuda_cdf.h.

◆ d_x

template<typename T>
T* Gpu::CDF< T >::d_x

Definition at line 37 of file cuda_cdf.h.

◆ multConst

template<typename T>
T Gpu::CDF< T >::multConst

Definition at line 42 of file cuda_cdf.h.

◆ Nsteps

template<typename T>
T Gpu::CDF< T >::Nsteps

Definition at line 43 of file cuda_cdf.h.

◆ probStep

template<typename T>
T Gpu::CDF< T >::probStep

Definition at line 41 of file cuda_cdf.h.

◆ type_size

template<typename T>
constexpr size_t Gpu::CDF< T >::type_size = sizeof(T)
static

Definition at line 34 of file cuda_cdf.h.

◆ volume_size

template<typename T>
size_t Gpu::CDF< T >::volume_size

Definition at line 40 of file cuda_cdf.h.


The documentation for this struct was generated from the following files: