8.3 TIDSVARIANT SPIKING FILTER OG INVERS Q-FILTRERING

Det er en del litteratur om fjerning av dempning fra seismogram tilgjengelig i dag. Clarke og Wang skrev som vi har nevnt artikler om emnet, der de anvendte et tidsvarierende spikingfilter som de knyttet til wienerfiltrering (Robinson (1966). Her vil vi utrede teorien for et tidsvarierende spikingfilter og anvende såkalt invers Q-filtrering (IQF). Dette er et meget aktuelt forskningsområde og Wang har nettopp skrevet en bok om emnet, Seismic Inverse Q Filtering, (2008)

Kapittel 3,4 og 7 har vist at vi kan inkludere dempning i Riccatiligningen ved å anvende forskjellige dempningmodeller på ligningen. Dempningsmodellene kunne tilknyttes filterfunksjoner som ble studert i forbindelse med generell viskoelastisk teori. Vi vil her vise at vi kan fjerne dempningen ved å anvende invers Q–filtrering som fjerner både viskoelastisk dempning og dispersjon.

I grove trekk kan vi si at våre filtre h demper og forsinker t’ en enhetspuls og det inverse filteret fjerner dempningen og flytter enhetspulsen tilbake til ankomsten. Robinson (1967) har gitt en regel som gjelder slike inverse filtere, og som vi vil komme tilbake til:

Enhver minimum fase operator(ensidig operator) har en entydig invers ensidig operator som også er minimum fase. Frekvensresponsene til operatorene har denne sammenheng med hverandre:

|H(iw]| |H*(iw)| = l

Der H*(iw) er den inverse til H(iw).

Det er svært enkelt i teorien å lage en transferfunksjon til den inverse generelle filterfunksjonen. Vi får ved å ta den inverse av ligning 7.2.1:

8.1

Der H*(t’,iw) er den kompleks konjungerte til H(t’,iw). Fra denne kan de enkelte inverse filtre til dempningsmodellene utledes som funksjoner i tiden. For eksempel vil det inverse filteret til dempningsmodell1 med et konstant logaritmisk dekrement bli:

H*(t’,iw) = exp( -iw(1 + ibft)

Vi kan også finne den inverse direkte i tidsdomenet. F.eks vil den inverse til H(t’,iw) være identisk med den inverse til fouriertransformen av H(t’,iw):

8.2

Der h(u) kan knyttes til vår diskusjon av dempningsdata i kapittel 6 med utrykket:

h(u) = ln((z+0,15)/0,15)

Det blir komplisert å anvende slike filtre i sin inverse analytiske form, men man har funnet enkle løsninger for det. En anvendelse er å utvikle det inverse filteret 8.1 i en Taylor-rekke. Jeg vil ikke gå mer inn på dette her.

En annen innfallsvinkel er å tolke det inverse filteret som en invers puls. Det inverse filteret kalles et spikingfilter og vil gjøre en hvilken som helst seismisk puls på et seismogram til en enhetspuls gjennom å ”spike” den. Grunnen til at dette kan være en nyttig innfallsvinkel er at den seismiske pulsen som er dempet ikke er en ideell enhetspuls som er påvirket av attenuasjon alene. Den er en skuddpuls som har gått gjennom jordfilteret, og blitt helt ugjenkjennelig i forhold til de enkle enhetspulsene vi har arbeidet med i den viskoelastiske teorien. Det vil derfor være svært vanskelig for oss å skille ut viskoelastisk dempning fra andre dempningseffekter, noe som er nødvendig dersom en korrekt fjerning av attenuasjon skal finne sted.

Fig.8.5 Et inverst (spiking) filter lages for hvert tidsområde. Vektingen øker med ankomsten.

Siden man innførte såkalt invers Q-filtrering har det vært mulig å separere den delen av ”spiking-prosessen” som gjelder fjerning av dempning i den seismiske pulsen til en egen prosess og dette skaper en god del forenklinger. Forutsetningen er at vi kjenner til hvilke dempningseffekter som har påvirket seismogrammet.

For å se dette bedre kan vi diskutere fig.8.5 i lys av den spiking-prosessen vi har begitt oss inn på. Vi kan begynne med den viskoelastiske filter-teorien fra kapittel 7. Vi så da at en isolert reflektor syntetisert ved Riccatiligningen kunne tilknyttes impulsresponsen til det tidsforsinkede filteret, der reflektorens ankomst var lik tidsforsinkelsen. Vi får en like enkel analyse for det inverse tilfellet (uten dispersjon) ut fra lign.4.7.5-8.

En reflektor i y(t) med en ankomst t’ syntetisert med Riccatiligningen når bare indre friksjon spiller inn, har en tosidet invers puls som representerer det skalerte inverse filteret til den generelle filterfunksjonen for samme tidsforsinkelse som ankomsten t'.

Hvordan er det dersom vi har flere reflektorer som løsning av Riccatiligningen når bare indre friksjon spiller inn? Det er klart vi kan benytte flere reflektorer i y(t) også i det inverse tilfellet og få representert alle mulige løsninger av det skalerte inverse filteret.

En rekke reflektorer med ankomster t’ syntetisert ved Riccatiligningen når bare indre friksjon spiller inn, har tosidede inverse pulser som representerer løsninger av den inverse til den generelle filterfunksjonen for samme tidsforsinkelser t’.

Fire slike reflektorer er representert på fig.8.5, og vi ønsker å definere inverse pulser innenfor hvert tidsområde på figuren. Når filteret er definert på en ønsket form kan det anvendes på Riccatiligningen. Selve anvendelsen blir på samme måte som når vi anvender den generelle filterfunksjonen, dvs. ved konvolvering.

I kap.4.7 viser vi hvordan filterfunksjonen anvendes i en dekonvolveringsprosess for å få den udempede refleksjonskoeffisientrekken tilbake. Dekonvolvering er derfor det samme som å konvolvere det dempede seismogrammet med den inverse til konvolveringsfunksjonen.

Fig.8.6 Det inverse filteret anvendes med sin filtervekt over inndelte områder på trasen.

Filtrering av seismogrammet ved konvolveringsprosessen er enkelt ved bruk av den generelle filterfunksjonen, fordi filteret er en kontinuerlig funksjon av sin tidsforsinkelse. Siden filteret er en funksjon av tiden vil også den inverse kunne være det. Men det er ikke gitt at det inverse filteret alltid vil kunne anvendes i en kontinuerlig konvolveringsproses – som et tidsvariant filter. Derfor må vi dele den seismiske trasen inn i områder (vinduer) og konvolvere hvert vindu med det inverse filteret som er definert for hvert vindu slik som på fig.8.5..

Er det stor variasjon i ankomstene pga. indre friksjon må vi ha mange intervaller i seismogrammet, som gir mange inverse pulser, dvs. vi optimaliserer filteret tett. Er det liten variasjon blir optimaliseringen spredt.

Anvendelsen av det definerte filteret blir tungvint da en inndeling av seismogrammet tilsvarende optimaliseringen må gjøres for hvert enkelt seismogram. De områder som seismogrammet er delt inn i kan tolkes som et tidsinvariant system. Det inverse filteret anvendes i sin tidsforsinkede form(den definerte form) og konvolveres som et tidsinvariant filter med reflektorer innenfor sitt optimaliseringsintervall. Dersom det er tett optimalisering, blir det mange områder i seismogrammet og mange numeriske beregninger. Jo tettere vi optimaliserer, jo mindre blir de tidsinvariante områdene i seismogrammet og jo mer nærmer inversjonsprosedyren seg en kontinuerlig prosess. Grensetilfellet har vi ved uendelig tett optimalisering og vi får et uendelig antall konvolveringer. Jo mer spredt vi optimaliserer jo enklere blir inversjonsprosessen. Grensetilfellet ved spredt optimalisering er at det inverse filteret er tidsinvariant og vi trenger bare en konvolvering med seismogrammet.

Tegningen på fig.8.6 viser gangen i prosedyren. Dette skjemaet er lignende det vi tegnet på fig.7.2 i det direkte tilfellet.

Det inverse filteret kan brukes til alle typer reflektorer og uansett ankomsttid innenfor det enkelte filters definerte område, men vi må regne med at jo mer komplisert seismogrammet er jo dårligere virker filteret.

Forsterkning av amplityden for pulsrekken og det inverse Q-filteret

Men vi hadde det problem at den skuddpulsen vi ser på seismogrammet har gått gjennom jordfilteret, og blitt helt ugjenkjennelig i forhold til de enkle enhetspulsene vi har arbeidet med i den viskoelastiske teorien, og som nevnt er det en forutsetning at vi kan skille ut den viskoelastiske dempningen før filteret anvendes. Derfor må vi ha en mer generell prosess i bakhånd når vi fjerner dempningen og se på invers Q-filtrering som et spesialtilfelle av dekonvolusjon.

Wiener-filtrering

Det har vært vanlig siden Robinson (1966) å lage et slikt generelt filter som er optimalisert som et inverst filter gjennom en prosedyre som kalles Wiener-filtrering. Da får man fjernet effekten av alle prosesser som har dempet ut den initiale pulsen. (som vi må kjenne til for at filteret skal fungere.) Når vi vet denne kan vi gjennom en numerisk prosess som kalles minste kvadraters metode finne den minste forskjellen vi kan regne oss til mellom en spike og hvor lik vi kan få skuddpulsen til å bli en spike gjennom vårt Wiener filter uten at vi trenger å vite hvilke prosesser som har virket på den. Vi lager det inverse filteret direkte ut fra seismogrammet ved å ta den inverse pulsen til ankomstene og finne ut hvor tett opp til en enhetspuls vi kan få den aktuelle pulsen ved å minimisere forskjellen på energien i enhetspulsen og den aktuelle pulsen ved minste kvadraters metode. Da får vi det tidsforsinkede inverse filteret direkte.

AGC-forsterkning

En langt enklere form for forsterkning får vi gjennom såkalt Automatic Gain Control (AGC). Metoden er basert på det samme prinsipp som vi har vist på figur 8.6. der man deler den seismiske trasen inn i områder som forsterkes gjennom ulik vekting. Vi definerer lengden på et inndelt område som 2L +1 og regner ut de vektede amplituder rundt dette samplingspunktet etter formelen:

8.4

der w_i er vekt-faktorer og x_i representerer amplityden til tidssamplet i.

Dersom S representerer et bruker-definert skaleringsnivå vil amplitydeverdien til midt-punktet av vinduet være:

x’ = x_i S/ x_i

Svakheten ved denne metoden er at man får forsterket opp alt innhold på seismogrammet uten at man får et skarpere bilde, og den må sees på som den groveste form for forsterkning vi kan tenke oss i seismikk. Styrken er at den er enkel i bruk og den brukes derfor i stor utstrekning.

Invers Q-filtrering

IQF er et bedre alternativ til forsterkning enn AGC da man får korrigert for den dempning på seismogrammet som skyldes attenuasjon og dette er en metode som er blitt tatt i bruk i stor utstrekning de siste årene. Jeg vil nå gå inn på det matematiske grunnlaget for denne metoden ved å knytte teorien for invers Q-filtrering til den viskoelastiske teorien i kapittel 7 og vår introduksjon av det inverse filteret i utrykket 8.1. Vi har diskutert teorien i kapittel 7 ved å se på attenuasjonskoeffisienten.

Det er enkelt å knytte denne til konstanten Q som er det mest brukte utrykk for dempning i dag. Det generelle filteret utviklet for Hortons modeller skrives på denne for oss velkjente formen fra 7.2.1:

Vår oppgave blir å finne ulike inverse funksjoner av denne.

Anvendelse av invers Q-filtrering i Riccatiligningen

Vår fremgangsmåte i anvendelsen vil følge skjema på fig. 8.8. Vi benytter en deltapuls, legger dempning på denne og anvender invers Q-filtrering for å forsterke tilbake den dempning som de viskoelastiske modellene har laget på pulsen. Dersom alt fungerer som det skal vil vi få tilbake den originale deltapulsen.Fig.8.9 utvider skjema i fig.8.8 til en full syntetiserings og inverterings-prosedyre med Riccatiligningen. Her trekker vi inn skuddpulsen som kan tolkes som en bandpassfiltrert utgave av refleksjonskoeffisientrekken og impedansen som vi har som modell for vår konstruksjon av det syntetiske seismogrammet, og som vi, dersom alt fungerer som det skal, vil komme tilbake til gjennom inversjonen.

Fig. 8.8 Skjema for testing av IQF-filter og Wienerfilter

8.9 Oversikt over syntetiserings og inverterings-prosedyren med Riccatiligningen

På fig.8.10. ser vi en pulsrekke som kan tolkes som primærer i seismogram. Vi har Q=136. Den er regnet ut etter formelen a=π/Q, der a er attenuasjonskoeffisienten vi har benyttet i kapittel 7 og som vi har satt til a=0.023.

8.4 Wiener-filtrering og viskoelastisk dempning

I kapittel 8 i Sørsdals reviderte hovedoppgave påsto vi – på bakgrunn av utrykket 4.7.1 i kapittel 4 - at vi kunne utføre inversjonen med Riccatiligningen med dekonvolusjon ved å definere den inverse til den seismiske pulsen over et avgrenset område i toveistid og konvolvere denne med seismogrammet i det valgte området. På fig. 8.5 illustrerte vi hvordan det kunne gjøres med å dele den 0,6 sek lange trasen inn i fire områder, og definere et inverst Q-filter for hvert område. Dersom dette filteret ble konvolvert med utrykket på venstre side i utrykket 4.7.1 over et ønsket tidsintervall ville vi få det dekonvolverte seismogrammet der effekten av viskoelastisk dempning ble fjernet.

Men hva om vi ikke kjenner den viskoelastiske dempningen i mediet? Da blir det ikke mulig å bruke Q-filteret fordi den seismiske pulsen består av mange andre effekter.

Riccatiligningen tar seg av noen av disse. Det har vi sett ut av iterasjonsprosessen. I utrykket på høyre side av 4.7.1. blir transmisjonstap og multippelrefleksjoner fjernet fra det opprinnelige seismogrammet. Til slutt blir impedansen beregnet. Den ferdig konvolverte sekvensen settes så inn i utrykket på høyre side, og slik forsetter vi inntil vi ikke har variasjoner i impedansen når vi itererer. Dermed har vi en svært omfattende og samtidig elegant inversjonsprosedyre for å prosessere våre data der også viskoelastisk dempning fjernes.

Men svakheten med prosessen er at Riccatiligningen krever en input-seismisk modell for γ som er mest mulig lik en refleksjonskoeffisientrekke (en spiket seismisk rekke) dersom vi skal få en nøyaktig inversjon. Og vi har ikke deltapulser som reflektorer i et reelt seismogram, men bandpassfiltrerte enhetspulser som kan defineres som Rickerpulser som på fig.3.7. Dette kan det gjøre det nødvendig å ikke bare bruke et Q-filter i dekonvolveringsprosessen, men også et generelt spiking-filter (Wiener-filter). Vi vil nå utlede teorien for dette og gjøre noen anvendelser. Vi vil på denne bakgrunn lage to definisjoner for den type filter som anvendes i inversjonen, og som vi har samlet i begrepet operatorer:

1). Det inverse filteret for dempede enhetspulser vil være et inverst Q-filter.

2). Det inverse filteret for dempede Rickerpulser vil være en invers puls.

Anvendelse av filteret definert ved 1) krever at vi kjenner den viskoelastiske dempningen og at vi kan fjerne multipler og transmisjonstap på en enkel måte i vår regneprosess.

Anvendelse av filteret definert ved 2) krever at vi først og fremst kjenner den initiale skuddpulsen.

På figur. 8.10.a ser vi grafen av dempede enhetspulser med dempning b=0.023 og figur 8.10 b viser en graf av en Rickerpulser tatt ut fra fig.8.9. På denne modellen kan vi anvende IQF på reflektorene i fig.8.10.a og få et bra resultat, men ikke på Rickerpulsen og forvente samme resultat.

Rickerpulsen motiverer oss til å bruke Wienerteori, og vi vil gå gjennom teorien for denne i den videre fremstilling. Da vil vi finne den inverse både til Rickerpulser av typen på fig. 8.10.b og til de dempede enhetspulsene i fig. 8.10 a. De ferdige data kan vi konvolvere med impulsrekken på fig.8.10 a. og på den måten få et dempet refleksjons-seismogram tilbake. Vi vil så kunne anvende det inverse Q-filteret på denne og få et ferdig invertert resultat der dempningen er fjernet.

Figur 8.10. a og b. a viser dempede enhetspulser og b er en Rickerpuls som er plukket ut på tidsaksen som tilsvarer en bandpassfiltret enhetspuls med samme ankomst

Jeg vil også nevne at det har vært nødvendig med et krafig regneverktøy for å få til Wienerfilteret og dekonvolusjon med Riccatiligningen. I appendix. 2 finner vi det FORTRAN program som utfører den numeriske beregning. Flytskjema viser gangen i prosedyren.

To bøker om dekonvolusjon ble utgitt som en Geophysics reprint series (1978) og gir et meget godt bilde av fagområdet, og jeg baserer min fremgangsmåte på disse to hefter.

Rice (1961) i sin artikkel ”inverse convolution filters” ble tatt med i reprint serien. Han har konsentrert seg om å skape et inverst ”spiking” filter som kan brukes til å rekonstruere de udempede refleksjonene fra målingene av de dempede pulser der han tar utgangspunkt i Robinson. Problemstillingen blir som på fig. 8.11. Rice bruker en negativ deltapuls som utgangspunkt, men det blir samme teori ved bruk av en positiv deltapuls som vi gjør.

Fig.8.11 Problemet er å bestemme et inverst filter som kan omforme en hvilken som helst puls til en deltapuls.

Vi tar utgangspunkt i en dempet enhetspuls i hvert enkelt område 1-4 på fig. 8.5. Vi vil da følge generell teori og plukker ut to sampler fra den dempede enhetspulsen som har z-transform:

a + bz 8.16

Så multipliserer vi dette utrykket med en responsfunksjon som kan gi den udempede enhetspulsen:

X(z) (a+bz) = 1 8.17

Da har vi :

X(z) = 1/(a+bz) som kan utvikles i potenser av z og da kan vi finne filterkoeffisientene til X(z) som dermed definerer et inverst filter.

[1/a, -b/a² , b²/a³, … ] 8.18

Dersom dette filteret konvolveres med samplene fra den dempede enhetspulsen får vi en udempet enhetspuls forutsatt at b<a og vi tar med nok koeffisienter. Vi sier at pulsen må være i minimum fase.

Dermed kan vi definere filterkoeffisienter for den inverse pulsen når vi har to samplinger av den halvdel av den dempede null-fase pulsen som har minimum-fase karakter (første sample er større enn den andre).

Vi kan få en lignende teori for den halvdelen der første sampel er mindre enn den andre og der vi får en maksimal-fase puls.

Vi kan ikke bruke utrykket over da det vil divergere, men vi kan omskrive det slik at vi får koeffisientene:

[1/b, -a/b² , a²/b³, … ] 8.19

Vi ser da at vi har to ulike filtre, men at vi kan definere et inverst filter til en dempet null-fasepuls ved å sette sammen alle filterkoeffisienter, skalere dem med maksimalkoeffisientene og benytte symmetri.

Her viser vi noen verdier for filterkoeffisientene for noen dempede enhetspulser med to samplepunkter (Område 1 a=1 b=0,5 | Område 2 a=1 b=0,5 | Område 3 a=1 b=0,5 |

	Område 1	Område 2	Område 3
-4	0,0625	0,00390625	0,0001
-3	-0,125	-0,015625	-0,001
-2	0,25	0,0625	0,01
-1	-0,5	-0,25	-0,1
0	1	1	1
1	-0,5	-0,25	-0,1
2	0,25	0,0625	0,01
3	-0,125	-0,015625	-0,001
4	0,0625	0,00390625	0,0001

Figur.8.12. Filterkoeffisienter for inverse dempede enhetspulser

Vi har tegnet et plot av hvert filter som en sammensatt invers puls.

Fig.8.13. Inverse filter for 3 dempede enhetspulser for verdier som er gitt på fig.8.12. Øverst til venstre område 1, Øverst til høyre område 2. Nederst til venstre område 3. Disse kan knyttes til pulsene fra fig.8.5.Vi får svakere oscillasjoner på den inverse pulsen, jo mer dempet den originelle pulsen er.

På grunnlag av avsnitt 8.3 i kapittel 8 kan vi sette opp en generell formel for prosedyren. Vi antar at vi opererer med et lineært tidsinvariant system slik at filtreringsprosessen defineres matematisk ved det kjente konvolveringsintegralet:

8.20

Rice påpeker at dersom det er gitt en impulsform f(t) vil det ikke finnes noen endelig realiserbar g(t) som lager en deltapuls ∂(t). Derfor er f(t) og ∂(t) inkompatible. Men vi kan finne en g(t) som lager den beste tilnærmelsen til deltapulsen.

Dersom vi antar at f(t) og g(t) er null overalt unntatt på et tidsintervall t, kan 8.20. skrives på formen:

8.21

Vi kan bruke Laplace-transformasjoner for å løse 8.20 og får:

H(s) = F(s) G(s) 8.22

8.22 kan presenteres som genererende Laplace-funksjoner (se Rice) og vi får:

der a_v er midtpunkt amplitude av f(t) 8.23

der x_v er midtpunkt amplitude av g(t)

der b_v er midtpunkt amplitude av h(t)

Siden vi har H(u) = F(u) G(u) , kan de ukjente amplituder finnes ved polynomisk divisjon av H(u) med F(u).

Vi får da en Laurant-rekke som konvergerer på enhetssirkelen:

Rice har gjort en grundig redegjørelse for konvergens av denne rekken og hva det betyr for løsningen og dannelsen av det inverse filteret. Vi skal ikke ta den diskusjonen opp her, men slå fast at dersom konvolusjon presenteres som polynomisk multiplikasjon og koeffisienter som er samsvarende samles på begge sider av likhetstegnet vil de følgende lineære ligninger i den ukjente amplitudene av x_i oppnåes:

a₁x₁ = b₁ 8.24

a₂x_{1 +} a₁x₂ = b₂

_.

_.

_.a_mx_{1 +} a _{m-1 x2} = b_n-m+1

a_mxn-m+1 _+….. = b_n

Dersom h(t) er lik en deltapuls blir b₁=1 og b_v = 0 kan vi løse ligningene som en iterasjonsprosess.

Minste kvadraters metode

Siden det ikke finnes noen eksakt løsning på problemet, må tilnærmede løsninger finnes. Den beste tilnærmelsen har vært den såkalte minste kvadraters metode. Da antar vi at alle x_v er null for v>N. Da spør vi etter det unike sett av x’er som vil mimimisere summen av kvadratene av avvikene fra de kalkulerte b_v ene fra de vi ønsker som er b₁=1 og b_v=0 for v>1.

Den beste måten å finne en løsning på er ved hjelp av en n x N - matrise. Da kan vi skrive 8.24 på formen:

AX – B = E og summen av kvadratet til avviket E blir:

X’A’AX – X’A’B – B’AX + B’B

Der merkingen ’ betyr transponerte.

Dette utrykket kan også skrives:

A’AX = A’B 8.25

Etter en del regninger kan avviket (least square errors) skrives på formen:

E’E = B’B – B’AX 8.26

Det er lett å se at X må tilfredstille AX=B som en nødvendig og tilskrekkelig betingelse for at E’E=0.

Rice har lansert et teorem:

Minste kvadraters avvik E’E avtar monotont etter som antall ikke-null termer i det inverse filteret er tillat å øke, Det betyr at man kan få en tilnærmet eksakt tilnærmelse til enhetspulsen ved å velge lengden på det inverse filteret stor nok.

Vi introduserer denne løsningen ved matrisen:

A’A = 8.27

Der der k = 0,1 …….., n-m 8.28

Med α_k = 0 for k>m-1

Nå er α_k amplitudene til autokorrelasjonen til den originale pulsen med utrykket:

8.29

Da blir høyre side av 10.1.6:

8.30

For det tilfellet at enhetspulsen plasseres på det m-te punktet (bm=1 bv=0 for v<> m) reduseres denne vektoren til utrykket:

-a_m

-a_m-1

.

.

.

-a_2m-m

Som er negative amplituder for den originale pulsen i revers orden ned til -a_2m-m

I det videre skal vi anvende Wienerfilteret på de dempede enhetspulsene fra fig, 8.10 a. og se om vi kan få inverse filtere som på fig.8.13.

Kryss og autokorrelasjon for dempet enhetspuls

Det kan være nyttig å se på løsningen av 8.28 (autokorrelasjonen) og 8.30 (krysskorrelasjonen) som en separat del av løsningen. Figur 8.14 viser dette.

8.5. Wiener-filter med minste kvadraters metode

Før vi kjører programmet vil vi forsøke å få til et spiking-filter som får likhet med de vi utledet analytisk på fig.8.13. Vi vil bruke wiener-filteret i sin opprinnelige form slik vi finner det hos Robinson (1966). Spiking-filteret er kjernen i programkoden der ”Least square” –subroutinen måler hvor nær en spike pulsen er kommet. Den best tilpassede filterkoeffisienten blir returnert fra programmet.

Fig.8.17. a.Beregning med Wienerfilter med en enkel enhetspuls av bredde 0,02 sek som filtreres til en spike. Filteret er sterkt oscillerende. Til venstre ser vi pulsen, i midten spiking-filteret og til høyre har vi grafen over feilutslag for minste kvadraters metode.

Fig.8.17. b.Beregning med Wienerfilter med en enkel enhetspuls av bredde 0,04 sek som filtreres til en spike og som representerer en puls som er mer dempet enn den over i a). Filteret er mindre oscillerende enn for den skarpere pulsen på fig.8.17.a. noe som er i overenstemmelse med våre analytiske data på fig.8.13. Til venstre ser vi pulsen, i midten spiking-filteret og til høyre har vi grafen over feilutslag for minste kvadraters metode.

Horisontalakse på alle grafer går fra 0 til 0,2 sek. (50 samplepunkter med en sampling ∆t=0,004). Effekten av dempning er introdusert med bredden av pulsen på denne aksen.

Fig.8.18. a.Beregning med Wienerfilter med en enkel enhetspuls av bredde 0,08 sek som filtreres til en spike og representerer en dempet enhetspuls. Filteret er enda mindre oscillerende enn for pulsen på fig.8.17.b. noe som er i overenstemmelse med våre analytiske data på fig.8.13. Til venstre ser vi pulsen, i midten spiking-filteret og til høyre har vi grafen over feilutslag for minste kvadraters metode.

Fig.8.18. b.(venstre) Konvolvering av dempet enhetspuls med spiking-filter fig. 8.17 a. (0,02 sek) Resultatet er pulsen til venstre som viser at vi har fått en spiking effekt med noe høyfrekvent støy. Pulsen i midten er spikingfilteret fig.8.17.b konvolvert med den 0,04 sek brede pulsen, og puls til høyre er puls fig.8.18.b 0,08 sek bred konvolvert med spikingfilteret på fig.8.17.a.

De ferdig spikede pulsene kan tyde på at Wienerfilteret ikke er god nok når vi har en null-fasepuls. Vi har fått introdusert støy sammen med pulsene. Men vi ser at spikingen er god. Dempningen er fjernet og vi har fått en skarpere puls enn inngang. Vi vil gjøre et forsøk med en puls som har en minimum-fase karakter.

Fig.8.19. a.Inversjon med spikingfilter på tre ensidige (minimum fase) enhetspulser. Bredden på pulsen på tidsaksen representerer dempningen som er 0,04 sek (fiolett), 0,08 sek (blå) og 0,12 sek (rød) I høyre graf b) ser vi de respektive spiking-filtre.(pulsbredden er oppgitt i antall sampler der samplingsintervallet=0,004 sek)

Fig.8.19 c.) Konvolvering av dempet enhetspuls med spikingfilter or inputpuls som er 0,04 sek (fiolett), 0,08 sek (blå) og 0,12 sek (rød) ) Resultatet er svært bra for alle pulser. Vi får litt støy i økende ankomst for hver puls

De ferdig spikede pulsene kan tyde på at Wienerfilteret er svært nøyaktig når vi har en ensidet dempet enhetspuls i minimum fase. Dempningen er fjernet og vi har fått en skarpere puls enn inngangen(input).

8.6.Inverse filtre for dempede enhetspulser laget med Wienerfilter

Vi har gjort noen beregninger av inverse filtre for dempede enhetspulser og hentet dem ut fra fig.8.10. På fig.8.20 a ser vi de aktuelle dempede pulsene for ulike ankomster med dempning a=0.023. Og på fig.8.20 b ser vi deres inverse filtre. Fig.8.20.c viser resultatet etter at den dempede pulsen er konvolvert med dets respektive spikingfilter.

8.7. Q-filtrering og Wiener-filter

Når vi skal anvende Wienerfilteret vil vi vanskelig komme utenom en form for forsterkning i tillegg til spikingen fra Wienerfilteret. AGC-forsterkning definert ved ligning 8.4 blir for unøyaktig - og unødvendig – for siden vi kjenner dempning fra indre friksjon kan vi anvende invers Q-filtrering på en enkel form. Vi vet at dempning kan defineres over et område med formelen fra 7.3.4.b:

8.31

Fig.8.21.a. Grafen viser utrykket 8.31 som en funksjon av t (rød strek). Den lineære grafen (blå) er dens inverse. Fig.8.21.b viser Wienerfilteret vi har definert for område 3 fra fig. 8.20 c.

Fig.8.21.c. viser en pulsrekke som er dempet med b=0.023 (som gir Q=136) Pulsrekken er identisk med 1. iterasjon i inversjonen av Riccatiligningen. Blå graf viser dempet pulsrekke. Rød graf er den dempede pulsrekken spiket med Wienerfilteret

Den enkleste form for invers Q-filtrering vil være å multiplisere seismogrammet med den inverse til dette utrykket. Et slikt inverst filter vil forsterke amplityden på reflektorene som en funksjon av tiden. Spiking-effekten vil Wiener-filteret ta seg av. Når vi kombinerer begge filtre samtidig kan vi tilpasse dem til hverandre og vi slipper å bruke så mange vinduer (områder) på seismogrammet. Vi vil for eksempel kunne få et resultat ved kun å bruke et vindu for spiking (her område 3) og la det inverse Q-filteret anvendes over resten av tidsområdet ved å skalere det riktig..

Fig.8.21.a. Q-dempning med inverst Q-filter b. Wiener-filter (puls 3 fra fig.8.20.b)

Fig.8.21.c. Vi har en pulsrekke dempet med b=0.023 (Q=136) som vi har kalt Modell (blå strek). På denne har vi anvendt Wienerfilteret optimalisert for område 3 fra fig.8.10. (spiket – Rød strek). Vi ser at vi får en skarpere puls som viser at Wienerfilteret har utført en spiking. Så har vi anvendt det inverse Q-filteret på pulsrekken som forsterker opp all energi. (brun strek).Vi plukker ut de tidlige ankomstene på grafen (sirklet inn) og har plottet disse på fig.8.22. Da ser vi bedre effekten både av invers Q-filtrering og spiking. Området tilsvarer område 1 og 2 fra fig.8.10. Vi burde ha forventet best spiking-resultat i område 3 (0,2-0,3 sek.), men i vårt tilfelle har vi fått best resultat i område 1 og 2 og har derfor ikke tatt med område 3 på den detaljerte grafen fig. 8.22. (0,12-0,24 sek)

Fig.8.22.Brun graf viser spikede og Q-filtrerte pulser mellom t=0,12 og 0,24 sek. Blå graf viser de dempede enhetspulsene, rød graf viser de dempede pulser spiket med Wienerfilteret og brun graf viser at vi har fått en god forsterkning av denne delen av seismogrammet. Det har blitt introdusert ringing og denne er forsterket kraftig opp av det inverse Q-filteret. Ringing kan enkelt fjernes ved glatting, eller ved å bruke en lagdeling som fjerner ringingen.

Nå har vi utviklet et filter som kan anvendes for Wiener og Q-filtrering for en dempet pulsrekke dersom vi ønsker å fjerne dempingen slik at vi kan få den originale pulsrekken tilbake.

I neste kapittel vil vi utvide den anvendelsen vi har gjort her i kapittel 8 til å gjelde i seismikk. Da vil transmisjonstap og multippelenergi komme inn som effekter i tillegg til de viskoelastiske som vi har studert til nå.