Clear Filters
Clear Filters

How to generate odd frequency sinusoid input using idinput

2 views (last 30 days)
Hi,
I am trying to generate a multisinusoid input signal with odd frequencies for non-linear system identification, as recommended in [1]. Going though the documentation of the System Identification Toolbox, i found the idinput function which allow you to generate a multisinusoid signal. According to the doccumentation [2], the function is designed to also generate even/odd frquences using the parameter "GridSkip", but I cannot figure out how to use the parameter. Can someone help me with some understanding of how this parameter works?
[1]: "Nonlinear system identification: a user-oriented road map" Johan Schoukens and Lennart Ljung, IEEE control system magazine 2019
[2]:

Answers (1)

Mathieu NOE
Mathieu NOE on 10 Jul 2023
I don't have anything against idinput , but why can't you do this directly with some basic code :
f = (1:2:7); % odd frequencies (array)
ampl = ones(size(f))./f; % amplitude (here decreasing in 1/f law)
dt = 1e-3;
samples = 5e3;
t = (0:samples-1)'*dt;
y = [];
for ci = 1:numel(f)
y = [y ampl(ci)*sin(2*pi*f(ci)*t)];
end
% signal sum of all frequencies and normalisation to +/- 1 amplitude
y = sum(y,2);
y = y./max(abs(y));
plot(t,y)
  3 Comments
Mathieu NOE
Mathieu NOE on 10 Jul 2023
fyi , some methods that can help you reduce the crest factor of the signal
clc
clearvars
f = (1:2:7); % odd frequencies (array)
ampl = ones(size(f))./f; % amplitude (here decreasing in 1/f law)
dt = 1e-3;
samples = 5e3;
t = (0:samples-1)'*dt;
% %multitone signals with low crest factor and as suggested in the Boyd's paper : quadratic phase distribution:
% % phas=pi*(f-1).^2/(numel(f)); % the original version
% phas=pi*(f-1).^2/(numel(f))^2; % I prefer this one
% % or Shapiro Rudin method
phas = shapiro(numel(f));
y = [];
for ci = 1:numel(f)
y = [y ampl(ci)*sin(2*pi*f(ci)*t+phas(ci))];
end
% signal sum of all frequencies and normalisation to +/- 1 amplitude
y = sum(y,2);
y = y./max(abs(y));
plot(t,y)
y_rms = sqrt(mean(y.^2));
CF = max(abs(y))/y_rms % compute crest factor of y
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function myphase = shapiro(n) %#ok<*DEFNU>
% phase according to Shapiro Rudin method (minimize crest factor)
if n == 1
myphase = 0;
else
% initialisation
str = [1 1];
k = ceil(log(n)/log(2));
for ci = 1:k-1
l =length(str);
str2 = str;
str2(l/2+1:l) = -str(l/2+1:l);
str = [str str2];
end
mystr = str(1:n);
myphase= zeros(1,n);
ind = mystr==-1;
myphase(ind) = pi;
end
end

Sign in to comment.

Categories

Find more on Linear Model Identification in Help Center and File Exchange

Products


Release

R2022a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!