MATLAB 8.0 (R2012b): создание, обработка и фильтрация сигналов, Signal Processing Toolbox

№ 11’2013
PDF версия
Новая матричная система компьютерной математики MATLAB 8.0 (R2012b) хорошо приспособлена для решения фундаментальных и прикладных задач создания и обработки сигналов. В этой статье впервые в отечественной литературе описаны возможности работы с сигналами как в среде самой системы MATLAB 8.0, так и с помощью Signal Processing Toolbox и нового пакета расширения этой системы — Simulink 8.0, реализующего блочное визуально-ориентированное имитационное моделирование. Автор благодарит корпорацию The MathWork [1] за предоставленную систему MATLAB 8.0 + Simulink 8.0, использованную для подготовки этой серии статей.

Введение

Окно справки по пакету расширения Signal Processing Toolbox системы MATLAB 8.0 представлено на рис. 1 на фоне рабочего окна самой системы с открытой вкладкой каталога пакетов расширения APPS. Одна из кнопок в панели инструментов Signal Analysis дает доступ к браузеру сигналов, фильтров и спектров, показанному в правой части окна справки. В левой части этого окна указаны наименования разделов пакета Signal Processing Toolbox.

Окно справки по пакету расширения Signal Processing Toolbox

Рис. 1. Окно справки по пакету расширения Signal Processing Toolbox

Как видно на рис. 1, пакет Signal Processing Toolbox состоит из следующих разделов:

  • Waveforms — создание сигналов с различной формой и разными законами модуляции;
  • Convolution and Correlation — свертка и корреляция сигналов;
  • Transform — преобразование сигналов;
  • Analog and Digital Filters — аналоговые и цифровые фильтры;
  • Spectral Analysis — спектральный анализ.

 

Создание сигналов

Многие сигналы представлены как функции времени s(t), параметры которой можно изменять с помощью модуляции того или иного вида. Модуляцией называют процесс изменения какого-либо параметра (амплитуды, частоты, фазы и т. д.) по определенному закону, в результате чего сигнал становится переносчиком информации.

MATLAB 8.0 со своими встроенными средствами позволяет создавать множество сигналов. Например, простейшим является синусоидальный сигнал:

s = Asin(2πft+j),

где A — амплитуда; f — частота; j — фаза сигнала. Задав эти параметры (t — как вектор отсчетов сигнала, число элементов которого определяет число отсчетов сигнала):

>> t = 1:0.01:10; A=1; f=1; phi=0;

>> s=A*sin(2*pi*f*t+phi);

можно легко построить график сигнала s(t):

>> plot(t, s)

Он показан на рис. 2 в графическом окне системы MATLAB. Такое построение осуществляется в командной строке, на что указывает приглашение к вводу >>. В командном окне также показано окно с краткими данными о системе MATLAB 8.0 (R2012b). Согласно этим данным система выпущена на рынок в 2012 году.

Окно командного режима

Рис. 2. Окно командного режима

В дальнейшем, учитывая множество сигналов, мы будем представлять их по четыре в каждом графическом окне, используя для этого функцию его разбиения на четыре под-окна subplot.

Окно с подокнами меандра, апериодического треугольного импульса, пилообразного периодического импульса и периодического треугольного импульса

Рис. 3. Окно с подокнами меандра, апериодического треугольного импульса, пилообразного периодического импульса и периодического треугольного импульса

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

subplot(2,2,1); t=0:0.01*pi:5*pi; plot(t, square(t,20))

xlabel(‘Время( с)’);ylabel(‘Уровень’); title(‘Square’)

subplot(2,2,2);t=–10:0.1:10; plot(t, tripuls(t,5,0.5))

xlabel(‘Время( с)’);ylabel(‘Уровень’); title(‘Tripuls’)

subplot(2,2,3);t=0:0.01:20; plot(t, sawtooth(t,1))

xlabel(‘Время( с)’);ylabel(‘Уровень’); tit(t,1’)

subplot(2,2,4);t=0:0.01:20; plot(t, sawtooth(t,1/2))

xlabel(‘Время( с)’);ylabel(‘Уровень’); title(‘Sawtooth(t,1/2’)
Следующая программа создает четыре непрерывных сигнала, относящихся к функциям Гаусса (рис. 4):



subplot(2,2,1); tc = gauspuls(‘cutoff’,50e3,0.6,[],–40);

t = –tc : 1e–6 : tc; yi = gauspuls(t,50e3,0.6); plot(t,yi)

xlabel(‘Время( с)’);ylabel(‘Уровень’); title(‘Gauspuls’)

subplot(2,2,2); fc = 2E9; fs=100E9; tc = gmonopuls(‘cutoff’,fc);

t = –2*tc : 1/fs : 2*tc; y = gmonopuls(t,fc); plot(t,y)

xlabel(‘Время( с)’);ylabel(‘Уровень’); title(‘Gmonopuls’)

subplot(2,2,3); fc = 2E9; fs=100E9;

D = [2.5 10 17.5]’ * 1e–9; tc = gmonopuls(‘cutoff’,fc);

t = 0:1/fs:160*tc; yp = pulstran(t,D,@gmonopuls,fc);

plot(t,yp): xlabel(‘Время( с)’);ylabel(‘Уровень’); title(‘Pulstrain@gmonopuls’)

subplot(2,2,4); x = linspace(0,4*pi,300);

plot(x,diric(x,4)); axis tight; hold on

plot(x,diric(x,8)); axis tight; hold off

xlabel(‘Время( с)’);ylabel(‘Уровень’); title(‘2 Diric’)

t = (1:10)’; x = randn(size(t)); ts = linspace(–5,15,600)’;

y = sinc(ts(:,ones(size(t))) – t(:,ones(size(ts)))’)*x; plot(t,x,’o’,ts,y)
Временные диаграммы четырех сигналов, относящихся к функциям Гаусса

Рис. 4. Временные диаграммы четырех сигналов, относящихся к функциям Гаусса 

Важное значение имеет функция sinc (или sin(πt)/πt при t ≠ 0 и 1 при t = 0). Функция sinc(t) представляет обратное преобразование Фурье для прямоугольного импульса с высотой 1 и шириной 2π:

Кроме того, эту функцию можно использовать как базисную для восстановления любого сигнала g(t) по его отсчетам, если спектр сигнала ограничен условием –p:

Получение непрерывной кривой, проходящей через точки (отсчеты) произвольного сигнала

Рис. 5. Получение непрерывной кривой, проходящей через точки (отсчеты) произвольного сигнала

Это положение, вытекающее из известной теоремы Котельникова, иллюстрирует приведенный ниже пример для десяти случайных точек сигнала (рис. 5):

t = (1:10)’; x = randn(size(t)); ts = linspace(–5,15,600)’;

y = sinc(ts(:,ones(size(t))) – t(:,ones(size(ts)))’)*x;

plot(t,x,’o’,ts,y)

Подобный метод восстановления аналогового сигнала из его цифровых отсчетов ныне применяется во всех цифровых осциллографах.

Сложные модулированные сигналы

Рис. 6. Сложные модулированные сигналы

Сложные модулированные сигналы показаны на рис. 6, здесь chirp — период частотно-модулированного сигнала, sawtooth — сигнал с частотной модуляцией по треугольному закону, pulstrain — скачок с линейным спадом, pulstran gauspuls — последовательность импульсов Гаусса с убывающей амплитудой. Для создания этих четырех модулированных сигналов служит программа:

subplot(2,2,1); t = 0:0.001:0.4; y = chirp(t,0,1,250); plot(t,y)

xlabel(‘Время( с)’);ylabel(‘Уровень’); title(‘Chirp’)

subplot(2,2,2); fs = 100; t = 0:1/fs:2;

y = vco(sawtooth(2*pi*t,0.75),[0.1 0.4]*fs,fs);plot(t,y);

xlabel(‘Время (с)’);ylabel(‘Уровень’); title(‘Sawtooth’)

subplot(2,2,3); t = 0 : 1/1e3 : 1; d = 0 : 1/3 : 1;

y = pulstran(t,d,‘tripuls’,0.1,–1); plot(t,y)

xlabel(‘Время (c)’);ylabel(‘Уровень’); title(‘Pulstrain’)

subplot(2,2,4); t = 0 : 1/50E3 : 10e–3;

d = [0 : 1/1E3 : 10e–3 ; 0.8.^(0:10)]’;

y = pulstran(t,d,‘gauspuls’,10e3,0.5); plot(t,y);

xlabel(‘Время (с)’);ylabel(‘Уровень’);

title(‘Pulstrain (gauspuls)’)

Здесь функция:

y = chirp(t,f0,t1,f1,[ ‘method’,phi])

формирует выборку (дискретные значения) косинусоидального сигнала с частотой от f0 в начальный момент времени t до f1 в конечный момент времени t1. Звук такого сигнала напоминает визг, откуда и его название (chirp). По умолчанию t = 0, f0 = 0 и f1 = 100. Необязательный параметр phi (по умолчанию 0) задает начальную фазу сигнала. Другой необязательный параметр — method — задает закон изменения частоты: linear — линейный закон (по умолчанию), quadratic — квадратичный и logarithmic — логарифмический.

Функция strips позволяет детально (по частям с длительностью 0,25 с) рассмотреть первые два сложных сигнала, например частотно-модулированного от функции vco. В разделе модуляция/демодуляция есть пример с GUI (рис. 7), иллюстрирующий такие сигналы с различными видами модуляции.

Визуализация модулированных сигналов

Рис. 7. Визуализация модулированных сигналов

 

Свертка и корреляция сигналов

Пусть имеется две последовательности, представленные векторами a и b. Сверткой называют одномерный массив, вычисляемый следующим образом:

При записи этого выражения учтено, что нумерация индексов массивов в MATLAB идет с единицы (1). Свертка реализуется функцией conv(a,b):

a = [2 4 6]; b = [3 5 7]; c = conv(a,b)

c =

           6            22           52          58          42

Операцию свертки часто используют для вычисления сигнала на выходе линейной системы y по сигналу на входе x при известной импульсной характеристике системы h:

Эту операцию можно использовать для осуществления простейшей фильтрации сигнала:

 

x = randn(5,1);      % Создание вектора из пяти случайных чисел

h = [1 1 1 1]/4;      % Усредняющий фильтр с длиной 4

y = conv(h,x)           % Результат свертки (фильтрации)

y =

           0.2977

           0.5950

           0.5856

           0.6674

           0.4134

           0.1161

           0.1255

           0.0437

Для двух векторов x и y с длиной m и n определена операция свертки:

Линейная и циклическая свертка

Рис. 8. Линейная и циклическая свертка

Обратная свертке функция — это [q,r] = deconv(z,x). Она фактически определяет импульсную характеристику фильтра (рис. 8):

x = [2 1 2 1]; y = [1 2 3]; clin = conv(x,y);

xpad = [x zeros(1,6–length(x))]; ypad = [y zeros(1,6–length(y))];

ccirc = ifft(fft(xpad).*fft(ypad)); subplot(211);

stem(clin,‘markerfacecolor’,[0 0 1]);

title(‘Linear Convolution of x and y’);

set(gca,‘ylim’,[0 11]); subplot(212)

stem(ccirc,‘markerfacecolor’,[0 0 1]); set(gca,‘ylim’,[0 11]);

title(‘Circular Convolution of xpad and ypad’);
Очистка от шума линейно-нарастающего сигнала

Рис. 9. Очистка от шума линейно-нарастающего сигнала

Свертку часто применяют для очистки сигналов от шума (рис. 9):

x = –4:0.01:4; rng default; y = x+randn(size(x));

coeffs = polyfit(x,y,1); yfit = coeffs(2)+coeffs(1)*x;

plot(x,y); hold on; plot(x,yfit,‘r’,‘linewidth’,2);hold off

Для двумерных массивов также существует функция свертки: Z = conv2(X,Y) и Z = conv2(X,Y,‘option’). Возможна и многомерная свертка — функция convn. Новая функция mscohere строит график зависимости квадрата модуля функции когерентности от частоты (рис. 10):

randn(‘state’,0); h = fir1(30,0.2,boxcar(31));

h1 = ones(1,10)/sqrt(10);

r = randn(16384,1); x = filter(h1,1,r); y = filter(h,1,x);

mscohere(x,y,512,[],[],256)
Зависимость квадрата модуля функции когерентности от частоты

Рис. 10. Зависимость квадрата модуля функции когерентности от частоты

Следующая программа строит отсчеты двух сигналов с задержкой на три отсчета — треугольного сигнала и шума (рис. 11):

rng default; x = triang(20);

y = [zeros(3,1); x]+0.3*randn(length(x)+3,1);

subplot(211);stem(x,‘markerfacecolor’,[0 0 1]);axis([0 22 –1 2]);

subplot(212);stem(y,‘markerfacecolor’,[0 0 1]);axis([0 22 –1 2]);
Два коррелированных сигнала — треугольный и шума

Рис. 11. Два коррелированных сигнала — треугольный и шума

Кросс-корреляцию этих сигналов обеспечивает программа (рис. 12):

[xc,lags] = xcorr(y,x); [~,I] = max(abs(xc));

figure; stem(lags,xc,‘markerfacecolor’,[0 0 1]);
Кросс-корреляция сигналов

Рис. 12. Кросс-корреляция сигналов

 

Преобразование сигналов

В MATLAB всегда много внимания уделялось различным преобразованиям сигнала [2–5], в частности спектральным. Спектр дискретного сигнала является периодическим, и прямое дискретное преобразование Фурье (ДПФ или Discrete Fourier Transform, DFT) определяется выражением:

Для предотвращения растекания (размазывания) спектра дискретных сигналов часто используются окна. Для этого достаточно в формуле прямого ДПФ под знаком суммы ввести еще один множитель — W(k). Соответственно, обратное дискретное преобразование Фурье задается выражением:

ДПФ легко обеспечивает восстановление непрерывных периодических сигналов с ограниченным спектром. Для этого нужно номер отсчета k поменять на нормированное время t/T. Тогда формула восстановления при четном числе отсчетов будет иметь вид:

Для получения полосы частот сигнала от 0 до π/T приходится смещать нумерацию отсчетов. При нечетном числе отсчетов суммирование ведется при n, меняющемся от –(N–1)/2 до (N–1)/2. Коэффициенты X·(n) с отрицательными номерами вычисляют из соотношения симметрии.

Частотным спектром случайного процесса является преобразование Фурье от корреляционной функции случайного процесса Rx:

В радиоэлектронике особый интерес представляет спектральная оценка сильно зашумленных сигналов. Для таких сигналов применяются два подхода: непараметрический — использующий только информацию, извлеченную из сигнала (реализован в методах периодограмм и Уэлча), и параметрический — предполагающий наличие некоторой статистической модели сигналов, параметры которой подлежат определению. Реализовано восемь классов алгоритмовспектрального анализа: Periodogram, Welch, MTM (Thomson multitaper method), Burg, Covariance, Modified Covariance, Yule-Walker, MUSIC (Multiple Signal Classification) и Eigenvector. Их подробное описание дано в [4–6].

 

Периодограммы, спектрограммы и их применение

Обычный спектр строится методом быстрого преобразования Фурье (БПФ, FFT) часто с применением временного окна, предотвращающего разрывы сигнала на концах интервала анализа спектра. Например, так реализованы периодограммы. Для построения спектрограмм используется разбивка интервала анализа на короткие окна. Короткое окно пробегает общий интервал анализа, и в каждом частичном интервале строится свой спектр. Их наложение дает спектрограмму, определенную в пространстве «уровень – частота – время» (на плоскости уровень представляется цветом), тогда как обычный спектр определен в плоскости «уровень – частота».

Построение спектрограмм сигналов, модулированных по различным законам: логарифмическому, линейному, треугольному и синусоидальному

Рис. 13. Построение спектрограмм сигналов, модулированных по различным законам: логарифмическому, линейному, треугольному и синусоидальному

На рис. 13 показано построение четырех спектрограмм для сигналов с различными законами частотной модуляции. Нетрудно увидеть, что во всех случаях закон модуляции отчетливо распознается и позволяет судить об области частот, в которой действует сигнал (chirp или vco). Этот рисунок строит следующая программа:

subplot(2,2,1); t = 0.5:0.001:2; y = chirp(log(t),0,1,500);

spectrogram(y,256,250,256,1E3,‘yaxis’)

subplot(2,2,2); t = –2:0.001:0; y = chirp(t,0,1,250);

spectrogram(y,256,250,256,1E3,‘yaxis’)

subplot(2,2,3);fs = 10000; t = 0:1/fs:2;

x = vco(sawtooth(2*pi*t,0.75),[0.1 0.4]*fs,fs);

spectrogram(x,kaiser(256,5),220,512,fs,‘yaxis’)

subplot(2,2,4);fs = 10000; t = 0:1/fs:2;

x = vco(sin(2*pi*(fs/5000)*t),[0.1 0.4]*fs,fs);

spectrogram(x,kaiser(256,5),220,512,fs,‘yaxis’)
Сравнение периодограммы и спектрограммы сигнала sinc(t)

Рис. 14. Сравнение периодограммы и спектрограммы сигнала sinc(t)

Из сказанного может сложиться неверное представление о явных преимуществах спектрограмм по сравнению с периодограммами (функция psd). То, что это далеко не так, показывает программа для сигнала sinc (рис. 14):

Fs = 10000; t = 0:0.1/Fs:0.0001; x = sinc(t);

Hs=spectrum.periodogram(‘Blackman’);

subplot(1,2,1); psd(Hs,x,‘Fs’,Fs)

subplot(1,2,2); spectrogram(x)

title(‘Sinc spectogram’)

Эта программа строит периодограмму и спектрограмму функции sinc(t). Теперь беспомощной оказывается спектрограмма, по которой ничего нельзя сказать о спектре сигнала и области его частот. А периодограмма более информативна: она указывает на вид спектра и область занимаемых им частот. В частности, хорошо видно постоянство спектра в начальной области частот, присущее этой функции. Правда, вид спектрограммы сильно зависит от типа короткого окна: в данном случае задано окно Блэкмана (Blackman). Вид периодограммы также зависит от выбора окна, но глобального.

Чем сложнее сигнал, тем более детальной и эффектной оказывается спектрограмма. На рис. 15 показан пример модуляции/демодуляции с построением спектрограммы сложного звукового сигнала. Кстати, оригинальный сигнал и сигнал, прошедший модуляцию/демодуляцию, можно воспроизвести на компьютере, оборудованном звуковой картой и акустической системой. Это можно сделать при различных видах модуляции.

Спектрограмма сложного акустического сигнала с амплитудной модуляцией

Рис. 15. Спектрограмма сложного акустического сигнала с амплитудной модуляцией

Похожие на спектрограммы картинки дают вейвлетограммы и скайлеграммы, получаемые при вейвлет-анализе сигналов [7]. Порою они более информативны. Но вейвлеты в Signal Ptocessing Toolbox не реализованы: они описаны и используются в отдельном пакете расширения Wavelet Toolbox.

Для выполнения дискретного преобразования Фурье различными методами служит GUI-окно Discrete Fourier Transform (рис. 16).

Спектр пилообразного сигнала в окне дискретного Фурье-преобразования

Рис. 16. Спектр пилообразного сигнала в окне дискретного Фурье-преобразования

Оно позволяет задать один из трех видов сигнала (синусоида, меандр и пилообразный) и его периодограмму (спектр) при одном из шести видов окон. Как сигнал, так и окно можно загружать извне из файла или рабочего пространства MATLAB. Нетрудно убедиться в большом влиянии на вид спектра выбранного окна.

 

Оконные функции и браузер окон

Учитывая важную роль окон, в Signal Processing Toolbox входит 21 N‑точечное окно (N — целое число), называемое по фамилии предложивших окно ученых, например Hamming, Blackman, Bartlett, Chebyshev, Taylor, Kaiser и т. д. Для просмотра временных и амплитудно-частотных характеристик всех окон есть соответствующие функции, но удобно пользоваться GUI-браузером окон wvtool.

Например, команда:

wvtool(rectwin(64),hamming(64),hann(64),gausswin(64))

позволяет строить четыре типа 64‑точечных окон — прямоугольное, Хемминга, Ханна и Гаусса. Можно задать окна и отдельными командами.

Сравнение в браузере трех окон Кайзера с параметром бетта 1,5, Блэкмана и Блэкмана-Харисса

Рис. 17. Сравнение в браузере трех окон Кайзера с параметром бетта 1,5, Блэкмана и Блэкмана-Харисса

Это, а также сравнительное построение трех других типов окон обеспечивает следующая программа (рис. 17):

H = sigwin.kaiser(128,1.5); H1 = blackman(128);

H2 = blackmanharris(128); wvtool(H,H1,H2)
Конструктор окон

Рис. 18. Конструктор окон

Существует также конструктор-анализатор окон с GUI-интерфейсом, окно которого (рис. 18) открывается командой:

wintool

или

wintool(obj1, obj2, …)

По умолчанию в нем открывается окно объекта sigwin, но можно открыть и другие окна (в том числе окно пользователя) из списка Current Window Information. Открытые окна появляются в третьем нижнем окне. Окна можно скопировать, добавить в список, установить в рабочее пространство MATLAB или стереть.

 

Изменение числа отсчетов и интерполяция сигналов

Изменение числа отсчетов широко используется в технике цифровой обработки сигналов. Для уменьшения числа отсчетов применяются операция децимации и функция decimate (рис. 19):

t=0:.00025:1;x=sin(2*pi*30*t)-sin(2*pi*70*t);

y=decimate(x,4);

subplot(121); stem(x(1:120)),

axis([0 120 -2 2]); title(‘Original Signal’);

subplot(122); stem(y(1:30)); title(‘Decimated Signal’)
Уменьшение частоты дискретизации сигнала (децимация)

Рис. 19. Уменьшение частоты дискретизации сигнала (децимация)

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

t = 0:0.001:1; x = sin(2*pi*30*t) - sin(2*pi*70*t);

y = interp(x,4); subplot(121); stem(x(1:20));axis([0 20 –2 2]);

title(‘Оригинальный сигнал’); subplot(122); stem(y(1:80));

title(‘Сигнал после интерполяции’); axis([0 80 –2 2]);
Интерполяция сигнала и увеличение частоты его дискретизации

Рис. 20. Интерполяция сигнала и увеличение частоты его дискретизации

 

Аналоговые и цифровые фильтры

Фильтры имеют особое значение при обработке сигналов. С их помощью осуществляется очистка сигналов от шума или реализуются избирательные свойства систем. Будем считать, что читатель знаком с теорией фильтров.

Фильтрующие цепи обычно задаются своей операторной передаточной характеристикой:

h(s) = a(s)/b(s).

Имея векторы коэффициентов полиномов a(s) и b(s), с помощью функции freqs можно построить АЧХ и ФЧХ фильтрующей цепи в логарифмическом масштабе (рис. 21):

a = [1 0.4 1]; b = [0.2 0.3 1];

w = logspace(–1,1); freqs(b,a)
Построение логарифмической АЧХ и ФЧХ системы по ее операторной передаточной характеристике

Рис. 21. Построение логарифмической АЧХ и ФЧХ системы по ее операторной передаточной характеристике

В Signal Processing Toolboox входит множество функций по расчету и проектированию различных фильтров — нижних, верхних частот и полосовых. Порою для получения важных характеристик фильтров достаточно задать нужную строку программного кода. Например, построение характеристик аналогового фильтра Бесселя (рис. 22) реализуется следующим программным фрагментом:

[b,a] = besself(5,10000); freqs(b,a);

Для просмотра практически всех характеристик фильтров можно использовать визуализатор фильтров fvtool. Приведенная ниже программа дает характеристики фильтра Бесселя нижних частот (рис. 23):

[z,p,k] = ellip(6,3,50,300/500); [sos,g] = zp2sos(z,p,k);

Hd = dfilt.df2tsos(sos,g); h = fvtool(Hd)

set(h,‘Analysis’,‘freq’)
Характеристики аналогового эллиптического фильтра нижних частот

Рис. 23. Характеристики аналогового эллиптического фильтра нижних частот

Открытое меню Analysis дает представление обо всех доступных видах анализа фильтров. Среди них построение АЧХ и ФЧХ фильтра, импульсная и переходная характеристики, групповая задержка и т. д.

Конструирование двух цифровых фильтров НЧ (FIR и Баттерворта) со сравнительным построением их характеристик в окне fvtool (рис. 24) обеспечивает следующая программа:

d=fdesign.lowpass(‘Fp,Fst,Ap,Ast’,1000,1200,1,60,20000);

Hd1=design(d,‘equiripple’); %FIR equiripple design

Hd2=design(d,‘butter’); %Butterworth design

fvtool([Hd1 Hd2],‘legend’,’on’); axis([0 2 –100 10])
АЧХ двух цифровых фильтров для сравнения

Рис. 24. АЧХ двух цифровых фильтров для сравнения

Программа строит и другие характеристики фильтров, например переходные и импульсные. Напомним, что переходная характеристика является реакцией системы (фильтра) на единичный скачок, а импульсная — на импульс единичной площади с длительностью, стремящейся к нулю (рис. 25).

Импульсные характеристики двух цифровых фильтров

Рис. 25. Импульсные характеристики двух цифровых фильтров

 

Интерактивный конструктор-анализатор фильтров fdatool

В Signal Proccesing Toolbox входит интерактивный конструктор-анализатор фильтров с GUI-интерфейсом. Он позволяет без какого-либо программирования анализировать и проектировать семь основных типов фильтров (таблица).

Таблица. Семь основных типов фильтров

Русскоязычное наименование

Метод конструирования

Функция

Баттерворта

Butterworth

butter

Чебышева I

Chebyshev Type I

cheby1

Чебышева II

Chebyshev Type II

cheby2

Эллиптический

Elliptic

ellip

С максимально плоской вершиной

Maximally Flat

maxflat

Наименьших квадратов

Least-squares

firls

Оконный

Window

fir1

Кроме того, поддерживается 13 расширенных методов анализа и проектирования фильтров: FIR (с конечной импульсной характеристикой), IIR (с бесконечной импульсной характеристикой), Найквиста, интерполяционного, FIR второго порядка и др.

При загрузке программы командой fdatool появляется окно, показанное на рис. 26. В нем достаточно установить тип фильтра и его параметры, в частности указанные в зоне возможного расположения АЧХ фильтра. Эта зона меняется в зависимости от вида выбранного типа. Открытое меню Analysis показывает, какие виды анализа возможны.

Чтобы выполнить конструирование фильтра, достаточно активизировать мышью кнопку Design Filter внизу окна. На экране появится окно с параметрами выбранного фильтра (рис. 27), обычно по умолчанию это АЧХ фильтра. С АЧХ можно вывести и ФЧХ фильтра. Меню Analysisпозволяет выбрать и другие параметры, но они будут отображены отдельно.

 

Браузер сигналов, фильтров и спектров SPTool

Многие из описанных возможностей объединяет в себе браузер сигналов, фильтров и спектров SPTool, GUI-окно которого показано на рис. 1 (окно справа). Он имеет три раздела:

  • Signals — просмотр сигналов;
  • Filters — просмотр фильтров;
  • Spectra — просмотр спектров.

В разделе Signals можно просмотреть три встроенных типа сигнала. Первый — mtlib — сложный сигнал (рис. 28). На этом рисунке также показано открытое меню Tools и списки со статистикой сигнала и данными измерений его параметров.

Просмотр параметров сигнала mtlib

Рис. 28. Просмотр параметров сигнала mtlib

Второй сигнал — chirp (визг). На рис. 29 показана пара вертикальных маркеров, предназначенных для точного измерения отмеченного ими пика спектра. Кроме того, выведено окно с данными нового пакета расширения Signal Processing Toolbox. Видно, что пакет имеет версию 6.18.

Просмотр параметров сигнала chirp

Рис. 29. Просмотр параметров сигнала chirp

Третий сигнал — train. На рис. 30 выведены два курсора и приведены данные измерений указанного ими уровня. Таким образом, о сигнале можно получить исчерпывающие сведения.

Просмотр параметров сигнала train

Рис. 30. Просмотр параметров сигнала train

В разделе Filters можно просматривать характеристики одного из трех фильтров или загруженного извне фильтра. На рис. 31 показаны АЧХ и ФЧХ одного из фильтров, увеличенные в размерах по оси X.

АЧХ и ФЧХ цифрового полосового фильтра FIRbp

Рис. 31. АЧХ и ФЧХ цифрового полосового фильтра FIRbp

Раздел Spectra служит для детального просмотра спектров как одного из трех встроенных сигналов (рис. 32 и 33), так и для сигналов, загружаемых извне.

Спектр сигнала chirp, полученный методом FFT

Рис. 32. Спектр сигнала chirp, полученный методом FFT

Возможен просмотр спектра сигналов различными методами и работа с курсорными измерениями.

Спектр сигнала, полученный методом Уэлча, и средства его анализа

Рис. 33. Спектр сигнала, полученный методом Уэлча, и средства его анализа

 

Signal Processing Toolbox в пакете блочного моделирования Simulink

В более ранних версиях MATLAB был набор инструментальных ящиков Blockset, расширяющих возможности пакета блочного имитационного моделирования Simulink. В этом наборе был пакет Signal Processing Blockset. В MATLAB 8.0 этот пакет исключен, но его функции включены как в сам Simulink, так и в пакет DSP System Toolbox. Они содержат все блоки для построения графических диаграмм моделей, создающих и обрабатывающих различные сигналы. Для составления диаграмм моделей знание программирования не требуется. Многие средства DSP System Toolbox аппаратно зависимы. Предусмотрена реализация фильтров с фиксированной запятой с применением встроенного в микропроцессор формата чисел. Это существенно ускоряет вычисления.

На рис. 34 дан пример диаграммы с блоками, реализующей обычное быстрое преобразование Фурье, для построения спектра трехчастотного синусоидального сигнала с частотами 15, 130 и 210 Гц. Спектр строит блок Vector Scope — виртуальный анализатор спектра. Показаны также окно установки анализируемого сигнала и окно с его спектром. На нем отчетливо видны три пика на указанных частотах.

Диаграмма модели проведения БПФ

Рис. 34. Диаграмма модели проведения БПФ

Сравнение трех методов спектрального анализа демонстрирует Simulink-диаграмма (рис. 35). Подобный пример сравнения в командном режиме приводился ранее. Здесь он реализован в Simulink на уровне блочного имитационного моделирования.

Сравнение спектрального анализа тремя методами

Рис. 35. Сравнение спектрального анализа тремя методами

Другой пример (рис. 36) показывает диаграмму модели LPC (кодирование с линейным предсказанием) анализа и синтеза звукового сигнала от системы MATLAB с возможностью прослушивания исходного и синтезированного сигнала (кнопки под диаграммой). В системе используются три цифровых фильтра и три блока Waterfall Scope, строящих развитие спектра во времени.

Диаграмма модели LPC-анализа и синтеза звукового сигнала

Рис. 36. Диаграмма модели LPC-анализа и синтеза звукового сигнала

Еще один пример построения диаграммы адаптивного LMS-фильтра приведен на рис. 37. На вход фильтра подан сложный сигнал — синусоида с шумом. Фильтр использует для очистки от шума алгоритм наименьших квадратов с временной задержкой [8]. Результаты моделирования представлены на рисунке.

Диаграмма адаптивного LMS-фильтра и результаты его моделирования

Рис. 37. Диаграмма адаптивного LMS-фильтра и результаты его моделирования

 

Заключение

При подготовке MATLAB 8.0 + Simulink 8.0 к выходу на рынок была проделана огромная работа по совершенствованию пакета расширения Signal Processing Toolbox, содержащего основные средства создания и обработки сигналов, включая методы спектрального анализа и фильтрации сигналов. В MATLAB 8.0 создание и обработка сигналов ведется в командном окне с применением объектно-ориентированного программирования, в окнах с GUI-интерфейсом без применения программирования и в пакетах Simulink 8.0 и DSP System Toolbox на основе блочного имитационного моделирования и визуально-ориентированного программирования. При этом везде используются новейшие методы обработки сигналов. Особое внимание уделено удобству и наглядности применения доступных пользователю средств. В частности, многие из них включены в каталоги, панели которых расположены в верхней части окна MATLAB 8.0

Литература
  1. www.mathworks.com
  2. Дьяконов  В. П. MATLAB R2006/2007/2008 + Simulink 5/6.7. Основы применения. М.: Солон-Пресс, 2008.
  3. Дьяконов  В. П. MATLAB. Полный самоучитель. М.: ДМК-Пресс, 2012.
  4. Дьяконов  В. П. MATLAB и Simulink для радиоинженеров. М.: ДМК-Пресс, 2011.
  5. Дьяконов  В. П., Пеньков А. А. MATLAB и Simulink в электроэнергетике. М.: Горячая линия – Телеком, 2009.
  6. Дьяконов  В. П . MATLAB — новые возможности в технологии спектроскопии и спектрометрии // Компоненты и технологии. 2010. № 11.
  7. Дьяконов  В. П . MATLAB — новые возможности в технологии осциллографии // Компоненты и технологии. 2009. № 10.
  8. Haykin  S. Adaptive Filter Theory. Prentice-Hall (USA), 1996.

Добавить комментарий

Ваш адрес email не будет опубликован. Обязательные поля помечены *