1 #ifndef NPSTAT_CONVOLUTIONENGINEND_HH_
2 #define NPSTAT_CONVOLUTIONENGINEND_HH_
21 #include "npstat/wrap/arrayNDToNumpy.hh"
22 #include "npstat/wrap/numpyArrayUtils.hh"
39 unsigned optimization = FFTW_ESTIMATE);
43 template <
typename Real>
44 void setFilter(
const Real* data,
const unsigned *dataShape,
51 template <
typename Real1,
typename Real2>
53 const unsigned *dataShape,
unsigned dataRank);
61 template <
typename Real>
63 const Real* data,
const unsigned *dataShape,
71 template <
typename Real1,
typename Real2>
73 const Real1* in, Real2* out,
74 const unsigned *dataShape,
unsigned dataRank);
80 template <
typename Real2>
82 const unsigned *dataShape,
89 template <
typename Real2>
91 const unsigned *dataShape,
108 inline unsigned rank()
const {
return shape_.size();}
111 inline const std::vector<unsigned>&
shape()
const {
return shape_;}
117 inline bool isFilterSet()
const {
return validFilter_;}
118 inline unsigned filterBankSize()
const {
return filterMap_.size();}
121 typedef std::map<unsigned long,FourierImage*> FilterMap;
127 template <
typename Real1,
typename Real2>
128 void convolveWithImage(
const fftw_complex* image,
129 const Real1* in, Real2* out,
130 const unsigned *dataShape,
unsigned dataRank);
136 fftw_complex* filterImage;
138 std::vector<unsigned> shape_;
139 unsigned long dataLen_;
140 unsigned long cmplLen_;
143 FilterMap filterMap_;
147 inline void setFilter2(PyObject* input)
149 if (!isNumpyShapeCompatible(input, &shape_[0], shape_.size()))
throw std::invalid_argument(
150 "In npstat::ConvolutionEngineND::setFilter2: incompatible array shape");
151 const int typenum = numpyArrayType(input);
152 if (typenum == NumpyTypecode<double>::code)
153 setFilter((
double*)PyArray_DATA((PyArrayObject*)input), &shape_[0], shape_.size());
154 else if (typenum == NumpyTypecode<float>::code)
155 setFilter((
float*)PyArray_DATA((PyArrayObject*)input), &shape_[0], shape_.size());
156 else if (typenum == NumpyTypecode<int>::code)
157 setFilter((
int*)PyArray_DATA((PyArrayObject*)input), &shape_[0], shape_.size());
158 else if (typenum == NumpyTypecode<long long>::code)
159 setFilter((
long long*)PyArray_DATA((PyArrayObject*)input), &shape_[0], shape_.size());
160 else if (typenum == NumpyTypecode<long int>::code)
161 setFilter((
long int*)PyArray_DATA((PyArrayObject*)input), &shape_[0], shape_.size());
162 else if (typenum == NumpyTypecode<unsigned char>::code)
163 setFilter((
unsigned char*)PyArray_DATA((PyArrayObject*)input), &shape_[0], shape_.size());
165 throw std::invalid_argument(
"In npstat::ConvolutionEngineND::setFilter2: "
166 "unsupported input array type");
169 inline void depositFilter2(
const unsigned long id,
172 const int typenum = numpyArrayType(input);
173 if (typenum != NumpyTypecode<double>::code)
throw std::invalid_argument(
174 "In npstat::ConvolutionEngineND::depositFilter2: input is not an array of doubles");
175 if (!isNumpyShapeCompatible(input, &shape_[0], shape_.size()))
throw std::invalid_argument(
176 "In npstat::ConvolutionEngineND::depositFilter2: incompatible array shape");
177 double* in = (
double*)PyArray_DATA((PyArrayObject*)input);
181 inline PyObject* convolveWithFilter2(PyObject* input,
182 const bool makeNonNegative =
false)
184 const int typenum = numpyArrayType(input);
185 if (typenum != NumpyTypecode<double>::code)
throw std::invalid_argument(
186 "In npstat::ConvolutionEngineND::convolveWithFilter2: input is not an array of doubles");
187 const unsigned dataRank = shape_.size();
188 assert(dataRank <= CHAR_BIT*
sizeof(
unsigned long));
189 if (!isNumpyShapeCompatible(input, &shape_[0], dataRank))
throw std::invalid_argument(
190 "In npstat::ConvolutionEngineND::convolveWithFilter2: incompatible array shape");
191 double* in = (
double*)PyArray_DATA((PyArrayObject*)input);
192 npy_intp sh[CHAR_BIT*
sizeof(
unsigned long)];
193 for (
unsigned i=0; i<dataRank; ++i)
195 PyObject* xarr = PyArray_SimpleNew(dataRank, sh, typenum);
198 double* out = (
double*)PyArray_DATA((PyArrayObject*)xarr);
203 catch (
const std::exception& e)
209 for (
unsigned long i=0; i<dataLen_; ++i)
216 inline PyObject* convolveWithDeposit2(
const unsigned long id,
218 const bool makeNonNegative =
false)
220 const int typenum = numpyArrayType(input);
221 if (typenum != NumpyTypecode<double>::code)
throw std::invalid_argument(
222 "In npstat::ConvolutionEngineND::convolveWithDeposit2: input is not an array of doubles");
223 const unsigned dataRank = shape_.size();
224 assert(dataRank <= CHAR_BIT*
sizeof(
unsigned long));
225 if (!isNumpyShapeCompatible(input, &shape_[0], dataRank))
throw std::invalid_argument(
226 "In npstat::ConvolutionEngineND::convolveWithDeposit2: incompatible array shape");
227 double* in = (
double*)PyArray_DATA((PyArrayObject*)input);
228 npy_intp sh[CHAR_BIT*
sizeof(
unsigned long)];
229 for (
unsigned i=0; i<dataRank; ++i)
231 PyObject* xarr = PyArray_SimpleNew(dataRank, sh, typenum);
234 double* out = (
double*)PyArray_DATA((PyArrayObject*)xarr);
239 catch (
const std::exception& e)
245 for (
unsigned long i=0; i<dataLen_; ++i)
252 inline PyObject* convolveFilterAndDeposit2(
const unsigned long id,
253 const bool makeNonNegative =
false)
255 const int typenum = NumpyTypecode<double>::code;
256 const unsigned dataRank = shape_.size();
257 assert(dataRank <= CHAR_BIT*
sizeof(
unsigned long));
258 npy_intp sh[CHAR_BIT*
sizeof(
unsigned long)];
259 for (
unsigned i=0; i<dataRank; ++i)
261 PyObject* xarr = PyArray_SimpleNew(dataRank, sh, typenum);
264 double* out = (
double*)PyArray_DATA((PyArrayObject*)xarr);
269 catch (
const std::exception& e)
275 for (
unsigned long i=0; i<dataLen_; ++i)
285 #include "npstat/nm/ConvolutionEngineND.icc"
Memory manager for Fourier transforms produced by FFTW.
Definition: ConvolutionEngineND.hh:31
bool discardFilter(unsigned long id)
void convolveFilterAndDeposit(unsigned long id, Real2 *out, const unsigned *dataShape, unsigned dataRank)
bool isShapeCompatible(const unsigned *shape, unsigned rank) const
unsigned long dataLength() const
Definition: ConvolutionEngineND.hh:114
void deconvolveFilterAndDeposit(unsigned long id, Real2 *out, const unsigned *dataShape, unsigned dataRank)
const std::vector< unsigned > & shape() const
Definition: ConvolutionEngineND.hh:111
void depositFilter(unsigned long id, const Real *data, const unsigned *dataShape, unsigned dataRank)
unsigned rank() const
Definition: ConvolutionEngineND.hh:108
void setFilter(const Real *data, const unsigned *dataShape, unsigned dataRank)
void convolveWithFilter(const Real1 *in, Real2 *out, const unsigned *dataShape, unsigned dataRank)
ConvolutionEngineND(const unsigned *shape, unsigned rank, unsigned optimization=FFTW_ESTIMATE)
void convolveWithDeposit(unsigned long id, const Real1 *in, Real2 *out, const unsigned *dataShape, unsigned dataRank)
Definition: AbsArrayProjector.hh:14