Внедряете avgpower Matlab в Октаве? - PullRequest
Новая
       40

Внедряете avgpower Matlab в Октаве?

1 голос
/ 19 марта 2009

Люди,

Matlab 2007b (7.5.0) имеет функцию avgpower. Смотрите здесь :

"Метод avgpower использует приближение прямоугольника к интегралу рассчитать среднюю мощность сигнала, используя данные PSD, хранящиеся в объект.

"Метод avgpower возвращает среднюю мощность сигнала, который площадь под кривой PSD. "

Пример вызова: 

    numSamples = 10000
    frequency = 20
    amplitude = 10
    Fs = 1000
    t = [0:1/numSamples:1];
    sig = amplitude * sin(2*pi*frequency*t);
    h = spectrum.periodogram('rectangular');
    hopts = psdopts(h, signal);
    set(hopts,'Fs',Fs);
    p = psd(h,signal,hopts);
    lower = 12
    upper = 30
    beta_power = p.avgpower([lower upper]);

Я хочу повторить такую ​​функциональность в Octave. Функция «pwelch» кажется возможной. Для остроумия: 

    ...
    sig = amplitude * sin(2*pi*frequency*t);
    pwelch('R12+');
    [spectra, freq]=pwelch(signal, [], [], [], Fs, plot_type='dB');

Теперь я думаю, что спектры имеют значения y PSD, а freq имеет x ценности. Таким образом, я мог найти образцы в частоте, которые попадают между "нижней" а "верхний" и ... эр, усреднить соответствующие значения в спектрах? Я довольно размышляю об этом.

Более того, значения в "freq" не обязательно соответствуют моим желаемый верх и низ, и я не уверен, что с этим делать. Что, если нижний или верхний попадают прямо в середину широкого диапазона частот? Например, взять ли половину ячейки (то есть линейно интерполировать)?

Также возможно получить одно значение из некоторого БПФ вместо использования pwelch.

Предложения

1 Ответ

2 голосов
/ 20 марта 2009

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

function[avgp] = oavgpower(signal, sampling_freq, lowfreq, highfreq, window)

[spectra, freq]=pwelch(signal, window, [], [], sampling_freq);

idx1=max(find(freq = highfreq));

% Index and actual frequency of lower and upper bins
%idx1
%freq(idx1)
%idx2
%freq(idx2)

% 0: don't include the last bin
width = [diff(freq); 0];

pvec = width(idx1:idx2).*spectra(idx1:idx2);
avgp = sum(pvec);
...