Xmipp  v3.23.11-Nereus
Classes | Functions | Variables
BSpline Transform Settings
Collaboration diagram for BSpline Transform Settings:

Classes

class  BSplineTransformSettings< T >
 
class  BSplineGeoTransformer< T >
 

Functions

void BSplineTransformSettings< T >::check () const override
 
 BSplineGeoTransformer< T >::BSplineGeoTransformer ()
 
 BSplineGeoTransformer< T >::BSplineGeoTransformer (const BSplineGeoTransformer &)=delete
 
 BSplineGeoTransformer< T >::BSplineGeoTransformer (const BSplineGeoTransformer &&)=delete
 
virtual BSplineGeoTransformer< T >::~BSplineGeoTransformer ()
 
BSplineGeoTransformerBSplineGeoTransformer< T >::operator= (const BSplineGeoTransformer &)=delete
 
BSplineGeoTransformerBSplineGeoTransformer< T >::operator= (const BSplineGeoTransformer &&)=delete
 
virtual void BSplineGeoTransformer< T >::setSrc (const T *data) override
 
virtual const T * BSplineGeoTransformer< T >::getSrc () const
 
virtual T * BSplineGeoTransformer< T >::getDest () const override
 
virtual void BSplineGeoTransformer< T >::copySrcToDest () override
 
virtual T * BSplineGeoTransformer< T >::interpolate (const std::vector< float > &matrices)
 
virtual void BSplineGeoTransformer< T >::sum (T *dest, const std::vector< float > &weights, size_t firstN, T norm)
 
virtual void BSplineGeoTransformer< T >::initialize (bool doAllocation) override
 
virtual void BSplineGeoTransformer< T >::release ()
 
virtual void BSplineGeoTransformer< T >::setDefault ()
 
virtual void BSplineGeoTransformer< T >::check () override
 
virtual bool BSplineGeoTransformer< T >::canBeReused (const BSplineTransformSettings< T > &s) const override
 

Variables

bool BSplineTransformSettings< T >::keepSrcCopy
 

Detailed Description

Function Documentation

◆ BSplineGeoTransformer() [1/3]

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

Definition at line 50 of file bspline_geo_transformer.h.

50  {
51  setDefault();
52  }

◆ BSplineGeoTransformer() [2/3]

template<typename T>
BSplineGeoTransformer< T >::BSplineGeoTransformer ( const BSplineGeoTransformer< T > &  )
delete

◆ BSplineGeoTransformer() [3/3]

template<typename T>
BSplineGeoTransformer< T >::BSplineGeoTransformer ( const BSplineGeoTransformer< T > &&  )
delete

◆ canBeReused()

template<typename T >
bool BSplineGeoTransformer< T >::canBeReused ( const BSplineTransformSettings< T > &  s) const
overrideprotectedvirtual

Implements AGeoTransformer< BSplineTransformSettings< T >, T >.

Definition at line 89 of file bspline_geo_transformer.cpp.

89  {
90  bool result = true;
91  if ( ! this->isInitialized()) {
92  return false;
93  }
94  auto &sOrig = this->getSettings();
95  result &= sOrig.dims.size() >= s.dims.size(); // previously, we needed more space
96  result &= !(( ! sOrig.keepSrcCopy) && s.keepSrcCopy); // we have a problem if now we need to make a copy and before not
97 
98  return result;
99 }
const BSplineTransformSettings< T > & getSettings() const
constexpr size_t size() const
Definition: dimensions.h:104

◆ check() [1/2]

template<typename T>
void BSplineTransformSettings< T >::check ( ) const
inlineoverridevirtual

Reimplemented from GeoTransformerSettings< T >.

Definition at line 41 of file bspline_geo_transformer.h.

41  {
43  }
virtual void check() const

◆ check() [2/2]

template<typename T >
void BSplineGeoTransformer< T >::check ( )
overrideprotectedvirtual

Implements AGeoTransformer< BSplineTransformSettings< T >, T >.

Definition at line 49 of file bspline_geo_transformer.cpp.

49  {
50  const auto &s = this->getSettings();
51  if (InterpolationDegree::Linear != s.degree) {
52  REPORT_ERROR(ERR_NOT_IMPLEMENTED, "Only linear interpolation is available");
53  }
54  if (s.doWrap) {
55  REPORT_ERROR(ERR_NOT_IMPLEMENTED, "Wrapping is not yet implemented");
56  }
57  if (InterpolationType::NToN != s.type) {
58  REPORT_ERROR(ERR_NOT_IMPLEMENTED, "Only NToN is currently implemented");
59  }
60 }
Case or algorithm not implemented yet.
Definition: xmipp_error.h:177
const BSplineTransformSettings< T > & getSettings() const
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211

◆ copySrcToDest()

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

Implements AGeoTransformer< BSplineTransformSettings< T >, T >.

Reimplemented in CudaBSplineGeoTransformer< T >.

Definition at line 64 of file bspline_geo_transformer.cpp.

64  {
65  bool isReady = this->isInitialized()
66  && this->isSrcSet();
67  if ( ! isReady) {
68  REPORT_ERROR(ERR_LOGIC_ERROR, "Instance is either not initialized or the 'src' has not been set.");
69  }
70  memcpy(m_dest.get(),
71  m_src,
72  this->getSettings().dims.size() * sizeof(T));
73 }
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>
virtual T* BSplineGeoTransformer< T >::getDest ( ) const
inlineoverridevirtual

Implements AGeoTransformer< BSplineTransformSettings< T >, T >.

Reimplemented in CudaBSplineGeoTransformer< T >.

Definition at line 71 of file bspline_geo_transformer.h.

71  {
72  return m_dest.get();
73  }

◆ getSrc()

template<typename T>
virtual const T* BSplineGeoTransformer< T >::getSrc ( ) const
inlinevirtual

Implements AGeoTransformer< BSplineTransformSettings< T >, T >.

Reimplemented in CudaBSplineGeoTransformer< T >.

Definition at line 67 of file bspline_geo_transformer.h.

67  {
68  return m_src;
69  }

◆ initialize()

template<typename T >
void BSplineGeoTransformer< T >::initialize ( bool  doAllocation)
overrideprotectedvirtual

Implements AGeoTransformer< BSplineTransformSettings< T >, T >.

Definition at line 31 of file bspline_geo_transformer.cpp.

31  {
32  const auto &s = this->getSettings();
33 
34  for (auto &hw : s.hw) {
35  if ( ! dynamic_cast<CPU*>(hw)) {
36  REPORT_ERROR(ERR_LOGIC_ERROR, "Instance of CPU is expected");
37  }
38  }
39 // m_threadPool.resize(s.hw.size()); // FIXME DS set to requested number of thread
40  m_threadPool.resize(CPU::findCores());
41 
42  if (doAllocation) {
43  release();
44  m_dest = std::unique_ptr<T[]>(new T[s.dims.size()]);
45  }
46 }
const BSplineTransformSettings< T > & getSettings() const
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211
static unsigned findCores()
Definition: cpu.h:41
void resize(int nThreads)
Definition: ctpl.h:70
Some logical error in the pipeline.
Definition: xmipp_error.h:147

◆ interpolate()

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

Reimplemented in CudaBSplineGeoTransformer< T >.

Definition at line 107 of file bspline_geo_transformer.cpp.

107  {
108  bool isReady = this->isInitialized()
109  && this->isSrcSet();
110  if ( ! isReady) {
111  REPORT_ERROR(ERR_LOGIC_ERROR, "Instance is either not initialized or the 'original' has not been set.");
112  }
113  const Dimensions dims = this->getSettings().dims;
114  const size_t n = dims.n();
115  const size_t z = dims.z();
116  const size_t y = dims.y();
117  const size_t x = dims.x();
118 
119  auto futures = std::vector<std::future<void>>();
120 
121  auto workload = [&](int id, size_t signalId){
122  size_t offset = signalId * dims.sizeSingle();
123  auto in = MultidimArray<T>(1, z, y, x, const_cast<T*>(m_src + offset)); // removing const, but data should not be changed
124  auto out = MultidimArray<T>(1, z, y, x, m_dest.get() + offset);
125  in.setXmippOrigin();
126  out.setXmippOrigin();
127  // compensate the movement
128  Matrix2D<double> m(3,3);
129  const float *f = matrices.data() + (9 * signalId);
130  for (int i = 0; i < 9; ++i) {
131  m.mdata[i] = f[i];
132  }
133  applyGeometry(xmipp_transformation::LINEAR, out, in, m, true, xmipp_transformation::DONT_WRAP);
134  };
135 
136  for (size_t i = 0; i < n; ++i) {
137  futures.emplace_back(m_threadPool.push(workload, i));
138  }
139  for (auto &f : futures) {
140  f.get();
141  }
142  return m_dest.get();
143 }
const BSplineTransformSettings< T > & getSettings() const
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211
auto push(F &&f, Rest &&... rest) -> std::future< decltype(f(0, rest...))>
Definition: ctpl.h:152
CUDA_HD constexpr size_t z() const
Definition: dimensions.h:69
static double * y
void applyGeometry(int SplineDegree, MultidimArray< std::complex< double > > &V2, const MultidimArray< std::complex< double > > &V1, const Matrix2D< double > &A, bool inv, bool wrap, std::complex< double > outside, MultidimArray< double > *BcoeffsPtr)
doublereal * x
#define i
CUDA_HD constexpr size_t x() const
Definition: dimensions.h:51
int in
double * f
CUDA_HD constexpr size_t sizeSingle() const
Definition: dimensions.h:100
CUDA_HD constexpr size_t y() const
Definition: dimensions.h:60
double z
CUDA_HD constexpr size_t n() const
Definition: dimensions.h:78
int m
int * n
Some logical error in the pipeline.
Definition: xmipp_error.h:147

◆ operator=() [1/2]

template<typename T>
BSplineGeoTransformer& BSplineGeoTransformer< T >::operator= ( const BSplineGeoTransformer< T > &  )
delete

◆ operator=() [2/2]

template<typename T>
BSplineGeoTransformer& BSplineGeoTransformer< T >::operator= ( const BSplineGeoTransformer< T > &&  )
delete

◆ release()

template<typename T >
void BSplineGeoTransformer< T >::release ( )
protectedvirtual

Definition at line 76 of file bspline_geo_transformer.cpp.

76  {
77  m_dest.release();
78  setDefault();
79 }

◆ setDefault()

template<typename T >
void BSplineGeoTransformer< T >::setDefault ( )
protectedvirtual

Definition at line 82 of file bspline_geo_transformer.cpp.

82  {
83  m_dest.reset();
84  m_src = nullptr;
85  m_threadPool.resize(1);
86 }
void resize(int nThreads)
Definition: ctpl.h:70

◆ setSrc()

template<typename T>
virtual void BSplineGeoTransformer< T >::setSrc ( const T *  data)
inlineoverridevirtual

Implements AGeoTransformer< BSplineTransformSettings< T >, T >.

Reimplemented in CudaBSplineGeoTransformer< T >.

Definition at line 62 of file bspline_geo_transformer.h.

62  {
63  m_src = data;
64  this->setIsSrcSet(nullptr != data);
65  }

◆ sum()

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

Reimplemented in CudaBSplineGeoTransformer< T >.

Definition at line 102 of file bspline_geo_transformer.cpp.

102  {
103  REPORT_ERROR(ERR_NOT_IMPLEMENTED, "This functionality is not yet available.");
104 }
Case or algorithm not implemented yet.
Definition: xmipp_error.h:177
#define REPORT_ERROR(nerr, ErrormMsg)
Definition: xmipp_error.h:211

◆ ~BSplineGeoTransformer()

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

Definition at line 56 of file bspline_geo_transformer.h.

56  {
57  release();
58  }

Variable Documentation

◆ keepSrcCopy

template<typename T>
bool BSplineTransformSettings< T >::keepSrcCopy

Definition at line 39 of file bspline_geo_transformer.h.