Astronomy Answers: Frequency Determination

Astronomy Answers
Frequency Determination


[Astronomy Answers] [Dictionary] [AnswerBook] [Universe Family Tree] [Science] [Starry Sky] [Planet Positions] [Calculate] [Colophon]

1. Uit tijdstippen van gebeurtenissen ... 1.1. De Methode ... 1.2. Voorbeelden ... 1.3. De Uitleg

\(\def\|{&}\DeclareMathOperator{\D}{\bigtriangleup\!} \DeclareMathOperator{\d}{\text{d}\!}\)

\( \DeclareMathOperator{\sinc}{sinc} \DeclareMathOperator{\sgn}{sgn} \)

This page explains some relatively quick methods for estimating characteristics (such as the frequency) of periodic components of a set of measurements.

1. Uit tijdstippen van gebeurtenissen

Stel, je hebt een reeks van \( n \) oplopende tijdstippen \( t_i \) (zonder dubbele), waarop een interessante gebeurtenis voorkwam. Hoe kun je perioden \( P \) in die reeks vinden?

Hieronder staat een methode beschreven die ik bedacht om dit probleem op te lossen. Het blijkt beter te passen om die beschrijving te doen aan de hand van frequentie \( ν \) in plaats van periode \( P \). Ze zijn verbonden door \( νP = 1 \), dus \( ν = 1/P \) en \( P = 1/ν \).

1.1. De Methode

  1. Bereken de tijdsduur \( T \) van de meetreeks, het verschil tussen de laatste en eerste meettijd:

    \begin{equation} T = t_{\text{max}} − t_{\text{min}} \end{equation}

  2. Bepaal het kleinste verschil \( \D t_{\text{min}} \) tussen de opvolgende tijdstippen:

    \begin{equation} \D t_{\text{min}} = \min_i(t_i − t_{i−1}) \end{equation}

  3. Doe nu voor alle frequenties \( f_m \) tussen \( 2/T \) en \( 1/\D t_{\text{min}} \) in stappen van \( 1/(2T) \) het volgende:

    1. Bereken de fase \( φ_i \) van elk tijdstip \( t_i \):

      \begin{equation} φ_i = f_mt_i \bmod 1 = f_mt_i − ⌊f_mt_i⌋ \end{equation}

      Hier betekent \( \bmod 1 \) dat je de rest bepaalt bij deling door 1, en betekent \( ⌊\cdot⌋ \) dat je afrond naar beneden naar het dichtstbijzijnde gehele getal. Er geldt \( 0 \le φ_i \lt 1 \).

    2. Sorteer de fases in oplopende volgorde. Die gesorteerde fases noemen we \( φ'_k \), met \( k \) lopend van 1 tot \( n \).

    3. De rest volgt later...

  4. Als de hoogste waarde van \( s \) duidelijk hoger is dan \( 1/n \) dan is dat een goede aanwijzing dat de bijbehorende frequentie \( ν_m \) in de buurt ligt van een frequentie \( ν_0 \) die bij de metingen past.

    Als \( s_m = 1 \) dan is \( ν_0 = ν_m \).

    Als \( s_m \lt 1 \) dan kan dat zijn omdat \( ν_m \) niet precies gelijk is aan \( ν_0 \) (omdat \( ν_0 \) niet in het rijtje uitgeprobeerde frequenties zit), of omdat er tijdstippen van meer dan één periode in de meetreeks zitten, of omdat de periode van het verschijnsel een beetje met de tijd varieert, of omdat er metingen ontbreken, of omdat er ruis op de metingen zit.

    Je kunt nu het frequentieinterval van \( ν_{m−1} \) tot \( ν_{m+1} \) verdelen in (bijvoorbeeld) 10 gelijke delen en dan bovenstaand recept herhalen voor die nieuwe lijst van frequenties en op die manier steeds dichter bij \( ν_0 \) te komen.

    Ik ben nog op zoek naar manieren om sneller bij \( ν_0 \) te komen die ook werken als de meetreeks op diverse manieren niet perfect is.

1.2. Voorbeelden

Als voorbeeld nemen we een meetreeks waarin twee frequenties zitten, \( ν_1 = 0.087 \) en \( ν_2 = 0.136 \). Dat betekent dat er in de eerste reeks elke \( 1/0.087 ≈ 11.49 \) seconden iets gebeurt, en in de tweede reeks elke \( 1/0.136 ≈ 7.35 \) seconden. De eerste reeks loopt van \( t = 7.3 \) tot \( t = 145.23 \) (13 punten), en de tweede reeks loopt van \( t = 2.3 \) tot \( t = 142.01 \) (20 punten).

We kijken hoe onze kwaliteitsmaten het doen apart voor de meetreeks van elk van de twee frequenties, en ook voor de meetreeks waar de metingen van beide frequenties vermengd zijn.

1.3. De Uitleg

Stel, je vermoedt dat alleen frequentie \( ν \) voorkomt in deze meetreeks, dus dat alle waarnemingen voldoen aan \( t_i = t_0 + iP = t_i = t_0 + i/ν \) dus \( νt_i = νt_0 + i \). Dan is \( νt_i = νt_0 \pmod{1} \) voor elke meting. Dit is de basis voor onze methode.

Voor onze meetreeks met alleen \( ν_1 \) vinden we dat \( ν_1t_i \bmod 1 = 0.6351 \) voor alle gebeurtenissen, en voor onze meetreeks met alleen \( ν_2 \) vinden we dat \( ν_2t_i \bmod 1 = 0.3125 \).

We definiëren de fase \( φ \):

\begin{equation} φ = νt \bmod 1 \end{equation}

De fase \( φ \) is tussen 0 en net onder 1, voor elke frequentie. Bereken de fase \( φ \) voor alle tijdstippen \( t_i \).

Fig. 1: Fase voor verschillende frequenties, één frequentie in de data
Fig. 1: Fase voor verschillende frequenties, één frequentie in de data
Figuur 1 toont de fase \( φ \) als een functie van tijd \( t \) voor onze voorbeeldmeetreeks met alleen \( ν_1 \). Voor frequentie 0.087 is de fase steeds hetzelfde. Voor de andere frequenties is de fase steeds anders.


Fig. 2: Fase voor één frequentie, twee frequenties in de data
Fig. 2: Fase voor één frequentie, twee frequenties in de data
Figuur 2 toont de fase \( φ \) als functie van tijd \( t \) voor de meetreeks met beide frequenties erin, met \( ν_1 \) als de testfrequentie waarmee de fase berekend is. We zien veel meetpunten met dezelfde fase 0.6351 − die horen bij \( ν_1 \). De andere meetpunten horen bij \( ν_2 \) en tonen fases die ongeveer gelijk verdeeld zijn tussen 0 en 1.


Sorteer de fase \( φ \) in oplopende volgorde. \( φ_k \) is de \( k \)-de fase in oplopende volgorde van fase (met \( k \) lopend van 0 tot en met \( n − 1 \)). Onthoud de \( i \) voor elke \( k \). Voor het gemak definiëren we

\begin{equation} x = \dfrac{k}{n} \end{equation}

dus \( x \) geeft aan welke fractie van de gesorteerde meetpunten we al gehad hebben. Dan lopen zowel \( φ \) als \( x \) van 0 tot net onder 1, ongeacht het aantal meetpunten.

Fig. 3: Gesorteerde fase voor verschillende perioden
Fig. 3: Gesorteerde fase voor verschillende perioden

Figuur 3 toont de gesorteerde fase \( φ \), voor dezelfde frequenties als in figuur 1. Voor frequentie 0.087 die echt in de meetreeks voorkomt is er een horizontaal plateau op fase 0.6351 − dat zijn de meetpunten die bij die frequentie horen. Voor frequentie 0.200 die ver van de juiste frequentie ligt zijn de fases redelijk gelijkmatig verdeeld tussen 0 en 1 en is de grafiek ongeveer een rechte lijn van linksonder naar rechtsboven.


Als er meetpunten bij frequentie \( ν_0 \) horen, dan kun je met op de fase gebaseerde detectiemethoden ook veelvouden en omgekeerde veelvouden van die frequentie detecteren. Als er meetpunten bij frequentie \( ν_0 \) horen dan is de bijbehorende fase \( ν_0t \bmod 1 = φ_ν \) constant voor die meetpunten. De fase voor een veelvoud \( mν_0 \) van die frequentie (met \( m \) een heel getal groter dan 0) is dan

\begin{equation} mν_0t_i = mν_0 \left( t_0 + \dfrac{i}{ν_0} \right) = mφ_ν + mi = mφ_ν \pmod{1} \end{equation}

dus ook dan geven al die meetpunten een constante fase, die een veelvoud (mod 1) is van \( φ_ν \). De fase van een omgekeerd veelvoud \( ν_0/m \) van de beginfrequentie is

\begin{equation} \dfrac{ν_0t_i}{m} = \dfrac{ν_0}{m} \left( t_0 + \dfrac{i}{ν_0} \right) = \dfrac{φ_ν + i}{m} \pmod{1} \end{equation}

dus zijn de meetpunten dan verdeeld over \( m \) klassen met in elke klasse een vaste fase die een veelvoud van \( 1/m \) verschoven is ten opzichte van \( φ_ν \).

Fig. 4: Gesorteerde fase voor (omgekeerde) veelvouden van de frequentie
Fig. 4: Gesorteerde fase voor (omgekeerde) veelvouden van de frequentie

Figuur 4 toont de gesorteerde fase \( φ \) voor veelvouden en omgekeerde veelvouden van frequentie \( ν = ν_0 = 0.087 \) waarbij een deel van de meetpunten horen.

Het horizontale plateau bij \( φ = φ_ν = 0.635 \) wordt gevormd door de meetpunten die bij \( ν_0 \) horen.

Voor een veelvoud van die basisfrequentie \( ν_0 \) verschuift het plateau naar hetzelfde veelvoud (modulus 1) van \( φ_ν \). Dus voor \( 2ν_0 \) verschuift het plateau naar \( 2φ_ν \bmod 1 = 2×0.635 \bmod 1 = 0.270 \), en voor \( 3ν_0 \) naar \( 3φ_ν \bmod 1 = 0.905 \).

Voor een omgekeerd veelvoud van \( ν_0 \) (dus \( ν_0 \) gedeeld door een heel getal) wordt het plateau opgesplitst in kortere plateau's, elk bij een andere fase. Voor frequentie 0.087/2 verschuift de helft van het plateau naar fase \( φ_ν/2 = 0.318 \) en de andere helft naar \( (φ_ν + 1)/2 = 0.818 \). Voor frequentie 0.087/3 wordt het plateau net zo verdeeld in drieën.


\( \D φ \) is het verschil tussen opeenvolgende gesorteerde fasen:

\begin{equation} \D φ_k = \left\{ \begin{split} φ_k − φ_{k−1} \| \qquad (1 \le k \lt n) \\ φ_0 + 1 − φ_{n−1} \| \qquad (k = 0) \end{split} \right. \end{equation}

De reden voor de aanwezigheid van \( \D φ_0 \) is dat zonder die regel de resultaten af kunnen hangen van de keuze van het nulpunt van de tijd (\( t = 0 \)), en dat is ongewenst. Nu is de som van de gesorteerdefaseverschillen altijd gelijk aan 1.

Het zou rekentechnisch fijn zijn als we als maat voor de kwaliteit van de periode konden gebruiken:

\begin{equation} s_1 = \sum_{k=0}^{n−1} f(\D φ_k) \end{equation}

waarin \( f(\D φ) \) een nog nader te bepalen functie is. Helaas blijkt dit in de praktijk alleen te kunnen werken voor meetreeksen waarin nagenoeg alle meetpunten bij één periode horen. We moeten dus niet alleen de losse faseverschillen gebruiken, maar ook naar hun volgorde kijken, of de fases zelf gebruiken.

Bekijk figuur 3 nog eens. Voor een verkeerde periode is de grafiek ongeveer een rechte lijn van linksonder naar rechtsboven. Die lijn noemen we hier de diagonaal. Voor de goede periode komen sommige punten van de grafiek ver van de diagonaal, en dat kunnen we gebruiken als maat voor de kwaliteit van de periode. Het is daarbij een probleem dat hoever de meetpunten voor de goede periode van de diagonaal zijn afhangt van het gekozen nulpunt van de tijd.

Dit is te zien in figuur 4 voor de grafieken die horen bij de basisfrequentie en veelvouden (2× en 3×) daarvan. De plateau's zijn dan even lang maar liggen bij verschillende fases. Merk op dat het plateau grotendeels rechts van de diagnonaal ligt als het plateau een lage fase heeft, en grotendeels links van de diagnonaal als het plateau een hoge fase heeft. Als het plateau een fase rond 0.5 heeft dan ligt het ongeveer evenveel links en rechts van de diagnonaal. Daardoor hangt een maat zoals de gemiddelde afstand van de meetpunten tot de diagonaal af van de fase van het plateau.

Als we een kwaliteitsmaat willen baseren op de afstand van de meetpunten van de diagonaal, dan is één manier om met die afhankelijkheid om te gaan om de kwaliteitsmaat uit te rekenen voor

\begin{equation} φ' = (φ + c) \bmod 1 \end{equation}

voor veel faseverschuivingen \( c \) regelmatig verdeeld tussen 0 en 1, en dan de uitkomsten op een toepasselijke manier samen te voegen. In de praktijk doen kwaliteitsmaten gebaseerd op de afstand tot de diagonaal het ook niet zo goed.

De beste kwaliteitsmaat die ik heb gevonden is gebaseerd op het ontdekken van opeenvolgende lage \( \D φ \), dus van de plateau's in grafieken van gesorteerde \( φ \) versus \( x \), zoals figuren 3 en 4.

Stel, de frequentie \( ν \) die we proberen ligt heel dicht bij frequentie \( ν_0 \) die in de data voorkomt, dus

\begin{equation} \D ν = ν − ν_0 \end{equation}

is klein. \( m \) van de \( n \) meetpunten horen bij frequentie \( ν_0 \), en de rest niet. Dan zal er in een grafiek zoals figuur 3 een plateau zijn dat een fractie \( m/n \) breed is in \( x \). De tijdstippen van de \( m \) meetpunten zijn

\begin{equation} t_i = t_0 + \dfrac{i}{ν_0} \end{equation}

voor hele getallen \( i \), en de overeenkomstige fases uitgerekend met frequentie \( ν \) zijn

\begin{equation} φ_i = νt_i = νt_0 + \dfrac{iν}{ν_0} = φ_0 + \dfrac{i\D ν}{ν_0} \pmod{1} \end{equation}

Als de \( i \) varieert van \( i_{\text{min}} \) tot \( i_{\text{max}} \) (met \( i_{\text{max}} − i_{\text{min}} + 1 \ge m \)) dan is het gesorteerde faseverloop

\begin{equation} Φ = |φ_{i_{\text{max}}} − φ_{i_{\text{min}}}| = (i_{\text{max}} − i_{\text{min}}) \dfrac{|\D ν|}{ν_0} ≈ T \left| \D ν \right|\end{equation}

Alleen als \( |Φ| \ll 1 \) kunnen we detecteren dat we in de buurt van de juiste frequentie zijn. We kiezen \( Φ_0 = 1/10 \ll 1 \) als toegestane bovengrens van \( Φ \). We willen de frequentie ook kunnen detecteren als alle \( n \) meetpunten bij de frequentie horen, dus de gemiddelde fasestap op het plateau moet voldoen aan \( \D φ \le Φ_0/n \).

  1. Begin bij \( k = −1 \).
  2. Zoek de eerstvolgende \( k \) waarvoor \( \D φ_k \le Φ_0/n \). Als die er niet is (vóór \( k = n \)) dan ben je klaar. Als die er wel is dan heb je het begin van een kandidaat-plateau gevonden. Die \( k \) noemen we \( k_1 \).
  3. Ga verder zolang de gemiddelde fasestap op het tot dan toe gevonden kandidaat-plateau niet groter is dan \( Φ_0/n \). Als je bij \( k = n − 1 \) komt en ook daarvoor aan die eis is voldaan, ga dan verder met \( k = 0 \). De laatste \( k \) waarvoor aan de eis is voldaan noemen we \( k_2 \). Als \( k_2 \lt k_1 \) en er is overlap met het eerstgevonden plateau, gooi dan dat eerstgevonden plateau weg.
  4. Als \( k_2 = k_1 \) dan is er maar 1 faseverschil in het kandidaat-plateau, en dus maar hooguit 2 bijbehorende meetpunten. Dan wijzen we het kandidaat-plateau af en gaan terug naar punt 2.
  5. De kwaliteitsmaat \( s_2 \) voor de frequentie \( ν \) is

    \begin{equation} s_2 = (k_2 − k_1 + 1) \bmod n \end{equation}

    Met deze kwaliteitsmaat kun je zien dat je in de buurt bent van een gezochte frequentie, en ongeveer hoeveel meetpunten er bij die frequentie horen, maar niet zo goed hoe ver je van de gezochte frequentie bent.

  6. Bekijk de verzameling van fasestappen van \( \D φ_{k_1} \) tot en met \( \D φ_{k_2} \). Als ze allemaal bij dezelfde frequentie horen dan zouden ze allemaal veelvouden van dezelfde waarde moeten zijn (namelijk van \( |\D ν|/ν_0 \)). Bepaal de deelverzamelingen die allemaal (praktisch) veelvouden zijn van dezelfde waarde. Als \( \D φ_a \bmod \D φ_b \lt \min(\D φ_a, \D φ_b)/100 \) dan vinden we dat \( \D φ_a \) en \( \D φ_b \) bij dezelfde deelverzameling horen. Alle deelverzamelingen met minder dan 3 leden wijzen we af. Als er geen deelverzamelingen overblijven, ga dan terug naar punt 2.
  7. Bereken voor elke overblijvende deelverzameling het gemiddelde \( μ \) van de fasestappen \( \D φ \) in die deelverzameling.

    Bekijk voor de eerste twee \( k \) in elke deelverzameling de bijbehorende \( i \). Als de \( i \) afnemen, dan nam de fase af met de tijd. Vermenigvuldig dan het gemiddelde \( μ \) met −1.

    Een schatting voor de gezochte frequentie \( ν_0 \) is dan \( ν/(1 + μ) \).

    Het aantal leden van de deelverzameling is een schatting voor hoeveel meetpunten er bij de gezochte frequentie horen.

    Als er meerdere deelverzamelingen naar dezelfde \( ν_0 \) wijzen dan tellen de meetpunten op.



[AA]

languages: [en] [nl]

//aa.quae.nl/en/reken/frequentie.html;
Last updated: 2020-07-18