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

La 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 periodograma modificado promediado de Welch de estimación espectral.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, 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 de columna por columna.

  • Si y son matrices con el mismo número de filas pero diferentes números de columnas, 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 en columna:.xycpsdpxy(:,n) = cpsd(x(:,n),y(:,n)) Para obtener una matriz de múltiples entradas/salidas múltiples, anexe a la lista de argumentos.'mimo'

De verdad y, devuelve un CPSD unilateral.xycpsd Para complejos o, devuelve un CPSD de dos lados.xycpsd

pxy = cpsd(x,y,window) 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 transformada de Fourier discreta.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 la 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.ffsfscpsd Para introducir una frecuencia de muestreo y seguir utilizar 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 en 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 en las frecuencias especificadas en.f

[___] = cpsd(x,y,___,freqrange) Devuelve la estimación de la densidad espectral de potencia cruzada sobre el rango de frecuencias 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 la densidad espectral de potencia cruzada en la ventana de la 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 solapamiento.

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 sinusoides de 2 2 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 varianzas unitarias.

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));

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

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

Por defecto, trabaja columna por columna en las entradas matriciales del mismo tamaño.cpsd Cada canal alcanza picos en las frecuencias de las 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 de y el primer canal de muestra un pico mejorado a la frecuencia común de 300 Hz.qr

Genere señales sinusoidales de 2 100 Hz muestreadas a 1 kHz para 296 ms. Una de las sinusoides se retrasa la otra por 2,5 ms, equivalente a un desfase de fase de/2.π Ambas señales están incrustadas en blanco Gaussiano ruido de varianza 1/4 ².

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));

Calcule y trace la magnitud de la densidad espectral de potencia cruzada. Utilice la configuración predeterminada para.cpsd Los picos de magnitud en la frecuencia donde hay una coherencia significativa entre las señales.

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

La magnitud de la trama-función de coherencia cuadrada y la fase del espectro transversal. La coordenada en la frecuencia de alta coherencia corresponde al desfase de fase entre las 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')

Genere dos

<math display="block">
<mrow>
<mi>N</mi>
</mrow>
</math>
-secuencias exponenciales de muestra,
<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 los efectos de tamaño finito.

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

Calcule y trace 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 sin solapamientos 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 grandes

<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 la ventana.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 están de acuerdo con 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 sonoras generadas cuando se marca 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. Muestrear 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

Trace el periodograma Welch de cada señal y anote las frecuencias de los componentes. Utilice una ventana de Hamming de 200 muestras para dividir las señales en segmentos no solapados y en la ventana de 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 corrompida 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)

Calcule 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.β Graficar 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 de 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 una sinusoide incrustada en el ruido Gaussiano blanco.cos(pi/4*(0:159))+randn(1,160)

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

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

  • Si es un número entero, a continuación, divide y en segmentos de longitud y ventanas de cada segmento con una ventana de Hamming de esa longitud.windowcpsdxywindow

  • Si es un vector, a continuación, divide y en segmentos de la misma longitud que el vector y ventanas de cada segmento utilizando.windowcpsdxywindow

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

Si especifica como vacío, a continuación, utiliza una ventana de Hamming tal que y se dividen en ocho segmentos con muestras superpuestas.windowcpsdxynoverlap

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

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

Tipos de datos: single | double

Número de muestras superpuestas, especificadas 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 los segmentos.noverlapcpsd Si la longitud del segmento no está especificada, la función se establece en ⌊/4,5 ⌋, donde es 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, a continuación, establece el parámetro paranfftcpsd max(256,2p)Dónde p = ⌈log2 N para las señales de entrada de longitud.N

Tipos de datos: single | double

Frecuencia de muestreo, especificada como un escalar positivo. La frecuencia de muestreo es el número de muestras por unidad de tiempo. Si la unidad de tiempo es de 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/sample.

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. La hora de la unidad 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 frecuencias para la estimación de la densidad espectral de potencia cruzada, especificada como, o.'onesided''twosided''centered' El valor predeterminado es para las señales con valores reales y para las señales con valores complejos.'onesided''twosided'

  • : Devuelve la estimación unilateral de la densidad espectral de potencia cruzada de dos señales de entrada con valores reales y.'onesided'xy Si es Even, tiene/2 + 1 filas y se calcula en el intervalonfftpxynfft [0,π] RAD/sample. Si es impar, tiene (+ 1)/2 filas y el intervalo esnfftpxynfft [0,π) RAD/sample. Si se especifica, los intervalos correspondientes son [0,/2] ciclos/unidad de tiempo para pares y [0,/2) ciclos/unidad de tiempo para impar.fsfsnfftfsnfft

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

  • : Devuelve la estimación centrada a dos caras de la densidad espectral de potencia cruzada de dos señales de entrada con valores complejos o de valor real, y.'centered'xy En este caso, tiene filas y se calcula a lo largo del intervalopxynfft (–π,π] RAD/sample para Even ynfft (–π,π) RAD/sample para Odd.nfft Si especifica, los intervalos correspondientes son (–/2,/2] ciclos/unidad de tiempo para pares y (–/2,/2) ciclos/unidad de tiempo para impar.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 es incluso ynfft [0,π) cuando es impar.nfft

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

  • Si está centrado en DC, abarca el intervalopxyw (–π,π] cuando es incluso ynfft (–π,π) cuando es impar.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 fijos de forma conjunta, –∞ < n < ∞, <n<Y E {· } es el operador de valor esperado.

Algoritmos

utiliza el método periodograma modificado 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