C# Moving Window Filter Example

using System;
using System.IO;

using CenterSpace.NMath.Core;

namespace CenterSpace.NMath.Examples.CSharp
{
  /// <summary>
  /// A .NET example in C# showing smoothing using the Savitzky-Golay moving filter.
  /// </summary>
  class MovingWindowFilterExample
  {

    static void Main( string[] args )
    {
      int signalLength = 2000;

      Console.WriteLine();

      DoubleVector signal = NoisySignal( signalLength );
      double noisySignalVariance = NMathFunctions.Variance( signal );
      Console.WriteLine( "Noisy signal variance = " + noisySignalVariance );

      // Set up a moving average filter with an asymmetric window. A moving window
      // filter replaces each input data point with a linear combination of 
      // its surrounding points.  A filter is thus described by the number of
      // of points to the left and right of the input point and the coefficients
      // of the linear combination.
      // This filter will replace each input data point with the average of 
      // its value with the 4 values to the left, and the 5 values to the right.
      // Thus the coefficients to use are all equal to one over the number of points 
      // in the window, 10 in this case.
      int numberLeft = 4;
      int numberRight = 5;
      DoubleVector filterCoefficients = MovingWindowFilter.MovingAverageCoefficients( numberLeft, numberRight );
      var movingAverageFilter = new MovingWindowFilter( numberLeft, numberRight, filterCoefficients );

      // Filter the signal. Note that we must specify a boundary option. When
      // replacing the first input value, we dont have any points to the left, similarly, when 
      // replacing the last input value, we dont have any point to the right. The "PadWithZeros"
      // boundary option prepends "numberLeft" zeros and appends "numberRight" zeros to the 
      // input vector to deal with this.
      DoubleVector filteredSignal = movingAverageFilter.Filter( signal, MovingWindowFilter.BoundaryOption.PadWithZeros );
      double filteredSignalVariance = NMathFunctions.Variance( filteredSignal );
      Console.WriteLine( "Moving Average: filtered signal variance = " + filteredSignalVariance );

      // Set up a Savitzky-Golay filter. This is a smoothing filter that replaces input values with
      // the value of a polynomial of specified degree fit through the input value and its
      // surrounding points. A least squares algorithm is used to determine the fitting polynomial.
      int degree = 4;
      filterCoefficients = MovingWindowFilter.SavitzkyGolayCoefficients( numberLeft, numberRight, degree );
      MovingWindowFilter savitzkyGolayFilter = new MovingWindowFilter( numberLeft, numberRight, filterCoefficients );

      // Filter the signal. Here we use a different boundary option: "DoNotFilterBoundaryPoints". 
      // This option will not filter or replace the first "numberLeft" or last "numberRight"
      // values of the input signal.
      filteredSignal = savitzkyGolayFilter.Filter( signal, MovingWindowFilter.BoundaryOption.DoNotFilterBoundaryPoints );

      filteredSignalVariance = NMathFunctions.Variance( filteredSignal );
      Console.WriteLine( "Savitzky-Golay: filtered signal variance = " + filteredSignalVariance );
      Console.WriteLine();

      // If you are filtering lots of signals with the same length, it is more
      // economic to use the "Filter" method which allows you to specify
      // the vector to place the output filtered signal in. This avoids having
      // to allocate potentially large vectors on every call to "Filter".
      int numSignals = 10;
      DoubleMatrix noisySignals = NoisySignals( signalLength, numSignals );
      var y = new DoubleVector( signal.Length );
      for ( int i = 0; i < numSignals; i++ )
      {
        DoubleVector s = noisySignals.Col( i );
        Console.WriteLine( string.Format( "Noisy signal {0} variance = {1}", i, NMathFunctions.Variance( s ) ) );
        savitzkyGolayFilter.Filter( s, MovingWindowFilter.BoundaryOption.PadWithZeros, ref y );
        Console.WriteLine( string.Format( "Savitzky-Golay filtered signal {0} variance = {1}", i, NMathFunctions.Variance( y ) ) );
        Console.WriteLine();
      }

      Console.WriteLine();
      Console.WriteLine( "Press Enter Key" );
      Console.Read();
    }

    static DoubleVector NoisySignal( int length )
    {
      var rng = new RandGenNormal();
      var signal = new DoubleVector( length );
      for ( int i = 0; i < length; i++ )
      {
        signal[i] = Math.Cos( .2 * i ) + rng.Next();
      }
      return signal;
    }

    static DoubleMatrix NoisySignals( int rows, int columns )
    {
      var rng = new RandGenNormal();
      var signals = new DoubleMatrix( rows, columns );
      for ( int i = 0; i < rows; i++ )
      {
        for ( int j = 0; j < columns; j++ )
        {
          signals[i, j] = Math.Cos( .2 * i * j ) + rng.Next();
        }
      }
      return signals;
    }
  }
}