Инструкция по выполнению работы

Лабораторная работа № 2

Исследование разомкнутой линейной системы при случайных возмущениях

Цели работы

· освоение методов анализа одномерной линейной непрерывной системы при случайных возмущениях с помощью среды Matlab

Задачи работы

· научиться вычислять среднеквадратическое отклонение и дисперсию ан выходе линейной системы, возбуждаемой единичным белым шумом

· научиться моделировать случайные процессы, используя в качестве источника сигнала белый шум (с ограниченной полосой)

· научиться оценивать СКВО и дисперсию случайного процесса, полученного при моделировании

· научиться вычислять автокорреляционную функцию случайного процесса

· научиться вычислять спектральную плотность случайного процесса по известной корреляционной функции

· научиться использовать быстрое преобразование Фурье (БПФ) для оценки спектральной плотности случайного процесса

· научиться использовать спектральные окна для сглаживания оценки спектральной плотности случайного процесса

Оформление отчета

Отчет по лабораторной работе выполняется в виде связного (читаемого) текста в файле формата Microsoft Word (шрифт основного текста Times New Roman, 12 пунктов, через 1,5 интервала, выравнивание по ширине). Он должен включать

· название предмета, номер и название лабораторной работы

· фамилию и инициалы авторов, номер группы

· фамилию и инициалы преподавателя

· номер варианта

· краткое описание исследуемой системы

· результаты выполнения всех пунктов инструкции, которые выделены серым фоном (см. ниже): результаты вычислений, графики, ответы на вопросы.

При составлении отчета рекомендуется копировать необходимую информацию через буфер обмена из рабочего окна среды Matlab. Для этих данных используйте шрифт Courier New, в котором ширина всех символов одинакова.

Инструкция по выполнению работы

Основная часть команд вводится в командном окне среды Matlab. Команды, которые надо применять в других окнах, обозначены иконками соответствующих программ.

Этап выполнения задания	Команды и иллюстрации
1. Очистите рабочее пространство Matlab (память).	clear all
2. Очистите окно Matlab.	clc
3. Введите передаточную функцию .	F = tf(1, [1 1])
4. Используя функцию norm, подсчитайте среднеквадратическое значение выхода этой системы при единичном белом шуме на входе.	norm ( F )
5. Подсчитайте дисперсию выхода системы при единичном белом шуме на входе.	norm ( F )^2
6. Найдите полосу пропускания этой системы (в рад/с).	bw = bandwidth ( F )
7. Найдите рекомендуемый максимальный интервал корреляции для моделирования по формуле	tau = 2*pi/100/bw
8. Запустите Simulink и создайте новую модель. Установите время моделирования 100 с (меню Simulation – Simulation Parameters – Stop Time).	на панели инструментов, в окне Simulink
9. Добавьте в модель блоки Band-Limited White Noise (белый шум с ограниченной полосой, группа Sources) и Scope (осциллограф, группа Sinks). Установите для белого шума параметр Noise Power (мощность) равный 1. Запустите модель и посмотрите, что представляет собой этот сигнал.
10. Так же, как и в предыдущей работе, подключите блоки Auto Correlator (автокорреляционная функция) и Power Spectral Density (спектральная плотность) из группы Simulink Extras – Additional Sinks). Посмотрите свойства этого сигнала.
11. Добавьте в схему звено с передаточной функцией так, как показано на схеме.
12. В параметрах блока Band-Limited White Noise уменьшите время корреляции (Sample Time) до значения, рассчитанного в п. 7. Для этого можно ввести в нужное поле имя переменной tau.
13. Откройте окно осциллографа, щелкните по кнопке и настройте параметры так, как показано на рисунке. На вкладке General в списке Sampling выберите вариант Sampling time (установить интервал вручную), а в соседнем поле введите имя переменной tau.На вкладке Data history нужно сбросить флажок Limit data points to last, установить флажок Save data to workspace, ввести имя переменной out и выбрать формат данных Array. Запустите моделирование.
14. Перейдите в окно Matlab, найдите среднеквадратическое отклонение (СКВО) и дисперсию сигнала на выходе звена. Сравните их со значениями, полученными в п. 4 и 5 по теоретическим формулам. Вычислите относительную ошибку при определении СКВО с помощью моделирования.	t = out(:,1); y = out(:,2); std ( y ) var ( y )
15. В окне блока Auto Correlator посмотрите, как выглядит корреляционная функция процесса, определенная по результатам моделирования.
16. Вычислите автокорреляционную функцию на выходе, используя функции Matlab, и постройте ее график. Сравните его с теоретической корреляционной функцией . Удобно создать новый m-файл и записать в него такие команды (без номеров строк): 1R = xcorr(y)/length(y); 2Rplus = R(floor(length(R)/2):end); 3M = 200; 4t = t(1:M); Rplus = Rplus(1:M); 5R_teor = 0.5*exp(-abs(t)); 6 figure(1); 7 plot(t, Rplus, t, R_teor) 8 xlim([0 max(t)]); Комментарий: 1 – вычисляем экспериментальную корреляционную функцию, используя стандартную функцию xcorr из пакета Signal Processing 2 – выделяем корреляционную функцию для положительных 3 – число точек для оценки корреляционной функции (меньше, чем длина сигнала) 4 – «обрезаем» массив отсчетов времени и корреляционную функцию 5 – для тех же значений времени находим теоретическую корреляционную функцию 6-7 – построение обоих графиков на одном поле 8 – установить точные границы по оси абсцисс
17. Подключите блок Averaging Power Spectral Density (усредненная спектральная плотность) из группы Simulink Extras – Additional Sinks). Выполните моделирование еще раз и сравните спектры, полученные с помощью двух разных блоков. Сделайте выводы.
18. Постройте спектральную плотность сигнала для частот от 0 до 5 рад/с. По полученным данным. Сравните ее с теоретической спектральной плотностью . 1T = t(2) - t(1); 2w = 0:0.02:5; 3Sw = []; Sw_teor =[]; 4for i=1:length(w) 5Sw(i) = sum(Rplus .* cos(w(i)*t)); 6Sw_teor(i) = 1 / (w(i)*w(i) + 1); 7end; 8Sw = 2*T*Sw; 9 figure(2); 10plot ( w, Sw, w, Sw_teor ); Комментарий: 1 – находим интервал дискретизации (он должен быть равен tau) 2 – задаем сетку частот, от 0 до 5 рад/с с шагом 0,02 рад/с 3 – опустошаем массивы 4-7 – цикл по всем выбранным частотам 5 – находим спектр как преобразование Фурье корреляционной функции 6 – теоретический спектр 8 – умножаем на 2T 9 – строим графики обоих спектров
19. Используя формулу , с помощью функции trapz (численное интегрирование методом трапеций), оцените дисперсию по экспериментальному и теоретическому спектрам. Объясните результаты, сравнив их со значением дисперсии, полученным в п. 5.	trapz(w,Sw)/pi trapz(w,Sw_teor)/pi
20. Постройте сглаженную оценку корреляционной функции с помощью окна Хэмминга. Для этого нужно добавить (в нужное место скрипта) команды hamm = 0.54 + 0.46*cos(pi*t/max(t)); Rhamm = Rplus .* hamm; и при построении корреляционных функций вывести третью линию plot(t, Rplus, t, R_teor, t, Rhamm) Перенесите в отчет полученный график.
21. Постройте оценку спектральной плотности, используя сглаженную корреляционную функцию. На графике должны быть три спектра (теоретический, оценка без сглаживания и оценка со сглаживанием). Скопируйте график в отчет.
22. Добавьте в скрипт команды для оценки спектральной плотности мощности с помощью быстрого преобразования Фурье (БПФ) и выполните его. 1N = 2*pi/0.5/T; 2N = 2^nextpow2(N); 3Fw = T * fft(y, N); 4Sw_fft = Fw .* conj(Fw) / N / T; 5Sw_fft = Sw_fft(1:N/2+1); 6w1 = 2*pi*[0:N/2] / N / T; 7plot ( w, Sw_teor, w1, Sw_fft ); 8xlim([0 max(w)]); Комментарий: 1 – считаем число точек для БПФ, чтобы шаг по частоте был равен 0,5 рад/с 2 – определяем ближайшую большую степень двойки 3 – выполняем БПФ 4 – считаем оценку спектра 5 – берем первую половину спектра до частоты Найквиста 6 – сетка угловых частот для построения графика 7 – строим теоретический спектр и оценку с помощью ДПФ 8 – устанавливаем пределы по оси абсцисс Сравните полученный результат с теоретической кривой. Сделайте выводы.
23. Повторите построение спектральной плотности, используя окно Хэмминга с масштабированием. Для этого добавьте в скрипт следующие команды: 1 scale = 1/sqrt(0.54^2 + 0.46^2/2); 2 hamm = hamming(N) * scale; 3 yHamm = y(1:N) .* hamm; 4 Fw = T * fft(yHamm, N); 5 Sw_fftHamm = Fw .* conj(Fw) / N / T; 6 Sw_fftHamm = Sw_fftHamm(1:N/2+1); 7 plot ( w, Sw_teor, w1, Sw_fft, w1, Sw_fftHamm ); 8 xlim([0 max(w)]); Комментарий: 1 – находим масштабирующий коэффициент для окна Хэмминга 2 – строим окно Хэмминга с масштабированием 3 – применяем окно к первым отсчетам сигнала 4 – выполняем БПФ 5 – считаем оценку спектра 6 – берем первую половину спектра до частоты Найквиста 7 – строим теоретический спектр и оценки 8 – устанавливаем пределы по оси абсцисс Запустите скрипт. Скопируйте полученный график в отчет. Сделайте выводы.