I am trying to use this c program to do the numerical modeling. It has been successfully compiled. However, when I am trying to input:
./forward 0.0 0.0 1.0 3.00 100.0
.
It keeps reminding me:
segmentation fault: 11
.
I do not know what goes wrong. My code is as following:
#include<stdio.h>
#include<stdlib.h>
#include<string.h>
#include<math.h>
#include"util_err.h"
#define pi 3.1415926
#define MAX_Array 16000
typedef struct source {
float vp;
float ev_radius;
} source;
source model_search(float evdp, float a_radius, char *file){
int i;
float *radius=NULL,*alpha_s=NULL,*beta_s=NULL;
float *z_s=NULL;
source sour;
FILE *f;
/*radius = vector(0,MAX_Array-1);
alpha_s = vector(0,MAX_Array-1);
beta_s = vector(0,MAX_Array-1);
z_s = vector(0,MAX_Array-1);*/
f = fopen(file,"r");
if (f==NULL)
{
nrerror("unable to open the file");
}
sour.vp = 0;
sour.ev_radius = 0;
for(i=0;!feof(f);i++){
fscanf(f,"%f %f %f",&radius[i],&alpha_s[i],&beta_s[i]);
z_s[i] = a_radius - radius[i];
if(i==0) continue;
if(evdp>=z_s[i-1]&&evdp<=z_s[i]){
sour.vp = alpha_s[i-1]+(evdp-z_s[i-1])/(z_s[i]-z_s[i-1])
*(alpha_s[i]-alpha_s[i-1]);
sour.ev_radius = a_radius-evdp;
printf("vp = %f\n",sour.vp);
printf("ev_radius = %f\n",sour.ev_radius);
break;
}
}
fclose(f);
return sour;
}
float* linespace(float fl, float fr, int n){
int i;
static float *u=NULL,interval;
interval = (fr-fl)/(n-1);
for(i=0;i<n-1;i++){
u[i] = fl + i * interval;
}
u[n-1] = fr;
return u;
}
/* Input: x, y, z, t0, evdp*/
int main(int argc,char **argv){
int i,j;
int sto,sto1,sa,sa1;
float x,y,z,origin_time;
float travel_time,evdp,a_radius=6371.1;
float *takeoff=NULL,*takeoff1=NULL,*azi=NULL,*azi1=NULL;
float t0,az;
char filename[160];
source sour;
FILE *fp=NULL,*fp1=NULL;
fp=fopen("data.dat","w");
if(fp==NULL){
nrerror("unable to allocate");
}
fp1=fopen("data1.dat","w");
if(fp1==NULL){
nrerror("unable to allocate");
}
if(argc<6){
nrerror("Too few arguments!");
}
x = atof(argv[1]);
y = atof(argv[2]);
z = atof(argv[3]);
origin_time = atof(argv[4]);
evdp = atof(argv[5]);
strcpy(filename,"/Users/a123/Desktop/Ear_relo/plot_seis/iasp91.dat");
sour = model_search(evdp,a_radius,filename);
printf("vp = %f, ev_radius = %f",sour.vp,sour.ev_radius);
sto = 91;
sto1 = 10;
sa = 361;
sa1 = 37;
takeoff = linespace(0.0,90.0,sto);
takeoff1 = linespace(0.0,90.0,sto1);
azi = linespace(0.0,360.0,sa);
azi1 = linespace(0.0,360.0,sa1);
/* write data.dat for the background values */
for(i=0;i<sto;i++){
t0 = takeoff[i];
for(j=0;j<sa;j++){
az = azi[j];
travel_time = -x*cos(az*pi/180)*sin(t0*pi/180)
-y*sin(az*pi/180)*sin(t0*pi/180)+z*cos(t0*pi/180);
travel_time = travel_time/sour.vp + origin_time;
fprintf(fp, "%5.1f %5.1f %13.8f\n", az,t0,travel_time);
}
}
/* write data1.dat */
for(i=0;i<sto1;i++){
t0 = takeoff1[i];
for(j=0;j<sa1;j++){
az = azi1[j];
travel_time = -x*cos(az*pi/180)*sin(t0*pi/180)
-y*sin(az*pi/180)*sin(t0*pi/180)+z*cos(t0*pi/180);
travel_time = travel_time/sour.vp + origin_time;
fprintf(fp1, "%5.1f %5.1f %13.8f\n", az,t0,travel_time);
}
}
free(takeoff);
free(takeoff1);
free(azi);
free(azi1);
fclose(fp);
fclose(fp1);
return 0;
}
This is my iasp91.dat file:
6371 5.8 3.36
6370 5.8 3.36
6369 5.8 3.36
6368 5.8 3.36
6367 5.8 3.36
6366 5.8 3.36
6365 5.8 3.36
6364 5.8 3.36
6363 5.8 3.36
6362 5.8 3.36
6361 5.8 3.36
6360 5.8 3.36
6359 5.8 3.36
6358 5.8 3.36
6357 5.8 3.36
6356 5.8 3.36
6355 5.8 3.36
6354 5.8 3.36
6353 5.8 3.36
6352 5.8 3.36
6351 5.8 3.36
6351 6.5 3.75
6350 6.5 3.75
6349 6.5 3.75
6348 6.5 3.75
6347 6.5 3.75
6346 6.5 3.75
6345 6.5 3.75
6344 6.5 3.75
6343 6.5 3.75
6342 6.5 3.75
6341 6.5 3.75
6340 6.5 3.75
6339 6.5 3.75
6338 6.5 3.75
6337 6.5 3.75
6336 6.5 3.75
6336 8.04 4.47
6331 8.0406 4.4718
6326 8.0412 4.4735
6321 8.0418 4.4753
6311 8.0429 4.4788
6301 8.0441 4.4824
6291 8.0453 4.4859
6281 8.0465 4.4894
6271 8.0476 4.4929
6261 8.0488 4.4965
6251 8.05 4.5
6251 8.05 4.5
6241 8.0778 4.502
6231 8.1056 4.504
6221 8.1333 4.506
And the "util_err.h" declares the function --- nrerror.
I do not know what the exact problem is but it seems to me that the problem is in the block linespace.