number for font all matlab signal-processing phase

matlab - for - identificar el cambio de fase entre señales



plot matlab (5)

He generado tres ondas idénticas con un cambio de fase en cada una. Por ejemplo:

t = 1:10800; % generate time vector fs = 1; % sampling frequency (seconds) A = 2; % amplitude P = 1000; % period (seconds), the time it takes for the signal to repeat itself f1 = 1/P; % number of cycles per second (i.e. how often the signal repeats itself every second). y1 = A*sin(2*pi*f1*t); % signal 1 phi = 10; % phase shift y2 = A*sin(2*pi*f1*t + phi); % signal 2 phi = 15; % phase shift y3 = A*sin(2*pi*f1*t + phi); % signal 3 YY = [y1'',y2'',y3'']; plot(t,YY)

Ahora me gustaría usar un método para detectar este cambio de fase entre las ondas. El punto de hacer esto es para que eventualmente pueda aplicar el método a datos reales e identificar los cambios de fase entre las señales.

Hasta ahora he estado pensando en calcular los espectros cruzados entre cada onda y la primera onda (es decir, sin el cambio de fase):

for i = 1:3; [Pxy,Freq] = cpsd(YY(:,1),YY(:,i)); coP = real(Pxy); quadP = imag(Pxy); phase(:,i) = atan2(coP,quadP); end

Pero no estoy seguro de si esto tiene algún sentido.

¿Alguien más ha hecho algo similar a esto? El resultado deseado debe mostrar un cambio de fase en 10 y 15 para las ondas 2 y 3, respectivamente.

Cualquier consejo sería apreciado.

Aquí está la pequeña modificación de su código: phi = 10 está realmente en grado, luego, en la función seno, la información de la fase se expresa principalmente en radianes, por lo que necesita cambiar deg2rad(phi) siguiente manera:

t = 1:10800; % generate time vector fs = 1; % sampling frequency (seconds) A = 2; % amplitude P = 1000; % period (seconds), the time it takes for the signal to repeat itself f1 = 1/P; % number of cycles per second (i.e. how often the signal repeats itself every second). y1 = A*sin(2*pi*f1*t); % signal 1 phi = deg2rad(10); % phase shift y2 = A*sin(2*pi*f1*t + phi); % signal 2 phi = deg2rad(15); % phase shift y3 = A*sin(2*pi*f1*t + phi); % signal 3 YY = [y1'',y2'',y3'']; plot(t,YY)

entonces utilizando el método de dominio de frecuencia como se menciona chipaudette

fft_y1 = fft(y1); fft_y2 = fft(y2); phase_rad = angle(fft_y1(1:end/2)/fft_y2(1:end/2)); phase_deg = rad2deg(angle(fft_y1(1:end/2)/fft_y2(1:end/2)));

ahora esto le dará una estimación de cambio de fase con error = +-0.2145

Con las señales correctas:

t = 1:10800; % generate time vector fs = 1; % sampling frequency (seconds) A = 2; % amplitude P = 1000; % period (seconds), the time it takes for the signal to repeat itself f1 = 1/P; % number of cycles per second (i.e. how often the signal repeats itself every second). y1 = A*sin(2*pi*f1*t); % signal 1 phi = 10*pi/180; % phase shift in radians y2 = A*sin(2*pi*f1*t + phi); % signal 2 phi = 15*pi/180; % phase shift in radians y3 = A*sin(2*pi*f1*t + phi); % signal 3

Lo siguiente debería funcionar:

>> acos(dot(y1,y2)/(norm(y1)*norm(y2))) >> ans*180/pi ans = 9.9332 >> acos(dot(y1,y3)/(norm(y1)*norm(y3))) ans = 0.25980 >> ans*180/pi ans = 14.885

Si eso es lo suficientemente bueno o no para sus señales "reales", solo usted puede decirlo.

El código de Matlab para encontrar la fase relativa utilizando la correlación cruzada:

fr = 20; % input signal freq timeStep = 1e-4; t = 0:timeStep:50; % time vector y1 = sin(2*pi*t); % reference signal ph = 0.5; % phase difference to be detected in radians y2 = 0.9 * sin(2*pi*t + ph); % signal, the phase of which, is to be measured relative to the reference signal [c,lag]=xcorr(y1,y2); % calc. cross-corel-n [maxC,I]=max(c); % find max PH = (lag(I) * timeStep) * 2 * pi; % calculated phase in radians >> PH PH = 0.4995

Hay varias maneras de medir el desplazamiento de fase entre señales. Entre su respuesta, los comentarios debajo de su respuesta y las otras respuestas, ha obtenido la mayoría de las opciones. La elección específica de la técnica se basa generalmente en temas como:

  • Ruidoso o limpio : ¿Hay ruido en su señal?
  • Componente múltiple o componente único : ¿hay más de un tipo de señal dentro de su grabación (múltiples tonos en múltiples frecuencias que se mueven en diferentes direcciones)? O, ¿hay una sola señal, como en tu ejemplo de onda sinusoidal?
  • Instantáneo o promediado : ¿Busca el promedio de desfase en toda su grabación o está buscando cómo cambia la fase a lo largo de la grabación?

Dependiendo de su respuesta a estas preguntas, podría considerar las siguientes técnicas:

  • Correlación cruzada : use el comando a como [c,lag]=xcorr(y1,y2); para obtener la correlación cruzada entre las dos señales. Esto funciona en las señales de dominio de tiempo originales. [maxC,I]=max(c); el índice donde c es máximo ( [maxC,I]=max(c); ) y luego obtienes tu valor de retraso en unidades de muestras lag = lag(I); . Este enfoque le da el retraso de fase promedio para toda la grabación. Requiere que su señal de interés en la grabación sea más fuerte que cualquier otra cosa en su grabación ... en otras palabras, es sensible al ruido y otras interferencias.

  • Dominio de frecuencia : aquí convierte sus señales en el dominio de frecuencia (usando fft o cpsd o lo que sea). Luego, encontrará el contenedor que corresponde a la frecuencia que le interesa y obtendrá el ángulo entre las dos señales. Entonces, por ejemplo, si bin # 18 corresponde a la frecuencia de su señal, obtendría el retardo de fase en radianes a través de phase_rad = angle(fft_y1(18)/fft_y2(18)); . Si sus señales tienen una frecuencia constante, este es un enfoque excelente porque, naturalmente, rechaza todo el ruido y la interferencia en otras frecuencias. Puede tener una interferencia realmente fuerte en una frecuencia, pero aún puede obtener limpiamente su señal en otra frecuencia. Esta técnica no es la mejor para señales que cambian de frecuencia durante la ventana de análisis FFT.

  • Transformada de Hilbert : una tercera técnica, a menudo pasada por alto, es convertir su señal de dominio de tiempo en una señal analítica a través de la transformada de Hilbert: y1_h = hilbert(y1); . Una vez que haces esto, tu señal es un vector de números complejos. Un vector que contenga una onda sinusoidal simple en el dominio del tiempo ahora será un vector de números complejos cuya magnitud es constante y cuya fase está cambiando en sincronía con la onda sinusoidal original. Esta técnica le permite obtener el retraso de fase instantáneo entre dos señales ... es potente: phase_rad = angle(y1_h ./ y2_h); o phase_rad = wrap(angle(y1_h) - angle(y2_h)); . La principal limitación de este enfoque es que su señal debe ser mono-componente, lo que significa que su señal de interés debe dominar su grabación. Por lo tanto, es posible que tenga que filtrar cualquier interferencia sustancial que pueda existir.

Para dos señales sinusoidales, la fase del coeficiente de correlación complejo le proporciona lo que desea. Solo puedo darte un ejemplo de python (usando scipy) ya que no tengo un matlab para probarlo.

x1 = sin( 0.1*arange(1024) ) x2 = sin( 0.1*arange(1024) + 0.456) x1h = hilbert(x1) x2h = hilbert(x2) c = inner( x1h, conj(x2h) ) / sqrt( inner(x1h,conj(x1h)) * inner(x2h,conj(x2h)) ) phase_diff = angle(c)

Hay una función corrcoeff en matlab, que también debería funcionar (el pitón uno descarta la parte imaginaria). Es decir, c = corrcoeff (x1h, x2h) debería funcionar en matlab.