c语言地震子波6,地震子波显示及合成地震记录
一、實驗目的
通過編制和運行相關C或matlab程序,理解地震子波的概念和特點,掌握地震合成記錄的制作流程。認識子波,對子波的波形有直觀的認識。(零相位子波,混合相位子波,最小相位子波;了解子波的分辨率與頻寬的關系;)
利用褶積公式合成一維地震記錄
二、主要內容和原理
1、理論子波的計算及顯示
(1)零相位子波、最小相位子波、鐘型子波的計算
零相位子波的表達式
最小相位子波的表達式為
;
n=m1/m2,為最大波峰m1和最大波谷m2之比,取2.0或1.5,本次試驗采用2.0.;
代表子波的中心頻率,t=i*dt,dt為時間采樣間隔,i為時間離散點序號;
本次試驗通過零相位子波和最小相位子波的公式,利用0.002s的采樣間隔,對零相位子波和最小相位子波進行離散采樣,得到一組數據。利用excel繪制圖像,得到中心頻率分別為10,25,40,100hz的8張圖像。
2、利用聲波測井資料計算層速度,密度和反射系數序列
本實驗通過C語言編程序處理獲得的聲波測井資料,把資料中的數據用于層速的計算,公式為
和密度計算,當速度與密度呈線性關系時,計算公式為:
本次試驗假設速度和密度為線性關系,采用這個公式參與計算。
3、 制作合成地震記錄
地震合成記錄是在一定的假設條件下進行的,如果地下實際介質有n層面都是反射界面的話,地面可以接收到每一層界面的反射波,于是一個實際地震記錄道上會記錄n個反射波。忽略介質吸收、透射損失影響,使用地震子波和反射系數離散卷積的方法得到。因為離散波形卷積需要采樣間隔相同,所以使用0.002s的間隔對第i個反射面到地表的垂直往返時間進行離散采樣.計算公式為:
得到反射系數的序號,然后求取每個序號求取反射系數,當然并不是所有的序號都有反射系數,對沒有反射系數的用0取代。得到反射系數的離散數據,在得到35hz和70hz的零相位子波的離散數據,根據離散卷積公式,
由此可得一系列隨時間變化的卷積結果即為對應時刻的波形振幅值,繪制散點圖即可觀察到所記錄的理論合成波記錄。
三、實驗結果分析(含圖)
1、理論子波的計算及顯示
(1)、零相位子波的波形
中心頻率為10HZ
中心頻率為25HZ
中心頻率為40HZ
中心頻率為100HZ
結果分析:根據零相位子波的計算結果和所繪制的圖像,我們可以發現當n一定時(n=m1/m2,為最大波峰m1和最大波谷m2之比),零相位子波圖像隨著子波的中心頻率fm的增加,而呈現有規律的變化。對于零相位子波來說,子波的中心頻率越大,最大幅度值越小。
(2)、最小相位子波的波形
中心頻率為10HZ
中心頻率為25HZ
中心頻率為40HZ
中心頻率為100HZ
結果分析:根據最小相位子波的計算結果和所繪制的圖像,我們可以發現當n一定時(n=m1/m2,為最大波峰m1和最大波谷m2之比),對于最小相位子波來說,隨著子波中心頻率的增加,振幅變化的速度加快了,但是最大幅值卻沒有明顯的變化。
2、利用聲波測井資料計算層速度,密度和反射系數序列
(1)、層速曲線
結果分析:通過編程計算得到層速度結果并繪制圖件,分析結果和圖件可以發現,層速度和深度并沒有明顯的線性關系,層速度應該與每層介質的密度有關。
(2)、密度曲線
結果分析:密度和深度并沒有明顯的線性關系,但是所有地層除了個別層密度較大外,基本變化不大。密度應該和層速度有密切的聯系。
(3)、反射系數
結果分析:反射系數是根據離散取樣而來,所以許多標號下,反射系數為零,因為相鄰兩層的的介質性質不同,所以,反射系數也不相同,有大有小,與相鄰底層的波阻抗性質有關。
3、 制作合成地震記錄
(1)、中心頻率為35HZ
結果分析:
(2)、中心頻率為70HZ
結果分析:觀察分析兩幅地震合成記錄圖件地震合成記錄70Hz相比35Hz其最大幅值要大并且相同水平距離70Hz的波數大于35Hz,原因在于同一介質頻率越高波長越短使得相同水平距離70Hz的波數大于35Hz,同樣頻率越高在相同時間間隔內疊加次數越多使得卷結果絕對值增大。
四、源程序代碼
1、理論子波的計算及顯示
#include
#include
#include
#define dz1 "file1.dat"
#define dz2 "file2.dat"
main ()
{ FILE *fp1,*fp2;//定義文件指針
double dt=0.002,t;int fm1=10,fm2=25,fm3=40,fm4=100,i,j;double r1,r2,r3,r4,w1,w2,w3,w4;
fp1=fopen(dz1,"w+");//文件新建
fp2=fopen(dz2,"w+");
fclose(fp1);
fclose(fp2);//關閉文件
fp1=fopen(dz1,"a+");
fp2=fopen(dz2,"a+");
for(i=1;i<=200;i++)
{ t=i*dt;
r1=exp(-pow(3.141592653*fm1*t,2))*(1-2*pow(3.141592653*fm1*t,2));//
r2=exp(-pow(3.141592653*fm2*t,2))*(1-2*pow(3.141592653*fm2*t,2));
r3=exp(-pow(3.141592653*fm3*t,2))*(1-2*pow(3.141592653*fm3*t,2));
r4=exp(-pow(3.141592653*fm4*t,2))*(1-2*pow(3.141592653*fm4*t,2));
w1=exp(-2*fm1*fm1*t*t*log(2))*sin(2*3.14*fm1*t);
w2=exp(-2*fm2*fm2*t*t*log(2))*sin(2*3.14*fm2*t);
w3=exp(-2*fm3*fm3*t*t*log(2))*sin(2*3.14*fm3*t);
w4=exp(-2*fm4*fm4*t*t*log(2))*sin(2*3.14*fm4*t);//計算子波
printf("%d %f %f %f\n",i,t,r1,w1);
fprintf(fp1,"%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\n",t,r1,r2,r3,r4,w1,w2,w3,w4);
}//輸出子波
fclose(fp1);
fclose(fp2);
}
2、利用聲波測井資料計算層速度,密度和反射系數序列
3、制作合成地震記錄
#include//速度密度,反射系數和制作地震記錄程序在一塊
#include//程序中有標明
#include
main(){
FILE *fp1,*fp2,*fp3;double vi[200],midu1[200],midu2[200],ti[200],rn[2000];int n[200];int s1[200];int s2[200];float s3[200];
double dt=0.002,t;int fm1=35,fm2=70,i,j;double r1[200],w1[200],r2[200],w2[200],xn1[2000],xn2[2000];
for(j=0;j<200;j++)
{ ???? t=(j+1)*dt;
r1[j]=exp(-pow(3.141592653*fm1*t,2))*(1-2*pow(3.141592653*fm1*t,2));
r2[j]=exp(-pow(3.141592653*fm2*t,2))*(1-2*pow(3.141592653*fm2*t,2));
w1[j]=exp(-2*fm2*fm2*t*t*log(2))*sin(2*3.141592653*fm1*t);
w2[j]=exp(-2*fm2*fm2*t*t*log(2))*sin(2*3.141592653*fm2*t); }//計算子波
if((fp1=fopen("file1","rb+"))==NULL){//第二小問層速,密度程序開始
printf("Cannot open file strike any key exit!"); }
if((fp2=fopen("file2","r+"))==NULL){
printf("Cannot open file strike any key exit!");
}
if((fp3=fopen("file3","r+"))==NULL){
printf("Cannot open file strike any key exit!"); }
for(i=0;i<200;i++)
fscanf(fp1,"%d %d %f \n",&s1[i],&s2[i],&s3[i]);/讀取文件中數據
for(i=0;i<200;i++)
{????vi[i]=1000000.0/s2[i];//層速計算
midu1[i]=1.62+0.00021*vi[i];
midu2[i]=0.31*pow(vi[i],0.25); }//密度計算,后邊結果采用線性密度公式
ti[0]=2*s3[0]/vi[0];//初值
for(i=1;i<200;i++)
{ ti[i]=ti[i-1]+2*(s3[i]-s3[i-1])/vi[i]; } //時間
for(i=0;i<200;i++)
{ n[i]=int(ti[i]/0.002); }
for(i=0;i<2000;i++)
{ ???????? rn[i]=0.0;???????? ???????? }
for(i=0;i<200;i++)
{rn[(n[i])]=(midu1[i+1]*vi[i+1]-midu1[i]*vi[i])/(midu1[i+1]*vi[i+1]+midu1[i]*vi[i] ;}//反射系數計算第二小問層速,密度,反射系數程序結束
for(i=0;i<2000;i++)//70hz子波卷積,地震記錄合成
{????xn1[i]=0.0;
for(j=0;j<200;j++)
{ if((i-j)>0)
xn1[i]=xn1[i]+r1[j]*rn[i-j] ;???? }
}
for(i=0;i<2000;i++)//70hz子波卷積,地震記錄合成
{????xn2[i]=0.0;
for(j=0;j<200;j++)
{ if((i-j)>0)
xn2[i]=xn2[i]+r2[j]*rn[i-j] ;}
}
for(i=0;i<200;i++)
{
fprintf(fp2,"%d\t%f\t%f\t%f\t%f\t%d\n",i+1,vi[i],midu1[i],midu2[i],ti[i],n[i]);
}//輸出層速,密度到文件
for(i=0;i<2000;i++)// 輸出反射系數,卷積到文件
{ fprintf(fp3,"%d\t%f\t%f\t%f\n",i+1,rn[i],xn1[i],xn2[i]); }
printf("處理完畢請打開文件尋找數據\n");
fclose(fp1);
fclose(fp2);
fclose(fp3); }
學生實驗 心得通過做這個實驗,我認識到了零相位子波和最小相位子波的波形,
認識了最大相位子波,強化了編程技能。學到了關于C語言關于文件
輸出的方式方法。練習了根據測井資料計算速度密度和反射系數,而且畫出曲線,并用卷積方法制作了地震合成記錄,掌握了分析的方法。
學習并掌握了了地震測井資料的分析處理和解釋方法。對勘探技能
的掌握也有了提高。
學生(簽名): 2015 年 10 月 26 日
指導
教師
評語成績評定:
指導教師(簽名):
年 月 日
總結
以上是生活随笔為你收集整理的c语言地震子波6,地震子波显示及合成地震记录的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 重磅!《中国迈向新一代人工智能》全文来了
- 下一篇: Ear Clipping算法简介