Реализуйте основной синтезатор цифровой музыки и используйте его, чтобы проигрывать традиционную песню в расположении с тремя речью. Задайте частоту дискретизации 2 кГц. Сохраните песню как расписание MATLAB®.

fs = 2e3;
t = 0:1/fs:0.3-1/fs;

l = [0 130.81 146.83 164.81 174.61 196.00 220 246.94];
m = [0 261.63 293.66 329.63 349.23 392.00 440 493.88];
h = [0 523.25 587.33 659.25 698.46 783.99 880 987.77];
note = @(f,g) [1 1 1]*sin(2*pi*[l(g) m(g) h(f)]'.*t);

mel = [3 2 1 2 3 3 3 0 2 2 2 0 3 5 5 0 3 2 1 2 3 3 3 3 2 2 3 2 1]+1;
acc = [3 0 5 0 3 0 3 3 2 0 2 2 3 0 5 5 3 0 5 0 3 3 3 0 2 2 3 0 1]+1;

song = [];
for kj = 1:length(mel)
    song = [song note(mel(kj),acc(kj)) zeros(1,0.01*fs)];
end
song = song'/(max(abs(song))+0.1);

% To hear, type sound(song,fs)

tune = timetable(seconds((0:length(song)-1)'/fs),song);

Откройте Signal Analyzer и перетащите расписание от браузера Рабочей области до таблицы Signal. Нажмите Display Grid ▼, чтобы создать two-two сетку отображений. Выберите лучшие два отображения и нижнее левое отображение и нажмите кнопку Spectrum, чтобы добавить представление спектра. Выберите нижнее правое отображение, нажмите Time-Frequency, чтобы добавить представление спектрограммы и нажать Time, чтобы удалить представление времени. Перетащите песню ко всем четырем отображениям. Выберите нижнее правое отображение, и во вкладке Spectrogram, задайте разрешение времени 0,31 вторых и 0%-х перекрытий (на 310 мс) между смежными сегментами. Установите Пределы Степени к $-50$ дБ и $-10$ дБ.

На вкладке Analyzer нажмите Duplicate три раза, чтобы создать три копии песни. Переименуйте копии как high, medium, и low путем двойного клика по столбцу Имени в таблице Signal. Переместите копии в лучшие два и нижние левые отображения.

Предварительно обработайте дублирующиеся сигналы с помощью фильтров.

На каждом из трех отображений, содержащих отфильтрованные сигналы:

Выберите три отфильтрованных сигнала путем нажатия на их столбец Имени в таблице Signal. На вкладке Analyzer нажмите Export и сохраните сигналы в MAT-файл под названием music.mat. В MATLAB загрузите файл к рабочей области. Постройте спектры трех сигналов.

load music

pspectrum(low)
hold on
pspectrum(medium)
pspectrum(high)
hold off

% To hear the different voices, type 
% sound(low.low,fs), pause(5), sound(medium.medium,fs), pause(5), sound(high.high,fs)

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

Всю информацию несут в сигнале с действительным знаком $\mathit{s}\left(\mathit{t}\right)$ может быть представлен соответствующим конвертом комплекса lowpass:

$\mathit{s}\left(\mathit{t}\right)=\mathrm{Ре}\left\{\mathit{g}\left(\mathit{t}\right){\mathit{e}}^{\mathit{j}2\pi {\mathit{f}}_{\mathrm{c}}\mathit{t}}\right\}=\mathit{i}\left(\mathit{t}\right)\mathrm{cos2}\pi {\mathit{f}}_{\mathit{c}}\mathit{t}+\mathit{q}\left(\mathit{t}\right)\mathrm{sin2}\pi {\mathit{f}}_{\mathit{c}}\mathit{t}$.

В этом уравнении:

Комплексный конверт модулируется к несущей частоте и отправляется через канал. В получателе шумная форма волны демодулируется с помощью несущей частоты.

Изменение фазы из-за несущей частоты предсказуемо и таким образом не передает информации. Комплексный конверт не включает изменение фазы и может быть произведен на более низком уровне.

Сгенерируйте сигнал, комплексный конверт которого состоит из синусоиды и щебета. Синфазный компонент является синусоидой с частотой 19 Гц. Квадратурный компонент является квадратичным щебетом чьи частотные диапазоны от 61 Гц до 603 Гц. Сигнал производится на уровне 2 кГц в течение 1 секунды.

fs = 2e3;
t = (0:1/fs:1-1/fs)';
inph = sin(2*pi*19*t);
quad = chirp(t-0.6,61,t(end),603,'quadratic');

Вычислите комплексный конверт и сохраните его как расписание MATLAB® частоты дискретизации fs.

env = inph + 1j*quad;
g = timetable(env,'SampleRate',fs);

Откройте Signal Analyzer и перетащите комплексный конверт от браузера Рабочей области до таблицы Signal. Отображение показывает синфазные компоненты и квадратурные компоненты конверта как линии того же оттенка и насыщения, но различной яркости. Первый цвет линии представляет синфазный компонент, и второй цвет линии представляет квадратурный компонент. Нажмите Spectrum ▼ на вкладке Display и выберите Spectrum. Отображения приложения набор осей со спектром сигнала. Комплексный конверт имеет двухсторонний спектр, отображенный как линия того же цвета синфазного компонента комплексного конверта.

На вкладке Display нажмите Panner, чтобы активировать регулятор панорамы. Используйте регулятор панорамы, чтобы создать окно изменения масштаба между 300 мс и 720 мс. Перетащите окно изменения масштаба так, чтобы оно было сосредоточено на уровне 0 Гц. Спектр имеет импульс на уровне 0,19 кГц и более широкий конический профиль на более высоких частотах. Область отрицательной частоты спектра является зеркальным отображением области положительной частоты.

Модулируйте сигнал с помощью несущей частоты 200 Гц. Умножьтесь $\sqrt{2}$ так, чтобы степень модулируемого сигнала равнялась степени исходного сигнала. Добавьте белый Гауссов шум, таким образом, что отношение сигнал-шум составляет 40 дБ.

fc = 200;
mod = sqrt(2)*real(env.*exp(2j*pi*fc*t));

SNR = 40;
mod = mod + randn(size(mod))*std(mod)/db2mag(SNR);
s = timetable(mod,'SampleRate',fs);

Нажмите Display Grid ▼, чтобы добавить второе отображение. Перетащите модулируемый сигнал к таблице Signal и введите информацию времени. Модуляция переместила спектр в положительные частоты, сосредоточенные на несущей частоте.

Вычислите аналитический сигнал и демодулируйте сигнал путем умножения аналитического сигнала с отрицательным экспоненциалом с комплексным знаком частоты 200 Гц.

dem = hilbert(mod).*exp(-2j*pi*fc*t)/sqrt(2);

Нажмите Display Grid ▼, чтобы создать три одним сетка отображений. Перетащите демодулируемый сигнал к таблице Signal. Добавьте время информация к комплексному конверту путем нажатия на Time Values во вкладке Analyzer. Двухсторонний спектр показывает восстановленное синфазное и квадратурные компоненты сгенерированного модулированного сигнала.

Нажмите Display Grid ▼, чтобы создать один за другим сетка отображений и построить демодулируемый сигнал. Нажмите Data Cursors ▼ и выберите Two. Установите курсоры временного интервала на уровне 300 мс и 900 мс, таким образом, они заключают спектральный peaks. Нажмите Extract Signals ▼ и выберите Between Time Cursors. Установите флажок Preserve Start Time. Очистите отображение и постройте извлеченный сигнал. Приложение извлекает и синфазные компоненты и квадратурные компоненты демодулируемого сигнала в необходимой области. Выберите извлеченный сигнал путем нажатия на его столбец Имени в таблице Signal. На вкладке Analyzer нажмите Export и сохраните сигнал в MAT-файл под названием dem_ROI.mat.

Загрузите dem_ROI файл к рабочей области. Вычислите демодулируемое синфазное и квадратурные компоненты путем взятия действительных и мнимых частей извлеченного сигнала. Храните информацию времени извлеченного сигнала в переменной t_dem времени.

load dem_ROI
inph_dem = real(dem_ROI);
quad_dem = imag(dem_ROI);
t_dem = 0.3+(0:length(dem_ROI)-1)/fs;

Сравните переданные формы волны и извлеченные необходимые области. Также сравните их спектры.

subplot(2,1,1)
plot(t,inph,t_dem,inph_dem,'--')
legend('Transmitted Inphase Signal','Received Inphase Signal')

subplot(2,1,2)
plot(t,quad,t_dem,quad_dem,'--')
legend('Transmitted Quadrature Signal','Received Quadrature Signal')

figure
subplot(2,1,1)
pspectrum(inph,fs)
hold on
pspectrum(inph_dem,fs)
legend('Transmitted Inphase Signal','Received Inphase Signal')
hold off

subplot(2,1,2)
pspectrum(quad,fs)
hold on
pspectrum(quad_dem,fs)
legend('Transmitted Quadrature Signal','Received Quadrature Signal')
hold off

Загрузите файл данных, содержащий импульс эхолотирования, испускаемый большой коричневой битой (Eptesicus fuscus) и измеренный с интервалом выборки 7 микросекунд. Создайте расписание MATLAB® с помощью сигнала и информации времени.

load batsignal

t = (0:length(batsignal)-1)*DT;
sg = timetable(seconds(t)',batsignal);

Откройте Signal Analyzer и перетащите расписание от браузера Рабочей области до таблицы Signal. Нажмите Display Grid ▼, чтобы создать два рядом друг с другом отображения. Выберите каждое отображение и нажмите кнопку Time-Frequency, чтобы добавить представление спектрограммы.

Перетащите расписание к обоим отображениям.

Выберите вкладку Spectrogram. На отображении в праве проверяйте Переприсвоение. Для каждого отображения:

Переприсвоенная спектрограмма ясно показывает три частотно-временных гребня. Чтобы отследить гребни, выберите отображение в праве. На вкладке Display нажмите Generate Script и выберите Spectrogram Script. Скрипт появляется в Редакторе.

% Compute spectrogram

% Generated by MATLAB(R) 9.7 and Signal Processing Toolbox 8.2.
% Generated on: 26-Dec-2018 17:21:44

% Parameters
timeLimits = seconds([0 0.002793]); % seconds
frequencyLimits = [0 71428.57]; % Hz
leakage = 0.9;
timeResolution = 0.00028; % seconds
overlapPercent = 85;
reassignFlag = true;

%%
% Index into signal time region of interest
sg_batsignal_ROI = sg(:,'batsignal');
sg_batsignal_ROI = sg_batsignal_ROI(timerange(timeLimits(1),timeLimits(2),'closed'),1);

% Compute spectral estimate
% Run the function call below without output arguments to plot the results
[P,F,T] = pspectrum(sg_batsignal_ROI, ...
    'spectrogram', ...
    'FrequencyLimits',frequencyLimits, ...
    'Leakage',leakage, ...
    'TimeResolution',timeResolution, ...
    'OverlapPercent',overlapPercent, ...
    'Reassign',reassignFlag);

Запустите скрипт. Постройте переприсвоенную спектрограмму.

mesh(seconds(T),F,P)
xlabel('Time')
ylabel('Frequency')
axis tight
view(2)
colormap pink

Используйте tfridge функционируйте, чтобы отследить гребни.

[fridge,~,lridge] = tfridge(P,F,0.01,'NumRidges',3,'NumFrequencyBins',10);

hold on
plot3(seconds(T),fridge,P(lridge),':','linewidth',3)
hold off

Благодаря Кертису Кондону, Кену Вайту и Аль Фэну из Центра Бекмана в Университете Иллинойса для bat данных и разрешения использовать его в этом примере.

Датчики могут возвратить отсеченные показания, если данные больше, чем данная точка насыщения. Чтобы восстановить показания, можно соответствовать полиному через точки, смежные с влажными интервалами. Запишите функцию, которая выполняет реконструкцию, и интегрируйте ее в Signal Analyzer.

Сгенерируйте сигнал с тремя каналами, произведенный на уровне 1 кГц в течение 14 секунд. Сигнал имеет несколько peaks различных размеров и форм. Датчик, который читает сигнал, насыщает на уровне 0,1 В.

fs = 1000;
t = 0:1/fs:14-1/fs;

sig = [chirp(t-1,0.1,17,2,'quadratic',1).*sin(2*pi*t/5);
    chirp(t-2,2,2,2.1,'quadratic',100).*exp(-(t-6.5).^2/20).*sin(2*pi*t*2);
    0.85*besselj(0,5*(sin(2*pi*(t+1.5).^2/20).^2)).*sin(2*pi*t/9)]';

sigsat = sig;
stv = 0.1;
sigsat(sigsat >= stv) = stv;

Откройте Signal Analyzer и перетащите исходный сигнал и влажный сигнал к таблице Signal. Перетащите каждый исходный и влажный канал к его собственному отображению.

Запишите функцию, которая использует полином, чтобы восстановить peaks сигнала:

function [y,tOut] = declip(x,tIn,varargin)
% Declip saturated signal by fitting a polynomial

    % Initialize the output signal

    y = x;

    % For signals with no time information, use sample numbers as abscissas
    
    if isempty(tIn)
        tOut = [];
        t = (1:length(x))';
    else
        t = tIn;
        tOut = t;
    end
    
    % Specify the degree of the polynomial as an optional input argument
    % and provide a default value of 4
    
    if nargin<3
        ndx = 4;
    else
        ndx = varargin{1};
    end

    % To implement your algorithm, you can use any MATLAB or Signal
    % Processing Toolbox function
    
    % Find the intervals where the signal is saturated and generate an 
    % array containing the interval endpoints
    idx = find(x==max(x)); 
    fir = [true;diff(idx)~=1];
    ide = [idx(fir) idx(fir([2:end 1]))];
    % For each interval, fit a polynomial of degree ndx over the ndx+1 points
    % before the interval and the ndx+1 points after the interval
    for k = 1:size(ide,1)
        bef = ide(k,1); aft = ide(k,2);
        intv = [bef-1+(-ndx:0) aft+1+(0:ndx)];
        [pp,~,mu] = polyfit(t(intv),x(intv),ndx);
        y(bef:aft) = polyval(pp,t(bef:aft),[],mu);
    end

end

Добавьте функцию в Signal Analyzer как пользовательская функция предварительной обработки. На вкладке Analyzer нажмите Preprocessing ▼ и выберите Add Custom Function. Введите имя функции и описание. Вставьте текст своей функции в окне редактора, которое появляется. Сохраните файл. Функция появляется в галерее предварительной обработки.

Продемонстрируйте, что функция, которую вы создали, восстанавливает влажные области.

Проверьте, что функция предварительной обработки работает, когда сигналы имеют информацию времени.

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

Используйте Signal Analyzer, чтобы вычислить огибающую спектра вибрации подшипника, сигнализируют и ищут дефекты. Сгенерируйте скрипты MATLAB® и функции, чтобы автоматизировать анализ.

[t,xBPFO,xBPFI,bpfi] = bearingdata;

envspec = preprocess(xBPFI,t);

% Compute power spectrum

% Generated by MATLAB(R) 9.6 and Signal Processing Toolbox 8.2.
% Generated on: 12-Nov-2018 15:13:34

% Parameters
timeLimits = [0 0.9999]; % seconds
frequencyLimits = [0 5000]; % Hz

%%
% Index into signal time region of interest
envspec_ROI = envspec(:);
sampleRate = 10000; % Hz
startTime = 0; % seconds
minIdx = ceil(max((timeLimits(1)-startTime)*sampleRate,0))+1;
maxIdx = floor(min((timeLimits(2)-startTime)*sampleRate,length(envspec_ROI)-1))+1;
envspec_ROI = envspec_ROI(minIdx:maxIdx);

% Compute spectral estimate
% Run the function call below without output arguments to plot the results
[Penvspec_ROI, Fenvspec_ROI] = pspectrum(envspec_ROI,sampleRate, ...
    'FrequencyLimits',frequencyLimits);

plot(Fenvspec_ROI,(Penvspec_ROI))
hold on
[X,Y] = meshgrid((1:10)*bpfi,ylim);
plot(X,Y,':k')
hold off
xlim([0 10*bpfi])

function y = preprocess(x,tx)
%  Preprocess input x
%    This function expects an input vector x and a vector of time values
%    tx. tx is a numeric vector in units of seconds.

% Generated by MATLAB(R) 9.6 and Signal Processing Toolbox 8.2.
% Generated on: 12-Nov-2018 15:09:44

y = detrend(x,'constant');
Fs = 1/mean(diff(tx)); % Average sample rate
y = bandpass(y,[2250 3750],Fs,'Steepness',0.85,'StopbandAttenuation',60);
[y,~] = envelope(y);
y = detrend(y,'constant');
end

function [t,xBPFO,xBPFI,bpfi] = bearingdata

p = 0.12;
d = 0.02;
n = 8;
th = 0;
f0 = 25;
fs = 10000;

t = 0:1/fs:1-1/fs;
z = [1 0.5 0.2 0.1 0.05]*sin(2*pi*f0*[1 2 3 4 5]'.*t);

xHealthy = z + randn(size(z))/10;

bpfo = n*f0/2*(1-d/p*cos(th));

tmp = 0:1/fs:5e-3-1/fs;
xmp = sin(2*pi*3000*tmp).*exp(-1000*tmp);

xBPFO = xHealthy + pulstran(t,0:1/bpfo:1,xmp,fs)/4;

bpfi = n*f0/2*(1+d/p*cos(th));

xBPFI = xHealthy + pulstran(t,0:1/bpfi:1,xmp,fs)/4;

end