Program to output Low and High pass filter impulse response coefficients for signal processing

Question

My issue that I was trying to solve was constructing a vector where the values were mirrored around the central element. The problem is here: https://stackoverflow.com/questions/58876916/construct-mirror-vector-around-the-centre-element-in-c

I solved the issue using the method in the accepted answer and am sharing my working code to hear any review suggestions or edit to make my code more efficient.

My code solution:

filters.cpp

#include <iostream>
#include "ErrorCode.h"
#include "variables.h"
#include <cmath>
#include <fstream>
#include <vector>
#include <algorithm>
#include <iterator>
#include <complex>

using namespace std;

/*
 *  Functions for calculating the coefficients of the selected filter.
 *  Options: Low Pass, High Pass, Band Pass, Band Reject.
 *
 *  (To do: put functions in separate file, use of headers and classes - constructors?)
 *
 */

vector <double> lowPass (vector<double> hL, double fcL, int M, int N, const double PI)
{
    std::vector<int> mArr;
    mArr.resize(N);
    int count;
    for(count = 0; count < N; count++)
        mArr[count] = count - M;
    int n;
    for(n = 0; n < M; n++)
    {
        hL[n] = (sin((fcL*mArr[n])-M))/((mArr[n]-M)*PI);
    }
    hL[M] = fcL/PI;

    for(n = 0; n < M; n++)
    {
        hL[N-n-1] = hL[n];
    }

    return hL;
}

vector <double> highPass (vector<double> hH, double fcH, int M, int N, const double PI)
{
    std::vector<int> mArr;
    mArr.resize(N);
    int count;
    for(count = 0; count < N; count++)
        mArr[count] = count - M;
    int n;
    for(n = 0; n < M; n++)
    {
        hH[n] = -((sin((fcH*mArr[n])-M))/(mArr[n]-(M*PI)));
    }
    hH[M] = 1+(fcH/PI);

    for(n = 0; n < M; n++)
    {
        hH[N-n-1] = hH[n];
    }

    return hH;
}

vector<double> bandPass (vector<double> bp, vector<double> hL, vector<double> hH, double fcL, double fcH, int M, const double PI)
{
    for(int n = 0; n < N; n++)
    {
        bp[n] = -(hH[n]) - hL[n];
        bp[M] = (fcH/PI) - (fcL/PI);
    }

    return bp;
}

vector<double> bandReject (vector<double> br, vector<double> hL, vector<double> hH, double fcL, double fcH, int M, const double PI)
{
    for(int n = 0; n < N; n++)
    {
        br[n] = hL[n] - (-(hH[n]));
        br[M] = 1 + (fcL/PI) - (fcH/PI);
    }

    return br;
}

/*
 *  Functions for calculating the coefficients of the selected window filter.
 *  Options: Rectangular, Hanning, Hamming, Blackman Harris
 *
 *  (To do: put functions in separate file, use of headers and classes - constructors?)
 *
 */

vector<double> hanning(vector<double> hann, int M, int N, const double PI)
{
    std::vector<int> mArr;
    mArr.resize(N);
    int count;
    for(count = 0; count < N; count++)
        mArr[count] = count - M;

    for(int n = 0; n < N; n++)
    {
        hann[n] = 0.5+(0.5*(cos((mArr[n]*PI)/M)));
    }
    return hann;
}

vector<double> hamming(vector<double> hamm, int N, const double PI)
{
    for(int n = 0; n < N; n++)
    {
        hamm[n] = 0.54 - (0.46*(cos(n*2*PI/N)));
    }
    return hamm;
}

vector<double> blackman (vector<double> bmhr, int N, const double PI)
{
    for(int n = 0; n < N; n++)
    {
        bmhr[n] = 0.42 - (0.5*(cos(2*PI*n/N))) + (0.08*cos(4*PI*n/N));
    }

    return bmhr;
}

/*
 *  Functions to calculate cut-off frequency for low and high
 */

double calcFcL (double f_cL, double f_s)
{
    double fcL = f_cL/f_s;
    return fcL;
}

double calcFcH (double f_cH, double f_s)
{
    double fcH = f_cH/f_s;
    return fcH;
}

/*
 * Check if cut off frequencies are lower than the sampling frequency
 * Check if filter length is odd
 */

ErrorCode checkFreq(double, double, double);
ErrorCode checkN(int);
ErrorCode checkNs(int);
ErrorCode checkF1(double);
ErrorCode checkF2(double);
ErrorCode checkA1(double);
ErrorCode checkA2(double);


/*
 *  User input for selecting filter type
 *  Options: Low Pass, High Pass, Band Pass, Band Reject
 *  Passes filter output of function by value
 *  Output coefficients to a csv file
 */

vector<double> filterType(vector<double> filter, double fcL, double fcH, int M, int N, const double PI)
{
    int choice;
    cout <<"Please select filter type:  \n";

    cout << "1: Low Pass\n"
            "2: High Pass\n"
            "3: Band Pass\n"
            "4: Band Reject\n";

    cout << "Enter your selection (1, 2, 3 or 4): ";
    cin >> choice;

    switch(choice)
    {
            case 1:
                {
                vector<double> hL;
                hL.resize(N);
                hL = lowPass(hL, fcL, M, N, PI);
                ofstream lpfile ("lowpass.csv");
                if (lpfile.is_open())
                {
                    for(int count = 0; count < N; count++){
                        lpfile << hL[count] << "," ;
                    }
                    lpfile.close();
                }
                else cout << "Unable to open file";
                filter = hL;
                break;
                }


            case 2:
                {
                vector<double> hH;
                hH.resize(N);
                hH = highPass(hH, fcH, M, N, PI);
                ofstream hpfile ("highpass.csv");
                if (hpfile.is_open())
                {
                for(int count = 0; count < N; count++){
                    hpfile << hH[count] << "," ;
                    }
                hpfile.close();
                }
                else cout << "Unable to open file";
                filter = hH;
                break;
                }


            case 3:
                {
                vector<double> hL;
                hL.resize(N);
                vector<double> hH;
                hH.resize(N);
                vector<double> bp;
                bp.resize(N);
                hL = lowPass(hL, fcL, M, N, PI);
                hH = highPass(hH, fcH, M, N, PI);

                bp = bandPass(bp, hL, hH, fcL, fcH, M, PI);
                ofstream bpfile ("bandpass.csv");
                if (bpfile.is_open())
                {
                for(int count = 0; count < N; count++){
                    bpfile << bp[count] << "," ;
                    }
                bpfile.close();
                }
                else cout << "Unable to open file";
                filter  = bp;
                break;
                }


            case 4:
                {
                vector<double> hL;
                hL.resize(N);
                vector<double> hH;
                hH.resize(N);
                vector<double> br;
                br.resize(N);
                hL = lowPass(hL, fcL, M, N, PI);
                hH = highPass(hH, fcH, M, N, PI);
                br = bandReject(br, hL, hH, fcL, fcH, M, PI);
                ofstream brfile ("bandreject.csv");
                if (brfile.is_open())
                {
                for(int count = 0; count < N; count++){
                    brfile << br[count] << "," ;
                    }
                brfile.close();
                }
                else cout << "Unable to open file";
                filter = br;
                break;
                }


            default:
                break;


    }
    return filter;
}

/*
 *  User input for selecting window type
 *  Options: Rectangular, Hanning, Hamming, Blackman Harris
 *  Passes window output of function by value
 *  Output coefficients to a csv file
 */

vector<double> windowType(vector<double> window, int M, int N, const double PI)
{
    int choice;
    cout <<"Please select Window type:  \n";

    cout << "1: Rectangular\n"
            "2: Hanning\n"
            "3: Hamming\n"
            "4: Blackman Harris\n";

    cout << "Enter your selection (1, 2, 3 or 4): ";
    cin >> choice;

    switch(choice)
    {
            case 1:
                {
                    vector<double> rect;
                    rect.assign (N, 1);
                    ofstream rectfile ("rectangular.csv");
                    if (rectfile.is_open())
                    {
                        for(int count = 0; count < N; count++){
                            rectfile << rect[count] << "," ;
                            }
                        rectfile.close();
                    }
                    else cout << "Unable to open file";
                    window = rect;
                    break;
                }

            case 2:
                {
                    vector<double> hann = {};
                    hann.resize(N);
                    hann = hanning(hann, M, N, PI);
                    ofstream hannfile ("hanning.csv");
                    if (hannfile.is_open())
                    {
                        for(int count = 0; count < N; count++){
                            hannfile << hann[count] << "," ;
                            }
                        hannfile.close();
                    }
                    else cout << "Unable to open file";
                    window = hann;
                    break;
                }

            case 3:
                {
                    vector<double> hamm = {};
                    hamm.resize(N);
                    hamm = hamming(hamm, N, PI);
                    ofstream hammfile ("hamming.csv");
                    if (hammfile.is_open())
                    {
                        for(int count = 0; count < N; count++){
                            hammfile << hamm[count] << "," ;
                            }
                        hammfile.close();
                    }
                    else cout << "Unable to open file";
                    window = hamm;
                    break;
                }

            case 4:
                {
                    vector<double> bmhr = {};
                    bmhr.resize(N);
                    bmhr = blackman(bmhr, N, PI);
                    ofstream bmhrfile ("blackman.csv");
                    if (bmhrfile.is_open())
                    {
                        for(int count = 0; count < N; count++){
                            bmhrfile << bmhr[count] << "," ;
                            }
                        bmhrfile.close();
                    }
                    else cout << "Unable to open file";
                    window = bmhr;
                    break;
                }

    }
    return window;
}

/*
 *  Combine window filter with filter to generate FIR coefficients
 */

vector<double> finiteIR(vector<double> fir, vector<double> filter, vector<double> window, int N)
{
    for(int n = 0; n < N; n++)
    {
        fir[n] = filter[n]*window[n];
    }
    ofstream firfile ("FIR.csv");
    if (firfile.is_open())
    {
        for(int count = 0; count < N; count++){
            firfile << fir[count] << "," ;
        }
        firfile.close();
    }
        else cout << "Unable to open file";
    return fir;
}

vector<complex<double> > genSignal(vector<complex<double> > signal, double f_s, int N_s, double f1, double f2, double a1, double a2, const double PI)
{
    complex<double> i = -1;
    i = sqrt(i);
    double tP = 1/f_s;
    vector<double> t;
    t.resize(N_s);
    for(int n = 0; n < N_s - 1; n++)
    {
        t[n] = n*tP;
    }
    vector<complex<double> > x = {};
    double x1Real;
    double x1Imag;
    double x2Real;
    double x2Imag;
    for (int n = 0; n < N_s; n++)
    {
      x1Real = real(a1*exp(i*(2*PI*f1*t[n])));
      x1Imag = imag(a1*exp(i*(2*PI*f1*t[n])));
      x2Real = real(a2*exp(i*(2*PI*f2*t[n])));
      x2Imag = imag(a2*exp(i*(2*PI*f2*t[n])));
      double xreal = x1Real + x2Real;
      double ximag = x1Imag + x2Imag;
      complex<double> iNum(xreal, ximag);
      x.push_back(iNum);
    }

    ofstream sigfile ("signal.csv");
    if (sigfile.is_open())
    {
        for(int count = 0; count < N_s; count++)
        {
            sigfile << real(x[count]) << "," ;
            sigfile << imag(x[count]) << ",\n" ;
        }
        sigfile.close();
    }
    else cout << "Unable to open file";

    signal = x;

    return signal;

}

// ----------------------------------------------------------------------------------------------------------------------------------

int main()
{
    double f_cL;
    double f_cH;
    cout << "Enter Low Cut-Off frequency (Hz) [Set to 0 for the highpass filter] : ";
    cin >> f_cL;
    cout << "Enter High Cut-Off frequency (Hz) [Set to 0 for the lowpass filter] : ";
    cin >> f_cH;

    // ------------------------------------------------------------------------------

    cout << "Enter sampling frequency (Hz): ";
    cin >> f_s;
    if (checkFreq(f_s, f_cL, f_cH) == ErrorCode::FAIL)
    {
        std::cout << "You entered a frequency less than cut-off!\n";
        cout << "Fail";
        return 0;
    }

    // ------------------------------------------------------------------------------

    cout << "Enter Order of the filter (filter length): ";
    cin >> N;
    if (checkN(N) == ErrorCode::FAIL)
    {
        cout << "You entered an invalid number!\n";
        cout << "Fail";
        return 0;
    }

    // ------------------------------------------------------------------------------

    double fcL;
    double fcH;
    fcL = calcFcL(f_cL, f_s);
    cout << "FcL is : " << fcL << endl;
    fcH = calcFcH(f_cH, f_s);
    cout << "FcH is: " << fcH << endl;

    // ------------------------------------------------------------------------------

    int M = (N-1)/2;
    vector<double> filter;
    filter.resize(N);
    filter = filterType(filter, fcL, fcH, M, N, PI);

    vector<double> window;
    window.resize(N);
    window = windowType(window, M, N, PI);

    vector<double> fir;
    fir.resize(N);
    fir = finiteIR(fir, filter, window, N);

    // ------------------------------------------------------------------------------

    cout << "\n\n" << "Generate Signal...\n";
    cout << "Enter signal sample length: ";
    cin >> N_s;
    if (checkNs(N_s) == ErrorCode::FAIL)
    {
        cout << "You entered an invalid number!\n";
        cout << "Fail";
        return 0;
    }

    cout << "Enter frequency 1 (Hz): ";
    cin >> f1;
    if (checkF1(f1) == ErrorCode::FAIL)
    {
        cout << "You entered an invalid frequency!\n";
        cout << "Fail";
        return 0;
    }

    cout << "Enter frequency 2 (Hz): ";
    cin >> f2;
    if (checkF2(f2) == ErrorCode::FAIL)
    {
        cout << "You entered an invalid frequency!\n";
        cout << "Fail";
        return 0;
    }

    cout << "Enter Amplitude 1: ";
    cin >> a1;
    if (checkA1(a1) == ErrorCode::FAIL)
    {
        cout << "You entered an invalid frequency!\n";
        cout << "Fail";
        return 0;
    }

    cout << "Enter Amplitude 2 (Hz): ";
    cin >> a2;
    if (checkA2(a2) == ErrorCode::FAIL)
    {
        cout << "You entered an invalid frequency!\n";
        cout << "Fail";
        return 0;
    }


    vector<complex<double> > signal;
    signal = genSignal(signal, f_s, N_s, f1, f2, a1, a2, PI);

    return 0;
}

error.cpp

#include <iostream>
#include "ErrorCode.h"

ErrorCode checkFreq(double value, double value2, double value3)
{
    // if value is less than cut-off
    if (value < value2 || value < value3)
       // early return an error code
        return ErrorCode::FAIL;

    // Success
    return ErrorCode::SUCCESS;
}

ErrorCode checkN(int value)
{
    // if value is even
    if(value % 2 == 0 || value < 0)
      // early return an error code
        return ErrorCode::FAIL;

    // Success
    return ErrorCode::SUCCESS;
}

ErrorCode checkNs(int value)
{
    // if value is negative
    if(value < 0)
      // early return an error code
        return ErrorCode::FAIL;

    // Success
    return ErrorCode::SUCCESS;
}

ErrorCode checkF1(double value)
{
    // if value is negative
    if(value < 0)
      // early return an error code
        return ErrorCode::FAIL;

    // Success
    return ErrorCode::SUCCESS;
}

ErrorCode checkF2(double value)
{
    // if value is negative
    if(value < 0)
      // early return an error code
        return ErrorCode::FAIL;

    // Success
    return ErrorCode::SUCCESS;
}

ErrorCode checkA1(double value)
{
    // if value is negative
    if(value < 0)
      // early return an error code
        return ErrorCode::FAIL;

    // Success
    return ErrorCode::SUCCESS;
}

ErrorCode checkA2(double value)
{
    // if value is negative
    if(value < 0)
      // early return an error code
        return ErrorCode::FAIL;

    // Success
    return ErrorCode::SUCCESS;
}

ErrorCode.h

#ifndef ERRORCODE_H_INCLUDED
#define ERRORCODE_H_INCLUDED

enum class ErrorCode
{
    SUCCESS = 0,
    FAIL = -1
};

#endif // ERRORCODE_H_INCLUDED

variables.h

#ifndef VARIABLES_H_INCLUDED
#define VARIABLES_H_INCLUDED


double f_c = 0.0, f_s = 0.0, fc = 0.0, w_c = 0.0;
double N_s = 0.0, f1 = 0.0, f2 = 0.0, a1 = 0.0, a2 = 0.0;
int N = 0;
const double PI  = 3.141592653589793238463;

#endif // VARIABLES_H_INCLUDED

Answer 1

Use good names

First impression is that the code is really hard to read, because the identifiers are terse and unhelpful. Take our first function:

vector <double> lowPass (vector<double> hL, double fcL, int M, int N, const double PI)

There's not even a comment to say what any of those values represent, or even why some are in all-caps like macros. (And PI is a terrible choice of name, as it looks a lot like the constant π).

And then we have uncommented lines such as

    hL[n] = (sin((fcL*mArr[n])-M))/((mArr[n]-M)*PI);

Even a little bit of whitespace would help readers parse this a bit more easily.

The whole thing reeks of something that's been hacked on until it produces the desired result, with no consideration given to future changes we might want to make.

Another one - without looking at the implementation, what does the function checkF2() do, and how is it different from checkF1()?

---

Don't import the whole standard namespace

using namespace std; is harmful, particularly so at file scope where it affects every function. The whole point of namespaces is to segregate identifiers so that upgrading your Standard Library (or just adding another include) won't suddenly change the meaning of your program.

---

Pass large input parameters by reference-to-const

bandPass() and bandReject both have arguments std::vector<double> hL, std::vector<double> hH which are never written to. We can avoid copying this if we accept std::vector<double> const& instead.

---

Validate external inputs

It's wrong to read from choice here, because we've not yet checked whether any value was assigned to it:

int choice;
cin >> choice;

switch(choice)

As a simplest possible fix, we could fail immediately:

    int choice;
    std::cin >> choice;
    if (!std::cin) {
        std::cerr << "Read failed.\n";
        return EXIT_FAILURE;
    }

    switch(choice)

A more user-friendly approach would use a function to prompt the user and retry when input cannot be parsed.

But it's probably easier for you (and for users, who might want to execute this from a script) to accept command-line arguments rather than pursue the interactive route.

---

Don't hard-code file names

The code writes to "lowpass.csv", "highpass.csv", "bandpass.csv" or "bandreject.csv" unconditionally, even if that overwrites the results from previous runs.

I'd argue that it's better to write to standard output instead, so that it can be directed wherever the user wants (perhaps piped into another program, for instance). Again, that makes scripted use easier.

---

Avoid global variables

The variables in variables.h are global, making it hard to reason about parts of the program in isolation. If they are related, consider creating one or more structures and/or classes to work with them, and pass them to the functions that need them.

---

Enable more warnings

It seems you're not enabling some important compiler warnings:

232446.cpp: In function ‘int main()’:
232446.cpp:533:17: warning: conversion from ‘double’ to ‘int’ may change value [-Wfloat-conversion]
  533 |     if (checkNs(N_s) == ErrorCode::FAIL)
      |                 ^~~
232446.cpp:578:37: warning: conversion from ‘double’ to ‘int’ may change value [-Wfloat-conversion]
  578 |     signal = genSignal(signal, f_s, N_s, f1, f2, a1, a2, PI);
      |                                     ^~~