Xmipp  v3.23.11-Nereus
Public Member Functions | List of all members
CudaBSplineGeoTransformer< T > Class Template Reference

#include <cuda_bspline_geo_transformer.h>

Inheritance diagram for CudaBSplineGeoTransformer< T >:
Inheritance graph
[legend]
Collaboration diagram for CudaBSplineGeoTransformer< T >:
Collaboration graph
[legend]

Public Member Functions

 CudaBSplineGeoTransformer ()
 
virtual ~CudaBSplineGeoTransformer ()
 
void setSrc (const T *data) override
 
const T * getSrc () const override
 
T * getDest () const override
 
void copySrcToDest () override
 
T * interpolate (const std::vector< float > &matrices) override
 
void sum (T *dest, const std::vector< float > &weights, size_t firstN, T norm) override
 
- Public Member Functions inherited from BSplineGeoTransformer< T >
 BSplineGeoTransformer ()
 
 BSplineGeoTransformer (const BSplineGeoTransformer &)=delete
 
 BSplineGeoTransformer (const BSplineGeoTransformer &&)=delete
 
virtual ~BSplineGeoTransformer ()
 
BSplineGeoTransformeroperator= (const BSplineGeoTransformer &)=delete
 
BSplineGeoTransformeroperator= (const BSplineGeoTransformer &&)=delete
 
- Public Member Functions inherited from AGeoTransformer< BSplineTransformSettings< T >, T >
 AGeoTransformer ()
 
virtual ~AGeoTransformer ()=default
 
void init (const BSplineTransformSettings< T > &s, bool reuse)
 
const BSplineTransformSettings< T > & getSettings () const
 

Additional Inherited Members

- Protected Member Functions inherited from BSplineGeoTransformer< T >
virtual bool canBeReused (const BSplineTransformSettings< T > &s) const override
 
- Protected Member Functions inherited from AGeoTransformer< BSplineTransformSettings< T >, T >
constexpr bool isInitialized () const
 
constexpr bool isSrcSet () const
 
void setIsSrcSet (bool status)
 

Detailed Description

template<typename T>
class CudaBSplineGeoTransformer< T >

Definition at line 36 of file cuda_bspline_geo_transformer.h.

Constructor & Destructor Documentation

◆ CudaBSplineGeoTransformer()

template<typename T>
CudaBSplineGeoTransformer< T >::CudaBSplineGeoTransformer ( )
inline

Definition at line 38 of file cuda_bspline_geo_transformer.h.

38  {
39  setDefault();
40  }

◆ ~CudaBSplineGeoTransformer()

template<typename T>
virtual CudaBSplineGeoTransformer< T >::~CudaBSplineGeoTransformer ( )
inlinevirtual

Definition at line 42 of file cuda_bspline_geo_transformer.h.

42  {
43  release();
44  }

Member Function Documentation

◆ copySrcToDest()

template<typename T >
void CudaBSplineGeoTransformer< T >::copySrcToDest ( )
overridevirtual

Reimplemented from BSplineGeoTransformer< T >.

Definition at line 78 of file cuda_bspline_geo_transformer.cpp.

78  {
79  bool isReady = this->isInitialized()
80  && this->isSrcSet();
81  if ( ! isReady) {
82  REPORT_ERROR(ERR_LOGIC_ERROR, "Instance is either not initialized or the 'src' has not been set.");
83  }
84  gpuErrchk(cudaMemcpy(
85  m_d_dest,
86  m_d_src,
87  this->getSettings().dims.size() * sizeof(T),
88  cudaMemcpyDeviceToDevice)); // FIXME DS use proper stream
89 }
#define gpuErrchk(code)
Definition: cuda_asserts.h:31
const BSplineTransformSettings< T > & getSettings() const
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211
Some logical error in the pipeline.
Definition: xmipp_error.h:147

◆ getDest()

template<typename T>
T* CudaBSplineGeoTransformer< T >::getDest ( ) const
inlineoverridevirtual

Reimplemented from BSplineGeoTransformer< T >.

Definition at line 52 of file cuda_bspline_geo_transformer.h.

52  {
53  return m_d_dest;
54  }

◆ getSrc()

template<typename T>
const T* CudaBSplineGeoTransformer< T >::getSrc ( ) const
inlineoverridevirtual

Reimplemented from BSplineGeoTransformer< T >.

Definition at line 48 of file cuda_bspline_geo_transformer.h.

48  {
49  return m_d_src;
50  }

◆ interpolate()

template<typename T >
T * CudaBSplineGeoTransformer< T >::interpolate ( const std::vector< float > &  matrices)
overridevirtual

Reimplemented from BSplineGeoTransformer< T >.

Definition at line 121 of file cuda_bspline_geo_transformer.cpp.

121  {
122  const auto& dims = this->getSettings().dims;
123  size_t bytes = dims.n() * 9 * sizeof(float);
124  auto stream = *(cudaStream_t*)m_stream->stream();
125  gpuErrchk(cudaMemcpyAsync(
126  m_d_matrices,
127  matrices.data(),
128  bytes,
129  cudaMemcpyHostToDevice, stream));
130  dim3 dimBlock(64, 1, 1);
131  dim3 dimGrid(
132  std::ceil(dims.x() / (float)dimBlock.x),
133  std::ceil(dims.y() / (float)dimBlock.y),
134  dims.n());
135 
136  interpolateKernel<<<dimGrid, dimBlock, 0, stream>>> (
137  m_d_src, m_d_dest, m_d_matrices,
138  dims.x(), dims.y());
139  m_stream->synch();
140 
141  return m_d_dest;
142 }
#define gpuErrchk(code)
Definition: cuda_asserts.h:31
void * stream() const
Definition: gpu.h:50
const BSplineTransformSettings< T > & getSettings() const
void synch() const
Definition: gpu.cpp:129

◆ setSrc()

template<typename T >
void CudaBSplineGeoTransformer< T >::setSrc ( const T *  data)
overridevirtual

Reimplemented from BSplineGeoTransformer< T >.

Definition at line 32 of file cuda_bspline_geo_transformer.cpp.

32  {
33  const auto& s = this->getSettings();
34  gpuErrchk(cudaMemcpy(
35  m_d_src,
36  h_data,
37  this->getSettings().dims.size() * sizeof(T),
38  cudaMemcpyHostToDevice));
39  this->setIsSrcSet(true);
40 }
#define gpuErrchk(code)
Definition: cuda_asserts.h:31
const BSplineTransformSettings< T > & getSettings() const

◆ sum()

template<typename T >
void CudaBSplineGeoTransformer< T >::sum ( T *  dest,
const std::vector< float > &  weights,
size_t  firstN,
norm 
)
overridevirtual

Reimplemented from BSplineGeoTransformer< T >.

Definition at line 145 of file cuda_bspline_geo_transformer.cpp.

147  {
148  const auto& dims = this->getSettings().dims.copyForN(firstN);
149  assert(weights.size() == dims.n());
150  assert(0 < dims.n());
151  dim3 dimBlock(32, 32, 1);
152  dim3 dimGrid(
153  std::ceil(dims.x() / (float)dimBlock.x),
154  std::ceil(dims.y() / (float)dimBlock.y),
155  1);
156  auto stream = *(cudaStream_t*)m_stream->stream();
157 
158  // copy weights
159  T *d_weights;
160  size_t bytesWeight = weights.size() * sizeof(float);
161  gpuErrchk(cudaMalloc(&d_weights, bytesWeight));
162  gpuErrchk(cudaMemcpyAsync(
163  d_weights,
164  weights.data(),
165  bytesWeight,
166  cudaMemcpyHostToDevice, stream));
167 
168  sumKernel<<<dimGrid, dimBlock, 0, stream>>> (
169  m_d_dest, m_d_src, d_weights,// FIXME DS this will damage src data
170  dims.x(), dims.y(), dims.n(), norm);
171 
172  size_t bytes = dims.sizeSingle() * sizeof(T);
173  gpuErrchk(cudaMemcpyAsync(
174  dest,
175  m_d_src,
176  bytes,
177  cudaMemcpyDeviceToHost, stream));
178 
179  gpuErrchk(cudaFree(d_weights));
180  m_stream->synch();
181 
182 }
#define gpuErrchk(code)
Definition: cuda_asserts.h:31
void * stream() const
Definition: gpu.h:50
const BSplineTransformSettings< T > & getSettings() const
T norm(const std::vector< T > &v)
Definition: vector_ops.h:399
void synch() const
Definition: gpu.cpp:129

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