signal_processing

#include <ipp.h>
#include <cmath>
#include <iostream>
#include <fstream>
#include <stdexcept>
#include <memory>

// Signal parameter.
Ipp32s signal_frequency [5] = { 1, 5, 10, 200, 500 };   // in Hz.
Ipp32s signal_amplitude = 20;                           // in levels.
Ipp32s noise_amplitude = 10;                            // in levels.
const Ipp32s sampling_frequency = 44100;                // in Hz.

// Free functor.
struct freeIpps {
    inline void operator()(void *_p)
    {
        if(!_p) return;

        ippsFree(_p);
    }
};

// Type references.
typedef std::unique_ptr<Ipp32f[], freeIpps> Ipp32f_ptr;
typedef std::unique_ptr<Ipp8u[], freeIpps> Ipp8u_ptr;
typedef std::unique_ptr<IppsFFTSpec_R_32f, freeIpps> IppsFFTSpec_R_32f_ptr;

// Signal data.
const Ipp32f time_step = 1;
const Ipp32f log_2 = std::log(2);
Ipp32f time_line[sampling_frequency];

// Check function.
inline bool ispwr_of_2(Ipp32s x)
{
    return (x && !(x & (x - 1)));
}

inline Ipp32s topwr_of_2(Ipp32s x)
{
    if(x > 0)
    {
        if(x & (x - 1))
        {
            --x;

            x |= x >> 1;
            x |= x >> 2;
            x |= x >> 4;
            x |= x >> 8;
            x |= x >> 16;

            ++x;
        }
    }

    return  x;
}

void init_time()
{
    time_line[0] = time_step;

    for(Ipp32u i = 1; i < sampling_frequency; ++i)
        time_line[i] = time_line[i-1] + time_step;
}

Ipp32f_ptr signal_gen(Ipp32s fi, Ipp32s &data_sz)
{
    // The length of the vector transformed by the FFT must be a power of 2.
    data_sz = topwr_of_2(sampling_frequency);
    Ipp32f_ptr data(ippsMalloc_32f(data_sz));

    Ipp32f y = IPP_2PI * signal_frequency[fi];
    y /= sampling_frequency;

    // Generate sine wave with sampling_frequency size.
    for(Ipp32u i = 0; i < sampling_frequency; ++i) {
        data[i] = y * time_line[i];
        data[i] = std::sin(data[i]);
        data[i] *= signal_amplitude;
    }

    Ipp32s zero_begin = sampling_frequency;
    Ipp32s zero_end = data_sz - zero_begin;

    // Zero-pad the input signals (add zeros to the end so that at least half of the wave is "blank").
    IppStatus status = ippsZero_32f(data.get() + zero_begin, zero_end);
    if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));

    std::cout << "The sine wave is generated. The length is equal to: " << sampling_frequency
              << "\nZeroed part of the signal is in the range: [ " << zero_begin << " ] [ " << zero_begin + zero_end << " ]\n";

    return data;
}

void do_processing()
{
    while(true)
    for(Ipp32u i = 0; i < 5;  i = (i+1 == 5)?(0):(i+1))
    {
        Ipp32s fftSz = 0; // Length of data with zeros.
        Ipp32f_ptr pSrc = signal_gen(i, fftSz);
        Ipp32f_ptr pDst(ippsMalloc_32f(fftSz));

        Ipp32s pSpecSize = 0;
        Ipp32s pSpecBuffSize = 0;
        Ipp32s pBuffSize = 0;
        Ipp32s order = std::log(fftSz)/log_2;

        IppStatus status = ippsFFTGetSize_R_32f(order, IPP_FFT_NODIV_BY_ANY, ippAlgHintFast, &pSpecSize, &pSpecBuffSize, &pBuffSize);
        if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));

        Ipp8u_ptr pBuffer(ippsMalloc_8u(pBuffSize));
        Ipp8u_ptr pSpec(ippsMalloc_8u(pSpecSize));
        Ipp8u_ptr pSpecBuffer(ippsMalloc_8u(pSpecBuffSize));

        IppsFFTSpec_R_32f *pFFTSpec = nullptr;

        status = ippsFFTInit_R_32f(&pFFTSpec, order, IPP_FFT_NODIV_BY_ANY, ippAlgHintFast, pSpec.get(), pSpecBuffer.get());
        if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));

        status = ippsFFTFwd_RToPack_32f(pSrc.get(), pDst.get(), pFFTSpec, pBuffer.get());
        if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));

        Ipp32fc *pDstC = (Ipp32fc*)(pDst.get());
        Ipp32s realSz = fftSz/2;          // Exclude zero part.
        Ipp32s magSz = realSz/2;          // Convert length to complex length.

        // Extract magniutde from input signal.
        status = ippsMagnitude_32fc(pDstC, pDst.get(), realSz);
        if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));

        Ipp32f lvl = 2.0/Ipp32f(realSz);
        Ipp32f B = 19.93156857;

        // Get amplitude in levels.
        status = ippsMulC_32f(pDst.get(), lvl, pDst.get(), magSz);
        if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));

        // Convert level to dB.
        status = ippsLn_32f(pDst.get(), pDst.get(), magSz);
        if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));

        status = ippsMulC_32f(pDst.get(), B, pDst.get(), magSz);
        if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));
    }
}

int main()
try{
    IppStatus status = ippInit();
    if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));

    std::cout << "The IPP library was successfully initialized.\n";

    status = ippSetAffinity(ippAffinityCompactFineHT, 0);
    if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));

    std::cout << "OpenMP succesfully bind all threads to OS processors. \n";

    init_time();
    std::cout << "Timeline was initialized.\n";

    do_processing();
    std::cout << "Signal processing has completed successfully.\n";

    return 0;
}
catch(std::runtime_error &e)
{
    std::cerr << e.what() << '\n';
    return 1;
}
catch(...)
{
    std::cerr << "Unhandled exception.\n";
    return 2;
}