1621.cpp

这是哈尔滨工业大学acmOJ的源代码

/*  This Code is Submitted by wywcgs for Problem 1621 on 2006-06-02 at 21:49:36 */ 
#include <cstdio>
#include <cmath>
#include <cassert>
#include <algorithm>
using namespace std;

const int SN = 32, RN = 11;
const double eps = 1e-4;

class Sphere {
public:
	int x, y, r;
	void make() { scanf("%d %d %d", &x, &y, &r); }
	bool out(double, double, double) const;
};
bool Sphere::out(double a, double b, double c) const {
	double d = (a*x+b*y+c)*(a*x+b*y+c) / (a*a+b*b);
	return (r*r-d < eps);
}

int main()
{
	Sphere sph[SN];
	int n, i, s, t;

	for(t = 1; scanf("%d", &n) != EOF && n != 0; t++) {
		for(i = 0; i < n; i++) sph[i].make();
		double x, y, dx, dy; scanf("%lf %lf %lf %lf", &x, &y, &dx, &dy);
		printf("Scene %d\n", t);
		bool end = true;
		for(s = 0; s < RN; s++) {
			int ro = -1; double k = 1e10, dc = dx*y-dy*x;
			assert(fabs(dy*x-dx*y+dc) < eps);
			for(i = 0; i < n; i++) {
				if(sph[i].out(dy, -dx, dc)) continue;
				double x1, x2, y1, y2, k1, k2;
				if(fabs(dx) < eps) {
					x1 = x2 = x;
					double dr = sqrt(sph[i].r*sph[i].r-(x-sph[i].x)*(x-sph[i].x));
					y1 = sph[i].y + dr; y2 = sph[i].y - dr;
					k1 = (y1-y)/dy; k2 = (y2-y)/dy;
				} else {
					double kl = dy/dx, bl = y-kl*x, db = bl-sph[i].y;
					double a = kl*kl+1, b = 2*(kl*db-sph[i].x), c = sph[i].x*sph[i].x+db*db-sph[i].r*sph[i].r,
						delta = sqrt(b*b-4*a*c);
					x1 = (-b+delta)/(2*a); x2 = (-b-delta)/(2*a);
					y1 = kl*x1+bl; y2 = kl*x2+bl;
					k1 = (x1-x)/dx; k2 = (x2-x)/dx;
				}
				if(k1 < eps && k2 < eps) continue;
				else if(k1 < eps && k > k2) { k = k2; ro = i; }
				else if(k2 < eps && k > k1) { k = k1; ro = i; }
				else if(k > min(k1, k2)) { k = min(k1, k2); ro = i; }
			}
			if(s == RN-1 && ro != -1) { end = false; break; }
			if(ro == -1) break;
			printf("%d ", ro+1);
			x += dx*k; y += dy*k;
			double dxm = x-sph[ro].x, dym = y-sph[ro].y, alpha = (dx*dym-dy*dxm)/(dx*dxm+dy*dym);
			double pdx = dx, pdy = dy;
			dx = alpha*dym-dxm; dy = -alpha*dxm-dym;
			double d1 = dym*pdx-pdy*dxm, d2 = dym*dx-dy*dxm;
			if(d1*d2 < 0) { dx = -dx; dy = -dy; d2 = -d2; }
			else if(fabs(d1*d2) < eps && dx*pdx > 0) { dx = -dx; dy = -dy; d2 = -d2; }
		}
		printf("%s\n\n", end ? "inf" : "...");
	}
	
	return 0;
}