5#ifndef BALL_MATHS_TFFT1D_H
6#define BALL_MATHS_TFFT1D_H
8#ifndef BALL_COMMON_EXCEPTION_H
12#ifndef BALL_DATATYPE_REGULARDATA1D_H
36 template <
typename ComplexTraits>
38 :
public TRegularData1D<std::complex<typename ComplexTraits::ComplexPrecision> >
42 typedef std::complex<typename ComplexTraits::ComplexPrecision>
Complex;
67 TFFT1D(
Size ldn,
double stepPhys = 1.,
double origin = 0.,
bool inFourierSpace =
false);
284 template <
typename ComplexTraits>
288 inFourierSpace_(false),
290 planCalculated_(false)
295 template <
typename ComplexTraits>
298 if (ComplexVector::size() == fft1d.
size() &&
311 double min = inFourierSpace_ ? minFourier_ : minPhys_;
312 double max = inFourierSpace_ ? maxFourier_ : maxPhys_;
313 double step = inFourierSpace_ ? stepFourier_ : stepPhys_;
315 for (
double pos=min; pos<=max; pos+=step)
317 if (getData(pos) != fft1d.
getData(pos))
329 template <
typename ComplexTraits>
332 Position internalOrigin = (int) rint(trans_origin/stepPhys_);
333 if (internalOrigin <= length_)
335 origin_ = trans_origin;
337 minPhys_ = ((-1.)*origin_);
338 maxPhys_ = (((length_-1)*stepPhys_)-origin_);
339 minFourier_ = ((-1.)*(length_/2.-1)*stepFourier_);
340 maxFourier_ = ((length_/2.)*stepFourier_);
350 template <
typename ComplexTraits>
359 stepPhys_ = new_width;
362 minPhys_ = ((-1.)*origin_);
363 maxPhys_ = (((length_-1)*stepPhys_)-origin_);
364 minFourier_ = ((-1.)*(length_/2.-1)*stepFourier_);
365 maxFourier_ = ((length_/2.)*stepFourier_);
371 template <
typename ComplexTraits>
377 template <
typename ComplexTraits>
383 template <
typename ComplexTraits>
389 template <
typename ComplexTraits>
395 template <
typename ComplexTraits>
401 template <
typename ComplexTraits>
407 template <
typename ComplexTraits>
410 return (length_ - 1);
413 template <
typename ComplexTraits>
416 return numFourierToPhys_;
420 template <
typename ComplexTraits>
423 if (!inFourierSpace_)
425 if (position >= ComplexVector::size())
432 r = -origin_ + (
float)position * stepPhys_;
438 if (position >= ComplexVector::size())
453 r = (
float)x * stepFourier_;
459 template <
typename ComplexTraits>
463 double normalization=1.;
465 if (!inFourierSpace_)
467 result = (*this)[pos];
468 normalization=1./pow((
double)length_,(
int)numFourierToPhys_);
472 result = (*this)[pos]*phase(pos);
474 normalization=1./(sqrt(2.*
Constants::PI))/pow((
double)length_,(int)numFourierToPhys_);
477 result *= normalization;
482 template <
typename ComplexTraits>
487 double min = inFourierSpace_ ? minFourier_ : minPhys_;
488 double max = inFourierSpace_ ? maxFourier_ : maxPhys_;
489 double step = inFourierSpace_ ? stepFourier_ : stepPhys_;
491 if ((pos < min) || (pos > max))
496 double h = pos - min;
497 double mod = fmod(h,step);
504 double before = floor(h/step)*step+ min;
505 double after = ceil(h/step)*step+ min;
507 double t = (pos - before)/step;
509 result = getData(before)*(
typename ComplexTraits::ComplexPrecision)(1.-t);
510 result += getData(after)*(
typename ComplexTraits::ComplexPrecision)t;
515 template <
typename ComplexTraits>
519 if (!inFourierSpace_)
521 dummy =
Complex(val.real()*((
typename ComplexTraits::ComplexPrecision)pow((
typename ComplexTraits::ComplexPrecision)(length_),(
int)numFourierToPhys_)),
522 val.imag()*((
typename ComplexTraits::ComplexPrecision)pow((
typename ComplexTraits::ComplexPrecision)(length_),(
int)numFourierToPhys_)));
528 val*=phase(pos)*(
typename ComplexTraits::ComplexPrecision)((sqrt(2*
Constants::PI)/stepPhys_))
529 *(
typename ComplexTraits::ComplexPrecision)pow((
typename ComplexTraits::ComplexPrecision)length_,(
int)numFourierToPhys_);
537 template <
typename ComplexTraits>
542 if (!inFourierSpace_)
544 internalPos = (
Index) rint((pos+origin_)/stepPhys_);
548 internalPos = (
Index) rint(pos/stepFourier_);
552 internalPos+=length_;
556 if ((internalPos < 0) || (internalPos>=(
Index) length_))
561 return operator [] ((
Position)internalPos);
564 template <
typename ComplexTraits>
569 if (!inFourierSpace_)
571 internalPos = (
Index) rint((pos+origin_)/stepPhys_);
575 internalPos = (
Index) rint(pos/stepFourier_);
579 internalPos+=length_;
583 if ((internalPos < 0) || (internalPos>=(
Index) length_))
588 return operator [] ((
Position)internalPos);
591 template <
typename ComplexTraits>
595 *(rint(origin_/stepPhys_))
602 template <
typename ComplexTraits>
605 return inFourierSpace_;
673 template <
typename ComplexTraits>
691 newGrid[i] = from[i].real()*normalization;
732 x-=(int)(lengthX/2.);
736 x+=(int)(lengthX/2.);
740 r = ((
float)xp * stepFourierX);
742 newGrid[i] = (from[i]*(
typename ComplexTraits::ComplexPrecision)normalization*from.
phase(r)).real();
#define BALL_CREATE(name)
BALL_EXPORT std::ostream & operator<<(std::ostream &os, const Exception::GeneralException &e)
TFFT1D< BALL_FFTW_DEFAULT_TRAITS > FFT1D
std::complex< BALL_COMPLEX_PRECISION > Complex
BALL_EXTERN_VARIABLE const double PI
PI.
BALL_INLINE size_type size() const
BALL_INLINE void setDimension(const CoordinateType &dimension)
BALL_INLINE void setOrigin(const CoordinateType &origin)
const ValueType & operator[](const IndexType &index) const
bool translate(double trans_origin)
double getFourierSpaceMax() const
double getGridCoordinates(Position position) const
TFFT1D(const TFFT1D &data)
Copy constructor.
void setNumberOfFFTTransforms(Size num)
ComplexTraits::FftwPlan planForward_
std::complex< typename ComplexTraits::ComplexPrecision > Complex
void setData(double pos, Complex val)
ComplexTraits::FftwPlan planBackward_
const Complex & operator[](const Position &pos) const
Complex & operator[](const Position &pos)
bool isInFourierSpace() const
Complex getInterpolatedValue(const double pos) const
virtual ~TFFT1D()
Destructor.
Complex phase(const double pos) const
double getPhysStepWidth() const
Size getNumberOfInverseTransforms() const
TRegularData1D< std::complex< typename ComplexTraits::ComplexPrecision > > ComplexVector
Complex getData(const double pos) const
bool setPhysStepWidth(double new_width)
void setNumberOfiFFTTransforms(Size num)
TFFT1D()
Default constructor.
bool operator==(const TFFT1D &fft1d) const
double getFourierSpaceMin() const
Complex & operator[](const double pos)
double getPhysSpaceMin() const
double getFourierStepWidth() const
double getPhysSpaceMax() const
const TFFT1D & operator=(const TFFT1D &fft1d)
Assignment operator.
TFFT1D(Size ldn, double stepPhys=1., double origin=0., bool inFourierSpace=false)