Cálculo de espectrogramas con Signal Processing Toolbox

Signal Processing Toolbox™ proporciona tres funciones que calculan el espectrograma de una señal no estacionaria. Cada una de las funciones tiene diferentes argumentos de entrada, valores predeterminados y salidas. La mejor opción para usted depende de su aplicación concreta.

Una señal no estacionaria es una señal cuyo contenido de frecuencia cambia con el tiempo. La transformada de Fourier de tiempo corto (STFT) se utiliza para analizar cómo cambia este contenido de frecuencia a medida que evoluciona la señal. La magnitud al cuadrado de la STFT se conoce como la representación de tiempo-frecuencia del espectrograma de la señal.

El espectrograma es solo una de las posibles representaciones de tiempo-frecuencia. Para obtener una visión general de otras representaciones de tiempo-frecuencia disponibles en Signal Processing Toolbox y Wavelet Toolbox™, consulte Time-Frequency Gallery. Para tratar señales estacionarias con la función periodogram, consulte Estimaciones de densidad del espectro de potencia utilizando FFT.

Funciones para el cálculo de espectrogramas

Signal Processing Toolbox tiene estas funciones que se pueden utilizar para calcular el espectrograma:

• spectrogram: diseñada para ofrecer la máxima flexibilidad. Calcula la STFT y las densidades espectrales de potencia o los espectros de potencia segmento a segmento. Admite reasignación.

• stft: diseñada para la invertibilidad y el máximo control. Calcula la STFT. La utilizan dlstft y stftLayer. Admite entradas multicanal. Su equivalente, la función istft, calcula la STFT inversa.

• pspectrum: diseñada para ofrecer facilidad de uso. Calcula los espectros de potencia. Se utiliza en los scripts de análisis generados por Signal Analyzer. Puede calcular espectros de señales estacionarias y espectros de persistencia. Admite reasignación.

La función spectrogram se utiliza como referencia en la siguiente descripción.

Categoría	Parámetros	Función
		spectrogram	stft	pspectrum
Entrada	Señal	Vector con Nx elementos	Vector con Nx elementos Matriz con Nx filas Horario con Nx filas	Vector con Nx elementos Horario con Nx filas
	Ventana, g(n)	Segundo argumento posicional (Predeterminado: Ventana de Hamming)	Argumento nombre-valor Window (Predeterminado: Ventana periódica de Hann)	Solo ventana de Kaiser
	Longitud de la ventana, M	Especificado como número de muestras (Predeterminado: ⌊Nx/4.5⌋)	Especificado como número de muestras (Predeterminado: 128)	Argumento nombre-valor TimeResolution
	Manchado, ℓ	Depende de la ventana Si se utiliza una ventana de Kaiser, debe ajustarse utilizando el factor de forma β	Depende de la ventana Si se utiliza una ventana de Kaiser, debe ajustarse utilizando el factor de forma β	Argumento nombre-valor Leakage relacionado con el factor de forma de la ventana de Kaiser β: ℓ = 1 – β/40
	Solapamiento, L	Número de muestras especificado como tercer argumento posicional (Predeterminado: 50% de la longitud de la ventana)	Número de muestras especificado con el argumento nombre-valor OverlapLength (Predeterminado: 75% de la longitud de la ventana)	Porcentaje de la longitud del segmento especificado con el argumento nombre-valor OverlapPercent (Predeterminado: $$\left(1-\frac{\text{1}}{2\times \text{ENBW}-1}\right)\times 100,$$ donde ENBW es el ancho de banda de ruido equivalente de la ventana)
	Número de puntos DFT, NDFT	Cuarto argumento posicional (Predeterminado: max(256,2⌈log2M⌉))	Argumento nombre-valor FFTLength (Predeterminado: 128)	Siempre 1024
	Información de tiempo	Tasa de muestreo especificada como quinto argumento posicional	Tasa de muestreo o vector temporal especificado como segundo argumento posicional	Tasa de muestreo o vector temporal especificado como segundo argumento posicional
Llamada a la función	fs = 100; x = exp(2j*pi*20* ... (0:1/fs:2-1/fs)); M = 200; lk = 0.5; g = kaiser(M,40*(1-lk)); L = 100; Ndft = 1024;	sps = abs( ... spectrogram(x,g,L, ... Ndft,fs,"centered") ... ).^2;	sts = abs( ... stft(x,fs,Window=g, ... OverlapLength=L, ... FFTLength=Ndft) ... ).^2;	pss = pspectrum(x,fs, ... "spectrogram", ... TimeResolution=M/fs, ... Leakage=lk, ... OverlapPercent=L/M*100 ... )*sum(g)^2;
Gráfica de conveniencia	fs = 6e2; ts = 0:1/fs:2.05; x = vco(sin(2*pi*ts).* ... exp(-ts),[0.1 0.4]*fs,fs); M = 32; lk = 0.9; g = kaiser(M,40*(1-lk)); L = 22; Ndft = 1024; Para spectrogram, añada "power","yaxis" Para stft, añada FrequencyRange="onesided"
Salida	Rango de frecuencia	Se controla con el argumento freqrange: "onesided": en los valores pares de NDFT, el intervalo de frecuencia se cierra en la frecuencia cero y en la frecuencia de Nyquist fs/2. En los valores impares de NDFT, el intervalo de frecuencia se cierra en la frecuencia cero y se abre en fs/2. (Predeterminado para señales de valor real). "twosided": en todos los valores de NDFT, el intervalo de frecuencia se cierra en la frecuencia cero y se abre en fs. "centered": en los valores pares de NDFT, el intervalo de frecuencia se abre en –fs/2 y se cierra en fs/2. En los valores impares de NDFT, el intervalo de frecuencia está abierto en ambos extremos. (Predeterminado para señales de valor complejo). El usuario puede especificar un vector de frecuencias en el que calcular la STFT y el espectrograma.	Se controla con el argumento nombre-valor FrequencyRange: "onesided": igual que en spectrogram. "twosided": igual que en spectrogram. "centered": igual que en spectrogram. (Predeterminado para señales de valores reales y complejos).	Se controla mediante el argumento lógico nombre-valor TwoSided: false: el intervalo de frecuencia se cierra en la frecuencia cero y en fs/2. (Predeterminado para señales de valor real). true: el intervalo de frecuencia se cierra en –fs/2 y en fs/2. (Predeterminado para señales de valor complejo).
	Intervalo de tiempo	Señal truncada tras el último segmento con ventana completa. Valores temporales en los centros de los segmentos.	Señal truncada tras el último segmento con ventana completa. Valores temporales en los centros de los segmentos.	Señal rellenada con ceros tras el último segmento con ventana completa. Valores temporales en los centros de los segmentos.
	Normalización	El primer argumento de salida es STFT. Su magnitud al cuadrado es el espectrograma. El cuarto argumento de salida es una magnitud al cuadrado. Para obtener el espectrograma, multiplique por (Σng(n))2 y especifique la opción "power".	El primer argumento de salida es STFT. Su magnitud al cuadrado es el espectrograma.	El primer argumento de salida es una magnitud al cuadrado. Para obtener el espectrograma, multiplique por (Σng(n))2.
	PSD y espectro de potencia	El cuarto argumento de salida de spectrogram contiene densidades espectrales de potencia de segmento o espectros de potencia de segmento. El espectrograma es igual al espectro de potencia multiplicado por el cuadrado de la suma de los elementos de la ventana. Argumento spectrumtype: "psd": multiplique por ENBW para obtener el espectro de potencia (Predeterminado) "power": espectro de potencia	El argumento de salida es STFT.	El argumento de salida es el espectro de potencia Para obtener la PSD, divida por ENBW
Ejemplos		Valores predeterminados del espectrograma Comparar la función spectrogram y la definición de STFT Seguimiento de chirps en una señal de audio	Transformada de Fourier de tiempo corto Comparar las funciones spectrogram y stft STFT de señales multicanal	Espectros de potencia de sinusoides Comparar las funciones spectrogram y pspectrum Filtrado eliminador de banda de señal musical

Definiciones de STFT y espectrograma

La STFT de una señal se calcula deslizando una ventana de análisis g(n) de longitud M sobre la señal y calculando la transformada discreta de Fourier (DFT) de cada segmento de los datos con ventana. La ventana salta sobre la señal original en intervalos de R muestras, lo que equivale a L = M – R muestras de solapamiento entre segmentos contiguos. La mayoría de las funciones de ventana se estrechan en los bordes para evitar el anillamiento espectral. La DFT de cada segmento con ventana se añade a una matriz de valor complejo que contiene la magnitud y la fase de cada punto en tiempo y frecuencia. La matriz STFT tiene

$$k=\lfloor \frac{N_x-L}{M-L}\rfloor$$

 columnas, donde N_x es la longitud de la señal x(n) y los símbolos ⌊⌋ denotan la función de suelo. En las transformadas centradas y bilaterales, el número de filas de la matriz es igual a N_DFT, el número de puntos de la DFT, y en las transformadas unilaterales de señales de valor real es un número impar cercano a N_DFT/2.

La columna m de la matriz de la STFT $$X(f)=\left[\begin{array}{ccccc}{X}_1(f)& X_2(f)& X_3(f)& \cdots & X_k(f)\end{array}\right]$$ contiene la DFT de los datos con ventana centrados en el tiempo mR:

$$X_m(f)={\displaystyle \sum _{n=-\infty }^{\infty }x(n)g(n-mR)e^{-j2\pi fn}}.$$

Comparar la función spectrogram y la definición de STFT

Genere una señal compuesta por un chirp cuadrático convexo de valor complejo muestreado a 600 Hz durante 2 segundos. El chirp tiene una frecuencia inicial de 250 Hz y una final de 50 Hz.

fs = 6e2;
ts = 0:1/fs:2;
x = chirp(ts,250,ts(end),50,"quadratic",0,"convex","complex");

Función spectrogram

Utilice la función spectrogram para calcular la STFT de la señal.

• Divida la señal en segmentos, cada uno de $\mathit{M}=49$ muestras.

• Especifique $\mathit{L}=11$ muestras de solapamiento entre segmentos contiguos.

• Descarte el último segmento, el más corto.

• Aplique una ventana de Bartlett a cada segmento.

• Evalúe la transformada discreta de Fourier de cada segmento en ${\mathit{N}}_{\mathrm{DFT}}=1024$ puntos. De forma predeterminada, spectrogram calcula transformadas bilaterales en señales de valor complejo.

M = 49;
L = 11;
g = bartlett(M);
Ndft = 1024;

[s,f,t] = spectrogram(x,g,L,Ndft,fs);

Calcule y muestre el espectrograma, definido como la magnitud al cuadrado de la STFT, con la función waterplot.

waterplot(s,f,t)

Definición de STFT

Calcule la STFT de la señal de ${\mathit{N}}_{\mathit{x}}$ muestras con la definición. Divida la señal en $$\lfloor \frac{N_x-L}{M-L}\rfloor$$ segmentos superpuestos. Aplique una ventana a cada segmento y evalúe la transformada discreta de Fourier en ${\mathit{N}}_{\mathrm{DFT}}$ puntos.

segs = framesig(1:length(x),M,OverlapLength=L);
X = fft(x(segs).*g,Ndft);

Calcule los intervalos de tiempo y frecuencia de la STFT.

• Para encontrar los valores de tiempo, divida el vector tiempo en segmentos superpuestos. Los valores de tiempo son los puntos medios de los segmentos, y cada segmento se trata como un intervalo abierto en el extremo inferior.

• Para encontrar los valores de frecuencia, especifique un intervalo de Nyquist cerrado en la frecuencia cero y abierto en el extremo inferior.

framedT = ts(segs);
tint = mean(framedT(2:end,:));

fint = 0:fs/Ndft:fs-fs/Ndft;

Compare la salida de spectrogram con la definición. Muestre el espectrograma.

maxdiff = max(max(abs(s-X)))

maxdiff = 
0

waterplot(X,fint,tint)

function waterplot(s,f,t)
% Waterfall plot of spectrogram
    waterfall(f,t,abs(s)'.^2)
    set(gca,XDir="reverse",View=[30 50])
    xlabel("Frequency (Hz)")
    ylabel("Time (s)")
end

Comparar las funciones spectrogram y stft

Genere una señal compuesta por un chirp muestreado a 1,4 kHz durante 2 segundos. La frecuencia del chirp disminuye linealmente de 600 Hz a 100 Hz durante el tiempo de medición.

fs = 1400;
x = chirp(0:1/fs:2,600,2,100);

Valores predeterminados de stft

Calcule la STFT de la señal con las funciones spectrogram y stft. Utilice los valores predeterminados de la función stft:

• Divida la señal en segmentos de 128 muestras y aplique una ventana periódica de Hann a cada segmento.

• Especifique 96 muestras de solapamiento entre segmentos contiguos. Esta longitud equivale al 75% de la longitud de la ventana.

• Especifique 128 puntos de la DFT y centre la STFT en la frecuencia cero, con la frecuencia expresada en hercios.

Compruebe que los dos resultados son iguales.

M = 128;
g = hann(M,"periodic");
L = 0.75*M;
Ndft = 128;

[sp,fp,tp] = spectrogram(x,g,L,Ndft,fs,"centered");

[s,f,t] = stft(x,fs);

dff = max(max(abs(sp-s)))

dff = 
0

Represente las dos salidas con la función mesh.

nexttile
mesh(tp,fp,abs(sp).^2)
title("spectrogram")
view(2), axis tight

nexttile
mesh(t,f,abs(s).^2)
title("stft")
view(2), axis tight

Valores predeterminados de spectrogram

Repita el cálculo utilizando los valores predeterminados de la función spectrogram:

• Divida la señal en segmentos de longitud $\mathit{M}=\lfloor {\mathit{N}}_{\mathit{x}}/4.5\rfloor$, donde ${\mathit{N}}_{\mathit{x}}$ es la longitud de la señal. Aplique una ventana de Hamming a cada segmento.

• Especifique el 50% de solapamiento entre segmentos.

• Para calcular la FFT, utilice $\max \left(256,2^{\lceil {\log }_2\mathit{M}\rceil }\right)$ puntos. Calcule el espectrograma solo para las frecuencias normalizadas positivas.

M = floor(length(x)/4.5);
g = hamming(M);
L = floor(M/2);
Ndft = max(256,2^nextpow2(M));

[sx,fx,tx] = spectrogram(x);

[st,ft,tt] = stft(x,Window=g,OverlapLength=L, ...
    FFTLength=Ndft,FrequencyRange="onesided");

dff = max(max(sx-st))

dff = 
0

Represente las dos salidas con la función waterplot. Divida el eje de la frecuencia por $\pi$ en ambos casos. Para la salida de stft, divida los números de muestra por la tasa de muestreo efectiva, $2\pi$.

figure
nexttile
waterplot(sx,fx/pi,tx)
title("spectrogram")

nexttile
waterplot(st,ft/pi,tt/(2*pi))
title("stft")

function waterplot(s,f,t)
% Waterfall plot of spectrogram
    waterfall(f,t,abs(s)'.^2)
    set(gca,XDir="reverse",View=[30 50])
    xlabel("Frequency/\pi")
    ylabel("Samples")
end

Comparar las funciones spectrogram y pspectrum

Genere una señal que conste de un oscilador controlado por tensión y tres átomos gaussianos. La señal se muestrea a ${\mathit{f}}_{\mathit{s}}=2$ kHz durante 2 segundos.

fs = 2000;
tx = 0:1/fs:2;
gaussFun = @(A,x,mu,f) exp(-(x-mu).^2/(2*0.03^2)).*sin(2*pi*f.*x)*A';
s = gaussFun([1 1 1],tx',[0.1 0.65 1],[2 6 2]*100)*1.5;
x = vco(chirp(tx+.1,0,tx(end),3).*exp(-2*(tx-1).^2),[0.1 0.4]*fs,fs);
x = s+x';

Transformadas de Fourier de tiempo corto

Utilice la función pspectrum para calcular la STFT.

• Divida la señal de ${\mathit{N}}_{\mathit{x}}$ muestras en segmentos con longitud $\mathit{M}=80$ muestras, correspondientes a una resolución temporal de $80/2000=40$ milisegundos.

• Especifique $\mathit{L}=16$ muestras o el 20% de solapamiento entre los segmentos contiguos.

• Añada una ventana de Kaiser a cada segmento y especifique un manchado de $\ell =0.7$.

M = 80;
L = 16;
lk = 0.7;

[S,F,T] = pspectrum(x,fs,"spectrogram", ...
    TimeResolution=M/fs,OverlapPercent=L/M*100, ...
    Leakage=lk);

Compárela con el resultado obtenido con la función spectrogram.

• Especifique la longitud de la ventana y el solapamiento directamente en las muestras.

• pspectrum siempre usa una ventana de Kaiser como $\mathit{g}\left(\mathit{n}\right)$. El manchado $\ell$ y el factor de forma $\beta$ de la ventana están relacionados por $\beta =40\times \left(1-\ell \right)$.

• pspectrum siempre usa ${\mathit{N}}_{\mathrm{DFT}}=1024$ puntos cuando calcula la transformada discreta de Fourier. Puede especificar este número si desea calcular la transformada en un rango de frecuencias bilateral o centrado. Sin embargo, en las transformadas unilaterales, que son las predeterminadas para las señales reales, spectrogram utiliza $1024/2+1=513$ puntos. Como alternativa, puede especificar el vector de frecuencias en el que desea calcular la transformada, como en este ejemplo.

• Si una señal no puede dividirse exactamente en $\mathit{k}=\lfloor \frac{{\mathit{N}}_{\mathit{x}}-\mathit{L}}{\mathit{M}-\mathit{L}}\rfloor$ segmentos, spectrogram la trunca, mientras que pspectrum la rellena con ceros para crear un segmento extra. Para que las salidas sean equivalentes, elimine el segmento final y el elemento final del vector de tiempo.

• spectrogram devuelve la STFT, cuya magnitud elevada al cuadrado es el espectrograma. pspectrum devuelve el espectro de potencia segmento a segmento, que ya está elevado al cuadrado pero se divide por un factor de ${\sum }_{\mathit{n}}\mathit{g}\left(\mathit{n}\right)$ antes de elevarlo al cuadrado.

• En las transformadas unilaterales, pspectrum añade un factor extra de 2 al espectrograma.

g = kaiser(M,40*(1-lk));

k = (length(x)-L)/(M-L);
if k~=floor(k)
    S = S(:,1:floor(k));
    T = T(1:floor(k));
end

[s,f,t] = spectrogram(x/sum(g)*sqrt(2),g,L,F,fs);

Para mostrar los espectrogramas calculados por las dos funciones, utilice la función waterplot.

subplot(2,1,1)
waterplot(sqrt(S),F,T)
title("pspectrum")

subplot(2,1,2)
waterplot(s,f,t)
title("spectrogram")

maxd = max(max(abs(abs(s).^2-S)))

maxd = 
2.4419e-08

Espectros de potencia y gráficas de conveniencia

La función spectrogram tiene un cuarto argumento que se corresponde con el espectro de potencia segmento a segmento o con la densidad espectral de potencia. Al igual que la salida de pspectrum, el argumento ps ya está elevado al cuadrado e incluye el factor de normalización ${\sum }_{\mathit{n}}\mathit{g}\left(\mathit{n}\right)$. En los espectrogramas unilaterales de señales reales, hay que incluir el factor adicional de 2. Establezca el argumento de escala de la función en "power".

[~,~,~,ps] = spectrogram(x*sqrt(2),g,L,F,fs,"power");

max(abs(S(:)-ps(:)))

ans = 
2.4419e-08

Cuando se llaman sin argumentos de salida, tanto pspectrum como spectrogram representan el espectrograma de la señal en decibelios. Incluya el factor de 2 en los espectrogramas unilaterales. Establezca que los mapas de colores sean los mismos en ambas gráficas. Fije los límites del eje x en los mismos valores para hacer visible el segmento extra al final de la gráfica de pspectrum. En la gráfica de spectrogram, muestre la frecuencia en el eje y.

subplot(2,1,1)
pspectrum(x,fs,"spectrogram", ...
    TimeResolution=M/fs,OverlapPercent=L/M*100, ...
    Leakage=lk)
title("pspectrum")
cc = clim;
xl = xlim;

subplot(2,1,2)
spectrogram(x*sqrt(2),g,L,F,fs,"power","yaxis")
title("spectrogram")
clim(cc)
xlim(xl)

function waterplot(s,f,t)
% Waterfall plot of spectrogram
    waterfall(f,t,abs(s)'.^2)
    set(gca,XDir="reverse",View=[30 50])
    xlabel("Frequency (Hz)")
    ylabel("Time (s)")
end

Cálculo de espectrogramas centrados y unilaterales

Genere una señal compuesta por un chirp de valor real muestreado a 2 kHz durante 2 segundos.

fs = 2000;
tx = 0:1/fs:2;
x = vco(-chirp(tx,0,tx(end),2).*exp(-3*(tx-1).^2), ...
    [0.1 0.4]*fs,fs).*hann(length(tx))';

Espectrograma bilateral

Calcule y represente la STFT bilateral de la señal.

• Divida la señal en segmentos, cada uno de $\mathit{M}=73$ muestras.

• Especifique $\mathit{L}=24$ muestras de solapamiento entre segmentos contiguos.

• Descarte el último segmento, el más corto.

• Aplique una ventana flat-top a cada segmento.

• Evalúe la transformada discreta de Fourier de cada segmento en ${\mathit{N}}_{\mathrm{DFT}}=895$ puntos, teniendo en cuenta que es un número impar.

M = 73;
L = 24;
g = flattopwin(M);
Ndft = 895;
neven = ~mod(Ndft,2);

[stwo,f,t] = spectrogram(x,g,L,Ndft,fs,"twosided");

Utilice la función spectrogram sin argumentos de salida para representar el espectrograma bilateral.

spectrogram(x,g,L,Ndft,fs,"twosided","power","yaxis")

Calcule el espectrograma bilateral con la definición. Divida la señal en segmentos de $\mathit{M}$ muestras con $\mathit{L}$ muestras de solapamiento entre segmentos contiguos. Aplique una ventana a cada segmento y calcule la transformada discreta de Fourier en ${\mathit{N}}_{\mathrm{DFT}}$ puntos.

y = framesig(x,M,Window=g,OverlapLength=L);
Xtwo = fft(y,Ndft);

Calcule los intervalos de tiempo y frecuencia.

• Para encontrar los valores de tiempo, divida el vector tiempo en segmentos superpuestos. Los valores de tiempo son los puntos medios de los segmentos, y cada segmento se trata como un intervalo abierto en el extremo inferior.

• Para encontrar los valores de frecuencia, especifique un intervalo de Nyquist cerrado en la frecuencia cero y abierto en el extremo superior.

framedT = framesig(tx,M,OverlapLength=L);
ttwo = mean(framedT(2:end,:));

ftwo = 0:fs/Ndft:fs*(1-1/Ndft);

Compare las salidas del spectrogram con las definiciones. Visualice los espectrogramas con la función waterplot.

diffs = [max(max(abs(stwo-Xtwo)));
    max(abs(f-ftwo'));
    max(abs(t-ttwo))]

diffs = 3×1
10-12 ×

         0
    0.2274
    0.0002

figure
nexttile
waterplot(Xtwo,ftwo,ttwo)
title("Two-Sided, Definition")

nexttile
waterplot(stwo,f,t)
title("Two-Sided, spectrogram Function")

Espectrograma centrado

Calcule el espectrograma centrado de la señal.

• Utilice los mismos valores de tiempo que utilizó en la STFT bilateral.

• Desplace el componente de frecuencia cero de la STFT hacia el centro del espectro con la función fftshift.

• En los valores impares de ${\mathit{N}}_{\mathrm{DFT}}$, el intervalo de frecuencia está abierto en ambos extremos. En los valores pares de ${\mathit{N}}_{\mathrm{DFT}}$, el intervalo de frecuencia está abierto en el extremo inferior y cerrado en el extremo superior.

Compare las salidas y muestre los espectrogramas.

tcen = ttwo;

if ~neven
    Xcen = fftshift(Xtwo,1);
    fcen = -fs/2*(1-1/Ndft):fs/Ndft:fs/2;
else
    Xcen = fftshift(circshift(Xtwo,-1),1);
    fcen = (-fs/2*(1-1/Ndft):fs/Ndft:fs/2)+fs/Ndft/2;
end

[scen,f,t] = spectrogram(x,g,L,Ndft,fs,"centered");

diffs = [max(max(abs(scen-Xcen)));
    max(abs(f-fcen'));
    max(abs(t-tcen))]

diffs = 3×1
10-12 ×

         0
    0.2274
    0.0002

figure
nexttile
waterplot(Xcen,fcen,tcen)
title("Centered, Definition")

nexttile
waterplot(scen,f,t)
title("Centered, spectrogram Function")

Espectrograma unilateral

Calcule el espectrograma unilateral de la señal.

• Utilice los mismos valores de tiempo que utilizó en la STFT bilateral.

• En los valores impares de ${\mathit{N}}_{\mathrm{DFT}}$, la STFT unilateral consiste en las primeras $\left({\mathit{N}}_{\mathrm{DFT}}+1\right)/2$ filas de la STFT bilateral. En los valores pares de ${\mathit{N}}_{\mathrm{DFT}}$, la STFT unilateral consiste en las primeras ${\mathit{N}}_{\mathrm{DFT}}/2+1$ filas de la STFT bilateral.

• En los valores impares de ${\mathit{N}}_{\mathrm{DFT}}$, el intervalo de frecuencia se cierra en la frecuencia cero y se abre en la frecuencia de Nyquist. En los valores pares de ${\mathit{N}}_{\mathrm{DFT}}$, el intervalo de frecuencia se cierra en ambos extremos.

Compare las salidas y muestre los espectrogramas. En señales de valor real, el argumento "onesided" es opcional.

tone = ttwo;

if ~neven
    Xone = Xtwo(1:(Ndft+1)/2,:);
else
    Xone = Xtwo(1:Ndft/2+1,:);
end

fone = 0:fs/Ndft:fs/2;

[sone,f,t] = spectrogram(x,g,L,Ndft,fs);

diffs = [max(max(abs(sone-Xone)));
    max(abs(f-fone'));
    max(abs(t-tone))]

diffs = 3×1
10-12 ×

         0
    0.1137
    0.0002

figure
nexttile
waterplot(Xone,fone,tone)
title("One-Sided, Definition")

nexttile
waterplot(sone,f,t)
title("One-Sided, spectrogram Function")

function waterplot(s,f,t)
% Waterfall plot of spectrogram
waterfall(f,t,abs(s)'.^2)
set(gca,XDir="reverse",View=[30 50])
xlabel("Frequency (Hz)")
ylabel("Time (s)")
end

Calcular las PSD de los segmentos y los espectros de potencia

La función spectrogram tiene como cuarto argumento de salida una matriz que contiene la densidad espectral de potencia (PSD) o el espectro de potencia de cada segmento. El espectro de potencia es igual a la PSD multiplicada por el ancho de banda de ruido equivalente (ENBW) de la ventana.

Genere una señal compuesta por un chirp logarítmico muestreado a 1 kHz durante 1 segundo. El chirp tiene una frecuencia inicial de 400 Hz que disminuye a 10 Hz al final de la medición.

fs = 1000;
tt = 0:1/fs:1-1/fs;
y = chirp(tt,400,tt(end),10,"logarithmic");

PSD de segmentos y espectros de potencia con tasa de muestreo

Divida la señal en segmentos de 102 muestras y aplique una ventana de Hann a cada uno. Especifique 12 muestras de solapamiento entre segmentos contiguos y 1024 puntos de la DFT.

M = 102;
g = hann(M);
L = 12;
Ndft = 1024;

Calcule el espectrograma de la señal con el tipo de espectro de la PSD predeterminado. Genere la STFT y el arreglo de densidades espectrales de potencia de los segmentos.

[s,f,t,p] = spectrogram(y,g,L,Ndft,fs);

Repita el cálculo con el tipo de espectro especificado como "power". Obtenga como salida la STFT y el arreglo de espectros de potencia de los segmentos.

[r,~,~,q] = spectrogram(y,g,L,Ndft,fs,"power");

Compruebe que el espectrograma es el mismo en ambos casos. Represente el espectrograma con una escala logarítmica para la frecuencia.

max(max(abs(s).^2-abs(r).^2))

ans = 
0

waterfall(f,t,abs(s)'.^2)
set(gca,XScale="log",...
    XDir="reverse",View=[30 50])

Compruebe que los espectros de potencia son iguales a las densidades espectrales de potencia multiplicadas por el ENBW de la ventana.

max(max(abs(q-p*enbw(g,fs))))

ans = 
1.1102e-16

Compruebe que la matriz de espectros de potencia del segmento es proporcional al espectrograma. El factor de proporcionalidad es el cuadrado de la suma de los elementos de la ventana.

max(max(abs(s).^2-q*sum(g)^2))

ans = 
0

PSD de segmentos y espectros de potencia con frecuencias normalizadas

Repita el cálculo, pero ahora trabaje en frecuencias normalizadas. Los resultados son los mismos cuando se especifica la tasa de muestreo como $2\pi$.

[~,~,~,pn] = spectrogram(y,g,L,Ndft);
[~,~,~,qn] = spectrogram(y,g,L,Ndft,"power");

max(max(abs(qn-pn*enbw(g,2*pi))))

ans = 
1.1102e-16

Consulte también

Apps

• Signal Analyzer

Funciones

• dlstft | dlistft | pspectrum | spectrogram | stft | istft | xspectrogram

Objetos

• stftLayer | istftLayer

Temas

• Time-Frequency Gallery
• Practical Introduction to Time-Frequency Analysis

Sitios web externos

• Análisis de Fourier (Recursos didácticos de MathWorks)