forked from geodynamics/hc
-
Notifications
You must be signed in to change notification settings - Fork 0
/
ggrd_test.c
72 lines (63 loc) · 1.65 KB
/
ggrd_test.c
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
#include "hc.h"
#include "hc_ggrd.h"
int main(int argc, char **argv)
{
struct ggrd_master *ggrd;
double time,vr[4],vphi[4],vtheta[4],dtrange;
double xloc[3];
static int order = 3;
hc_boolean calc_derivatives ;
double lon,lat,age;
int mode;
mode = 1;
switch(mode){
case 1:
if(argc>1)
ggrd_grdinfo(argv[1]);
break;
case 2:
/*
initialize velocity structure
*/
ggrd = (struct ggrd_master *)calloc(1,sizeof(struct ggrd_master));
ggrd_init_master(ggrd);
ggrd->v.history = TRUE; /* expect history */
ggrd->age_control = TRUE; /* expect seafloor age files */
ggrd->age_bandlim = 900;
/*
read in velocities
*/
if(ggrd_read_vel_grids(ggrd,1.0,FALSE,TRUE,
"/home/walter/becker/data/plates/past/clb/hall/",FALSE)){
fprintf(stderr,"error opening grids\n");
exit(-1);
}
if(argc>1)
sscanf(argv[1],"%lf",&time);
else
time = 0.0;
dtrange = 1.0; /* transition width, in Ma */
fprintf(stderr,"%s: using time %g\n",argv[0],time);
calc_derivatives = FALSE;
xloc[HC_R] = HC_ND_RADIUS(0.0);
for(lat=-89;lat<=89;lat+=2)
for(lon=0;lon<=358;lon+=2){
//lon=270;lat=-15;{
xloc[HC_THETA] = LAT2THETA(lat);
xloc[HC_PHI] = LON2PHI(lon);
/*
interpolate
*/
if(ggrd_find_vel_and_der(xloc,time,dtrange,ggrd,
order,calc_derivatives,
TRUE,vr,vtheta,vphi))
exit(-1);
if(interpolate_seafloor_ages(xloc[HC_THETA],
xloc[HC_PHI],time,ggrd, &age))
exit(-1);
fprintf(stdout,"%11g %11g\t%11g %11g %11g\t%11g\n",lon,lat,vphi[0],-vtheta[0],vr[0],age);
}
break; /* end mode 2 */
}
return 0;
}