Analyse av Treghetsnavigasjons Systemer

UNIVERSITETET I OSLO

Fysisk institutt

Analyse av Treghets- navigasjons Systemer

Masteroppgave

(30 studiepoeng)

Petter Moagi Lefoka

21. Mai 2013

Forord

Denne oppgaven er en del av min mastergrad elektronikk og datateknologi ved fysisk institutt på universitetet i Oslo. Jeg vil gjerne takke veilederen min, Oddvar Hallingstad for all støtte og innspill til denne oppgaven, og jeg vil takke moren min for å ha lest korrektur.

Petter Moagi Lefoka

Kjeller 21. Mai 2013

Sammendrag

I denne oppgaven er en analyse av treghetsnavigasjons systemer. Målet har hvert å finne minimums kvalitet akselerometer og gyrometer sensorer må opprettholde for treghetsnavigasjons systemer til forskjellige bruksområder.

Dette har blitt gjort ved å lage en mattematiskmodell av systemet som inneholder først èn banegenerator og et støymodellerings program der det blir generert målinger fra en kjent bane så genereres støy som inneholder hvitstøy og fargetstøy som adderes på målingene. Deretter et TNS system som gjør navigasjonssimuleringer på de stokastiske målingene, med et Kalmanfilter som minimerer driften gitt av støyen og integreringsfeil. Tilslutt er det et kovariansanalyseprogram for og se hva som må bli gjort av endringer i det teoretiske systemet for at det skal kunne opprettholde minimumsdriftkravene.

Det har blitt simulert på to forskjell forskjellige senarioer, der det første er et senario der målet er og kunne måle konturene til et lite objekt med en mobiltelefon, der bevegelses område er lite og det blir brukt nullhastighets måleoppdateringer for å korrigere for driften. I det andre senarioet er det ett modell fly der bevegelses området er mye større og måleoppdateringer blir gjort ved hjelp av GPS, i dette senarioet må sensorene være gode nokk til at driften i TNS’et blir mindre en 1 meter når det går minst 60 sekunder uten måleoppdateringer.

Notasjon 10 A INNLEDNING. . . .11 B MATEMATISKGRUNNLAG. . . .13

B.1 Navigasjonsrammer 13

B.1.1 Jordramme { e} 13

B.1.2 Geogra…skramme, { n} 13

B.1.3 Legemeramme, { b} . 13

B.1.4 Akselerometerramme og Gyrometerramme, { a} og { g} . 13

B.2 Vektorer 13

B.3 Skjevsymmetriskmatrise 13

B.4 Koordinattransformasjonsmatrise (KTM) 14

B.5 Støy 14

B.5.1 Hvitstøy 14

B.5.2 Fargetstøy 15

B.6 Diskretisering 15

C BANEGENERATOR. . . 17

C.1 Støymodellering 18

C.1.1 Mattematiskgrunnlag 19

C.1.2 Modellering 19

Akselerometerstøy 19

C.1.3 Diskretisering 21

Gyrometerstøy 21

C.1.4 Pesudokode 22

C.2 Resultat 24

D TREGHETSNAVIGASJONS SYSTEMER . . . .27

D.1 Eulers metode 27

D.2 Eulers-Heuns metode 27

E KALMANFILTER. . . .29

E.0.1 KF 29

E.0.2 LKF 29

E.1 Filtermodell 30

E.2 Nav Prosessmodell 31

INNHOLD

E.2.1 Sann-prosessmmodell 31

E.2.2 Filtermodell 32

E.2.3 Målemodell 34

E.2.4 Spektraltettheter og kovarians 34

E.3 Lineærisert diskretisert …ltermodell 36

E.3.1 LKF likningene 36

Måleoppdatering 36

E.4 Resultater av TNS med et LKF 37

E.4.1 Euler mot Euler-Heuns metoden 38

F FEILBUDSJETT . . . .41

F.1 Feilgrupper 41

F.2 Resultat av kovarians analysen 42

F.2.1 Pseudo-kode Feilbudsjett 49

G MODELFLY . . . .51

G.1 Banegenerator for modell ‡y. 51

G.1.1 Generering av deterministiske målinger 52

Koordinattransformasjonsmatrisen 53

G.1.2 Psudokode 53

G.1.3 Resultat fra den deterministiske banegenereringen for modell‡y 54

G.1.4 Deterministisk simuleringer av TNS 55

Resultat av deterministiske simuleringer 55

G.1.5 Stokastiske simuleringer av TNS 56

G.1.6 Resultat fra stokastiske simuleringer 57

H KOVARIANSANALYSE AV TNS FOR MODELFLY. . . .61 I KONKLUSJON OG VIDERE ARBEID . . . .63

I.1 Konklusjon 63

I.2 Videre arbeid 64

Referanser 65

Appendiks A: BKG3.M: 67

Appendiks B: DBG.M: 75

Appendiks C: STOY.M: 81

Appendiks D: ERRBUD.M: 85

Appendiks E: PLOT_ BANEGEN.M: 89

Appendiks F: PLOT_ KALMAN.M: 91

Appendiks G: PLOT_ NAV_ EULER.M: 95

Appendiks H: PLOT_ NAV_ HEUN.M: 97

Notasjon

Banegenerator:

R_n^b Transformasjonsmatrise fra n til b pⁿ Sann posisjonsvektor

vⁿ Sann hastighetsvektor aⁿ Sann akselerasjonsvektor fⁿ akselerasjonsvektor sett fra n

!^nb_b Vinkelhastighetsvektor sett fra n Kalman:

x Sann tilstandsvektor x feilvektor

x Prediktert feilvektor b

x Estimert feilvektor z Målevektor

H Målematrise Kontinuerlig:

x(t) Tilstandsvektor u(t) Pådragsvektor v(t) Hvitstøyvektor F(t) Prosessmatrise L(t) Pådragsmatrise G(t) Støymatrise P(t) Kovariansmatrise Q~ Spektraltetthetsmatrise Diskret:

x_k Tilstandsvektor u_k Pådragsvektor v_k Hvitstøyvektor

k Prosessmatrise

k Pådragsmatrise

k Støymatrise P_k Kovariansmatrise Q Spektraltetthetsmatrise NAV Prosessmodell:

x Sann tilstandsvektor

~

x Estimert tilstandsvektor

~

pⁿ Estimert posisjonvektor sett fra n

~

vⁿ Estimert hastighetvektor sett fra n

R~_bⁿ Estimert transformasjonsmatrise fra b til n f^b Akselerasjon sett fra b

!^bb_b Vinkelhastighet sett fra b Støyvektor gamma Støyvektor beta

v Hvitstøyvektor 10

Del A:

Innledning

Denne oppgave er basert på ‡ere tidligere masteroppgaver der temaet har hvert treghetsnavigasjons systemer. Der det blant annet tidligere har blitt lagd en Android app. Som kan logge data fra akselerome- ter og gyrometer sensorer, og det blitt lagd et program som estimerer gyroparametre fra stokastiske målinger. Det har også blitt gjort en simu- leringsoppgave der det ble generert deterministiske måledata til sensor- modeller for en enkel bane og et Kalman…lter, men Kalman…lteret er enda ikke testet på stokastisk måledata.

Målet med denne oppgaven er og hovedsakelig gjøres en analyse av treghetsnavigasjons systemer der det må undersøkes hvilken e¤ekt forskjel- lige typer støy har på et treghetsnavigasjons system og det skal un- dersøkes hvilke spesi…kasjoner de forskjellige sensorene må ha for at TNS’et skal fungere optimalt. Da må det bygges videre på den sist- nevnte oppgaven, der det må lages et støymodelleringsprogram for få testet Kalman…lteret på stokastiske målinger. Det må også lages et kovariansanalyseprogram som kan isolere e¤ekten de forskjellige støyk- ildene har på et TNS.

Videre skal det også lages en mattematiskmodell for et TNS beregnet modell‡y med en banegenerator som genererer stokastiske målinger, der sensorene er en klasse bedre enn mobiltelefon sensorer. Og det skal gjøres en kovariansanalyse på dette systemet for og se hvilken e¤ekt de samme støykildene har på et system med høyere nøyaktighetskrav.

11

B.3 Skjevsymmetriskmatrise

Del B:

Matematiskgrunnlag

Treghetsnavigasjon er navigasjon der man beregner posisjon, hastighet og orientering i rommet, ved å integrere opp akselerasjon og vinkelhastighets målinger over tid. Denne metoden trenger i teorien ikke noen eksterne referanse kilder til å beregne posisjonen, og er derfor mye brukt i romfart, luftfart og maritime sammenhenger. I praksis så er driften i et TNS en funksjon av tiden, og over tid vil derfor være nødvendig med en ekstern korrigerings kilde. I denne oppgaven skal det modelleres et TNS basert på MEMS (MicroElectroMechanical Systems) sensorer, også kjent som lavkost sensorer.

B.1 Navigasjonsrammer

B.1.1 Jordramme {e}

Jordrammen er viser posisjonen i forhold til jordkloden, den er ikke i bruk i oppgaven men er som regel en ECEF (Earth centered, earth …xed) ramme, med origo i jordens massesenter, og akser som sammenfaller med IRP (International reference pole) og IRM (International reference merdian).

B.1.2 Geogra…skramme, {n}

Den geogra…ske rammen, er en lokal representasjon av område det navigeres i. Den geogra…ske rammen kan være en geogra…skrepresentasjon av område legemet oppholder seg i der aksene til rammen sammenfaller med nord, øst og opp. Eller den kan være et enkelt kartesisk koordinatsys- tem, der origo ligger i et kjent punkt ofte p₀.

B.1.3 Legemeramme, {b}.

Legemerammen har origo i massesenteret til legemet som det navigeres fra, aksene sammenfaller med roll, pitch og yaw aksene til legemet.

B.1.4 Akselerometerramme og Gyrometerramme, {a} og {g}.

Akselerometer og gyrometer rammene ligger i hvert sitt kartesiske koordinatsystem, de har sitt origo i massesenteret til henholdsvis akselerometeret og gyrometeret. Det er her antatt at de sammenfaller med legemerammen.

B.2 Vektorer

Normen til en vektor xer de…nert ved jjxjj=p

x²₁+x²₂+x²₃

Enhetsvektoren er en retningsvektor med lengde en og er de…nert ved b= _jj^x_xjj

B.3 Skjevsymmetriskmatrise

En skjevsymmetriskmatrise A er en nxn matrise der den transponerte er lik den negative av

B Matematiskgrunnlag

matrisen, slik atA^T = A. Den skjevsymmetriskmatrisen av en vektorx;er gitt av S(x) =

2

4 0 x3 x2

x₃ 0 x₁ x2 x1 0

3

5 (B- 1)

Når den skjevsymmetriskematrisen er de…nert slik kan kryssproduktet av to vektorer skrives slik:

~ax~b <=> S(a) b (B- 2)

B.4 Koordinattransformasjonsmatrise (KTM)

Koordinattransformasjonsmatrisen transformerer koordinater mellom koordinatsystemer. I denne oppgaven blir det mellom geogra…skramme {n} og legemerammen {b}. Matrisen Rⁿ_b beskriver en transformasjon fra legemerammen {b} til geogra…skramme {n}. Hvis vektoren er bygd opp av eulervinklene ' ^T som står for roll, pitch og yaw. Når eulervinklene beskriver orienterin- gen til {b} i forhold til {n}, kan koordinattransformasjonsmatrisen de…ners som Rⁿ_b = R321( ).

Der rotasjonsmatrisenR₃₂₁( ) beskriver en rotasjon rundt z-y-x aksene.

R321( ) er di…nert som:

R321( ) =R3( ) R2(') R1( ) (B- 3)

Der R3 en rotasjon rundt z-aksen; R2 en rotasjon rundt y-aksen og R1 en rotasjon rundt x-aksen er gitt av

R₃( ₃) = 2

4 cos( ₃) sin( ₃) 0 sin( 3) cos( 3) 0

0 0 1

3

5 (B- 4)

R2( 2) = 2

4 cos( ₂) 0 sin( ₂)

0 1 0

sin( ₂) 0 cos( ₂) 3

5 (B- 5)

R1( 1) = 2

4 1 0 0

0 cos( 1) sin( 1) 0 sin( ₁) cos( ₁)

3

5 (B- 6)

Det gjør at R₃₂₁( )kan skrives slik, når c_x = cos(x) og s_x= sin(x):

R₃₂₁( ) = 2 4

c ₃c ₂ c ₃s₂s ₁ s ₃c ₁ c ₃s ₂c ₁+s ₂s ₁ s ₃c ₂ s ₃s ₂s ₁+c ₃c ₁ s ₃s₂c ₁ c ₃s ₁

s ₂ c ₂s ₁ c ₂c ₁

3

5 (B- 7)

R^b_n= (Rⁿ_b) ¹ =R^{n T}_b (B- 8)

B.5 Støy

B.5.1 Hvitstøy

Hvitstøy er de…nert som et signal med null som forventningsverdi og en uniform sannsynlighet- stetthetsfordeling med uendelig båndbredde. Et hvitstøy signal med en båndbredde som dekker hele frekvensspekteret er vanskelig å modellere, og det blir derfor ofte brukt gausiskhvitstøy som en tilnærming.

Browns bevegelse er integrert hvitstøy som fører til en tilfeldig vandring fra forventingsverdien også kalt virrevandring (random walk). Denne typen støy er hvit støy gjeldende på!_ og f_nivå gitt av

B.6 Diskretisering

_

x = w (B- 9)

der w v N(0; ²) (B- 10)

B.5.2 Fargetstøy

Fargetstøy forekommer når hvitstøy blir …lteret. Den fargetstøy prosessen som beskrives i oppgaven er en første ordens Gaus-Markov prosess. En Gaus-Markov prosess er beskrevet av

_

x = 1

Tx+w (B- 11)

der w v N(0; ²) (B- 12)

Der _T¹ er korelasjons tiden, og w er et hvitstøy element. En prosess er en første ordens Markov prosess hvis

p[x(t_k)jx(t_{k 1}); :::; x(t₁)] =p[x(t_k)jx(t_{k 1})] (B- 13) Dvs. hvis sannsynlighetstetthetsfunksjonen kun er avhengig av verdiene fra forrige tidsskritt. Hvis da hvitstøyelementet, w, har en gausisk sannsynlighetstetthetsfunksjon, kan prosessen kalles en Gauss-Markov prosess.

B.6 Diskretisering

Hvis man har et kontinuerlig-diskret system på formen _

x(t) = F x(t) +Lu(t) +Gv(t) (B- 14)

P_(t) = F P(t) +P(t)F^T +GQGe ^T (B- 15)

z_k = Hx_k+w_k (B- 16)

Skrives diskretiseringslikningene slik

(t_k+1;t_k) = e^{F(tk+1 tk)} =e^{F t} (B- 17)

u_k =

tk+1

Z

tk

(t_k+1; )Lu( )d (B- 18)

Q ^T = S =

tk+1

Z

tk

(t_k+1; )GQG~ ^{T T}(t_k+1; )d (B- 19)

Hvis man antar at u(t_k) =u_k kan det lages et nytt system på formen _dt^dxe=Fexe _

x _

u = F L

0 0 x

u der Fe= F L 0 0 e =e^{Fe t}= e₁₁ e₁₂

0 I (B- 20)

Slik at nå er x_k+1

u_k+1 = e₁₁ e₁₂

0 I

x_k

u_k (B- 21)

Nå kan systemet skrives slik x_k+1 = x_k+ u_k = e₁₁x_k+ e₁₂u_k man har da funnet og

B Matematiskgrunnlag matrisene.

For å …nne lambda kan det konstrueers en Fee matrise på formen Fee = F GQG~ ^T

0 F^T (B- 22)

e^{Fee t} = " ee₁₁ ee₁₂ 0 ee₂₂

#

(B- 23)

Slik at man kan …nne Q ^T slik Q ^T = ee₁₂ee₂₂¹:Ut ifra Q ^T kan bli funnet ved hjelp av UD faktorisering. DerS =U DU^T = Q ^T

Hvis Q = I kan stykket skrives om til I ^T = U D¹² I (U D¹²)^T ved hjelp av Cholesky faktorisering, slik at =U D¹².

Der S= ^T =U D¹² (U D¹²)^T

T =

2

4 ¹¹ 0 0

21 22 0

31 32 33

3 5

2

4 0^{11 12}22 ¹³23

0 0 ₃₃

3

5 (B- 24)

ii = vu utA_ii

i 1

X

k=1

L²_i;k (B- 25)

ji = 1

jj

(Aji i 1

X

k=1

L_j;kL_i;k) (B- 26)

Dermed har man funnet , og får det diskrete systemet

x_k+1 = x_k+ u_k+ v_k (B- 27)

P_k+1 = P_k ^T + Q ^T der Q=I og P₀= t P(t₀) (B- 28)

z_k = Hx_k+w_k (B- 29)

Del C:

Banegenerator

Banegeneratoren gir ut realistiske målinger for en forhåndsbestemt bane, de målingene brukes til og teste nøyaktigheten til navigasjons systemet. I første omgang modelleres banen til et kvadrat, der legemet (mobiltelefonen) som er i bevegelse stopper i to sekunder i hvert hjørne. Der banen sett rett oven i fra ser slik ut som i …guren under.

Banen generert av banegeneratoren

Figuren, prosessmodell banegenerator, viser blokkskjemat til banegeneratoren. Iden der er at det først blir generert støyfrie målinger som etter og ha kjørt igjennom treghetsnavigasjonslikningene gir sannposisjon, sannhastighet og sannorientering. Deretter blir støyen modellert og lagt til målin- gen slik at man har realistisk støyfylte målinger, som brukes til å teste treghetsnavigasjonssystemet.

Prosessmodell, Banegenerator

For simuleringen av mobiltelefon MEMS-sensorer så ser målingene som blir generert fra banegen- eratoren slik ut, …guren under. Der ses det at legemet roteres ikke. Langs z-aksen dvs. ned, måles det en konstant akselerasjon på 9.81 m/s^2. Og det kan ses at legemet står stille fra 0-2 sek, 4-6sek, 8-10 sek og 12-14 sek. (Legemet står også stille fra 16 til 18 sekunder, det er noe som ble lagt til senere og jeg kommer tilbake til det.) Der det kan bli gjort en nullhastighets måleoppdatering av Kalman…lteret.

C Banegenerator

Simulert målt-akselerasjon og vinkelhastighet

Banegeneratoren gir også ut simulert akselerasjon, simulert hastighet og simulert posisjon, for og kunne sammenligne avviket fra beregnet hastighet og beregnet posisjon fra navigasjonssystemet med den sanne tilstanden. Den simulerte akselerasjonen, hastigheten og posisjonen ser slik ut, illustrert i …guren ove.

Sann akselerasjon, hastighet og posisjon.

Den simulerte hastigheten og posisjonen som …guren over illustrerer er det optimale resultatet man har som mål og sitte igjen med når alt er sagt og gjort og systemet er ferdig. Resultatet videre kommer til og bli vist i avvik fra denne banen.

C.1 Støymodellering

C.1 Støymodellering

C.1.1 Mattematiskgrunnlag

Støyligningen som må modelleres er f og ! de er gitt av

f = ^a+ ^a+ ^a+v^f (C- 30)

! = ^g+ ^g+ ^g+v^! (C- 31)

De inneholder et fargetstøy element, et integrert hvitstøy element og et hvitstøy element. De inneholder også en tilfeldig konstant, men det ses bort ifra siden det leddet kan inkluderes i det integrerte hvitstøy elementet. Slik at likningene som modelleres blir

f = ^a+ ^a+v^f (C- 32)

! = ^g+ ^g+v^! (C- 33)

Der selve gamma, beta og de tilhørende hvitstøyelementene er gitt av følgende likninger, som er like for både akselerometre og gyrometre.

_ (t) = 1

T (t) +v (t) (C- 34)

der (t0) = 0 og v (t)~N(0; ²)

_ (t) =v (t) (C- 35)

der v (t)~N(0; ²)og v^x(t)~N(0; ²_x)

C.1.2 Modellering

Størrelsen til variansen i hvitstøyelementene er proposjonal med driften over tid i systemet, i første omgang er det ønsket å ha en sammenlagt drift på cirka 1 deg/s og 1 mG/s.

Ved å bruke disse likningene kan det lages to prosessmodeller som beskriver støyen. Modellene er på formen:

_

x(t) = F x(t) +Gv(t) (C- 36)

P_(t) = F P(t) +P(t)F^T +GQG~ ^T (C- 37)

Akselerometerstøy

For akselerometerstøyen kan likningene skrives slik.

_^a(t) = v^a(t) (C- 38)

_^a(t) = 1 T^a

a(t) +v^a(t) (C- 39)

v^f(t)~N(0; ²_f) (C- 40)

der v^a(t)~N(0; ²); ^a(t0) = 0 og v^a(t)~N(0; ²a) Der variansen til hvitstøyelementene settes til

C Banegenerator

2a = (1mG=s)² = (0:00981m=s³)² (C- 41)

2a = (1mG=s)² = (0:00981m=s³)² (C- 42)

2

f = q~

f = (1mG=s)² = (0:00981m=s³)² (C- 43)

og korrelasjonstiden blir satt til T^a= 2

Hvis tilstandsvektoren og støyvektoren blir de…nert slik, x^a= ^{a a T} og v^a= v^a v^{a T}: Kan systemmodellen skrives slik,

_

x^a(t) = F^ax^a(t) +G^av^a(t) der x^a(t0) er kjent: (C- 44) P_^a(t) = F^aP^a(t) +P^a(t)F^aT +G^aQ~^aG^aT der P^a(t0) ogQ~^a er kjent (C- 45) så blirF^a og G^a lik

F^a = 0 0

0 _T^Ia

(C- 46) G^a = I 0

0 I (C- 47)

Systemmodellen sier at initialverdiene x^a(t₀) og P^a(t₀); og spektraltetthetmatrisen Q~^a er kjent.

Initialverdien til tilstandsvektoren settes til null,x^a(t₀) = 0:For og …nne spektraltetthet og initiel kovorians må man ses på hvordan kovariansen utvikler seg. Utviklingen til kovariansen …nner man ved og sette opp likningene for utviklingen til varians for beta og gamma elementene, det gir:

Beta : p_^a(t) =qe^a (C- 48)

Gamma : p_^a(t) = 2

T^ap^a(t) +qe^a (C- 49)

Beta: Variansen til betaelementene øker linært siden den deriverte er lik en konstant, som er spektraltettheten til gitt varibel. Hvis den likningen integreres opp ser man at den initiele kovariansen er lik 0. p_^a(t) =qe^a=> p^a(t) =qe^a tSom gir at initiel variansp^a(0)må være lik 0 og spektraltettheten må være lik variansen ved 1 sekund.

p^a(0) = 0 (C- 50)

e

q^a = p^a(1) (C- 51)

Gamma: Kovariansen til gamma går mot en stabil tilstand når tiden går mot uendelig slik at _

p^a(1) = 0 = _T²ap^a(1) +eq^a=>

2

T^ap^a(1) =qe^a

Der vi setterp^a(1) =p^a(0)slik at

p^a(0) = ²a (C- 52)

e

q^a = 2

T^ap^a(0) (C- 53)

C.1 Støymodellering Spektraltettheten kan da skrives slik

Q~^a=

"

I qe^a 0 0 I qe^a

#

= P^a(1) 0 0 _T²aP^a(0) Q~^a= I ²a 0

0 I _T²a

2a

(C- 54) og initiel kovarians kan skrives

P^a(t₀) = 0 0 0 I ²a

(C- 55)

C.1.3 Diskretisering

Diskretiseringen av denne prosessen gir

a=e^{Fa t} (C- 56)

Spektraltettheten kan diskretiseres slik

Q^a =

I qe^aT^a

2 (1 e ^{2T at}) (C- 57)

Q^a = t I qe^a (C- 58)

q_f = t qe_f = t ²_f (C- 59)

Q^a = Q^a 0

0 Q^a (C- 60)

De diskrete hvitstøyelementene bli da beskrevet slik v^a_k~N(0; Q^a)

v^f_k~N(0; q_f)

Når spektratetthetmatrisen blir diskretisert slik blir støymatrisen lik G slik at ^a=G^a

x^a_k+1 = ^ax^a_k+ ^av^a(t_k) (C- 61)

P_k+1^a = ^aP_k^{a aT} + ^aQ^{a aT} der P_k^a= t P^a(t0) (C- 62) Der gamma og beta er gitt ved

a

k+1 = ^a_k+v^a_;k og v^a_;k~N(0; q^a) (C- 63)

a

k+1 = e ^{T a2 t a}

k+v^a_;k og v^a_;k~N(0; q^a) (C- 64)

Deretter adderes alle elementene samme til

f_k= ^a_k+ ^a_k+v^f_k (C- 65)

Gyrometerstøy

For gyrometerstøyen gjelder samme framgangsmåte. Men ved og bruke diskretiserings metoden beskrevet i innledingen (denne metoden blir brukt for alle diskretiseringer), blir koden mindre komplisert.

_^g(t) =v^g(t) der v^g(t)~N(0; ²)

C Banegenerator _^g(t) = _T¹g g(t) +v^g(t) der ^g(t0) = 0 og v^g(t)~N(0; ²g) v^!(t)~N(0; ²_!)

Der korrelasjons tidenT^g er likT^g = 1. Og variansen til hvitsøyelementene blir her satt til

2g = P^g(1) = (1 deg=s)² = (0:017453rad=s)² (C- 66)

2

g = P^g(0) = (1 deg=s)² = (0:017453rad=s)² (C- 67)

2

! = qe

f = (1 deg=s)² = (0:017453rad=s)² (C- 68)

P^g(t₀) = 0 0 0 ^2g Qe^g =

"

e q^g 0

0 qe^g

#

= P^g(1) 0 0 _T²gP^g(0)

Hvis tilstandsvektoren og hvistøyvektoren skrivesx^g = ^{g g T} og v^g = v^g v^{g T} : blir prosessmatrisenF^g og støymatrisenG^g lik

F^g = 0 0

0 _T^Ig

(C- 69) G^g = I 0

0 I (C- 70)

Systemmodellen kan da settes opp slik _

x^g(t) = F^gx^a(t) +G^gv^a(t) (C- 71)

P_^g(t) = F^gP^g(t) +P^g(t)F^gT +G^gQ~^gG^gT (C- 72) Diskretisert

x^g_k+1 = ^gx^g_k+ ^gv^g(t_k) (C- 73)

P_k+1^g = ^gP_k^{g gT} + ^gQ^{g gT} der Q^g =I og P_k^g =P^g(t0) (C- 74) diskretiseringen av v^!(t) blir her gjort på samme måte som forrige gang som var og multiplisere variansen med samplingstiden t, slik at

q! = t qe

! = t ²_! (C- 75)

v^!_k~N(0; q_!) (C- 76)

Deretter adderes alle elementene samme til

!_k= ^g_k+ ^g_k+v^!_k (C- 77)

C.1.4 Pesudokode

I =eye(3);

N = 9 200;

C.1 Støymodellering T^a= 2;

T^g = 1;

t= ₁₀₀¹ ;

2a= (0:00981^m=s_s ²)²

2g = (0:017453^rad_s )²

P^a(0) =P^a(1) =diag ²_a ²_a ²_a ; P^g(0) =P^g(1) =diag ²_g ²_g ²_g ; Qe^g = P^g(1) 0

0 _T²gP^g(0) ; Qe^a= P^a(1) 0

0 _T²aP^a(0) ; P^g(0) =P^a(0) = 0;

F^g = 0 0 0 _T^Ig

;

F^a= 0 0 0 _T^Ia

; G=G^g =G^a=I;

x^g(0) = ^g(0) ^g(0) ^T = 0 chol(P^a(0); lower) rand(3;1) ^T ; x^a(0) = ^a(0) ^a(0) ^T = 0 chol(P^a(0); lower) rand(3;1) ^T ;

g g =diskretisering(F^g; G; t;Qe^g);

a a =diskretisering(F^a; G; t;Qe^a);

P₁^g = t P^g(0) 0 0 P^g(0) ; P₁^a= t P^a(0) 0

0 P^a(0) ; f or k= 1 :N 1

f

x^g_k+1 = ^gx^g_k+ ^g rand(6;1);

P_k+1^g = ^gP_k^{g gT} + ^{g gT}; x^a_k+1 = ^ax^a_k+ ^a rand(6;1);

P_k+1^a = ^aP_k^{a aT} + ^{a aT}; g;

fk=x^a_k(4 : 6) +x^a_k(1 : 3) + ²_a t rand(3;1);

!_k=x^g_k(4 : 6) +x^g_k(1 : 3) + ²_g t rand(3;1);

C Banegenerator

C.2 Resultat

Forventet resultat fra støymodelleringen er å ha gamma og beta verdier innen for estimert standard avvik ca 1/3 av tiden. Figuren under viser støyen akselerometerstøy elementene gamma og beta i x-akse sensoren.

Modellering av gamma og beta for akselerometerstøyen.

Det som kan ses her er at over tid kan potensielt random walk prosessen (beta) vandre et stykke unna forventningsverdien. Mens Gauss-Markov prosessen (gamma) potensielt ikke vil være like dominerende, vil den fortsatt gi et solid bidrag til usikkerheten rundt forventningsverdien.

For gyrometerstøyen forventes det samme, men her burde støyen ligge rundt 1 deg/s. Resultatet er cirka det samme som for akselerometerstøyen, det kan ses i …guren under. Der ser også at støyen ligger rundt 1 deg/s.

C.2 Resultat

Modellering av gamma og beta for gyrometerstøyen.

Selve akselerometer og gyrometer målingene er gitt av

f~^b = f^b f (C- 78)

~

!^nb_n = !^nb_n ! (C- 79)

Når de blir addert sammen, blir resultatet som i …gurene under. Der ser man at gjennomsnittet av støyen faktisk ligger litt over 1 deg/s og 1 mG/s. Noe som vil føre til litt større drift i systemet en antatt.

C Banegenerator

Gyrometermålinger Akselerometermålinger

D.2 Eulers-Heuns metode

Del D:

Treghetsnavigasjons Systemer

Treghetsnavigasjon er en posisjonerings metode der man over tid integrerer opp akselerasjon og vinkelhastighet for å bestemme posisjon, hastighet og orientering (Roll, Pitch, Yaw). Et TNS er brukt i systemer der man ikke kan/vil bruke eksterne kilder til posisjonering, som for eksempel undervannsfartøyer eller satteliter eller det blir brukt til og forbedre nøyaktigheten til eksisterende målinger som for eksempel GPS målinger. Problemet med TNS er at driften er en funksjon av t så hvis ikke sensorene er veldig nøyaktige vil driften bli betydelig over kort tid.

D.1 Eulers metode

Eulers metode for treghets navigasjon er som følger.

d

dtpeⁿ(t) = evⁿ(t) (D- 80)

d

dtevⁿ(t) = R^b_n(t) fe^b(t) gⁿ (D- 81)

d

dtRe^b_n(t) = Re^b_n(t) S(e!^nb_b ) (D- 82)

diskretisert gir dette e

pⁿ_k+1 = peⁿ_k+ev_kⁿ t (D- 83)

e

vⁿ_k+1 = ev_kⁿ+R^b_n;k fe^b(t) t gⁿ(t) t (D- 84)

Re^b_n;k+1 = Re^b_n;k+Re^b_n;k S(!e^nb_b;k t) (D- 85)

D.2 Eulers-Heuns metode

Eulers-Heuns metode for treghets navigasjon gir et mer nøyaktig resultat, Euler-Heuns kan utrykkes ved å først skrive Euler likningene på standardform, som blir gjort på følgende måte:

_

y = f(u(t); y(t);Re^b_n(t)) (D- 86)

d

dtRe^b_n(t) = f^R(u(t);Re_n^b(t)) (D- 87)

Nåru(t) = fe^b(t) e

!^bb_b (t) og y(t) = peⁿ(t) e

vⁿ(t) vektorene blir de…nert slik, kan systemmatrisen skrives slik.

f(u(t); y(t);Re_n^b(t)) = evⁿ(t)

R^b_n(t) fe^b(t) gⁿ (D- 88)

f^R(u(t);Re_n^b(t)) = Re^b_n(t) S(!e^nb_b ) (D- 89) Euler-Heuns kan da utrykkes diskret slik

D Treghetsnavigasjons Systemer

Re^m_k+1 = Re_n;k^b +f^R(u_k;Re^b_n;k) t (D- 90)

y^m_k+1 = y_k+f(u_k; y_k;Re^b_n;k) t (D- 91)

Re^b_n;k+1 = Re_n;k^b + (f^R(u_k;Re^b_n;k) +f^R(u_k+1;Re^m_k+1)) t

2 (D- 92)

y_k+1 = y_k+ (f(u_k; y_k;Re^b_n;k) +f(u_k+1; y^m_k+1;Re^m_k+1)) t

2 (D- 93)

som gir

R^m_k+1 = Re^b_n;k+Re^b_n;k S(!e^nb_b;k t) (D- 94)

v_k+1^m = peⁿ_k+ev_kⁿ t (D- 95)

Re^b_n;k+1 = Re_n;k^b + (Re^b_n;k S(!e^nb_b;k) +R^m_k+1 S(!e^nb_b;k+1)) t

2 (D- 96)

e

pⁿ_k+1 = peⁿ_k+ (ev_kⁿ+v_k+1^m ) t

2 (D- 97)

e

vⁿ_k+1 = ev_kⁿ+ (R^b_n;k fe_k^b+R^m_k+1fe_k+1^b ) t

2 gⁿ t (D- 98)

Deterministiske simuleringer av Euler og Euler-Heuns likningene viser at de gir en integreringsfeil over tid, der Euler metodens integreringfeil synker linært mens samplingsfrekvensen øker, synker Euler-Heuns metodens integreringsfeilen kvadratisk mens samplingsfrekvensen øker. I praksis blir integreringsfeilen til Euler likningene er ca. 1000 ganger større en Euler-Heuns likningene, der det simuleres på en bane som gir konstant bevegelse over 180 sekunder og en samplingsfrekvens på 40Hz. Det er da klart at simuleringer med Euler-Heuns metoden er og foretrekke.

Deterministisk Euler Nav Deterministisk Euler-Heuns Nav

Del E:

Kalman…lter

Et Kalman…lter er en rekursiv estimator som på bakgrunn av kunnskap om prosessmodell og initialverdier gir et forbedret estimat av tilstand og feilkovarians på en måling med tilfeldig støy(, en stokastisk måling). Den estimerte tilstanden Kalman…lteret gir er ofte nærmer sannverdi fordi den har en bedre estimert usikkerhet enn målt tilstand. Filteret er rekursivt og er bare avhengig av forrige tidskrit, det gjør det til et veldig nyttig verktøy til tilstandsestimering og signalbehandling.

For å gjøre dette er første skritt og prediktere tilstanden og feilkovariansen, deretter beregnes en minimums gjennomsnitts sum til alle lineære kombinasjoner av tilsdandsfeilene (”least mean square estimation”) også korrigeres den prediktert tilstanden for gitt minimum.

E.0.1 KF

Den mattematiskemodellen for et lineær kontinuerlig-diskret system uttrykt slik _

x(t) = F(t)x(t) +L(t)u(t) +G(t)v(t) (E- 99)

z_k = Hkx_k+w_k (E- 100)

De kontinuerlig-diskrete Kalman…lterlikningen har formen

Tidsoppdatering, TO, er prediksjons fasen der likningene er gitt av:

d

dtx(t) = F(t)x(t) +L(t)u(t) (E- 101)

d

dtP(t) = F(t)P(t) +P(t)F^T(t) +G(t)QG^T(t) (E- 102) der P(t_k) = ^P_k og P(t₀) er gitt:

Måleoppdatering, MO, er korrigeringsfasen der eventuelle predikterings feil blir korrigert:

~

z_k = H_kex_k+w_k (E- 103)

K_k = P_kH^T(HP_kH^T Rd) ¹ ; P_k=P(t_k) (E- 104)

^

x_k = x_k+Kk(~z_k Hx_k) (E- 105)

P^_k = (I K_kH)P_k (E- 106)

E.0.2 LKF

Lineærisert Kalman…lter er en lineærisert variant av Kalman…lteret. Det beregner feilen mellom sann og målt verdi.

De mattematiskelikningene for en lineærisert kontinuerlig-diskret systemmodell uttrykt slik _

x(t) = f(x(t); u(t)) +G(t)v(t) (E- 107)

_

x(t) = x(t) ex(t) (E- 108)

z_k = H_kx_k (E- 109)

e

z_k = H_kex_k+w_k (E- 110)

De kontinuerlig-diskrete Kalman…lterlikningen for da formen

E Kalman…lter Tidsoppdatering TO

_

x(t) = F(t) x(t) (E- 111)

P(t) = F(t)P(t) +P(t)F^T(t) +G(t)QG^T(t) (E- 112) der P(t_k) = ^P_k og P(t₀) er gitt:

Måleoppdatering MO

z_k = z_k ze_k (E- 113)

K_k = P_kH^T(HP_kH^T Rd) ¹; P_k =P(t_k) (E- 114)

^

x_k = x_k+K_k( ~z_k H x_k) (E- 115)

^

x_k = xe_k+ ^x_k (E- 116)

P^_k = (I K_kH)P_k (E- 117)

E.1 Filtermodell

Blokkskjemaet til navigasjonssystemet er gitt ved.

Bloksjema, systemmodell uten tilbakekobling

Lineriseringen av navigasjonssystemetmodellen blir gjort fordi et Kalman…lter er egnetlig en linær opperator, og et navigasjonssystem er per de…nisjon ulineært. Dette blir gjort som illustrert i blokkskjema for en systemmodell uten tilbakekobling Der det estimeres avvik fra den drift i sys- temet og det ukjente er sann tilstand. Det løses ved å ta måleoppdatering når man vet hva sann hastigheten eller posisjonen er. Sann posisjon er kjent ved start eller når man får en ekstern po- sisjons avlesning. Sann hastighet er kjent, når man vet at systemet står stille, eller når det blir gitt en ekstern hastighets måling. For å få til det trengs en sann prosessmodell, en nav-prosessmmodell og tilhørende målemodeller.

E.2 Nav Prosessmodell

E.2 Nav Prosessmodell

Nav prosessmodellen er Euler navigasjonslikningene skrevet på standardform.

d

dtpeⁿ(t) = veⁿ(t) (E- 118)

d

dtevⁿ(t) = Rbⁿ_b(t)(fe^b(t)) gⁿ (E- 119)

d

dtReⁿ_b(t) = Rbⁿ_b(t)S(!e^nb_b (t)) (E- 120)

Der det er ønskelig og utrykke det slik F_KP^{U S} _dt^dex(t) =f(x;e u;e Reⁿ_b)

d

dtReⁿ_b(t) =fR(u;e Reⁿ_b) F_DP^{U S} ez=Hxe

Tilstandsvektoren og pådragsvektoren blir her e

x= peⁿ ev^{n T} ue=h

fe^b !e^nb_b iT

Dermed kan man utrykke systemet på standardform slik

f(x;e eu;Reⁿ_b) =

"

e vⁿ Rⁿ_b(fe^b) gⁿ

#

(E- 121) fR(eu;Reⁿ_b) =

hRbⁿ_b(t)S(!e^nb_b (t)) i

(E- 122)

E.2.1 Sann-prosessmmodell

Likningene til den sanne prosessmodell beskriver hvordan èn matematisk kan …nne sann posisjon, hastighet og orientering, uti fra Euler-navigasjons likningene hvis de sanne hvitstøyverdiene hadde vært kjent.

_

pⁿ(t) = vⁿ (E- 123)

_

vⁿ(t) = Rⁿ_b(fe^b+ f) gⁿ=Rⁿ_b(fe^{b a a} v^f(t)) gⁿ (E- 124) _^a(t) = 1

T^a

a(t) +v^a(t) (E- 125)

_^g(t) = 1 T^g

g(t) +v^g(t) (E- 126)

_^a(t) = v^a(t) (E- 127)

_^g(t) = v^g(t) (E- 128)

R_ⁿ_b(t) = Rⁿ_b S(!e^nb_b + !) =Rⁿ_b S(!e^nb_b ^{g g} v^!(t)) (E- 129) Der det er ønskelig og de…nere systemmodellen slik

S_KP^{U S} x_ =f (x; u; R_bⁿ) +G (Rⁿ_b)v R__bⁿ(t) =f_R(x; u; Rⁿ_b; v) S_DP^{U S} z_k =H x_k

E Kalman…lter

hvis tilstands vektoren, pådragsvektoren og støyvektoren de…neres slik

x = pⁿ v^{n a g a g T} (E- 130)

v = h

v^f v^! v^a v^g v^a v^g iT

(E- 131) e

u =

h fe^b !e^nb_b iT

(E- 132) vil systemet utryket på standardform bli slik

f (x; u; Rⁿ_b) = 2 66 66 66 64

vⁿ

Rⁿ_b(fe^{b a a}) gⁿ

1 T^a

a(t)

1 T^g g(t)

0 0

3 77 77 77 75

G (Rⁿ_b) = 2 66 66 66 4

0 0 0 0 0 0

0 Rⁿ_b 0 0 0 0

0 0 I 0 0 0

0 0 0 I 0 0

0 0 0 0 I 0

0 0 0 0 0 I

3 77 77 77 5 f_R(x; u; R_bⁿ; v) =

h

Rⁿ_b S(!e^nb_b ^{g g} v^!(t)) i

E.2.2 Filtermodell

Filtermodellen utledes ved å trekke sann tilstand fra beregnet tilstand. x=x x:~ Blokkskjemaet blir da slik som i …gur sett fra Kalman…lteret.

Blokkskjema av …ltermodell _

p= _dt^dpⁿ _dt^depⁿ=vⁿ evⁿ= v _

v= _dt^dvⁿ _dt^devⁿ=Rⁿ_b(fe^b+ f) gⁿ (Re_bⁿfe^b+gⁿ)

=Rⁿ_bfe^b+Rⁿ_b f Reⁿ_bfe^b

=R(")Reⁿ_bfe^b+R(")Reⁿ_b f Reⁿ_bfe^b

= (I+S("))Reⁿ_bfe^b Reⁿ_bfe^b+ (I+S("))Reⁿ_b f

E.2 Nav Prosessmodell

= (I+S(") I)Re_bⁿfe^b+ (I+S("))Reⁿ_b f

=S(")Reⁿ_bfe^b+Reⁿ_b f +S(")Reⁿ_b f

= S(Re_bⁿfe^b)"+Re_bⁿ f S(Re_bⁿ f)"

Det siste leddet i S(Reⁿ_b f)"likningen v, blir feilen kvadrert, som når feilen er liten blir ett veldig_ lite tall. Det antas derfor at S(Reⁿ_b f)"går mot null og det kan ses bort i fra. Slik at

_

v= S(Reⁿ_bfe^b)"+Reⁿ_b f, der f = ^{a a} v^f

d

dtR_bⁿ _dt^dReⁿ_b =Rⁿ_bS(!e^nb_b + !) Reⁿ_bS(!e^nb_b )

d

dt(R(")Reⁿ_b) _dt^dReⁿ_b =R(")Reⁿ_b(S(!e^nb_b ) +S( !)) Reⁿ_bS(!e^nb_b )

S(_dt^d")Reⁿ_b + (I+S("))_dt^dReⁿ_{b dt}^dReⁿ_b = (I+S("))Reⁿ_b(S(!e^nb_b ) +S( !)) Reⁿ_bS(!e^nb_b ) S(_dt^d")Reⁿ_b +_dt^dRe_bⁿ+S(")_dt^dReⁿ_{b dt}^dReⁿ_b = (Re_bⁿ+S(")Reⁿ_b)(S(!e^nb_b ) +S( !)) Reⁿ_bS(e!^nb_b )

S(_dt^d")Reⁿ_b +S(")Reⁿ_bS(!e^nb_b ) =Reⁿ_bS(!e^nb_b ) +Reⁿ_bS( !) +S(")Re_bⁿS(!e^nb_b ) +S(")Reⁿ_bS( !) Re_bⁿS(!e^nb_b ) S(_dt^d")Reⁿ_b =Reⁿ_bS( !) +S(")Reⁿ_bS( !)

I det siste leddet i S(")Reⁿ_bS( !) likningen, blir feilen også kvadrert, når feilen er liten blir det ett veldig lite tall. Det antas derfor at S(")Reⁿ_bS( !) går mot null og det kan ses bort i fra. Slik at S(_dt^d")Reⁿ_b =Reⁿ_bS( !)

S(_dt^d") =Reⁿ_bS( !)Re^b_n

d

dt"=Reⁿ_b !¹, der != ^{g g} v^!

Dermed kan feil likningene utrykes slik d

dt p(t) = v(t) (E- 133)

d

dt v(t) = S(Rⁿ_b(t) fe^b(t)) "(t) +Rⁿ_b(t) ( ^a(t) ^a(t) v^f(t)) (E- 134) d

dt"(t) = R_bⁿ(t)( ^g(t) ^g(t) v^!(t)) (E- 135) d

dt

a(t) = 1 T^a

a(t) +v^a(t) (E- 136)

d dt

g(t) = 1 T^g

g(t) +v^g(t) (E- 137)

d dt

a(t) = v^a(t) (E- 138)

d dt

g(t) = v^g(t) (E- 139)

Det er ønskelig og ha …ltermodellen uttrykt på denne formen K_KP^{U S} x_ =F(u;Re_bⁿ) x+G(Reⁿ_b)v;

K_DP^{U S} ez=z ez+w Siden z alltid er 0

¯ blir ez= ez+w

1 TNS med lavkostsensorer, Diego Mugisha, 2012

E Kalman…lter Med tilstandsvektor, pådragsvektor og støyvektor de¢ nert slik

x = pⁿ vⁿ " ^{a g a g T} (E- 140)

e u =

h fe^b !e^nb_b iT

(E- 141) v =

h

v^f v^! v^a v^g v^a v^g iT

(E- 142) vil de tilhørende matrisene bli

F(u;R~ⁿ_b) = 2 66 66 66 66 4

0 I 0 0 0 0 0

0 0 S( ~R_bⁿf~^b) R~_bⁿ 0 R~ⁿ_b 0 0 0 0 0 R~ⁿ_b 0 R~ⁿ_b

0 0 0 _T¹aI 0 0 0

0 0 0 0 _T¹gI 0 0

0 0 0 0 0 0 0

0 0 0 0 0 0 0

3 77 77 77 77 5

G^K( ~Rⁿ_b) = 2 66 66 66 66 4

0 0 0 0 0 0

R~ⁿ_b 0 0 0 0 0 0 R~ⁿ_b 0 0 0 0

0 0 I 0 0 0

0 0 0 I 0 0

0 0 0 0 I 0

0 0 0 0 0 I

3 77 77 77 77 5

E.2.3 Målemodell

Det blir kun gjort måleoppdateringer når en av tilstandene til systemmodellen er kjent. Som når posisjonen eller hastigheten blir gitt av en eksternkilde eller man vet at systemet står stille.

Målematrisene H og H^{*} er blir da gitt etter hva det gjøres en måleoppdatering på.

Der målemodellene skal være på formenS_DP^{U S} z_k=H x_k; F_DP^{U S} ez=Hxe Måleoppdateringer på posisjonen, gir, H= I 0 ,H = I 0 0 0 0 0 0 Måleoppdateringer på Hastigheten, gir,H = 0 I ,H = 0 I 0 0 0 0 0

E.2.4 Spektraltettheter og kovarians

Spektraltetthetmatrisen er gitt av

Qe = 2 66 66 66 66 4

e

q^a_fI 0 0 0 0 0

0 eq^g_!I 0 0 0 0 0 0 eq^aI 0 0 0 0 0 0 eq^gI 0 0 0 0 0 0 eq^aI 0

0 0 0 0 0 eq^gI

3 77 77 77 77 5 Kovariansmatrisen er gitt av

E.2 Nav Prosessmodell

P(t0) = 2 66 66 66 66 4

p _p(0)I 0 0 0 0 0 0

0 p _v(0)I 0 0 0 0 0

0 0 p_"(0)I 0 0 0 0

0 0 0 P^a(0) 0 0 0

0 0 0 0 P^a(0) 0 0

0 0 0 0 0 P^g(0) 0

0 0 0 0 0 0 P^g(0)

3 77 77 77 77 5

Standaravviket til initielposisjonsfeil settes til lite, og kan for mobiltelefonsensor programmet egentlig settes til null: Da initiel posisjon og orientering ikke er interessant siden det er konturen til et objekt som måles og da er banen i forhold til start interessant..

p _p(0) = (0:01m)² p v(0) = (0:01 m=s)²

p"(0) = (0:1 deg)²= (0:00175rad)²

Standardavviket til hvitstøyelementen blir i første omgang i banegeneratoren satt til 1 deg/s og 1mG/s, for akselerometer og gyrometerelementene. Spektraltetthet verdiene kan da utledes på ved hjelp av kovarianslikningene.

P^a(0) =I (1mG=s)² P^a(1) =I (1mG=s)² P^g(0) =I (1 deg=s)² P^g(1) =I (1 deg=s)²

Ut i fra det kan initialkovarians utledes på følgende måte:

d

dtP^g(t) = 2

T^gP^g(t) +eq^g (E- 143)

Kovariansen til ^g går mot en stasjonær tilstand når t går mot uendelig, _dt^dP^g(1) = 0 , slik at

2

T^gP^g(1) =qe^g derP^g(0) =P^g(1), det girP^g(0) = ^T₂^gqe^g d

dtP^g(t) =eq^g => P^g(t) =qe^g t (E- 144)

Når vi vet at kovariansen til ^g utvikler seg slik P^g(t) = qe^g t , det gir atP^g(0) = 0 De samme formlene gjelder for akslerometerstøyen slik at.

P^a(0) = T^a

2 qe^a (E- 145)

P^g(0) = T^g

2 eq^g (E- 146)

P^a(0) = 0 I (E- 147)

P^g(0) = 0 I (E- 148)

Spektraltettheten blir da satt til.