meanfreq

Generate 1024 samples of a chirp sampled at 1024 kHz. The chirp has an initial frequency of 50 kHz and reaches 100 kHz at the end of the sampling. Add white Gaussian noise such that the signal-to-noise ratio is 40 dB.

nSamp = 1024;
Fs = 1024e3;
SNR = 40;

t = (0:nSamp-1)'/Fs;
x0 = chirp(t,50e3,nSamp/Fs,100e3);

x1 = x0 + randn(size(x0))*std(x0)/db2mag(SNR);

Estimate the mean frequency of the chirp. Plot the power spectral density (PSD) and annotate the mean frequency.

meanfreq(x1,Fs);

Generate another chirp. Specify an initial frequency of 200 kHz, a final frequency of 300 kHz, and an amplitude that is twice that of the first signal. Add white Gaussian noise.

x2 = 2*chirp(t,200e3,nSamp/Fs,300e3);
x2 = x2 + randn(size(x2))*std(x2)/db2mag(SNR);

Concatenate the chirps to produce a two-channel signal. Estimate the mean frequency of each channel.

y = meanfreq([x1 x2],Fs)

y = 1×2
105 ×

    0.7503    2.4999

Plot the PSDs of the two channels and annotate their mean frequencies.

meanfreq([x1 x2],Fs);

Add the two channels to form a new signal. Plot the PSD and annotate the mean frequency.

meanfreq(x1+x2,Fs);

Generate 1024 samples of a 100.123 kHz sinusoid sampled at 1024 kHz. Add white Gaussian noise such that the signal-to-noise ratio is 40 dB.

nSamp = 1024;
Fs = 1024e3;
SNR = 40;

t = (0:nSamp-1)'/Fs;
x0 = sin(2*pi*t*100.123e3);

x1 = x0 + randn(size(x0))*std(x0)/db2mag(SNR);

Use the periodogram function to compute the power spectral density (PSD) of the signal. Specify a Kaiser window with the same length as the signal and a shape factor of 38. Estimate the mean frequency of the signal and annotate it on a plot of the PSD.

[Pxx,f] = periodogram(x1,kaiser(nSamp,38),[],Fs);

meanfreq(Pxx,f);

Generate another sinusoid, this one with a frequency of 257.321 kHz and an amplitude that is twice that of the first sinusoid. Add white noise.

x2 = 2*sin(2*pi*t*257.321e3);
x2 = x2 + randn(size(x2))*std(x2)/db2mag(SNR);

Concatenate the sinusoids to produce a two-channel signal. Estimate the PSD of each channel and use the result to determine the mean frequency.

[Pyy,f] = periodogram([x1 x2],kaiser(nSamp,38),[],Fs);

y = meanfreq(Pyy,f)

y = 1×2
105 ×

    1.0013    2.5732

Annotate the mean frequencies of the two channels on a plot of the PSDs.

meanfreq(Pyy,f);

Add the two channels to form a new signal. Estimate the PSD and annotate the mean frequency.

[Pzz,f] = periodogram(x1+x2,kaiser(nSamp,38),[],Fs);

meanfreq(Pzz,f);

Generate a signal whose PSD resembles the frequency response of an 88th-order bandpass FIR filter with normalized cutoff frequencies $$0.25\pi$$ rad/sample and $$0.45\pi$$ rad/sample.

d = fir1(88,[0.25 0.45]);

Compute the mean frequency of the signal between $$0.3\pi$$ rad/sample and $$0.6\pi$$ rad/sample. Plot the PSD and annotate the mean frequency and measurement interval.

meanfreq(d,[],[0.3 0.6]*pi);

Output the mean frequency and the band power of the measurement interval. Specifying a sample rate of $$2\pi$$ is equivalent to leaving the rate unset.

[mnf,power] = meanfreq(d,2*pi,[0.3 0.6]*pi);

fprintf("Mean = %.3f*pi, power = %.1f%% of total \n", ...
    mnf/pi,power/bandpower(d)*100)

Mean = 0.373*pi, power = 75.6% of total 

Add a second channel with normalized cutoff frequencies $$0.5\pi$$ rad/sample and $$0.8\pi$$ rad/sample and an amplitude that is one-tenth that of the first channel.

d = [d;fir1(88,[0.5 0.8])/10]';

Compute the mean frequency of the signal between $$0.3\pi$$ rad/sample and $$0.9\pi$$ rad/sample. Plot the PSD and annotate the mean frequency of each channel and the measurement interval.

meanfreq(d,[],[0.3 0.9]*pi);

Output the mean frequency of each channel. Divide by $$\pi$$.

mnf = meanfreq(d,[],[0.3 0.9]*pi)/pi

mnf = 1×2

    0.3730    0.6500

Since R2026b

Generate a timetable that comprises a three-second quadratic swept-frequency cosine signal with a sample rate of 1 kHz. The frequency sweeping ranges from 150 Hz to 350 Hz.

fs = 1000;
t = (1/fs:1/fs:3)';
f1 = 150;
f2 = 350;
sg = chirp(t,f1,t(end),f2,"quadratic") + 1e-3*randn(size(t));

sgTT = timetable(sg,SampleRate=fs);

Estimate the mean frequency of the signal timetable. Plot the PSD and annotate the mean frequency.

meanfreq(sgTT);

You can estimate the mean frequency of the input signal timetable. Remove one sample from the signal timetable every 0.1 seconds.

sgTT(mod(seconds(sgTT.Time),0.1)<1/(2*fs), :) = [];

Estimate the mean frequency of the modified timetable. Plot the PSD and annotate the mean frequency. The estimated mean-frequency value approximates well with the mean-frequency estimation from the unmodified timetable.

meanfreq(sgTT);