合成记录制作示意程序_油香飘飘

http://blog.sina.com.cn/u/6985345888

首页博文目录关于我

个人资料

微博

加好友发纸条

写留言加关注

博客等级：
博客积分：

博客访问：
关注人气：
获赠金笔：0支
赠出金笔：0支
荣誉徽章：

正文字体大小：大中小

合成记录制作示意程序

(2022-11-11 16:28:21)

标签：

上机作业（2）

分类：《信号分析与处理》课

#include "stdio.h" // standard input/output head file.

#include "math.h" // math. head file.

#include "malloc.h" // memory allocation head file.

#define fm 35 // The main frequency of Ricker wavelet, you can change it.

#define dt 0.002 // The sample interval, you can also set it to 0.001 or 0.004.

#define XL 0.060 // The length of wavelet(s), you can change it.

#define HL 0.300 // The length of h(t), namely reflection coefficient.

// Here the length is 0.3 second, you can change it.

#define PI 3.1415926

main()

{

FILE *fp1,*fp2,*fp3;

float *x,*h,*y; // define the dynamic array.

float t;

int i,M,N;

void conv();

M=XL/dt+1; // The number of x(t). For this example, M=31.

N=HL/dt+1; // The number of h(t). For this example, N=151.

// allocate dynamic array.

x=(float *)calloc(M,sizeof(float));

h=(float *)calloc(N,sizeof(float));

y=(float *)calloc(M+N-1,sizeof(float));

// Generate Ricker wavelet.

fp1=fopen("Ricker.xls","w");

for(i=-M/2;i<=M/2;i++)

{ t=i*dt;

x[i+M/2]=(1.0-2.0*pow(PI*fm*t,2.0))*exp(-pow(PI*fm*t,2.0));

fprintf(fp1,"%8.3f \t .4f\n",t*1000,x[i+M/2]);

}

fclose(fp1);

// Generate reflection coefficient.

fp2=fopen("Ksai.xls","w");

// set initial value to zero.

for(i=0;i<=N-1;i++)h[i]=0.0;

// Generate different ksai[i].

// You can design different reflection coefficient sequence.

h[20]=0.4;

h[80]=-0.3;

h[90]=-0.5;

h[130]=0.5;

for(i=0;i<=N-1;i++)

{ t=i*dt;

fprintf(fp2,"%8.3f \t .4f\n",t*1000,h[i]);

}

fclose(fp2);

// Do convolution x(n)*h(n).

conv(x,M,h,N,y,M+N-1);

// Output the results. Let the length of y(n) be equal to h(t), cut y(n).

fp3=fopen("Seis.xls","w");

for(i=M/2;i<=N+M/2-1;i++)

{

t=(i-M/2)*dt;

fprintf(fp3,"%8.3f \t .4f\n",t*1000,y[i]);

}

fclose(fp3);

}

// The function of convolution.

void conv(float x[],int m,float h[],int n,float y[],int l)

{ int k,i;

for(k=0;k

{ y[k]=0.0;

for(i=0;i

if(k-i>=0&&k-i<=n-1)y[k]=y[k]+x[i]*h[k-i];

}

阅读┊ 收藏 ┊ 喜欢 ▼ ┊打印┊举报/Report

前一篇：一般卷积实现程序

后一篇：王岩、李志远、郭文静在课堂完成演算、上机

新浪BLOG意见反馈留言板　欢迎批评指正