Generate C Code from Spectral Analysis Algorithm and Verify Using Parameterized Equivalence Tests

Define MATLAB Entry-Point Function

Define an entry-point function powerSpectrumPeaks that calculates the Fourier transform of a time-domain signal sampled at frequency fs. The function then squares the absolute value of the Fourier transform to obtain the power spectrum of the signal in the frequency domain. Finally, it returns the frequencies at which the peaks of the power spectrum, which have a pre-specified minimum prominence, are located.

type powerSpectrumPeaks.m

function peaks = powerSpectrumPeaks(signal,fs,prominence) %#codegen
arguments
    signal (1,:) double {coder.mustBeComplex}
    fs (1,1) double
    prominence (1,1) double
end

n = length(signal);
signalIndices = 1:n;

y = fft(signal);
power = abs(y).^2/n;

peakIndices = signalIndices(islocalmax(power,MinProminence=prominence));
peaks = (peakIndices - 1)*fs/n;
end

Define a complex time-domain signal x that has two frequency components at 50 Hz and 300 Hz. Here fs is the sampling frequency and L is the length of the signal.

fs = 1000;                                
T = 1/fs;                    
L = 1500;             
t = (0:L-1)*T;        
x = 0.7*exp(2*1i*pi*50*t) + 1.3*exp(2*1i*pi*300*t);  

To filter out the noise from the power spectrum, define a minimum peak prominence that the islocalmax function uses.

prominence = 0.001;

Call the entry-point function powerSpectrumPeaks with these inputs. The function returns a 2-by-1 array that contains the two component frequencies of the input signal.

powerSpectrumPeaks(x,fs,prominence)

ans = 1×2

    50   300

Check For Run-Time Issues by Generating and Running MEX

Use the codegen command to generate a MEX for the powerSpectrumPeaks function. The arguments block in the function body specifies to the code generator that the input signal is an unbounded variable-size row vector of complex double values. Similarly, the sampling frequency and prominence are scalar doubles. Specify the name of the code generation folder as mexcode using the -d option.

codegen powerSpectrumPeaks -d mexcode

Code generation successful.

Run the generated MEX with the same inputs. The output of the MEX matches that of MATLAB execution.

powerSpectrumPeaks_mex(x,fs,prominence)

ans = 1×2

    50   300

Author Parameterized Equivalence Test

Author a class tEquivalence that contains:

This class provides a flexible interface that you can use to test both the generated MEX and the generated standalone code across a range of inputs.

type tEquivalence.m

classdef tEquivalence < matlabtest.coder.TestCase
    properties
        Target {mustBeMember(Target,["mex","sil","pil"])} = "mex"
    end

    properties (TestParameter)
        SignalFunction = {@(t) 0.7*exp(2*1i*pi*50*t) + 1.3*exp(2*1i*pi*300*t), ...
            @(t) 0.5*exp(2*1i*pi*100*t) + 2*exp(2*1i*pi*500*t)}
        LengthOfSignal = {1500,2000}
        SamplingFrequency = {1000,3000}
        Prominence = {0.001,0.01}
    end

    methods
        function testCase = tEquivalence(target)
            testCase.Target = target;
        end
    end

    methods (Test)
        function tMexForPowerSpectrumPeaks(testCase,SignalFunction,LengthOfSignal,SamplingFrequency,Prominence)
            T = 1/SamplingFrequency;
            t = (0:LengthOfSignal-1)*T;
            signal = SignalFunction(t);

            mexPath = fullfile(pwd,strcat("powerSpectrumPeaks_", ...
                testCase.Target,".",mexext));
            executionResults = execute(testCase,mexPath,Inputs= ...
                {signal,SamplingFrequency,Prominence});
            verifyExecutionMatchesMATLAB(testCase,executionResults)
        end
    end
end

Validate Generated MEX by Using Equivalence Testing

Create a test case for the MEX target.

testcase = tEquivalence("mex");

Run this test. The test passes, indicating that the generated MEX is functionally equivalent to your source MATLAB code for an exhaustive combinations of the inputs you provided in the class definition.

run(testcase)

Running tEquivalence
.......... ......
Done tEquivalence
__________

ans = 
  1×16 TestResult array with properties:

    Name
    Passed
    Failed
    Incomplete
    Duration
    Details

Totals:
   16 Passed, 0 Failed, 0 Incomplete.
   0.38024 seconds testing time.

Generate Static Library and SIL Interface

Use the codegen command to generate a static library and a SIL interface for the powerSpectrumPeaks function. Specify the name of the code generation folder as silcode using the -d option.

cfg = coder.config('lib', 'ecoder', true);
cfg.VerificationMode = 'SIL';
codegen -config cfg powerSpectrumPeaks -d silcode -report

Code generation successful: View report

Run the generated SIL MEX with sample inputs. The output of the SIL MEX matches that of MATLAB execution.

powerSpectrumPeaks_sil(x,fs,prominence)

### Starting SIL execution for 'powerSpectrumPeaks'
    To terminate execution: clear powerSpectrumPeaks_sil

ans = 1×2

    50   300

Validate Generated Static Library By Using Equivalence Testing

Create a test case for the SIL target using the input parameters you used when running the SIL MEX.

testcaseSil = tEquivalence("sil");

Run this test. The test passes, indicating that the generated static library is functionally equivalent to your source MATLAB code for an exhaustive combinations of the inputs you provided in the class definition.

run(testcaseSil)

Running tEquivalence
.......... ......
Done tEquivalence
__________

ans = 
  1×16 TestResult array with properties:

    Name
    Passed
    Failed
    Incomplete
    Duration
    Details

Totals:
   16 Passed, 0 Failed, 0 Incomplete.
   15.9465 seconds testing time.

Finally, clear the SIL MEX function from memory.

clear powerSpectrumPeaks_sil