17 |
|
* each of NI2DIR surrounding directions. To speed this |
18 |
|
* calculation, we sort the data into half-planes and apply |
19 |
|
* simple tests to see which neighbor is closest in each |
20 |
< |
* direction. Once we have our approximate neighborhood |
20 |
> |
* angular slice. Once we have our approximate neighborhood |
21 |
|
* for a sample, we can use it in a modified Gaussian weighting |
22 |
|
* with allowing local anisotropy. Harmonic weighting is added |
23 |
|
* to reduce the influence of distant neighbors. This yields a |
90 |
|
{ |
91 |
|
if (mind <= 0) |
92 |
|
return; |
93 |
< |
if ((.998*ip->dmin <= mind) && (mind <= 1.002*ip->dmin)) |
93 |
> |
if ((.998*ip->dmin <= mind) & (mind <= 1.002*ip->dmin)) |
94 |
|
return; |
95 |
|
if (ip->da != NULL) { /* will need to recompute distribution */ |
96 |
|
free(ip->da); |
238 |
|
static double |
239 |
|
get_wt(const INTERP2 *ip, const int i, double x, double y) |
240 |
|
{ |
241 |
< |
double dir, rd, d2; |
241 |
> |
double dir, rd, r2, d2; |
242 |
|
int ri; |
243 |
|
/* get relative direction */ |
244 |
|
x -= ip->spt[i][0]; |
251 |
|
rd -= (double)ri; |
252 |
|
rd = (1.-rd)*ip->da[i][ri] + rd*ip->da[i][(ri+1)%NI2DIR]; |
253 |
|
rd = ip->smf * DECODE_DIA(ip, rd); |
254 |
+ |
r2 = 2.*rd*rd; |
255 |
|
d2 = x*x + y*y; |
256 |
+ |
if (d2 > 21.*r2) /* result would be < 1e-9 */ |
257 |
+ |
return(.0); |
258 |
|
/* Gaussian times harmonic weighting */ |
259 |
< |
return( exp(d2/(-2.*rd*rd)) * ip->dmin/(ip->dmin + sqrt(d2)) ); |
259 |
> |
return( exp(-d2/r2) * ip->dmin/(ip->dmin + sqrt(d2)) ); |
260 |
|
} |
261 |
|
|
262 |
|
/* Assign full set of normalized weights to interpolate the given position */ |