CodeArchive CircleBoundaryTraversal
From Ubcacm
Circle Boundary Traversal
Don't mess up with comments.
#include <iostream> #include <cmath> #include <complex> using namespace std; typedef long double ld; typedef complex < ld > point; const ld eps = 1e-9; const ld PI = 2*acos((ld)0); int vis[128], done[128][128], ind[128][128], sz[128], cnt, n; point c[128]; ld out[128][128], in[128][128], r[128], x, y; void dfs(int now){ if(vis[now]) return; vis[now] = 1; for(int i = 0; i < n; i++) if(!vis[i] && abs(c[now] - c[i]) < r[now] + r[i]) dfs(i); } int main(){ int t; cin >> t; while(t-- && cin >> n){ for(int i = 0; i < n && cin >> x >> y >> r[i]; i++) c[i] = point(x, y); memset(vis, 0, sizeof(vis)), memset(sz, 0, sizeof(sz)), memset(done, 0, sizeof(done)); for(int i = 0; i < n; i++) for(int j = 0; j < n; j++){ ld d = abs(c[i] - c[j]); if(i != j && abs(r[i] - r[j]) < d && d < r[i] + r[j]){ point v = c[j] - c[i]; ld cosa = (r[i]*r[i] + d*d - r[j]*r[j])/(2*r[i]*d); ld sina = sqrt(1 - cosa*cosa); v *= point(cosa, -sina); // v is direction to intersection v /= abs(v), v *= r[i]; point inter = v + c[i]; // now inter is the intersection point; bool good = true; for(int k = 0; k < n; k++) good &= (abs(inter - c[k]) + eps > r[k]); if(good){ // only if not inside a circle out[i][sz[i]] = atan2(v.imag(), v.real()); if(out[i][sz[i]] < 0) out[i][sz[i]] += 2*PI; ind[i][sz[i]] = j; v = c[i] - c[j]; cosa = (r[j]*r[j] + d*d - r[i]*r[i])/(2*r[j]*d); sina = sqrt(1 - cosa*cosa); v *= point(cosa, sina); // the other turn in[i][sz[i]++] = atan2(v.imag(), v.real()); if(in[i][sz[i] - 1] < 0) in[i][sz[i] - 1] += 2*PI; } } } cnt = 1; for(int i = 0; i < n; i++){ if(sz[i] && !vis[i]) cnt--, dfs(i); for(int j = 0; j < sz[i]; j++) if(!done[i][j]){ cnt++; ld angle = in[i][j] - eps; int now = ind[i][j]; while(true){ ld best = 1e99; int p = -1; for(int k = 0; k < sz[now]; k++){ ld dif = out[now][k] - angle; if(dif < 0) dif += PI*2; if(dif < best) best = dif, p = k; } if(done[now][p]) break; done[now][p] = 1; angle = in[now][p] - eps; now = ind[now][p]; } } } cout << cnt - 1 << endl; } return 0; }
-- Main.MatthewChan - 29 Mar 2006