/
ch03_ex07_time_varying_oscillator.m
84 lines (57 loc) · 1.44 KB
/
ch03_ex07_time_varying_oscillator.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
%% Example 3.7: Stochastic resonators with time-varying frequency
%
% Copyright:
% 2018 - Simo Särkkä and Arno Solin
%
% License:
% This software is provided under the MIT License. See the accompanying
% LICENSE file for details.
%%
% Lock random seed
if exist('rng') % Octave doesn't have rng
rng(0,'twister')
else
randn('state',1);
end
% Time
t = linspace(0,10,512);
% Define the sinc function
sinc = @(x) sin(x)./x;
% Set time-varying frequency trajectory
f = 1+.5*abs(t-5);
% Parameters
H = [1 0];
figure(1); clf; hold on
for j=1:4
% Allocate space
x = zeros(2,numel(t));
z = x;
% For harmonics
for i=1:2
% Diffusion
Qc = 10^-i;
% Initial state
x(:,1) = [0; 1];
% Loop through time points
for k=2:numel(t)
% Time delta
dt = t(k)-t(k-1);
% Define F and L
F = [0 2*pi*i*f(k); -2*pi*i*f(k) 0];
L = [0; 1];
% Solve A and Q
[A,Q] = lti_disc(F,L,Qc,dt);
% The next state
x(:,k) = A*x(:,k-1) + chol(Q,'lower')*randn(size(Q,1),1);
end
z = z + H*x;
end
plot(t,5*(j-1)+z,'-k','LineWidth',.5)
end
xlabel('Time, $t$')
xlim([0 10.5])
ylim([-3 18])
figure(2); clf
plot(t,f)
xlabel('Time, $t$')
ylabel('Frequency [Hz]')