% M03-de+corr.diary % Running an LCCDE, plus cross-correlation % 2007-09-18 % Evaluating an LCCDE algorithmically % y[n] + y[n-1] - 6y[n-2] = x[n] % with y[-2] = -1, y[-1] = 1, x[n] = 8.u[n] % Calculate for n = 0..4, but include -2 and -1 nn = -2:4; % So setup initial y values; zero is index 3 zz = 3; y(zz+ -2) = -1; y(zz+ -1) = 1; % Now evaluate the LCCDE step by step x = zeros(1,7); x(zz+ [0:4]) = 8; for n = 0:4; y(zz+ n) = x(zz+ n) - y(zz+ n-1) + 6*y(zz+ n-2); end y % y = % -1 1 1 13 1 85 -71 % Compare to algebraic solution l1 = -3; l2 = 2; beta = -2; a1 = -1.8; a2 = 4.8; yt = a1*(l1.^nn) + a2*(l2.^nn) + beta % yt = % -1.0000 1.0000 1.0000 13.0000 1.0000 85.0000 -71.0000 nn % nn = % -2 -1 0 1 2 3 4 % In this case, yt[-2, -1] *do* match ICs, but not always true % Use built-in filter() command to calculate impulse response % (you can't (easily) use filter() if there are initial conditions) % x DE coefficients {p_k} = {1}, y DE coefficients {d_k} = {1,1,-6} ir = filter([1], [1 1 -6], [1 0 0 0 0 0 0 0]) % ir = % Columns 1 through 4 % 1 -1 7 -13 % Columns 5 through 8 % 55 -133 463 -1261 % Compare to closed form h = 0.6*((-3).^[0:7]) + 0.4*(2.^[0:7]) % h = % 1.0e+03 * % Columns 1 through 5 % 0.0010 -0.0010 0.0070 -0.0130 0.0550 % Columns 6 through 8 % -0.1330 0.4630 -1.2610 % different formatting, but the same values... % Cross-correlation visualization, to find a pulse in noise % Basic pulse p = [-1 3 2 1]; % Noise sequence n = randn(1,20); % Make the signal be the noise plus pulse at time 10 x = n; x(10+[1:4]) = x(10+[1:4]) + p; % Look at them subplot(311) stem(10:13,p); axis([-2 20 -5 5]); grid on; title('Pulse p'); subplot(312) stem([0:19],n); axis([-2 20 -5 5]); grid on; title('Noise n'); subplot(313) stem([0:19],x); axis([-2 20 -5 5]); grid on; title('x = noise + pulse'); % Calculate the cross-correlation step by step; max should be at 10 corrviz(p,x)