Main Content

Esta página aún no se ha traducido para esta versión. Puede ver la versión más reciente de esta página en inglés.

cpsd

Densidad espectral de potencia cruzada

Descripción

ejemplo

pxy = cpsd(x,y) estima la densidad espectral de potencia cruzada (CPSD) de dos señales de tiempo discreto y, utilizando el método de estimación espectral de periodograma modificado y promediado de Welch.xy

  • Si y son ambos vectores, deben tener la misma longitud.xy

  • Si una de las señales es una matriz y la otra es un vector, entonces la longitud del vector debe ser igual al número de filas de la matriz. La función expande el vector y devuelve una matriz de estimaciones de densidad espectral de potencia cruzada columna por columna.

  • Si y son matrices con el mismo número de filas pero diferentes números de columnas, a continuación, devuelve una matriz tridimensional, , que contiene estimaciones de densidad espectral de potencia cruzada para todas las combinaciones de columnas de entrada.xycpsdpxy Cada columna de corresponde a una columna de , y cada página corresponde a una columna de : .pxyxypxy(:,m,n) = cpsd(x(:,m),y(:,n))

  • Si y son matrices de igual tamaño, entonces opera columna: .xycpsdpxy(:,n) = cpsd(x(:,n),y(:,n)) Para obtener una matriz de entrada múltiple/salida múltiple, anexe a la lista de argumentos.'mimo'

Para real y , devuelve un CPSD unilateral.xycpsd Para complejo o , devuelve un CPSD de dos lados.xycpsd

pxy = cpsd(x,y,window) se utiliza para dividir y en segmentos y realizar ventanas.windowxy

pxy = cpsd(x,y,window,noverlap) utiliza muestras de superposición entre segmentos adyacentes.noverlap

pxy = cpsd(x,y,window,noverlap,nfft) utiliza puntos de muestreo para calcular la transformación discreta de Fourier.nfft

ejemplo

pxy = cpsd(___,'mimo') calcula una matriz de múltiples entradas/salidas múltiples de estimaciones de densidad espectral de potencia cruzada. Esta sintaxis puede incluir cualquier combinación de argumentos de entrada de sintaxis anteriores.

[pxy,w] = cpsd(___) devuelve un vector de frecuencias normalizadas, en el que se estima la densidad espectral de potencia cruzada.w

[pxy,f] = cpsd(___,fs) devuelve un vector de frecuencias, expresado en términos de la frecuencia de muestreo, en el que se estima la densidad espectral de potencia cruzada. debe ser la sexta entrada numérica para .ffsfscpsd Para introducir una frecuencia de muestreo y seguir utilizando los valores predeterminados de los argumentos opcionales anteriores, especifique estos argumentos como vacíos, .[]

ejemplo

[pxy,w] = cpsd(x,y,window,noverlap,w) devuelve las estimaciones de densidad espectral de potencia cruzada a las frecuencias normalizadas especificadas en .w

ejemplo

[pxy,f] = cpsd(x,y,window,noverlap,f,fs) devuelve las estimaciones de densidad espectral de potencia cruzada a las frecuencias especificadas en .f

[___] = cpsd(x,y,___,freqrange) devuelve la estimación de densidad espectral de potencia cruzada sobre el rango de frecuencia especificado por .freqrange Las opciones válidas para son , , y .freqrange'onesided''twosided''centered'

ejemplo

cpsd(___) sin argumentos de salida traza la estimación de densidad espectral de potencia cruzada en la ventana de figura actual.

Ejemplos

contraer todo

Genere dos señales de ruido de color y trace su densidad espectral de potencia cruzada. Especifique una longitud-1024 FFT y una ventana triangular de 500 puntos sin superposición.

r = randn(16384,1);  hx = fir1(30,0.2,rectwin(31)); x = filter(hx,1,r);  hy = ones(1,10)/sqrt(10); y = filter(hy,1,r);  cpsd(x,y,triang(500),250,1024)

Genere dos sinusoides de dos canales muestreados a 1 kHz durante 1 segundo. Los canales de la primera señal tienen frecuencias de 200 Hz y 300 Hz. Los canales de la segunda señal tienen frecuencias de 300 Hz y 400 Hz. Ambas señales están incrustadas en el ruido gaussiano blanco de varianza unitaria.

fs = 1e3; t = (0:1/fs:1-1/fs)';  q = 2*sin(2*pi*[200 300].*t); q = q+randn(size(q));  r = 2*sin(2*pi*[300 400].*t); r = r+randn(size(r));

Calcular la densidad espectral de potencia cruzada de las dos señales. Utilice una ventana Bartlett de 256 muestras para dividir las señales en segmentos y ventanar los segmentos. Especifique 128 muestras de superposición entre segmentos adyacentes y 2048 puntos DFT. Utilice la funcionalidad integrada de para trazar el resultado.cpsd

cpsd(q,r,bartlett(256),128,2048,fs)

De forma predeterminada, funciona columna por columna en entradas de matriz del mismo tamaño.cpsd Cada canal alcanza los picos de las frecuencias de los sinusoides originales.

Repita el cálculo, pero ahora anexe a la lista de argumentos.'mimo'

cpsd(q,r,bartlett(256),128,2048,fs,'mimo')

Cuando se llama con la opción, devuelve una matriz tridimensional que contiene estimaciones de densidad espectral de potencia cruzada para todas las combinaciones de columnas de entrada.'mimo'cpsd La estimación del segundo canal y el primer canal de muestra un pico mejorado en la frecuencia común de 300 Hz.qr

Generar dos señales sinusoidales de 100 Hz muestreadas a 1 kHz durante 296 ms. Uno de los sinusoides se relaga el otro en 2,5 ms, equivalente a un desfase de fase de /2.π Ambas señales están incrustadas en ruido gaussiano blanco de varianza 1/42.

Fs = 1000; t = 0:1/Fs:0.296;  x = cos(2*pi*t*100)+0.25*randn(size(t)); tau = 1/400; y = cos(2*pi*100*(t-tau))+0.25*randn(size(t));

Calcular y trazar la magnitud de la densidad espectral de potencia cruzada. Utilice la configuración predeterminada para .cpsd La magnitud alcanza la frecuencia en la que hay una coherencia significativa entre las señales.

cpsd(x,y,[],[],[],Fs)

Función de coherencia cuadrada de magnitud de la gráfica y la fase del espectro cruzado. La coordenada en la frecuencia de alta coherencia corresponde al desfase de fase entre los sinusoides.

[Cxy,F] = mscohere(x,y,[],[],[],Fs);  [Pxy,F] = cpsd(x,y,[],[],[],Fs);  subplot(2,1,1) plot(F,Cxy) title('Magnitude-Squared Coherence')  subplot(2,1,2) plot(F,angle(Pxy))  hold on plot(F,2*pi*100*tau*ones(size(F)),'--') hold off  xlabel('Hz') ylabel('\Theta(f)') title('Cross Spectrum Phase')

Generar dos

<math display="block">
<mrow>
<mi>N</mi>
</mrow>
</math>
-muestra secuencias exponenciales,
<math display="block">
<mrow>
<msub>
<mrow>
<mi>x</mi>
</mrow>
<mrow>
<mi>a</mi>
</mrow>
</msub>
<mo>=</mo>
<msup>
<mrow>
<mi>a</mi>
</mrow>
<mrow>
<mi>n</mi>
</mrow>
</msup>
</mrow>
</math>
Y
<math display="block">
<mrow>
<msub>
<mrow>
<mi>x</mi>
</mrow>
<mrow>
<mi>b</mi>
</mrow>
</msub>
<mo>=</mo>
<msup>
<mrow>
<mi>b</mi>
</mrow>
<mrow>
<mi>n</mi>
</mrow>
</msup>
</mrow>
</math>
Con
<math display="block">
<mrow>
<mi>n</mi>
<mo></mo>
<mn>0</mn>
</mrow>
</math>
. Especificar
<math display="block">
<mrow>
<mi>a</mi>
<mo>=</mo>
<mn>0</mn>
<mo>.</mo>
<mn>8</mn>
</mrow>
</math>
,
<math display="block">
<mrow>
<mi>b</mi>
<mo>=</mo>
<mn>0</mn>
<mo>.</mo>
<mn>9</mn>
</mrow>
</math>
, y una pequeña
<math display="block">
<mrow>
<mi>N</mi>
</mrow>
</math>
para ver efectos finitos.

N = 10; n = 0:N-1;  a = 0.8; b = 0.9;  xa = a.^n; xb = b.^n;

Calcular y trazar la densidad espectral de potencia cruzada de las secuencias en el intervalo completo de frecuencias normalizadas,

<math display="block">
<mrow>
<mo stretchy="false">[</mo>
<mo>-</mo>
<mi>π</mi>
<mo>,</mo>
<mi>π</mi>
<mo stretchy="false">]</mo>
</mrow>
</math>
. Especifique una ventana rectangular de longitud
<math display="block">
<mrow>
<mi>N</mi>
</mrow>
</math>
y no hay superposición entre segmentos.

w = -pi:1/1000:pi; wind = rectwin(N); nove = 0;  [pxx,f] = cpsd(xa,xb,wind,nove,w);

El espectro de potencia cruzada de las dos secuencias tiene una expresión analítica para

<math display="block">
<mrow>
<mi>N</mi>
</mrow>
</math>
:

<math display="block">
<mrow>
<mi>R</mi>
<mo stretchy="false">(</mo>
<mi>ω</mi>
<mo stretchy="false">)</mo>
<mo>=</mo>
<mrow>
<mfrac>
<mrow>
<mn>1</mn>
</mrow>
<mrow>
<mrow>
<mn>1</mn>
<mo>-</mo>
<mi>a</mi>
<msup>
<mrow>
<mi>e</mi>
</mrow>
<mrow>
<mo>-</mo>
<mi>j</mi>
<mi>ω</mi>
</mrow>
</msup>
</mrow>
</mrow>
</mfrac>
</mrow>
<mspace width="0.16666666666666666em"></mspace>
<mrow>
<mfrac>
<mrow>
<mn>1</mn>
</mrow>
<mrow>
<mrow>
<mn>1</mn>
<mo>-</mo>
<mi>b</mi>
<msup>
<mrow>
<mi>e</mi>
</mrow>
<mrow>
<mi>j</mi>
<mi>ω</mi>
</mrow>
</msup>
</mrow>
</mrow>
</mfrac>
</mrow>
<mo>.</mo>
</mrow>
</math>

Convierta esta expresión en una densidad espectral de potencia cruzada dividiéndola por

<math display="block">
<mrow>
<mn>2</mn>
<mi>π</mi>
<mi>N</mi>
</mrow>
</math>
. Compare los resultados. La ondulación en el resultado es una consecuencia de las ventanas.cpsd

nfac = 2*pi*N;  X = 1./(1-a*exp(-1j*w)); Y = 1./(1-b*exp( 1j*w)); R = X.*Y/nfac;  semilogy(f/pi,abs(pxx)) hold on semilogy(w/pi,abs(R)) hold off legend('cpsd','Analytic')

Repita el cálculo con

<math display="block">
<mrow>
<mi>N</mi>
<mo>=</mo>
<mn>2</mn>
<mn>5</mn>
</mrow>
</math>
. Las curvas acuerdan seis cifras para
<math display="block">
<mrow>
<mi>N</mi>
</mrow>
</math>
tan pequeño como 100.

N = 25; n = 0:N-1; xa = a.^n; xb = b.^n;  wind = rectwin(N);  [pxx,f] = cpsd(xa,xb,wind,nove,w); R = X.*Y/(2*pi*N);  semilogy(f/pi,abs(pxx)) hold on semilogy(w/pi,abs(R)) hold off legend('cpsd','Analytic')

Utilice la densidad espectral de potencia cruzada para identificar un tono altamente dañado.

Las señales de sonido generadas al marcar un número o símbolo en un teléfono digital son sumas de sinusoides con frecuencias tomadas de dos grupos diferentes. Cada par de tonos contiene una frecuencia del grupo bajo (697 Hz, 770 Hz, 852 Hz o 941 Hz) y una frecuencia del grupo alto (1209 Hz, 1336 Hz o 1477 Hz).

Generar señales correspondientes a todos los símbolos. Muestre cada tono a 4 kHz durante medio segundo. Prepare una tabla de referencia.

fs = 4e3; t = 0:1/fs:0.5-1/fs;  nms = ['1';'2';'3';'4';'5';'6';'7';'8';'9';'*';'0';'#'];  ver = [697 770 852 941]; hor = [1209 1336 1477];  v = length(ver); h = length(hor);  for k = 1:v     for l = 1:h         idx = h*(k-1)+l;         tone = sum(sin(2*pi*[ver(k);hor(l)].*t))';         tones(:,idx) = tone;     end end

Trazar el periodograma Welch de cada señal y anotar las frecuencias de los componentes. Utilice una ventana Hamming de 200 muestras para dividir las señales en segmentos no superpuestos y avistar los segmentos.

[pxx,f] = pwelch(tones,hamming(200),0,[],fs);  for k = 1:v     for l = 1:h         idx = h*(k-1)+l;         ax = subplot(v,h,idx);         plot(f,10*log10(pxx(:,idx)))         ylim([-80 0])         title(nms(idx))         tx = [ver(k);hor(l)];         ax.XTick = tx;         ax.XTickLabel = int2str(tx);     end end

Una señal producida marcando el número 8 se envía a través de un canal ruidoso. La señal recibida está tan dañada que el número no puede ser identificado por la inspección.

mys = sum(sin(2*pi*[ver(3);hor(2)].*t))'+5*randn(size(t'));  % To hear, type soundsc(mys,fs)

Calcular la densidad espectral de potencia cruzada de la señal dañada y las señales de referencia. Ventana las señales utilizando una ventana Kaiser de 512 muestras con factor de forma 5.β Trazar la magnitud de cada espectro.

[pxy,f] = cpsd(mys,tones,kaiser(512,5),100,[],fs);  for k = 1:v     for l = 1:h         idx = h*(k-1)+l;         ax = subplot(v,h,idx);         plot(f,10*log10(abs(pxy(:,idx))))         ylim([-80 0])         title(nms(idx))         tx = [ver(k);hor(l)];         ax.XTick = tx;         ax.XTickLabel = int2str(tx);     end end

El dígito en la señal dañada tiene el espectro con los picos más altos y el valor RMS más alto.

[~,loc] = max(rms(abs(pxy)));  digit = nms(loc)
digit =  '8' 

Argumentos de entrada

contraer todo

Señales de entrada, especificadas como vectores o matrices.

Ejemplo: especifica un sinusoides incrustado en el ruido gaussiano blanco.cos(pi/4*(0:159))+randn(1,160)

Tipos de datos: single | double
Soporte de números complejos:

Ventana, especificada como un entero o como un vector de fila o columna. Se utiliza para dividir la señal en segmentos.window

  • Si es un entero, divide y en segmentos de longitud y ventanas cada segmento con una ventana Hamming de esa longitud.windowcpsdxywindow

  • Si es un vector, se divide y se divide en segmentos de la misma longitud que el vector y windows cada segmento utilizando .windowcpsdxywindow

Si la longitud de y no se puede dividir exactamente en un número entero de segmentos con muestras superpuestas, entonces las señales se truncan en consecuencia.xynoverlap

Si especifica como vacío, utilice una ventana Hamming de forma que y se dividan en ocho segmentos con muestras superpuestas.windowcpsdxynoverlap

Para obtener una lista de las ventanas disponibles, consulte .Windows

Ejemplo: y ambos especifican una ventana De Naden de longitud + 1.hann(N+1)(1-cos(2*pi*(0:N)'/N))/2N

Tipos de datos: single | double

Número de muestras superpuestas, especificada como un entero positivo.

  • Si es escalar, debe ser menor que .windownoverlapwindow

  • Si es un vector, debe ser menor que la longitud de .windownoverlapwindow

Si especifica como vacío, utiliza un número que produce una superposición del 50% entre segmentos.noverlapcpsd Si la longitud del segmento no está especificada, la función se establece en /4.5 , donde está la longitud de las señales de entrada y salida.noverlapNN

Tipos de datos: double | single

Número de puntos DFT, especificado como un entero positivo. Si especifica como vacío, establece el parámetro ennfftcpsd max(256,2p)Dónde p = ⌈log2 N para señales de entrada de longitud .N

Tipos de datos: single | double

Frecuencia de muestreo, especificada como escalar positiva. La frecuencia de muestreo es el número de muestras por unidad de tiempo. Si la unidad de tiempo es segundos, entonces la frecuencia de muestreo tiene unidades de Hz.

Frecuencias normalizadas, especificadas como un vector de fila o columna con al menos dos elementos. Las frecuencias normalizadas están en rad/muestra.

Ejemplo: w = [pi/4 pi/2]

Tipos de datos: double

Frecuencias, especificadas como un vector de fila o columna con al menos dos elementos. Las frecuencias están en ciclos por unidad de tiempo. El tiempo unitario se especifica mediante la frecuencia de muestreo, .fs Si tiene unidades de muestras/segundo, entonces tiene unidades de Hz.fsf

Ejemplo: fs = 1000; f = [100 200]

Tipos de datos: double

Rango de frecuencia para la estimación de densidad espectral de potencia cruzada, especificado como , , o .'onesided''twosided''centered' El valor predeterminado es para las señales de valor real y para las señales de valor complejo.'onesided''twosided'

  • — Devuelve la estimación unilateral de la densidad espectral de potencia cruzada de dos señales de entrada de valor real, y .'onesided'xy Si es par, tiene /2 + 1 filas y se calcula en el intervalonfftpxynfft [0,π] rad/muestra. Si es impar, tiene ( + 1)/2 filas y el intervalo esnfftpxynfft [0,π) rad/muestra. Si especifica , los intervalos correspondientes son [0, /2] ciclos/tiempo unitario para ciclos pares y [0, /2) ciclos/tiempo unitario para impares .fsfsnfftfsnfft

  • — Devuelve la estimación de dos lados de la densidad espectral de potencia cruzada de dos señales de entrada de valor real o de valores complejos, y .'twosided'xy En este caso, tiene filas y se calcula a lo largo del intervalopxynfft [0,2π) rad/muestra. Si especifica , el intervalo es [0, ) ciclos/tiempo unitario.fsfs

  • — Devuelve la estimación centrada a dos lados de la densidad espectral de potencia cruzada de dos señales de entrada de valor real o de valores complejos, y .'centered'xy En este caso, tiene filas y se calcula a lo largo del intervalopxynfft (–π,π] rad/muestra para par ynfft (–π,π) rad/muestra para odd .nfft Si especifica , los intervalos correspondientes son (– /2, /2] ciclos/tiempo unitario para ciclos pares y (– /2, /2) ciclos/tiempo de unidad para impares .fsfsfsnfftfsfsnfft

Argumentos de salida

contraer todo

Densidad espectral de potencia cruzada, devuelta como vector, matriz o matriz tridimensional.

Frecuencias normalizadas, devueltas como un vector de columna de valor real.

  • Si es unilateral, abarca el intervalopxyw [0,π] cuando está parejo ynfft [0,π) cuando es extraño.nfft

  • Si es de dos lados, abarca el intervalopxyw [0,2π).

  • Si está centrado en CC, abarca el intervalopxyw (–π,π] cuando está parejo ynfft (–π,π) cuando es extraño.nfft

Tipos de datos: double | single

Frecuencias, devueltas como un vector de columna de valor real.

Tipos de datos: double | single

Más acerca de

contraer todo

Densidad espectral de potencia cruzada

La densidad espectral de potencia cruzada es la distribución de potencia por unidad de frecuencia y se define como

Pxy(ω)=m=Rxy(m)ejωm.

La secuencia de correlación cruzada se define como

Rxy(m)=E{xn+myn}=E{xnynm},

Dónde xn Y yn son procesos aleatorios estacionarios, –∞ < n < ∞, <n<Y E {· } es el operador de valor esperado.

Algoritmos

utiliza el método de periodograma modificado y promediado de Welch de estimación espectral.cpsd

Referencias

[1] Rabiner, Lawrence R., and B. Gold. Theory and Application of Digital Signal Processing. Englewood Cliffs, NJ: Prentice-Hall, 1975, pp. 414–419.

[2] Welch, Peter D. “The Use of the Fast Fourier Transform for the Estimation of Power Spectra: A Method Based on Time Averaging Over Short, Modified Periodograms.” IEEE® Transactions on Audio and Electroacoustics, Vol. AU-15, June 1967, pp. 70–73.

[3] Oppenheim, Alan V., Ronald W. Schafer, and John R. Buck. Discrete-Time Signal Processing. 2nd Ed. Upper Saddle River, NJ: Prentice Hall, 1999.

Capacidades ampliadas

Generación de código C/C++
Genere código C y C++ mediante MATLAB® Coder™.

Introducido antes de R2006a