2013-03-03 1 views
0

Обычно я вычислить спектр сигнала с использованием pmtm:Как можно векторизовать вычисление спектра с несколькими конусами?

signal = rand(1000,1); 
NW = 4; 
Fr = 1:50; 
Fs = 200; 
[p, fr] = pmtm(signal, NW, Fr, Fs); 

Однако я ищу способ векторизации это так, я могу вычислить несколько спектров одновременно. Я пробовал:

signal = rand(1000,10); %<--- notice I have 10 columns instead of 1 
NW = 4; 
Fr = 1:50; 
Fs = 200; 
[p, fr] = pmtm(signal, NW, Fr, Fs); 

но он производит ошибку, которая на самом деле не говорит мне, что я сделал неправильно. Я знаю, что могу заключить звонок в pmtm в цикле.

Здесь ошибка:

Error using .* Matrix dimensions must agree.

Error in pmtm>mtm_spectrum (line 231) [Xx,w] = computeDFT(E(:,1:k).*x(:,ones(1,k)),nfft,Fs);

Error in pmtm (line 142) [S,k,w] = mtm_spectrum(x,params);

Это заставляет меня подозревать, что не существует Векторизованный способ добиться того, чего я хочу. Я надеялся, что кто-то здесь будет знать, как это сделать.

+0

'mtm_spectrum' заставляет' x' быть вектором в строке 'x = x (:); 'Если' x' является матрицей, она будет иметь большее значение размера n_col, и именно поэтому вы получаете ошибку размера матрицы. Вы пробовали arrayfun? Это псевдо-векторизованное решение, но может немного ускорить время выполнения. Что-то вроде ' arrayfun (@ (x) pmtm (сигнал (:, x), NW, Fr, Fs), 1: размер (сигнал, 2)) '. Но вы можете получить только одну выходную переменную с этим. – yuk

+0

Похоже, что E 'и' x' не соответствуют точно таким же размерам ... –

+1

Вы, вероятно, не можете избежать цикла for здесь, что довольно дорого, потому что 'pmtm' медленный. Если производительность ваша, то вы можете значительно ускорить ее предварительно вычисляя слепые последовательности с помощью 'dpss' и передавая их на' pmtm'. –

ответ

0

Что заставляет вас думать, что pmtm предназначен для работы с матрицами? Меню помощи специально ожидает вектор:

>> help pmtm 
pmtm Power Spectral Density (PSD) estimate via the Thomson multitaper 
    method (MTM). 
    Pxx = pmtm(X) returns the PSD of a discrete-time signal vector X in 
    the vector Pxx. Pxx is the distribution of power per unit frequency. 
    The frequency is expressed in units of radians/sample. pmtm uses a 
    default FFT length equal to the greater of 256 and the next power of 
    2 greater than the length of X. The FFT length determines the length 
    of Pxx. 
... 

Если вы после выполнения, вы можете получить спектральное представление вашего сигнала с 2D FFT (fft2), а затем вычислить распределение мощности по своему усмотрению. Это позволит вам работать с 2D-массивами без каких-либо циклов, но это будет означать для вас более точное кодирование :-(

Смежные вопросы