Low pass filter in Visual c++

[sdherzo]

Hi
How do i implement a low pass filter function in Visual c++ with the following function :

Y(n+1) = µ ( x- y(n) ) + y(n)

with µ = t/ ( t + constant )

x being the input
yn being the output
y(n+1) being the next output which has to be found

Thanks !
Regards,
S.D

[salman108]

U will do this by a recursive function call.

Since the function F is eating it's own output.

U create a function CalculateY(int y)

and in the function u replace the y(n) with a call to the function it self.

Enjoy

[cilu]

Well, I disagree. I don't know what a low passing filter is, but I'm pretty sure the recursive approach would have a low performance. You can use an iterative approach.

double µ = t / (t + constant);

// initialize y[0]
double yn = 0;
// initialize y[n+1] with y[n]
double yn1 = yn;
// run a while until a contition is satisfied
while(condition)
{
  yn1 = µ * (x - yn) + yn;
  yn = yn1;
}

// here you have the yn1 as y[n+1]

[sdherzo]

Thanks ! let me try it out ...

[TheCPUWizard]

A low pass filter is used to smooth out a value which varies over time. The BASS control on a stereo is one example.

Cilu's samply is on the money, but can also be looked at:
[code]
double µ = t / (t + constant);

double input[10000]; // ASSUME this is initialized somewhere or is read from an exernal source inside the loop.

double output[10000];

double output[0] = input[0];
// run a while until a contition is satisfied
while(for t=1; t<10000])
{
output[t] = µ * (input[t] - output[t-1]) + output[t-1];

OR

output[t] = µ * input[t] + (1-µ) * output[t-1];
}

[code]

[sdherzo]

Thanks guys for your replies !! That code works ... but unfortunately i have to do a bit more with the low pass... I am getting values from a A/D converter through a USB Port into my comp , these values have to be sent into the low pass ( into the code above - the variable (input ) & i have to plot the values i get from the output ( always - "output t ) in real time ... any suggestions there ??

I will be greatful for ur help !

Regards,

S.D
This is my code for normally giving out the values ... i real time! how do i get the output values of this prog ( "rawADvalue" - CurValue" ) get into my low pass & then plot the output ??

/
 *
 * Functions:  WinMain
 *             MainMessageHandler
 *             CleanUp
 *
 **********************************************************************/

/* Include files */

#include <windows.h>                    /* Compiler's include files's */
#include <string.h>                   
#include <stdio.h>
#include "..\cbw.h"                        /* Universal Library's include file */

#define BOARD_NUM      1                /* Number of A/D board as defined with InstaCal */
#define BUF_SIZE       30000            /* Size of buffer */
#define NUM_SECS       50               /* Number of secs to collect data */
#define ADRANGE        BIP5VOLTS        /* A/D voltage range */
#define TIMER_NUM      1                /* Windows timer used by this program */

HWND hWndMain;                          /* handle for main window */

LONG FAR PASCAL MainMessageHandler (HWND, UINT, WPARAM, LPARAM);



/************************************************************************
*
* Name:      WinMain
*
* Arguments: hInstance - the instance handle for the program
*            hPrevInstance - the class name of the application (not used)
*            CmndLine - command line was called with (not used)
*            CmndShow - indicates how to display window
*                                                                      
* This is the entry point to the program. It gets called by Windows
* to start the program up.  This routine creates the main window,
* initializes a timer, and then falls into the main Windows Get
* message/Dispatch message loop.
*                                                                      
************************************************************************/

int PASCAL
WinMain(HANDLE hInstance, HANDLE hPrevInstance, LPSTR CmdLine, int nCmdShow)
{
    MSG       msg;                      /* MSG structure to pass to windows proc */
    WNDCLASS  wndclass;
    char      *AppName;                 /* Name for the window */

    cbErrHandling (PRINTALL, STOPALL);  /* Set library's error handling */

    CmdLine = NULL;                     /* Not used */
    AppName = "WINCDEMO";               /* The name of this application */
    if(!hPrevInstance)
        {
        wndclass.style      = CS_HREDRAW | CS_VREDRAW;
        wndclass.lpfnWndProc= MainMessageHandler;
        wndclass.cbClsExtra = 0;
        wndclass.cbWndExtra = 0;
        wndclass.hInstance  = hInstance;
        wndclass.hIcon      = LoadIcon (hInstance, AppName);
        wndclass.hCursor    = LoadCursor (NULL, IDC_ARROW);
        wndclass.hbrBackground  = GetStockObject (WHITE_BRUSH);
        wndclass.lpszMenuName   = AppName;
        wndclass.lpszClassName  = AppName;
        RegisterClass (&wndclass);
        }

    /* create application's Main window                                    */
    hWndMain = CreateWindow (AppName,                  /* Window class name          */
                             "AInScan Background",
                             WS_OVERLAPPEDWINDOW,
                             CW_USEDEFAULT,           /* Use default X, Y            */
                             CW_USEDEFAULT,           /* Use default X, Y            */
                             GetSystemMetrics(SM_CXSIZE) * 12,   /* x - fit text      */
                             GetSystemMetrics(SM_CYSIZE) * 12,  /* y - fit text      */
                             NULL,                    /* Parent window's handle      */
                             NULL,                    /* Default to Class Menu       */
                             hInstance,               /* Instance of window          */
                             NULL);                   /* Create struct for WM_CREATE */


    if (hWndMain == NULL)
        {
        MessageBox(NULL, "Could not create window in WinMain", NULL, MB_ICONEXCLAMATION);
        return (1);
        }

    ShowWindow(hWndMain, nCmdShow);     /* Display main window      */
    UpdateWindow(hWndMain);

    /* Start a 500ms timer to update display */
    if(!SetTimer(hWndMain, TIMER_NUM, 500, NULL))
        {
        MessageBox(NULL, "Error starting Windows timer", NULL, MB_OK |
                   MB_ICONEXCLAMATION);
        return (1);
        }                          
        
    while(GetMessage(&msg, NULL, 0, 0)) /* Main message loop */
        {
        TranslateMessage(&msg);
        DispatchMessage(&msg);
        }

    UnregisterClass (AppName, hInstance);
    return (msg.wParam);
} 
                                                                        
                                                                        

/************************************************************************
*
* Name:      MainMessageHandler
*
* Arguments: hWnd - handle to the window
*            Message - message that needs to be handled
*            hParam - message parameter
*            lParam - message parameter
*                                                                      
* This is the message dispatcher that gets called whenever Windows
* sends a message to this window.  WinMain started up a timer that
* sends a message every 1/2 sec.  When the message (WM_TIMER)is received
* by this routine, it reads the A/D and stores the values in ADValues[].
* It also causes a screen update which eill automatically generate a 
* WM_PAINT message.  The WM_PAINT handler takes care of converting the
* raw A/D values to voltages and printing them in the Window.
*                                                                      
************************************************************************/

LONG FAR PASCAL
MainMessageHandler(HWND hWnd, UINT Message, WPARAM wParam, LPARAM lParam)
{
    HDC        hDC;                     /* handle for the display device  */
    PAINTSTRUCT ps;                     /* holds PAINT information        */
    TEXTMETRIC tm;                      /* holds TEXT information         */

    static HRGN hRgn;                   /* Rectangle region Handles       */
    static int  CharWidth, CharHeight;
    
    static unsigned short *ADValues;    /* Win32 pointer                  */
    static HGLOBAL  MemHandle32;        
    static short  CurStatus;
    static HGLOBAL  MemHandle; 
    static long CurCount, CurIndex;    
    int        x, y;
    long       Count, Rate;
    char       OutString[80], *str, *StatusStr;
    float      Voltage;
    unsigned short CurValue;

    switch (Message)                    /* Windows Message Loop           */
        {
        case WM_CREATE:
            CurStatus = RUNNING;
            CurIndex = 0l;
            CurCount = 0l;

            hDC = GetDC(hWndMain);      /* Get the device context for window */
            GetTextMetrics(hDC, &tm);   /* From context, get size of font */
            CharWidth = tm.tmAveCharWidth;
            CharHeight = tm.tmHeight + tm.tmExternalLeading;
            ReleaseDC(hWndMain, hDC);
            hRgn = CreateRectRgn(0,0, CharWidth*30, CharHeight*12);    

            /* Allocate a Windows buffer */
            MemHandle = cbWinBufAlloc ((long)BUF_SIZE);
            
            
            /* WIN32: This code is added to be able access more than  */
            /* 32K words of data.                                     */
            
            /* WIN32: allocate a local WIN32 buffer to hold A/D Data  */
             MemHandle32 = GlobalAlloc(GMEM_FIXED | GMEM_DDESHARE, BUF_SIZE *sizeof(short));
				
            /* Get a 32-bit pointer to the WIN32 buffer */
            ADValues = (unsigned short *)GlobalLock(MemHandle32);
                
            /* Start up background A/D scan */
            Count = BUF_SIZE;
            Rate = BUF_SIZE / NUM_SECS;
            cbAInScan (BOARD_NUM, 0, 3, Count, &Rate, ADRANGE, MemHandle,
                       BACKGROUND);
            break;

        case WM_TIMER:                  /* All Timer Events Processed Here */
            /* Get status of background A/D scan */
            cbGetStatus (BOARD_NUM, &CurStatus, &CurCount, &CurIndex, AIFUNCTION);
            InvalidateRgn(hWndMain, hRgn, FALSE);  /* Force screen update  */
            break;

        case WM_PAINT:                  /* Repaint client area of window */
            hDC = BeginPaint(hWndMain, &ps);     
            
            x = CharWidth * 2;          /* Position cursor within window */     
            y = CharHeight;             /* One line down and 2 chars in */
            str = "         A/D Info";  /* Print title */
            TextOut(hDC, x, y, str, strlen (str));

            if (CurStatus == RUNNING)   /* Print status of background scan */
                StatusStr = "RUNNING";
            else
                StatusStr = "COMPLETE";
            y += CharHeight;
            sprintf (OutString,"Status = %s    ", StatusStr);
            TextOut(hDC, x, y, OutString, strlen (OutString));
                                              
            y += CharHeight;            /* Print current count */
            sprintf (OutString,"Count = %ld    ", CurCount);
            TextOut(hDC, x, y, OutString, strlen (OutString));

            y += CharHeight;            /* Print current index */
            sprintf (OutString,"Index = %ld    ", CurIndex);
            TextOut(hDC, x, y, OutString, strlen (OutString));

            if (CurIndex >= 0)
               {
               
               /* WIN32: This code is added to be able access more than  */
               /* 32K words of data.                                     */
               
               /* Copy Data from memory to 32 bit memory buffer */
               cbWinBufToArray(MemHandle, &ADValues[CurIndex], CurIndex ,1);
               
               /* use cbAConvertData to convert the 16-bit values to 12-bit values */
               //cbAConvertData(BOARD_NUM, 1, &ADValues[CurIndex], NULL);
            
               CurValue = ADValues[CurIndex];
               }
            else
                CurValue = 0;                
                                    
            y += CharHeight*2;            /* Print raw data value */
            sprintf (OutString,"Raw A/D value = %u    ", CurValue);
            TextOut(hDC, x, y, OutString, strlen (OutString));
                                
            y += CharHeight;    /* Convert raw A/D to volts and print */
            cbToEngUnits(BOARD_NUM,ADRANGE,CurValue,&Voltage);
            sprintf (OutString,"Voltage = %.2f    ", Voltage);
            TextOut(hDC, x, y, OutString, strlen (OutString));

            SetTextAlign(hDC, TA_LEFT | TA_TOP);
            EndPaint(hWndMain, &ps);
            break;

        case WM_CLOSE:                      /* Close the window */
            KillTimer (hWnd, TIMER_NUM);    /* Stop the timer */
            cbStopBackground (BOARD_NUM, AIFUNCTION);   /* Stop background scan */
            DestroyWindow(hWnd);                          
            if (hWnd == hWndMain)
                PostQuitMessage(0);         /* Send message to Quit application */
            if (MemHandle)                  
                cbWinBufFree (MemHandle);   /* Free allocated memory */

            if (MemHandle32)
               GlobalFree (MemHandle32);

            break;

        default:
            return (DefWindowProc(hWnd, Message, wParam, lParam));
        }
    return (0l);
}

[djai]

Go to following link and you will find every thigh you want.
http://www.nr.com/

[djai]

There are many ways to do this if you use simple lo pass filter you have to maintain a buffer other wise you can do this very rapidly by FFT. If you want realtime simulations use FFT.

[sdherzo]

but fft only give ma a fouier transformation. How do i apply it to a low pass filter ??

Or cld anyone suggest me how to edit the code above in order to put the ow pass filter in it :-)

S.D:

[TheCPUWizard]

There are many ways to do this if you use simple lo pass filter you have to maintain a buffer other wise you can do this very rapidly by FFT.

For a single pole low -pass filter, all you need is the previous iterations filtered value as shown in the loop I posted earlier. The only change (which I thought would be obvious, is switching away from using an input array and reading the sensor (A/D)

[Martin O]

Heres a music-dsp site incase you didn't know of it:

http://www.musicdsp.org

[TheCPUWizard]

but fft only give ma a fouier transformation. How do i apply it to a low pass filter ??

You don't need it, but here is houw you would do it.

Take Samples[0..n], run FFT, clear all bins above the cutoff frequence, do inverse FFT, save just sample1.

Take Samples[1..n+1], run FFT, clear all bins above the cutoff frequence, do inverse FFT, save just sample1.

Take Samples[2..n+2], run FFT, clear all bins above the cutoff frequence, do inverse FFT, save just sample1.

Take Samples[3..n+3], run FFT, clear all bins above the cutoff frequence, do inverse FFT, save just sample1.

class Filter
{
public:
   filter(double t)
   { 
      µ = t / (t + constant);  
     last = 0;
    }
   double µ;   double last;

   double Update(double input)
   {
        double last = µ * (input - last)+ last;
        return last;
    }
};


int main()
{
   Filter filt(10.0);
   while (true)
    {
        double ADvalue = ReadAD();
        double filteredValue=  filter.Update(ADValue);
        cout << filteredValue;
        // Some Delay to match desired sampling rate;
     }
}