#ifndef myfft_h
#define myfft_h
/** \file myfft.h
 * Header file for \ref myfft.c
 * 
 * Contains declarations of the functions implemented in \ref myfft.c and
 * structures around the fftransformations containg plans fft_plan_3d, plans
 * on a level LevelPlan and plan_weight.
 * 
  Project: ParallelCarParrinello
  Jan Hamaekers
  2000

  File: myfft.h
  $Id: myfft.h,v 1.15 2006/08/21 09:16:16 foo Exp $
*/

// use double precision fft when we have it
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#ifdef HAVE_DFFTW_H
#include "dfftw.h"
#else
#include "fftw.h"
#endif

#ifdef HAVE_DRFFTW_H
#include "drfftw.h"
#else
#include "rfftw.h"
#endif


#define fftw_flag FFTW_ESTIMATE   //!< flag used when calling fftw routine to create plans

struct Problem;

/** Plan weight structure
 * \todo what is plan weight?
 */
struct plan_weight {
  double weight;  //!< the acual weight
  int Index;      //!< ???
  int PE;         //!< resides on which process
}; 

/** Structure around the fft plan for one level
 */
struct LevelPlan {
  int LevelFactor;        //!< ratio of level size between maximum level and this level
  int nyz;                //!< mesh points in yz-plane (nx times ny)
  int local_nyz;          //!< mesh points in yz-plane per process (local_nx times local_ny)
  int nx;                 //!< mesh points in x direction
  int local_nx;           //!< mesh points in x direction per process (nx/MaxPE)
  int start_nx;           //!< start index for mesh point set of individual process
  int ny;                 //!< mesh points in y direction
  int local_ny;           //!< mesh points in y direction per process (ny/MaxPE) 
  int nz;                 //!< mesh points in z direction (halved plus 1)
  int LevelBlockSize;     //!< cubed plan::LevelSize
  int SendRecvBlockSize;  //!< size of block 'pen' that is sent and received between processes in FFT
  int PermBlockSize;      //!< size of block for permutation (2x2x2 generally)
  int LocalSizeC;         //!< number of complex mesh points (nx times ny times nz)
  int LocalSizeR;         //!< number of complex mesh points (nx times ny times 2*nz-1)
  int LevelNo;            //!< continuous cardinal number of the level
  int N[NDIM];            //!< number of nodes on the grid for each dimension
  fftw_plan *xplanCR;
  fftw_plan *xplanRC;
  rfftwnd_plan *yzplanCR;
  rfftwnd_plan *yzplanRC;
};

/** Internal fftw plan structure.
 */
struct fft_plan_3d {
  MPI_Comm comm;            //!< involved MPI communicator group
  double E_cut;             //!< energy cutoff
  double GBasisSQ[NDIM];    //!< squared norm of reciprocal basis vectors
  int MaxPE;                //!< Size of group in the communicator
  int myPE;	                //!< Rank of process in the communicator
  int LevelSize;            //!< factor of number of mesh points between 0th and MaxLevel's level
  int Maxpw;                //!< maximum number of plan weights
  int LocalMaxpw;           //!< number of plan weights stored on this process
  int MaxLevel;             //!< maximum number of LatticeLevel's
  int LocalDataSize;
  int LocalWorkSize;
  int LevelBlockSize;
  int *MaxNoOfnFields;      //!< Maximum number of NField's per level
  int **NFields;            //!< NField per level per field
  int N[NDIM];		          //!< number of mesh points per dimension
  int NLSR[NDIM];	          //!< ratio of mesh points N over LevelSize
  int *PENo;
  int *LevelSizes;          //!< factor of mesh points between this and upper level per level
  int *Lev_CR_RC;
  struct plan_weight *pw;
  struct plan_weight *pwFS;
  struct plan_weight *Localpw;
  struct LevelPlan *Levplan;
  fftw_complex *local_data;
  fftw_complex *work;
  struct Problem *P;         //!< Problem at hand
};

struct fft_plan_3d *fft_3d_create_plan(struct Problem *P, MPI_Comm comm, const double E_cut, const int MaxLevel, const int *LevelSizes, const double GBasisSQ[NDIM], const int N[NDIM], const int *MaxNoOfnFields, int **NFields);
void fft_3d_complex_to_real(const struct fft_plan_3d *plan, const int Level, const int nFields, fftw_complex *local_data, fftw_complex *work); /* local data -> local_data, work destroyed */
void fft_3d_real_to_complex(const struct fft_plan_3d *plan, const int Level, const int nFields, fftw_complex *local_data, fftw_complex *work); /* local data -> local_data, work destroyed */
void fft_3d_destroy_plan(const struct Problem *P, struct fft_plan_3d *plan);
void InitDataTestR(const struct fft_plan_3d *plan, const int Level, const int nFields,fftw_real *local_data);
/*void fft_3d_c2r_r2c_Test(struct Problem *P, const struct fft_plan_3d *plan, const int nFields);*/
#endif