ParallelFFT.cs

Go to the documentation of this file.

using System;
using System.Threading;
using System.Diagnostics;

namespace FFT
{
   ///
   /// Pararelized version of IterativeFFT. 
   ///
   /// Allows to take advantage of multiprocessor computing enviroment. 
   ///
   public class ParallelFFT : IterativeFFT
   {
      ///
      /// Desired number of seperate threads.
      ///
      /// Can be assumed to be a power of 2. 
      ///
      protected int threadNo;

      ///
      /// Batch of FFT work. 
      ///
      protected struct FFTBatch
      {
         public Complex[] a;     ///< the working table
         public Direction d;     ///< the direction 
         public int fromLevel;   ///< strarting level
         public int toLevel;     ///< ending level
         public int fromRow;     ///< starting row
         public int toRow;       ///< ending row

         /// 
         /// Simple constructor.
         ///
         public FFTBatch(Complex[] a, Direction d,
                         int fromLevel, int toLevel, int fromRow, int toRow) 
         {
            this.a = a;
            this.d = d;
            this.fromLevel = fromLevel;
            this.toLevel = toLevel;
            this.fromRow = fromRow;
            this.toRow = toRow;
         }

         ///
         /// Performing the designated batch of work. 
         /// 
         public void Do()
         {
            int sgn = (d == Direction.Forward) ? (1) : (-1);

            for (int s = fromLevel; s < toLevel; s++) {
               int m = Utils.Exp2(s);
               
               for (int j = 0; j < m / 2; j++) {
                  Complex twiddler = Complex.FromPolar(1, sgn * 2 * j * Math.PI / m);
                  for (int k = fromRow+j; k < toRow; k += m) {
                     Complex t = twiddler * a[k + m / 2];
                     Complex u = a[k];
                     a[k] = u + t;
                     a[k + m / 2] = u - t;
                  }
               }
            }
         }
      } // FFTBatch


      ///
      /// Sets the desired number of seperate threads. 
      ///
      /// This number will be normalized, i.e. the closest 
      /// power of 2 greater than <tt>noThreads</tt> will be used.
      ///
      public ParallelFFT(int noThreads)
      {
         this.threadNo = Utils.Exp2(Utils.Log2(noThreads));
      }
      
      ///
      /// \copydoc FFT::IterativeFFT::FFT(Complex[], Direction)
      ///
      protected override Complex[] FFT(Complex[] a, Direction d)
      {
         Debug.Assert(a != null, "a is null");
         Debug.Assert(Utils.IsPowerOf2(a.Length),"length of a is not a power of 2!");

         int n = a.Length;
         int logN = Utils.Log2(n);
         int m = this.threadNo;
         int logM = Utils.Log2(m);
         
         a = BitReverseCopy(a);

         Thread[] thread = new Thread[m];

         //Dispatching the job to threads
         for (int i = 0; i < m; i++) {
            thread[i] = new Thread(new FFTBatch(a, d, 0, logN - logM, i*n/m, (i+1)*n/m).Do);
            thread[i].Start();
         }
         

         //Waiting for all threads to finish.
         for (int i = 0; i < m; i++) {
            thread[i].Join();
         }

         //Finishing the remaining job (in the main thread)
         new FFTBatch(a, d, logN - logM, logN + 1, 0, n).Do();
         
         if (d == Direction.Backward) {
            for (int i = 0; i < n; i++) {
               a[i] /= n;
            }
         }
         
         return a;

      } //FFT



      //
      // Executes a batch of FFT work. 
      //
      // \param param a FFTParams object specifying the work to be done. 
      //
      //protected void FFTBatch(object param)
      //{
      //   FFTParams p = (FFTParams) param;
      //   p.Do();
      //} //FFTBatch

      ///
      /// \copybrief FFT::IDFT::ToString()
      ///
      /// Returs string <tt>"Parallel"</tt> followed by the number 
      /// of threads (surrounded by regular paranthesis). 
      ///
      public override string ToString()
      {
         return "Parallel(" + threadNo + ")";
      }
      
   }//ParallelDFT
}

---