地震記錄數(shù)值模擬的褶積模型法

1、本科生實(shí)驗(yàn)報(bào)告實(shí)驗(yàn)課程 數(shù)值模型模擬 學(xué)院名稱 地球物理學(xué)院 專業(yè)名稱 勘測(cè)技術(shù)與工程 學(xué)生姓名 學(xué)生學(xué)號(hào) 指導(dǎo)教師 熊高君 實(shí)驗(yàn)地點(diǎn) 5417 實(shí)驗(yàn)成績(jī) 2015年5月 成都理工大學(xué)地震數(shù)值模擬實(shí)驗(yàn)報(bào)告實(shí)驗(yàn)時(shí)間2015年 5月開課單位地球物理學(xué)院指導(dǎo)教師熊高君實(shí)驗(yàn)題目：地震記錄數(shù)值模擬的褶積模型法姓名學(xué)號(hào)班級(jí)專業(yè)勘測(cè)技術(shù)與工程（石油物探）院（系）地球物理學(xué)院 地球探測(cè)與信息技術(shù)系單項(xiàng)成績(jī)內(nèi)容理解寫作結(jié)構(gòu)程序設(shè)計(jì)模型設(shè)計(jì)計(jì)算結(jié)果結(jié)果分析總成績(jī)實(shí)驗(yàn)報(bào)告一、 實(shí)驗(yàn)題目：地震記錄數(shù)值模擬的褶積模型法二、 實(shí)驗(yàn)?zāi)康模赫莆振薹e模型基本理論、實(shí)現(xiàn)方法與程序編制，由褶積模型初步分析地震信號(hào)的分辨率問(wèn)題。三、

2、原理公式1、褶積原理地震勘探的震源往往是帶寬很寬的脈沖，在地下傳播、反射、繞射到測(cè)線，傳播經(jīng)過(guò)中 高頻衰減，能量被吸收。吸收過(guò)程可以看成濾波的過(guò)程，濾波可以用褶積完成。在濾波中， 反射系數(shù)與震源強(qiáng)弱關(guān)聯(lián)，吸收作用與子波關(guān)聯(lián)。最簡(jiǎn)單的地震記錄數(shù)值模擬，可以看成反射系數(shù)與子波的褶積。通常，反射系數(shù)是脈沖，子波取雷克子波。（1） 雷克子波：wave(t)=cos(2ft)*e-22f2t（2） 反射系數(shù)：（3） 褶積公式：數(shù)值模擬地震記錄trace(t): trace(t) =rflct(t)*wave(t)；反射系數(shù)的參數(shù)由 z 變成了 t，怎么實(shí)現(xiàn)？在簡(jiǎn)單水平層介質(zhì)，分垂直和非垂直入射兩種實(shí)現(xiàn)，

3、分別如圖 1 和圖 2 所示。圖1 圖21） 垂直入射：t=2h/v；2）非垂直入射：t=2h2+x2v2、褶積方法（1） 離散化（數(shù)值化）計(jì)算機(jī)數(shù)值模擬要求首先必須針對(duì)連續(xù)信號(hào)離散化處理。反射系數(shù)在空間模型中存在，不同深度反射系數(shù)不同，是深度的函數(shù)。子波是在時(shí)間記錄上一延續(xù)定時(shí)間的信號(hào)，是時(shí)間的概念。在離散化時(shí)，通過(guò)深度采樣完成反射系數(shù)的離散化，通過(guò)時(shí)間采樣完成子波的離散化。如果記錄是 Trace(t)，則記錄是時(shí)間的函數(shù)，以時(shí)間采樣離散化。時(shí)間采樣間距以t 表示，深度采樣間距以z 表示。在做多道的數(shù)值模擬時(shí)，還有橫向x 的概念，橫向采樣間隔以x 表示。離散化的實(shí)現(xiàn)：t=It×t；

4、x=Ix×x；z=Iz×z；或：It=t/t; Ix=x/x; Iz=z/z（2） 離散序列的褶積trace(It)= Itao-+rflct Itao* wave(It -Itao)四、實(shí)驗(yàn)內(nèi)容1、 垂直入射地震記錄數(shù)值模擬的褶積模型；2、 非垂直入射地震記錄數(shù)值模擬的褶積模型。五、 方法路線1、 根據(jù)垂直入射褶積模型理論算法，填充程序（附后）的下劃線部分，使程序完整， 調(diào)試程序，算出結(jié)果，用“Fimage”顯示軟件顯示褶積結(jié)果；2、 根據(jù)非零偏移距算法，編制非零偏移距褶積模型程序，算出結(jié)果，用“Fimage”顯示軟件顯示褶積結(jié)果。（參考垂直入射褶積模型理論算法和程序，子

5、波與反射層不變）；3、 變換子波的主頻：fm（10hz到300hz范圍），重復(fù)1和2；4、 變換子波的長(zhǎng)度：Nw（80ms到160ms范圍），重復(fù)1和2；5、 改變反射層深度：h（800m到1600m范圍），重復(fù)1和2；6、 改變介質(zhì)速度：v（2000m/s到7000m/s 范圍），重復(fù)1和2。六、實(shí)驗(yàn)結(jié)果1、結(jié)果顯示1）垂直入射圖31 Nw=32,h=1000,v=3000,fm=100地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖32 Nw=32,h=1000,v=3000,fm=200地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖33 Nw=32,h=1000,v=3000,fm=300

6、地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖34 h=1000,v=3000,fm=25，Nw=20地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖35 h=1000,v=3000,fm=25，Nw=30地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖36 h=1000,v=3000,fm=25，Nw=40地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖37 Nw=32,v=3000,fm=25，h=1000地震記錄數(shù)值模擬的褶積模型圖38 Nw=32,v=3000,fm=25，h=1200地震記錄數(shù)值模擬的褶積模型圖39 Nw=32,v=3000,fm=25，h=1400地震記錄數(shù)值模擬的褶積模

7、型圖310 Nw=32, h=1000,fm=25，v=2000地震記錄數(shù)值模擬的褶積模型圖311 Nw=32, h=1000,fm=25，v=4000地震記錄數(shù)值模擬的褶積模型圖312 Nw=32, h=1000,fm=25，v=6000地震記錄數(shù)值模擬的褶積模型2）非垂直入射圖41 Nw=32,h=1000,v=3000,fm=100地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖42 Nw=32,h=1000,v=3000,fm=200地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖43 Nw=32,h=1000,v=3000,fm=300地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖44

8、h=1000,v=3000,fm=25，Nw=20地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖45 h=1000,v=3000,fm=25，Nw=30地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖46 h=1000,v=3000,fm=25，Nw=40地震記錄數(shù)值模擬的褶積模型（左）和子波（右）圖47 Nw=32,v=3000,fm=25，h=1000地震記錄數(shù)值模擬的褶積模型圖48 Nw=32,v=3000,fm=25，h=1200地震記錄數(shù)值模擬的褶積模型圖49 Nw=32,v=3000,fm=25，h=1400地震記錄數(shù)值模擬的褶積模型圖410 Nw=32, h=1000,fm=25，

9、v=2000地震記錄數(shù)值模擬的褶積模型圖411 Nw=32, h=1000,fm=25，v=4000地震記錄數(shù)值模擬的褶積模型圖412 Nw=32, h=1000,fm=25，v=6000地震記錄數(shù)值模擬的褶積模型子波振幅譜：圖51 fm=100的子波振幅譜圖52 fm=100的子波振幅譜2、對(duì)比分析a) 由圖31、圖32、圖33（或圖41、圖42、圖43）可知，當(dāng)子波長(zhǎng)度（Nw=32）、深度（h=1000）、速度(v=3000)不變，子波頻率變化時(shí)，褶積模型不變，且均在0.5到1s之間；b) 由圖34、圖35、圖36可知，當(dāng)深度（h=1000）、速度（v=3000）、子波頻率（fm=25）不

10、變, 子波長(zhǎng)度變化，且垂直入射時(shí)，褶積模型為直線模型，是因?yàn)榇怪比肷鋾r(shí)，時(shí)間與深度為線性關(guān)系；c) 由圖44、圖45、圖46可知，當(dāng)深度（h=1000）、速度（v=3000）、子波頻率（fm=25）不變, 子波長(zhǎng)度變化，且非垂直入射時(shí)，褶積模型前半部分為雙曲線模型，后半部分為直線模型，是由于計(jì)算的褶積結(jié)果的實(shí)際長(zhǎng)度小于所取的長(zhǎng)度，計(jì)算機(jī)賦的隨機(jī)數(shù)所致，所以，非垂直入射時(shí)的褶積模型只有圖件上顯示的前半部分的雙曲線，是由于非垂直入射時(shí)，時(shí)間與深度為雙曲線關(guān)系；d) 由圖37、圖38、圖39，當(dāng)子波長(zhǎng)度（Nw=32）、速度（v=3000）、子波頻率（fm=25）不變,深度變大時(shí)，垂直入射時(shí)，圖件上顯

11、示的褶積模型的位置逐漸向下移，即地震波的旅行路程變大，旅行時(shí)變長(zhǎng)；e) 由圖47、圖48、圖49可知，當(dāng)子波長(zhǎng)度（Nw=32）、速度（v=3000）、子波頻率（fm=25）不變,深度變大時(shí)，非垂直入射時(shí)，圖件上顯示的褶積模型的位置不變，這是由于非垂直入射時(shí)，相當(dāng)于地下有一半圓形界面，在圓心處自激自收；f) 由圖310、圖311、圖312（或圖410、圖411、圖412）可知，當(dāng)子波長(zhǎng)度（Nw=32）、深度（h=1000）、子波頻率（fm=25）不變,速度變大時(shí)，圖件上顯示的褶積模型的位置逐漸向上移，這是由于隨著速度變大，地震波的旅行時(shí)變??；g) 由圖51與圖52可知，當(dāng)子波頻率變化時(shí)，其振幅譜

12、不變。七、討論建議1、實(shí)驗(yàn)收獲通過(guò)此次試驗(yàn)，初步掌握了褶積模型基本理論、實(shí)現(xiàn)方法與程序編制，由褶積模型初步分析地震信號(hào)的分辨率問(wèn)題，掌握了褶積模型與子波主頻、子波長(zhǎng)度、界面深度、介質(zhì)中地震波速度的關(guān)系。2、存在問(wèn)題對(duì)褶積模型的理論實(shí)現(xiàn)過(guò)程不是十分清楚，對(duì)結(jié)果的物理含義理解不夠深入。3、其他問(wèn)題由于不同的計(jì)算機(jī)，計(jì)算精度不一樣，計(jì)算得到的數(shù)據(jù)結(jié)果可能會(huì)有部分差異，但總體趨勢(shì)不變。4、心得體會(huì)在此次試驗(yàn)中，應(yīng)特別注意褶積模型的理論實(shí)現(xiàn)過(guò)程的理解，以及程序調(diào)試時(shí)，要特別仔細(xì)地去檢查每一個(gè)錯(cuò)誤，每修改一處錯(cuò)誤，就重新運(yùn)行一次程序。附程序代碼：/=1. 預(yù)處理部分=/#include<math.h

13、>#include<stdio.h>#include<string.h>float Cnltn(float,float);float Rflct(float,float,float);float Wave(float,float);#define Nx 128#define Nt 256#define Nw 32#define PI 3.1415926/=2. 主程序波分=/void main()float dt=0.004,dx=20,fm=25,h=1000,v=3000; int iflag_Co,iflag_Re,iflag_Wv; if(iflag_Wv

14、=Wave(fm,dt)!=1)printf("Wave is error"); if(iflag_Re=Rflct(dt,h,v)!=1)printf("Reflection is error"); if(iflag_Co=Cnltn(dt,dx)!=1)printf("Convosion is error");/ =3.函數(shù)實(shí)現(xiàn)部分=/ =3.1 Wave Formaing function=/float Wave(float fm,float dt)FILE *fpw; int It; float WaNw,t; if(fpw=

15、fopen("wave.dat","wb")=NULL)printf("Connot open file ""wave"""); for(It=0;It<Nw;It+)t=It*dt;WaIt=cos(2*PI*fm*t)*exp(-2*PI*PI*fm*fm*t);/形成子波fwrite(&WaIt,sizeof(WaIt),1,fpw); fclose(fpw); return(1);/ =3.2 Reflect Formaing function=/float Rflct(f

16、loat dt,float h,float v)FILE *fpr; int It,Ix,J,Ltdpth; float t,dx=20,x; float ReNt;printf("請(qǐng)輸入J：n");scanf("%d",&J); if(fpr=fopen("Reflect.dat","wb")=NULL)printf("Connot open file ""Reflect"""); for(Ix=0;Ix<Nx;Ix+)for(It=0;I

17、t<Nt;It+)ReIt=0.;if(J=1)t=2*h/v;/垂直入射反射界面由深度轉(zhuǎn)換為自激自收時(shí)間if(J=2)x=Ix*dx;t=2*sqrt(h*h+x*x)/v;/非垂直入射反射界面由深度轉(zhuǎn)換為自激自收時(shí)間if(J!=1)&&(J!=2)printf("輸入錯(cuò)誤n"); Ltdpth=(int)(t/dt); ReLtdpth=1; for(It=0;It<Nt;It+)fwrite(&ReIt,sizeof(ReIt),1,fpr); fclose(fpr); return(1);/ =3.3 Convolution fu

18、nction=/float Cnltn(float dt,float dx)FILE *fpc,*fpw,*fpr;int It,Ix,Itao; float Wa1Nw,WaNw,ReNt+Nw+Nw,Re1Nt,t; float ConNt+Nw; if(fpc=fopen("Convosion.dat","wb")=NULL)printf("Connot open file ""Convosion"""); if(fpw =fopen("wave.dat","rb")=NULL)printf("Connot open file ""wave"""); if(fpr =fopen("Reflect.dat"