Как найти частоту сигнала матлаб

SampPerFrame = 1000; Fs = 1000; SW = dsp.SineWave('Frequency',100,. 'SampleRate',Fs,. 'SamplesPerFrame',SampPerFrame); TS = timescope('SampleRate',Fs,. 'TimeSpanSource','property',. 'TimeSpan',0.1,. 'YLimits',[-2, 2],. 'ShowGrid',true); SA = dsp.SpectrumAnalyzer('SampleRate',Fs); tic; while toc < 5 sigData = SW() + 0.05*randn(SampPerFrame,1); SA(sigData); TS(sigData); end release(TS) release(SA)

x = 3.5*ecg(2700).'; y = repmat(sgolayfilt(x,0,21),[1 13]); sigData = y((1:30000) + round(2700*rand(1))).';

TS_ECG = timescope('SampleRate',4000,. 'TimeSpanSource','Auto',. 'ShowGrid', true); TS_ECG(sigData); TS_ECG.YLimits = [-4, 4];

Измерения курсора

Включите Измерения Курсора из вкладки Measurements путем нажатия соответствующей кнопки панели инструментов. Курсоры появляются на осциллографе времени с полем, показывающим изменение вовремя и значение между этими двумя курсорами. Можно перетащить курсоры и использовать их, чтобы измерить время между событиями в форме волны. Когда вы перетаскиваете курсор, время, значение в курсоре появляется. Этот рисунок показывает, как использовать курсоры, чтобы измерить временной интервал между peaks в форме волны ECG. Измерение в поле курсора демонстрирует, что временной интервал между двумя peaks составляет 0,675 секунды, соответствуя сердечному ритму 1,482 Гц или 88,9 ударам/min.

Статистика сигнала и двухуровневые измерения

Можно также выбрать Signal Statistics и различные двухуровневые измерения от вкладки Measurements. Статистика сигнала может использоваться, чтобы определить минимальные и максимальные значения сигнала, а также другие метрики как от пика к пику, среднее значение, медиана и значения RMS. Двухуровневые измерения могут использоваться, чтобы определить информацию о повышении и падающих переходах, аберрациях перехода, перерегулировании и информации об отклонении от номинала, времени урегулирования, ширине импульса и рабочем цикле. Чтобы читать больше об этих измерениях, смотрите, Конфигурируют Time Scope Объект MATLAB.

Измерения частотного диапазона

Этот раздел объясняет, как сделать измерения частотного диапазона со спектром анализатор.

Спектр анализатор обеспечивает следующие измерения:

• Измерения курсора — устанавливают курсоры на отображение спектра.
• Пиковое Средство поиска — отображает максимумы и частоты, на которых они происходят.
• Измерения канала — отображают занимаемую полосу и измерения канала ACPR.
• Измерения искажения — отображают измерения искажения гармонической и интермодуляции.
• Измерения CCDF — отображают дополнительные измерения кумулятивной функции распределения.

Можно включить и отключить эти измерения от спектра панель инструментов анализатора или в меню Tools> Measurements.

Измерения искажения

Чтобы проиллюстрировать использование измерений со спектром анализатор, создайте синусоиду на 2,5 кГц, произведенную на уровне 48 кГц с аддитивным белым Гауссовым шумом. Оцените старший полином (9-я степень) в каждом значении сигналов, чтобы смоделировать нелинейное искажение. Отобразите сигнал в спектре анализатор.

Fs = 48e3; SW = dsp.SineWave('Frequency',2500,. 'SampleRate',Fs,. 'SamplesPerFrame',SampPerFrame); SA_Distortion = dsp.SpectrumAnalyzer('SampleRate',Fs,. 'PlotAsTwoSidedSpectrum',false); y = [1e-6 1e-9 1e-5 1e-9 1e-6 5e-8 0.5e-3 1e-6 1 3e-3]; tic; while toc < 5 x = SW() + 1e-8*randn(SampPerFrame,1); sigData = polyval(y, x); SA_Distortion(sigData); end clear SA_Distortion;

Включите гармонические измерения искажения путем нажатия на соответствующий значок на панели инструментов или путем нажатия на пункт меню Tools> Measurements> Distortion Measurements. В Измерениях Искажения измените значение для Цифры. Гармоники к 9 и проверка флажок Label Harmonics. В панели вы видите значение основного принципа близко к 2 500 Гц и 8 гармоник, а также их ОСШ, SINAD, THD и значений SFDR, на которые ссылаются относительно основной выходной мощности.

Пиковое средство поиска

Можно отследить изменяющиеся во времени спектральные компоненты при помощи Пикового диалогового окна измерения Средства поиска. Можно показать и опционально пометить до 100 peaks. Можно вызвать Пиковое диалоговое окно Средства поиска от пункта меню Tools> Measurements> Peak Finder, или путем нажатия на соответствующий значок на панели инструментов.

. 'SampleRate',Fs,. 'SamplesPerFrame',SampPerFrame); SW2 = dsp.SineWave(1e-1,15e3,0,. 'SampleRate',Fs,. 'SamplesPerFrame',SampPerFrame); SW3 = dsp.SineWave(1e-2,25e3,0,. 'SampleRate',Fs,. 'SamplesPerFrame',SampPerFrame); SA_Peak = dsp.SpectrumAnalyzer('SampleRate',Fs,. 'PlotAsTwoSidedSpectrum',false); tic; while toc < 5 sigData = SW1() + SW2() + SW3() + 1e-4*randn(SampPerFrame,1); SA_Peak(sigData); end clear SA_Peak;

Позвольте Пиковому Средству поиска пометить три частоты синусоиды. Значения частоты и степени в dBm отображены в Пиковой панели Средства поиска. Можно увеличить или сократить максимальное число peaks, задать минимальное пиковое расстояние и изменить другие настройки от панели Настроек в Пиковой панели Измерения Средства поиска.

Чтобы узнать больше об использовании измерений со спектром анализатор, смотрите Спектр пример Измерений Анализатора.

Открытый пример

У вас есть модифицированная версия этого примера. Вы хотите открыть этот пример со своими редактированиями?

Нахождения частоты сигнала

Извлечение сигнала заданной частоты из более сложного сигнала
Нужен совет. Предположим, у меня есть некий сложный сигнал. Я знаю, что в него входит синусоида.

Перенос частоты сигнала
Доброго времени суток! Возникла необходимость смоделировать в симулинке перенос частоты скажем.

Определение частоты сигнала из ДПФ
День добрый! В задаче определения частоты по отклонению фазы изображающего вектора получаю.

Анализ сигнала и разложение его на частоты
Доброго времени суток. Есть сигнал, записанный в wav файл. После Фурье преобразования получаю.

3390 / 1913 / 571
Регистрация: 09.04.2015
Сообщений: 5,365

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

Регистрация: 02.04.2018
Сообщений: 21
А как тогда точно найти частоту?

3390 / 1913 / 571
Регистрация: 09.04.2015
Сообщений: 5,365

Многое зависит от сигнала.
Если у Вас сигнал почти чистая синусоида, то надо определить период и по нему частоту.
Для повышения точности если есть возможность определяют длительность сразу нескольких периодов и затем для определения периода делят на число периодов.
При определении периодов можно использовать расстояние между максимумами / минимумами или пересечение средней линии синусоиды.

Регистрация: 02.04.2018
Сообщений: 21

Сигнал сгенерирован по принципу DTMF , там 2 синусоиды высокой и низкой частоты, мне нужно определить каждую из этих частот

Добавлено через 15 минут
вот у меня уже есть такой код, в котором я разбиваю код на части и достаточно точно вижу частоты, которые мне нужны, как я могу из этого графика найти частоту ?. Я понимаю, что мне не хватает каких-то элементарных знаний в языке матлаб, но к сожалению так устроена учеба в моем универе(

Спасибо за помощь!

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24

[sig1,Fs]=audioread('9234567.wav'); NFFT=512; WINDOW=NFFT; NOVERLAP=NFFT/2; t1=0:Fs/NFFT:Fs/2-Fs/NFFT; subplot(2,1,1) X1 = fft(sig1,NFFT); stem(t1,1/(NFFT/2)*abs(X1(1:NFFT/2)),'.'); subplot(2,1,2) nbWindows = floor(length(sig1)/ NFFT); for i=1:nbWindows coeffs = fft(sig1((i-1)*NFFT+1:i*NFFT)); plot(abs(coeffs(1:NOVERLAP))); end; f = (Fs/NFFT)*(0:NFFT-1); plot(f,abs(X1)); peaks=findpeaks(abs(X1)); axis([500 1500 0 100]) xlabel('f(Hz)') title('Power')

Measurements with MATLAB

В этом вебинаре будет сделан обзор основных возможностей MATLAB в области цифровой обработки сигналов.

Первая часть рассказа будет посвещена спектральному анализу и проектированию фильтров. Особое внимание при этом будет уделено возможностям среды, способствующим скорейшему ее освоению пользователями, плохо знакомыми с MATLAB.

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

Будут освещены следующие аспекты моделирования и реализации систем ЦОС:

• Разработка цифровых фильтров и спектральный анализ в MATLAB;
• Работа с измерительным оборудованием и платами сбора данных;
• Потоковая обработка данных.

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

Recorded: 15 May 2014

Пример расчета реакции сигнала с применением преобразования Фурье в среде МАТЛАБ

При решении задач передачи данных через линии, представленные частотными характеристиками, применяются преобразования Фурье – перевод сигналов из временной области в частотную область и обратно. Среда МАТЛАБ имеет полный набор функций для решения подобных задач. В этой работе разобран пример вычисления в МАТЛАБ реакции сигнала прошедшего через линию, характеристика которой измерена на частотах, не совпадающих с частотой передачи данных. Надеюсь, что этот пример позволит легче разобраться с особенностями технологии преобразования сигналов в среде МАТЛАБ.

Условие задачи

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

Среда, используемая для вычисления и отображения данных – MATLAB R2015а.
В качестве примера исходных данных взяты следующие отношения, опубликованные на сайте www.StatEye.org для версии метода StatEye 3.0 GUI [1, 2, 3].

Скорость передачи данных bps = 10,3125 Гбит/с. Постоянные времени нормированного фильтра второго порядка совпадают, их обратная величина составляет ¾ частоты передачи данных. Сигнальная линия представлена частотной характеристикой. Измерение характеристики выполнено на частотах channel.f = 0,006495:0,0012475:20 ГГц. Заданное число точек дискретизации преобразования Фурье: points = 2^13.

На Рисунок 1 показаны передача данных, последовательность и результаты обработки данных, которые рассмотрены в этой работе. Переход из временной области в частотную область и обратно выполняется при помощи алгоритма Быстрого Преобразования Фурье (БПФ, FFT).

Рисунок 1. Канал передачи данных. Входной сигнал iSignal.Tx, выходной сигнал фильтра iSignal.Filter_out, выход сигнальной линии iSignal.Rx. Представленные на диаграмме характеристики рассмотрены ниже.

Последовательность вычислений

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

Скорость передачи данных в два раза выше частоты, на которой происходит передача данных. Максимальная частота измеренной сигнальной линии max(channel.f) = 20 ГГц. На этой частоте можно передавать данные со скоростью 40 Гбит/с (как 2*max(channel.f)).

Максимальная скорость передачи данных, которая не превышает максимальную скорость передачи по сигнальной линии 40 Гбит/с и кратная скорости передачи bps = 10,3125 Гбит/с, равна fmax = 30,9375 Гбит/с, кратность N = 3 (N = fmax/bps). Далее, fmax используется как предельная частота для расчета реакции сигнала с применением преобразования Фурье.

Перевод входного сигнала в частотную область

Дискретность по времени для построения входного сигнала (бита данных) во временной области Ts = 1/fmax; Ts = 3,232e-11 с. Нормированная по отношению к длительности сигнала, временная шкала time состоит из 2^13 точек (points), шкала включает следующий массив точек time = bps/Ts .* (1:points). Дискретный единичный сигнал при скорости передачи bps = 10,3125 Гбит/с и квантовании с периодом Ts = 1/fmax состоит из трех точек в диапазоне от 10 до 11 единиц нормализованного времени. Сигнал единичной амплитуды можно построить и в любом другом месте временной шкалы, но лучше отступить с краев, чтобы полностью видеть предысторию и переходный процесс выходного сигнала. Импульсный сигнал (бит данных), построенный с использованием следующих команд МАТЛАБ, показан на Рисунок 2.

iSignal.Tx(1:size(time,2)) = 0; t0 = max(find(time<=10)); t1 = max(find(time<11)); iSignal.Tx(t0:t1) = 1.0;

Рисунок 2. Входной импульсный сигнал iSignal.Tx, бит данных.

Перевод сигнала iSignal.Tx в частотную область выполняют следующие БПФ функции.

iSignal.shiftedPSD = fft(iSignal.Tx); iSignal.PSD = fftshift(iSignal.shiftedPSD);

Функция преобразования Фурье fft строит симметричный спектр сигнала в областях положительных и отрицательных частот, максимальная частота которого находится в центре спектра (см. Рисунок 3). Функция fftshift восстанавливает спектр, сдвигая в центр нулевую частоту сигнала как показано на Рисунок 4.

Разрешение частоты спектра равно fs = fmax/points; Частоты спектра лежат в диапазоне от -fmax/2 до fmax/2-fs и равны f = -fmax/2:fs:fmax/2-fs;

Рисунок 3. Амплитудная характеристика сдвинутого спектра сигнала iSignal.Tx полученного с использованием БПФ.

Рисунок 4. Амплитудная характеристика восстановленного спектра сигнала iSignal.Tx показанного на Рисунок 3. Представлено 2^13 отсчетов. Средний отсчет в точке 4097 соответствует нулевой частоте. В левой части (от 1 до 4096 точки) располагаются отрицательные частоты, в правой части (от 4098 до 8192 точки) – область положительных частот.

Передаточная функция нормализованного фильтра нижних частот

В этом примере передаточная функция фильтра второго порядка имеет вид

где Т1 и Т2 – постоянные времени фильтра. Значения частот 1/T1 равны и 1/T2 заданы относительно частоты, на которой передаются данные: 1/T1 = 1/T2 = 0,75*bps (bps = 10,3125 Гбит/с).

Полоса частот нормализованного фильтра

f_nrm =fmax/bps/points.*(-points/2:points/2-1).

s = f_nrm .* j;

Амплитудно-фазовая характеристика нормализованного фильтра для положительных и отрицательных частот нормализованных относительно частоты передачи сигнала показана на Рисунок 5. Логарифмическая амплитудно-частотная характеристика фильтра показана на Рисунок 6.

Рисунок 5. Амплитудно-фазовая характеристика нормализованного фильтра

Рисунок 6. Логарифмическая амплитудно-фазовая частотная характеристика нормализованного фильтра. Синяя штриховая линия показывает положение частоты фильтра со значением 0,75 от частоты, на которой идет передача данных. На этой частоте (1/T1 = 1/T2) коэффициент передачи фильтра второго порядка равен -6 децибел. Красная штриховая линия показывает единичную частоту, на которой идет передача данных.

Перевод результатов измерения сигнальной линии к виду передаточной функции

Измеренная амплитудно-фазовая характеристика сигнальной линии включает 1599 отсчетов в полосе до 20 ГГц с фиксированным шагом 12,475 МГц. Она содержит следующие значения частот: channel.f = 0,006495:0,0012475:20 ГГц. Изначально, сигнальная линия была представлена характеристикой четырехполюсника. Эта характеристика была преобразована и в примере используется в виде одномерной комплексной функции.

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

1. Вычисление амплитуды характеристики линии на нулевой частоте путем ее экстраполяции. Для этого по десяти точкам амплитудной характеристики, ближайших к нулевой частоте, найдены коэффициенты линейного полинома, аппроксимирующего амплитудную характеристику:

[a] = polyfit(channel.f(1:10), channel.abs(1:10), 1);

Найденный второй коэффициент полинома равен амплитуде характеристики на нулевой частоте:

channel.dc = a(2);

2. Фазовая характеристика на нулевой частоте принята равной нулю.

channel.dcPhase = 0.00;

3. Пересчет амплитудной channel.abs и фазовой channel.phase характеристик сигнальной линии со значениями на нулевой частоте выполняется на частоты спектра входного сигнала (f = -fmax/2:fmax/points:fmax/2-fmax/points) с экстраполяцией характеристик в область нулевой и отрицательных частот:

ichannel.abs = interp1([0 channel.f], [channel.dc channel.abs], abs(f), 'linear', 'extrap'); ichannel.phase = interp1([0 channel.f], [channel.dcPhase unwrap(channel.phase)], abs(f), 'linear', 'extrap'); ichannel.s = ichannel.abs .* exp(+j.*ichannel.phase); ichannel.tf = real(ichannel.s) + j*imag(ichannel.s) .* sign(f);

Полученная передаточная функция — амплитудно-фазовая частотная характеристика канала в области низких частот показана на Рисунок 7. Амплитудно-частотные характеристики измеренной сигнальной линии и расчетной передаточной функции в полных частотных диапазонах показаны на Рисунок 8. Эти же характеристики в фазовом пространстве показаны на Рисунок 9.

Рисунок 7. Передаточная функция сигнальной линии в области низких частот. Красными и синими точками обозначены дискретные амплитудная и фазовая характеристики соответственно. Амплитудная характеристика показана в децибелах, фазовая — в радианах. Розовой линией отмечена самая низкая частота измеренной характеристики сигнальной линии. Коэффициент передачи на нулевой частоте равен 0,992.

Рисунок 8. Амплитудно-частотные характеристики сигнальной линии. Синими точками обозначены комплексные данные измеренной линии. Расчетная симметричная зависимость усиления сигнальной линии на частотах спектра входного сигнала выделена красным. В области нулевых частот эта характеристика показана на Рисунок 7.

Рисунок 9. Амплитудно-фазовые частотные характеристики измеренной линии передачи данных и ее нормированного спектра.

Вычисление реакции сигнала

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

Выходной сигнал фильтра во временной области iSignal.Filter_out вычисляется как

TransFunction.PSD = iSignal.PSD .* Filter.PSD_Tx; TransFunction.shiftedPSD = ifftshift(TransFunction.PSD); iSignal.Filter_out = real(ifft(TransFunction.shiftedPSD));

Выходной сигнал линии iSignal.Rx равен произведению спектра входного сигнала на передаточные функции фильтра и сигнальной линии с последующим переводом полученного сигнала из частотной области во временную область.

TransFunction.PSD = TransFunction.PSD .* ichannel.tf; TransFunction.shiftedPSD = ifftshift(TransFunction.PSD); iSignal.Rx = real(ifft(TransFunction.shiftedPSD)); 

Реакция фильтра на входной идеальный импульс и реакция канала показаны на Рисунок 10.

Рисунок 10. Выходной сигнал фильтра (красный график) и выходной сигнал линии передачи данных (зеленый график). Входной сигнал фильтра – единичный импульс показан на Рисунок 2. Входом сигнальной линии является выходной сигнал фильтра.

Приложение. Используемый m-код MATLAB

Листинг

clear all %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Ini data %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% bps = 1.03125e+10; FilterParam = [0.75 0.75]; points = 2^13; load('channel'); N = floor(max(channel.f)*2/bps); fmax = N*bps; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Signal %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % normalise all the scales for the bit rate time = bps/fmax .* (1:points); iSignal.Tx(1:size(time,2)) = 0; t0 = max(find(time<=10)); t1 = max(find(time<11)); iSignal.Tx(t0:t1) = 1.0; figure plot(time(1:t1+10), iSignal.Tx(1:t1+10),'b'); hold on plot(time(1:t1+10), iSignal.Tx(1:t1+10),'xb'); grid on xlabel('Normalised Time, tick Ts = 1/fmax'); ylabel('Normalised Amplitude'); title(['Pulse, data bit']); iSignal.shiftedPSD = fft(iSignal.Tx); figure plot(abs(iSignal.shiftedPSD),'c'); grid on xlabel('points, num'); ylabel('Amplitude'); title(['abs(fft(iSignal.Tx))']); iSignal.PSD = fftshift(iSignal.shiftedPSD); figure plot(abs(iSignal.PSD),'r'); grid on xlabel('points, num'); ylabel('Amplitude'); title(['abs(fftshift(fft(iSignal.Tx)))']); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Filter %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% f_nrm =fmax/bps/points.*(-points/2:points/2-1); s = f_nrm .* j; Filter_PSD = 1 ./(1 + s/FilterParam(1)) ./ (1 + s/FilterParam(2)); figure [AX,H1,H2] = plotyy (f_nrm, abs(Filter_PSD), f_nrm, phase(Filter_PSD)); hold(AX(1)); hold(AX(2)); set(H1,'LineWidth',2); grid(AX(2),'on'); xlabel('Normalised Frequency (Hz)'); set(get(AX(1),'Ylabel'),'String','Gain'); set(get(AX(2),'Ylabel'),'String','Phase, rad'); title(['Twopole filter [' sprintf(' %3.2f ', FilterParam) '] normalised to baud rate frequency']); figure plot_handles_Filter = plot(f_nrm(points/2 + 1:points), 20*log10(abs(Filter_PSD(points/2 + 1:points))), 'r', 'linewidth', 2); hold on stem_handles_br = stem(1, 20*log10(abs(Filter_PSD(max(find(f_nrm < 1))))), '-.ro'); hold on stem_handles_c = stem(FilterParam, [20*log10(abs(Filter_PSD(max(find(f_nrm < FilterParam(1)))))) 20*log10(abs(Filter_PSD(max(find(f_nrm < FilterParam(2))))))], '-.bo'); grid legend_handles = [plot_handles_Filter, stem_handles_br(1), stem_handles_c(1)]; legend(legend_handles, 'transfer function', 'filter attenuation at normalised baud rate', 'filter attenuation at normalised cutoff frequency', 3); xlabel('Normalised Frequency (Hz)'); ylabel('Magnitude (dB)'); title(['Twopole filter [' sprintf(' %3.2f ', FilterParam) '] normalised to baud rate frequency']); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Channel %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % create negative frequencies, convert data to complex value, taking care about negative frequency channel.abs = abs(channel.s); channel.phase = angle(channel.s); %channel.s = channel.abs .* exp(+j.*channel.phase); [a] = polyfit(channel.f(1:10), channel.abs(1:10), 1); channel.dc = a(2); channel.dcPhase = 0.00; fs = fmax/points; % frequency step f = -fmax/2:fs:fmax/2-fs; % frequency matrix % create new data structure with linearly interpolated data ichannel.abs = interp1([0 channel.f], [channel.dc channel.abs], abs(f), 'linear', 'extrap'); ichannel.phase = interp1([0 channel.f], [channel.dcPhase unwrap(channel.phase)], abs(f), 'linear', 'extrap'); % correct for negative frequencies ichannel.s = ichannel.abs .* exp(+j.*ichannel.phase); ichannel.tf = real(ichannel.s) + j*imag(ichannel.s) .* sign(f); figure disp_points = 2*round(channel.f(1)/fs); stem_handles_br = stem(channel.f(1), angle(ichannel.tf(max(find(f < channel.f(1))))), '-.mo'); hold on plot_abs = plot(f(points/2-disp_points:points/2+disp_points), 20*log10(abs(ichannel.tf(points/2-disp_points:points/2+disp_points))), '.r', 'linewidth', 3); hold on plot_phase = plot(f(points/2-disp_points:points/2+disp_points), angle(ichannel.tf(points/2-disp_points:points/2+disp_points)), '.b', 'linewidth', 3); grid legend_handles = [plot_abs, plot_phase, stem_handles_br(1)]; legend(legend_handles, 'absolute value (dB)', 'phase (rad)', 'min data freq', 3); xlabel('Relative Frequency (Hz)'); ylabel('Magnitude'); title(sprintf('dc extrapolation. dc trans function=%4.3f, dc phase=%4.3f rad', abs(ichannel.tf(points/2+1)), angle(ichannel.tf(points/2+1)))); figure plot(channel.f, 20*log10(channel.abs), '.r', 'linewidth', 3); hold on plot(f, 20*log10(ichannel.abs), 'g'); grid on legend('Measured Data', 'Interpolated Data', 3); xlabel('Frequency (Hz)'); ylabel('Magnitude (dB)'); title(['Chаnnel interpolated Data : ']); figure plot3(channel.f, real(channel.s), imag(channel.s),'r'); hold on plot3(f, real(ichannel.tf), imag(ichannel.tf),'g'); grid on legend('Measured Data', 'Interpolated Data'); xlabel('Frequency in Hz'); ylabel('Re(fwd transfer)'); zlabel('Im(fwd transfer)'); title(['Chаnnel interpolated Data : ']); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Response %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % filter Output TransFunction.PSD = iSignal.PSD .* Filter_PSD; TransFunction.shiftedPSD = ifftshift(TransFunction.PSD); iSignal.Filter_out = real(ifft(TransFunction.shiftedPSD)); % pass through channel TransFunction.PSD = TransFunction.PSD .* ichannel.tf; TransFunction.shiftedPSD = ifftshift(TransFunction.PSD); iSignal.Rx = real(ifft(TransFunction.shiftedPSD)); figure plot(time, iSignal.Filter_out,'r'); hold on [max_Tx, time_maxTx] = max(iSignal.Filter_out); [min_Tx, time_minTx] = min(iSignal.Filter_out); [max_Rx, time_maxRx] = max(iSignal.Rx); dtime_p5= round((time_maxRx - time_maxTx)*time(1) -1); plot(time - dtime_p5, iSignal.Rx,'g'); hold on plot(time, iSignal.Filter_out,'rx'); axis([(time_maxTx*time(1) - 3) (time_maxTx*time(1) + 5) (min_Tx-0.15) (max_Tx+0.1)]) grid on legend('Filter out','Rx', 2); xlabel('Normalised Time'); ylabel('Normalised Amplitude'); title(sprintf('Transmit pulse (Tx) max= %4.3f; Response (Rx) max (h0)= %4.3f', max(iSignal.Filter_out), max(iSignal.Rx))); 

Библиографический список

1. IEEE802.3ap. 10.3125Gbps NRZ Simulation Results Using “StatEye” and “Signal to Interference Model” on Cascaded Channel Components. Shannon Sawyer and Charles Moore / Agilent Technologies. January 24, 2005 www.ieee802.org/3/ap/public/jan05/sawyer_01_0105.pdf

• Передача данных
• обработка сигнала
• БПФ