--- a
+++ b/Fourier.def
@@ -0,0 +1,167 @@
+DEFINITION MODULE Fourier;
+
+  (*------------------------------------------------------------------------*)
+  (* Module provides some routines for Fourier transformations              *)
+  (*                                                                        *)
+  (* It might be a good idea to have a look on Peter Moylans Fourier module *)
+  (* as well.                                                               *)
+  (*------------------------------------------------------------------------*)
+  (* Implementation : Michael Riedl                                         *)
+  (* Licence        : GNU Lesser General Public License (LGPL)              *)
+  (*------------------------------------------------------------------------*)
+
+  (* $Id: Fourier.def,v 1.3 2018/03/05 11:21:42 mriedl Exp mriedl $ *)
+
+PROCEDURE RealDFT(VAR A,B : ARRAY OF LONGREAL; (* Koeff.       *)
+                  VAR F   : ARRAY OF LONGREAL; (* Werte        *)
+                      N   : CARDINAL;          (* Anzahl Werte *)
+                      rev : BOOLEAN);
+
+          (*----------------------------------------------------------------*)
+          (* Berechnet die diskreten Fourierkoeffizenten der in F "uber-    *)
+          (* gebenen Funktionswerte in den Felder A und B. Wenn rev = wahr  *)
+          (* wird aud A,B der Vektor F zurueckberechnet                     *)
+          (*                                                                *)
+          (* Calculates in A & B the discret fourier coefficients of the    *)
+          (* function values given in array F, if rev = true the reverse    *)
+          (* calculation is performed.                                      *)
+          (*                                                                *)
+          (* This routine is mainly intended for test purpose               *)
+          (*----------------------------------------------------------------*)
+
+PROCEDURE CmplxFT(VAR X   : ARRAY OF LONGCOMPLEX;
+                  VAR Y   : ARRAY OF LONGCOMPLEX;
+                      N   : CARDINAL;
+                      rev : BOOLEAN);
+
+          (*----------------------------------------------------------------*)
+          (* If rev = false perform a forward complex Fourier transform     *)
+          (* as given by:                                                   *)
+          (*                                                                *)
+          (*   Y_i = \sum_{j=1}^N  X_j \exp{{2 \pi I J} \over N}            *)
+          (*                                                                *)
+          (*   \forall i \in [1..N]                                         *)
+          (*                                                                *)
+          (* If rev = true perform a backward complex Fourier transform     *)
+          (* as given by:                                                   *)
+          (*                                                                *)
+          (*   X_i = {1 \over N} \sum_{j=1}^N Y_j \exp{{2 \pi I J} \over N} *)
+          (*                                                                *)
+          (*   \forall i \in [1..N]                                         *)
+          (*                                                                *)
+          (* This routine is mainly intended for test purpose               *)
+          (*----------------------------------------------------------------*)
+
+PROCEDURE CmplxFFT(VAR X   : ARRAY OF LONGCOMPLEX;
+                       N   : CARDINAL;
+                       rev : BOOLEAN);
+
+          (*----------------------------------------------------------------*)
+          (* Komplex nach komplex FFT                                       *)
+          (*                                                                *)
+          (* X    : Array of data points (real parts) for analyse as input, *)
+          (*        overwritten by complex fourier cofficents on output.    *)
+          (* X    : Array of complex fourier coefficiens for synthesis as   *)
+          (*        input, overwritten by re-formed data points (real part) *)
+          (*        on outpt.                                               *)
+          (* n    : Number of data points in X, must be a power of two      *)
+          (* rev  : true  : Analyse (analysis)                              *)
+          (*        false : Synthese (synthesis)                            *)
+          (* For more explanaions on analysis and synthesis see remarks for *)
+          (* procedure CfstFT below as well.                                *)
+          (*                                                                *)
+          (* This routine seems to be the most precise one in this library  *)
+          (*----------------------------------------------------------------*)
+
+PROCEDURE CfstFT(VAR X   : ARRAY OF LONGCOMPLEX;
+                     N   : CARDINAL;
+                     rev : BOOLEAN);
+
+          (*----------------------------------------------------------------*)
+          (*                                                                *)
+          (* Calculates the finite Fourier transform of a complex periodic  *)
+          (* sequence X[0:n-1] whose period n must be a power of two.       *)
+          (*                                                                *)
+          (* Either the direct transform (eq. 1)                            *)
+          (*                                                                *)
+          (*   $X_j = \sum_{k=0}^{n-1} X_k exp({{-i 2\pi jk} \over n})$     *)
+          (*                                                                *)
+          (* with                                                           *)
+          (*                                                                *)
+          (*   $ j=0,1,\dots,n-1 $                                          *)
+          (*                                                                *)
+          (* or the unscaled inverse transform (eq. 2)                      *)
+          (*                                                                *)
+          (*   $\alpha_k = \sum_{j=0}^{n-1} X_j exp({{i2\pi jk} \over n})$  *)
+          (*                                                                *)
+          (* with                                                           *)
+          (*                                                                *)
+          (*   $ k=0,1,\dots,n-1 $                                          *)
+          (*                                                                *)
+          (* where X_k, \alpha_k and X_j are complex numbers, may be        *)
+          (* calculated. To ensure optimum use of storage, the same array   *)
+          (* is used for input and output and all intermediate calculations *)
+          (* are carried out in this array.                                 *)
+          (*                                                                *)
+          (*   N   : the number of values in array X, must be a power of 2  *)
+          (*   X   : One dimensional array of dimension (0:d),              *)
+          (*         where d >= n-1                                         *)
+          (*         For rev = true:                                        *)
+          (*           On entry, X(k) = X_k (k=0,n-1)                       *)
+          (*           On exit,  X(j) = X_j (j=0,n-1), as defined by (1).   *)
+          (*         For rev = false:                                       *)
+          (*           On entry, X(j) = X_j (j=0,n-1), as defined by (2).   *)
+          (*           On exit , X(k) = X_k (k=0,n-1)                       *)
+          (*   rev : true  : The direct  transform is calculated.           *)
+          (*         false : The inverse transform is calculated.           *)
+          (*                                                                *)
+          (* This procedure is a Modula-2 translation of CERN Lib CFSTFT    *)
+          (* provided by K.S. K\"olbig and H.-H. Umst\"atter                *)
+          (*                                                                *)
+          (* This routine is about a factor two faster than CmplxFFT        *)
+          (*----------------------------------------------------------------*)
+
+PROCEDURE GlFFT(VAR U     : ARRAY OF LONGCOMPLEX;
+                    N     : CARDINAL;
+                    Invrs : BOOLEAN);
+
+          (*----------------------------------------------------------------*)
+          (* Glassmans General N FFT                                        *)
+          (*                                                                *)
+          (* Berechne die Fouriertransformation des Feldes U oder die       *)
+          (* inverse Fouriertransformierte (Invrs=wahr). Die Anzahl der     *)
+          (* Datenpunkte muss keine Potenz von 2 sein.                      *)
+          (*                                                                *)
+          (* Calculates the Fourier transform of U or the inverse Fourier   *)
+          (* transformaiton if Invrs = true. Note that the length N does    *)
+          (* not need to be a power of 2 !                                  *)
+          (*                                                                *)
+          (* Input                                                          *)
+          (*                                                                *)
+          (*  U     : N-vector to be transformed                            *)
+          (*  N     : length of vector to be transformed                    *)
+          (*  Invrs : false if the forward transform is desired             *)
+          (*          true  if the inverse transform is desired             *)
+          (*                                                                *)
+          (* Output                                                         *)
+          (*                                                                *)
+          (*  U     : the DFT of U if Invrs is false or N times             *)
+          (*          the inverse dft of U if Invrs is true                 *)
+          (*                                                                *)
+          (* This routine is nearly as fast as CfstFM but less exact, but   *)
+          (* the difference in precison should be absolutly OK for normal   *)
+          (* usage. It differs from the "usual" Cooley und Tukey algorithm  *)
+          (* based routines.                                                *)
+          (*                                                                *)
+          (* It might be worth the exercise to improve the code !           *)
+          (*                                                                *)
+          (* This routine is directly derived form the Fortran code         *)
+          (* provided in the reference. A "lineraized" Fortran 90 version   *)
+          (* can also be provided by the author of this library             *)
+          (*----------------------------------------------------------------*)
+          (* Ref: Ferguson, Warren E. jr.; "A simple derivation of          *)
+          (*      Glassman's general n fast Fourier transform",             *)
+          (*      Comp. & Math. with Applc., 8/6, pp 401-311 (1982)         *)
+          (*----------------------------------------------------------------*)
+
+END Fourier.