2013-06-16 15 views
8

sto il porting del codice MATLAB per Python usando SciPy e rimasto bloccato con la seguente riga:Periodogramma in Octave/Matlab vs SciPy

Matlab/Octave codice

[Pxx, f] = periodogram(x, [], 512, 5) 

codice Python

f, Pxx = signal.periodogram(x, 5, nfft=512) 

Il problema è che ottengo output diversi sugli stessi dati. Più specificamente, i vettori Pxx sono diversi. Ho provato diverse finestre per signal.periodogram, eppure non ho avuto fortuna (e sembra che la finestra box car di default scypy sia la stessa della finestra rettangolare di matlab di default) Un altro strano comportamento è che in python, il primo elemento di Pxx è sempre 0, non importa quali dati l'input è.

Mi manca qualcosa? Qualsiasi consiglio sarebbe molto apprezzato!


semplice codice/Octave Matlab con dati reali: http://pastebin.com/czNeyUjs
Python Semplice + codice SciPy con i dati reali: http://pastebin.com/zPLGBTpn

+0

Una trama sarebbe molto utile. –

risposta

5

Dopo la ricerca di un'ottava e il codice sorgente periodogramma di SciPy ho scoperto che usano algoritmo diverso per calcolare stima della densità spettrale di potenza. Octave (e MATLAB) use FFT, mentre il periodogramma di scipy usa lo Welch method.

Come ha già indicato @georgesl, l'output sembra abbastanza simile, ma tuttavia differisce. E per porting la ragione era critica. Alla fine, ho semplicemente scritto una piccola funzione per calcolare la stima PSD usando FFT, e ora l'output è lo stesso. Secondo il test timeit, funziona ~ 50% più veloce (1.9006s vs 2.9176s su un ciclo con 10.000 iterazioni). Penso che sia dovuto al fatto che FFT sia più veloce di Welch nell'implementazione di scipy, di essere semplicemente più veloce.

Grazie a tutti coloro che hanno mostrato interesse.

0

Hai guardano i risultati?

Outputs

Le lievi differenze tra i due risultati possono proviene da ottimizzazioni/finestre di default/implementazioni/qualunque ecc

1

Ho affrontato lo stesso problema, ma poi mi sono imbattuto la documentazione di scipy's periodogram

Come si vedrebbe lì che detrend = 'costante' è l'argomento di default. Ciò significa che python sottrae automaticamente la media dei dati di input da ciascun punto. (Read here). Mentre Matlab/Octave non fanno nulla del genere. Credo che questo sia il motivo per cui le uscite sono diverse. Prova a specificare detrend = False, mentre chiami il periodogramma di scipy dovresti ottenere lo stesso risultato di Matlab.

2

Dopo aver letto la documentazione Matlab e Scipy, un altro contributo ai diversi valori potrebbe essere che utilizzano diverse funzioni della finestra predefinita. Matlab usa una finestra di Hamming e Scipy usa un Hanning. Le due funzioni della finestra e simili ma non identiche.