Saltar al contenido

Trazar el espectro de frecuencia de una serie de datos usando Fourier

Contamos con el arreglo a esta inconveniente, al menos eso creemos. Si continuas con interrogantes dínoslo y sin tardanza

Solución:

Como mencionó @Rahul, no ha muestreado su onda sinusoidal con la suficiente frecuencia y ha introducido artefactos debido al alias. La frecuencia de Sin[500 x]=Sin[2 Pi f x] es $f=500/(2pi)$, que es aproximadamente 80 Hz. Se requieren al menos dos muestras por ciclo para evitar la formación de alias, por lo tanto, el intervalo predeterminado $x$ de 1 en x,0,100 debe reducirse a menos de aproximadamente $1/(2*80)=0,006$.

Además, la transformada de Fourier rápida discreta asume periodicidad. Por lo tanto, se debe tener cuidado para hacer coincidir los puntos finales con precisión. Un intervalo sin un múltiplo integral exacto de las longitudes de onda del seno devolverá funciones delta de Dirac borrosas.

Establecí el intervalo de muestreo en $(1/f)/4$, que es lo suficientemente pequeño para evitar el alias. Hay un número entero de longitudes de onda, 100*(2 Pi/500)/4, y los extremos coinciden.

testData = Table[[email protected][500 x], 
              x, 0, 100*(2 Pi/500)/4 - (2 Pi/500)/4, (2 Pi/500)/4];
ListLinePlot[Abs[Fourier[testData]], PlotRange->Full]

Este es un código que escribí mucho antes de que se agregaran funciones como DataRange a MMA. Es parte de un paquete más grande que escribí hace mucho tiempo y está lejos (De Verdad lejos) de elegante. Quizás el OP pueda encontrarlos útiles y, con un poco de reescritura, incluso podría hacer que se vean mejor en la versión más nueva de MMA.

El propósito de estos fragmentos es automatizar la escala de los ejes de tiempo y frecuencia en función de uno de los parámetros proporcionados como una opción para los procedimientos AddXXXRange. Se entiende que la hora de inicio y la frecuencia de inicio son cero.

Options[Domains] =
          DeltaT -> None, TotalT -> None, DeltaF -> None, TotalF -> None,
      Toffset -> 0, Centered -> False,OnlyParameters -> False;

Domains[n_, opts___] :=
  Module[
    tos, Ts, Ttot, Fs, Ftot,
      centered, param, timedomain, frequencydomain,

  Ts = DeltaT /. opts /. Options[Domains];
  Ttot = TotalT /. opts /. Options[Domains];
  Fs = DeltaF /. opts /. Options[Domains];
  Ftot = TotalF /. opts /. Options[Domains];
  tos = Toffset /. opts /. Options[Domains];
  centered = Centered /. opts /. Options[Domains];
  param = OnlyParameters /. opts /. Options[Domains];

  Switch[
      First[Flatten[Position[N[Ts, Ttot, Fs, Ftot, tos], _?NumberQ]]],
      1, Ttot = n Ts; Fs = 1/Ttot; Ftot = 1/Ts,
      2, Ts = Ttot/n; Fs = 1/Ttot; Ftot = 1/Ts,
      3, Ts = 1/(n Fs); Ttot = 1/Fs; Ftot = n Fs,
      4, Ts = 1/Ftot; Ttot = n Ts; Fs = Ftot/n,
      _, Ts = 1; Ttot = n; Fs = 1/n; Ftot = 1
  ];

  If[param,
      Return[Ts, Ttot, Fs, Ftot],
      timedomain = tos + (Range[n] - 1)Ts;
      frequencydomain =
          If[centered,
              (Range[n] - Ceiling[n/2])Fs,
              (Range[n] - 1)Fs
          ];
      Return[timedomain, frequencydomain]
  ]
]

Options[AddSignalRange] = Toffset -> 0;

AddSignalRange[data_List, opts___] :=
  Transpose[Domains[Length[data], opts][[1]], data]

Options[AddSpectrumRange] = 
      Centered -> False,
      HalfSpectrum -> False
      ;

AddSpectrumRange[data_List, opts___] :=
  Module[
      n, lista, half, centered,

      half = HalfSpectrum /. opts /. Options[AddSpectrumRange];
      If[half,
          centered = False,
          centered = Centered /. opts /. Options[AddSpectrumRange];
      ];

      n = Length[data];
      lista = If[centered,

              Transpose[
                          Domains[n, Centered -> True, opts][[2]],
                              RotateRight[data, Ceiling[n/2] - 1]
              ],

              Transpose[
                          Domains[n, Centered -> False, opts][[2]],
                              data
              ]   
      ];
      If[half, (*then centered is forced to False*)
                      Take[lista, 1, Ceiling[(n + 1)/2]],
                      lista (*can be centered or not *)
  ]
]

Ahora, tomemos una señal de ejemplo con componentes a 4 y 12 Hz, y probemos con una frecuencia de muestreo de 30 Hz. Si queremos 150 puntos de muestra, abarcaremos un intervalo de tiempo que va de 0 a 149 1/Fs. Siendo Fs la frecuencia de muestreo.

n = 150; Fs = 30;
data = Table[.7Sin[2Pi 4 t] - .4Cos[2Pi 12 t] + .2(Random[] - .5),
        t, 0, (n - 1)/Fs, 1/Fs];

Cuando graficamos los datos desnudos, tenemos n en los ejes x

ListLinePlot[data, PlotRange -> Full]

Podemos agregar la escala de tiempo correcta con AddSignalRange, especificando la frecuencia de muestreo utilizada para recopilar esos datos

signal = AddSignalRange[data, TotalF -> Fs];
ListLinePlot[signal, PlotRange -> Full]

Lo mismo ocurre con la FFT de los datos de muestra. Podemos realizar nuestra operación en los datos desnudos.

mag = Abs[Fourier[data]];

y luego aplicamos el rango de frecuencia correcto usando AddSpectumRange. También podemos especificar que queremos un espectro desplegado con frecuencias negativas y positivas

spectrum = AddSpectrumRange[mag, TotalF -> Fs, Centered -> True];
ListLinePlot[spectrum, PlotRange -> Full]

o simplemente obtenga la mitad positiva

spectrum = AddSpectrumRange[mag, TotalF -> Fs, HalfSpectrum -> True];
ListLinePlot[spectrum, PlotRange -> Full]

Los picos están a 4 y 12 hercios.

Si posees algún inconveniente o capacidad de refinar nuestro escrito eres capaz de ejecutar una nota y con gusto lo observaremos.

¡Haz clic para puntuar esta entrada!
(Votos: 0 Promedio: 0)



Utiliza Nuestro Buscador

Deja una respuesta

Tu dirección de correo electrónico no será publicada. Los campos obligatorios están marcados con *