Déterminer la fréquence et la période d'une onde

11

Je collecte des données de température dans un réfrigérateur. Les données ressemblent à une vague. Je voudrais déterminer la période et la fréquence de l'onde (afin de pouvoir mesurer si les modifications apportées au réfrigérateur ont un effet).

J'utilise R, et je pense que je dois utiliser une FFT sur les données, mais je ne sais pas où aller à partir de là. Je suis très nouveau dans l'analyse R et du signal, donc toute aide serait grandement appréciée!

Voici la vague que je produis:

Ma vague

Voici mon code R jusqu'à présent:

require(graphics)
library(DBI)
library(RSQLite)

drv <- dbDriver("SQLite")
conn <- dbConnect(drv, dbname = "s.sqlite3")

query <- function(con, query) {
  rs <- dbSendQuery(con, query)
  data <- fetch(rs, n = -1)
  dbClearResult(rs)
  data
}

box <- query(conn, "
SELECT id,
       humidity / 10.0 as humidity,
       temp / 10.0 as temp,
       ambient_temp / 10.0 as ambient_temp,
       ambient_humidity / 10.0 as ambient_humidity,
       created_at
FROM measurements ORDER BY id DESC LIMIT 3600
")

box$x <- as.POSIXct(box$created_at, tz = "UTC")

box$x_n <- box$temp - mean(box$temp)
png(filename = "normalized.png", height = 750, width = 1000, bg = "white")
plot(box$x, box$x_n, type="l")

# Pad the de-meaned signal so the length is 10 * 3600
N_fft  <- 3600 * 10
padded <- c(box$x_n, seq(0, 0, length= (N_fft - length(box$x_n))))
X_f    <- fft(padded)
PSD    <- 10 * log10(abs(X_f) ** 2)

png(filename = "PSD.png", height = 750, width = 1000, bg = "white")
plot(PSD, type="line")

zoom <- PSD[1:300]

png(filename = "zoom.png", height = 750, width = 1000, bg = "white")
plot(zoom, type="l")

# Find the index with the highest point on the left half
index <- which(PSD == max(PSD[1:length(PSD) / 2]))

# Mark it in green on the zoomed in graph
abline(v = index, col="green")

f_s     <- 0.5 # sample rate in Hz
wave_hz <- index * (f_s / N_fft)
print(1 / (wave_hz * 60))

J'ai posté le code R avec la base de données SQLite ici .

Voici un tracé du graphique normalisé (avec la moyenne supprimée):

graphique normalisé

Jusqu'ici tout va bien. Voici le tracé de densité spectrale:

densité spectrale

Ensuite, nous zoomons sur le côté gauche de l'intrigue et marquons l'indice le plus élevé (qui est 70) avec une ligne verte:

zoomer sur le tracé spectral

Enfin, nous calculons la fréquence de l'onde. Cette vague est très lente, nous la convertissons donc en minutes par cycle et imprimons cette valeur qui est 17,14286.

Voici mes données au format délimité par des tabulations si quelqu'un d'autre veut essayer.

Merci pour l'aide! Ce problème était difficile (pour moi) mais j'ai passé un bon moment!

Aaron Patterson
la source
Aaron, je pense que la meilleure chose ici est que vous mettiez un lien vers votre fichier de données (sous forme de texte ou quelque chose) sur une boîte de dépôt, afin que je puisse le télécharger et vous donner la réponse. Sinon, il y aura beaucoup de va-et-vient. Je ne peux pas distinguer les chiffres à l'extrême gauche. :-) (Donnez-vous également la fréquence d'échantillonnage - c'est-à-dire la fréquence à laquelle vous effectuez une lecture de température).
Spacey
Ah désolé. Les données contiennent la température en degrés C, j'ai converti en degrés F pour le graphique. Ce sont les données correctes (c'est la colonne "temp").
Aaron Patterson
Le problème avec la mesure de la fréquence de cette façon est que si vous obtenez une variabilité considérable d'un cycle à l'autre, il sera beaucoup plus difficile de déterminer la fréquence moyenne - les pics se saliront ensemble - alors que le simple comptage du temps entre les excursions vous permettra de bien faire la moyenne des choses (et aussi calculer std dev, etc.). L'utilisation de l'approche FFT serait plus nécessaire s'il y avait beaucoup de bruit, mais cela ne semble pas être le cas ici.
Daniel R Hicks
+1 pour la publication, le code, les données, les tracés et un lien vers github.
nibot
@DanielRHicks Dans ce cas particulier, je ne pense pas que cela soit important, mais oui, la FFT vous donnera la moyenne de tous, alors que si nous faisions quelque chose comme un passage à zéro, nous mesurerions la durée (fréquence) de chaque cycle, puis nous pouvons déterminer si nous voulons prendre la moyenne, la médiane, le mode, etc. Bon point!
Spacey

Réponses:

7

Projet intéressant que vous avez en cours! :-)

D'un POV d'analyse de signal, c'est en fait une question simple - et oui, vous avez raison d'utiliser la FFT pour ce problème d'estimation de fréquence.

real2+imag2

Ensuite, très simplement, trouvez le maximum de l'emplacement de votre PSD. L'abscisse de ce max correspondra à votre fréquence.

Caveat Emptor, je vous donne un aperçu général, et je soupçonne que le résultat de la FFT en R sera une fréquence normalisée, auquel cas vous devrez connaître votre taux d'échantillonnage, (ce que vous faites), afin de le reconvertir en Hz. Il y a de nombreux autres détails importants que je laisse de côté, tels que votre résolution de fréquence, la taille de la FFT et le fait que vous souhaitiez probablement réduire le signal en premier, mais ce sera bien de voir d'abord un tracé.

ÉDITER:

Prenons en compte votre signal. Après que je déteste, il ressemble à ceci:

entrez la description de l'image ici

x[n]

Nfft=103600=36000.

fs=0.5Hz

x[n]X(f)10log10(|X(f)|2)

entrez la description de l'image ici

Vous pouvez voir comment il est symétrique. Si vous ignorez la dernière moitié, et regardez simplement la première moitié et zoomez sur vous, vous pouvez voir ceci:

entrez la description de l'image ici

FsNfft=1.3889e005Hz01.3889e005=0Hz11.3889e005=1.3889e005Hz701.3889e005=9.7222e004Hz

1(9.7222e004)60=17.14

Spacey
la source
@AaronPatterson J'ai édité le message, veuillez voir. Vous pouvez également ajouter vos photos directement à votre message d'origine. :-). Veuillez ajouter une image du résultat PSD que vous obtenez.
Spacey
1
Pas exactement correct si la fréquence se situe entre les cases de résultats FFT.
hotpaw2
@ hotpaw2 C'est pourquoi j'ai averti OP que je donne des perspectives générales et pourquoi j'ai besoin de voir l'intrigue. Tout de même, j'ai modifié pour ajouter les mises en garde supplémentaires.
Spacey
1
@AaronPatterson Pas de problème, heureux de vous aider. En ce qui concerne les livres, regardez Richard Lyons "Understanding DSP" - c'est un livre rapide pour commencer.
Spacey
1
1.3x105
4

Pour une forme d'onde aussi lisse et stationnaire, le comptage des points d'échantillonnage entre les croisements positifs d'une certaine valeur de seuil moyenne vous donnera une estimation de période. Regardez plusieurs périodes de franchissement de seuil pour obtenir une estimation plus moyenne ou détecter toute tendance.

hotpaw2
la source
3

Il n'est pas nécessaire de faire quoi que ce soit de compliqué: il suffit de mesurer la durée entre les pics de la forme d'onde. C'est la période. La fréquence est juste 1 divisée par la période.

Avec environ 8 cycles sur 2 heures, la fréquence est de 4 cycles par heure, soit environ 1 mHz.

nibot
la source
3
Comment puis-je faire cela par programme?
Aaron Patterson