| 126 |
|
return(sgn2 == sgn3); |
| 127 |
|
} |
| 128 |
|
|
| 129 |
< |
/* Test and set for edge */ |
| 129 |
> |
/* Test (and set) bitmap for edge */ |
| 130 |
|
static int |
| 131 |
|
check_edge(unsigned char *emap, int nedges, const MIGRATION *mig, int mark) |
| 132 |
|
{ |
| 305 |
|
rbf->rbfa[n].crad = ANG2R(sqrt(rad0*rad0*(1.-t) + |
| 306 |
|
rad1*rad1*t)); |
| 307 |
|
ovec_from_pos(v, rbf1j->gx, rbf1j->gy); |
| 308 |
< |
geodesic(v, v0, v, t, GEOD_REL); |
| 308 |
> |
geodesic(v, v0, v, t*full_dist, GEOD_RAD); |
| 309 |
|
pos_from_vec(pos, v); |
| 310 |
|
rbf->rbfa[n].gx = pos[0]; |
| 311 |
|
rbf->rbfa[n].gy = pos[1]; |
| 331 |
|
float mbfact, mcfact; |
| 332 |
|
int n, i, j, k; |
| 333 |
|
FVECT v0, v1, v2; |
| 334 |
< |
double s, t; |
| 334 |
> |
double s, t, s_full, t_full; |
| 335 |
|
|
| 336 |
|
VCOPY(sivec, invec); /* find triangle/edge */ |
| 337 |
|
sym = get_interp(miga, sivec); |
| 360 |
|
normalize(v2); |
| 361 |
|
fcross(v1, sivec, miga[1]->rbfv[1]->invec); |
| 362 |
|
normalize(v1); |
| 363 |
< |
s = acos(DOT(v0,v1)) / acos(DOT(v0,v2)); |
| 363 |
> |
s = acos(DOT(v0,v1)); |
| 364 |
|
geodesic(v1, miga[0]->rbfv[0]->invec, miga[0]->rbfv[1]->invec, |
| 365 |
< |
s, GEOD_REL); |
| 366 |
< |
t = acos(DOT(v1,sivec)) / acos(DOT(v1,miga[1]->rbfv[1]->invec)); |
| 365 |
> |
s, GEOD_RAD); |
| 366 |
> |
s /= s_full = acos(DOT(v0,v2)); |
| 367 |
> |
t = acos(DOT(v1,sivec)) / |
| 368 |
> |
(t_full = acos(DOT(v1,miga[1]->rbfv[1]->invec))); |
| 369 |
|
n = 0; /* count migrating particles */ |
| 370 |
|
for (i = 0; i < mtx_nrows(miga[0]); i++) |
| 371 |
|
for (j = 0; j < mtx_ncols(miga[0]); j++) |
| 406 |
|
rad1j = R2ANG(rbf1j->crad); |
| 407 |
|
srad2 = (1.-s)*(1.-t)*rad0i*rad0i + s*(1.-t)*rad1j*rad1j; |
| 408 |
|
ovec_from_pos(v1, rbf1j->gx, rbf1j->gy); |
| 409 |
< |
geodesic(v1, v0, v1, s, GEOD_REL); |
| 409 |
> |
geodesic(v1, v0, v1, s*s_full, GEOD_RAD); |
| 410 |
|
for (k = 0; k < mtx_ncols(miga[2]); k++) { |
| 411 |
|
float mb = mtx_coef(miga[1],j,k); |
| 412 |
|
float mc = mtx_coef(miga[2],i,k); |
| 421 |
|
rad2k = R2ANG(rbf2k->crad); |
| 422 |
|
rbf->rbfa[n].crad = ANG2R(sqrt(srad2 + t*rad2k*rad2k)); |
| 423 |
|
ovec_from_pos(v2, rbf2k->gx, rbf2k->gy); |
| 424 |
< |
geodesic(vout, v1, v2, t, GEOD_REL); |
| 424 |
> |
geodesic(vout, v1, v2, t*t_full, GEOD_RAD); |
| 425 |
|
pos_from_vec(pos, vout); |
| 426 |
|
rbf->rbfa[n].gx = pos[0]; |
| 427 |
|
rbf->rbfa[n].gy = pos[1]; |