forked from geodynamics/hc
-
Notifications
You must be signed in to change notification settings - Fork 0
/
grdinttester.c
64 lines (49 loc) · 1.6 KB
/
grdinttester.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
/*
read in two grid files and compute the correlation of
points within a certain radius around each location
*/
#include "ggrd_grdtrack_util.h"
int main(int argc, char **argv)
{
struct ggrd_gt g1[1],g2[1];
char char_dummy;
float x,y,x1,y1,dx,dy,rad_km;
int n;
double *z1,*z2;
z1=(double *)malloc(sizeof(double));
z2=(double *)malloc(sizeof(double));
rad_km = 500;
if(argc < 3){
fprintf(stderr,"%s: usage:\n\t%s file1.grd file2.grd [radius, %g km]\n\n",
argv[0],argv[0],rad_km);
exit(-1);
}
if(argc > 3 )
sscanf(argv[3],"%f",&rad_km);
fprintf(stderr,"%s: assuming grids %s and %s are global, computing for radius %g km\n",
argv[0],argv[1],argv[2],rad_km);
if(ggrd_grdtrack_init_general(FALSE,argv[1],&char_dummy,
"-Lx",g1,TRUE,FALSE,FALSE)){
fprintf(stderr,"%s: error reading %s\n",argv[0],argv[1]);
exit(-1);
}
if(ggrd_grdtrack_init_general(FALSE,argv[2],&char_dummy,
"-Lx",g2,TRUE,FALSE,FALSE)){
fprintf(stderr,"%s: error reading %s\n",argv[0],argv[2]);
exit(-1);
}
/*
compute total correlation and best-fit slopes
*/
dy = 2;
for(n=0,y1 = 90-dy/2+1e-7,y = -90+dy/2;y <= y1;y += dy){
dx = dy/cos(y/GGRD_PIF); /* adjust for sphericity */
for(x=0,x1=360-dx+1e-7;x<=x1;x+=dx,n++){
ggrd_grdtrack_interpolate_xy((double)x,(double)y,g1,(z1+n),FALSE);
ggrd_grdtrack_interpolate_xy((double)x,(double)y,g2,(z2+n),FALSE);
fprintf(stderr,"%g %g %g %g %.6e\n",x,y,z1[n],z2[n],fabs(z1[n]-z2[n]));
z1=(double *)realloc(z1,(n+2)*sizeof(double));
z2=(double *)realloc(z2,(n+2)*sizeof(double));
}
}
}