标签：Shadows 1146 sum th theta C2 C1 rot AIZU

传送门：Secrets in Shadows

题意：
给定\(n(n\le 100)\)个半径为\(1\)的圆柱体，假设太阳光是从无限远处来的平行于地面直线，圆柱会在地面上投射出无限长的矩形阴影。
宽度定义为这些矩形的宽度并，问宽度最大的太阳角度。

题解：
一开始以为是枚举两个圆心，然后垂直连线啥的，发现错了，但是网上找不到好的题解。。于是看std理解了思路，遂写此题解。

显然，我们有了太阳的角度，求宽度是很简单的，只需要逆时针旋转\(\theta\)度，然后计算\(y\)轴上的宽度并，排序+扫一遍就行。
但是如何确定角度呢。
注意到，假如圆的相交关系不变，宽度的变化跟角度的变化是连续的。
先给出这个夹角的一个特解，然后找到每一段连续相交的圆的\(y\)最小和最大的两个，那么

\[w = \sum (rot(p_a, \theta + t).y - rot(p_b, \theta + t).y) + 2 \]

\[w = \sum (A.xsin\theta + A.ycos\theta - B.xsin\theta - B.ycos\theta) + 2 \]

求个导

\[dw = \sum (A.x cos\theta - A.ysin\theta - B.x cos \theta + B.y sin \theta) d\theta \]

\[dw = \sum (A.x - B.x) cos\theta - \sum (A.y - B.y)sin \theta) d\theta \]

令导数为\(0\)，得到\(tan\theta = \frac{\sum (A.y - B.y)}{\sum (A.x - B.x)}\)，于是可以得到\(\theta\)。

那么如何得到所有的特解呢，我们发现，两个圆投影是否相交当且仅当在两圆公切线处发生改变，于是我们枚举所有公切线角度，排序后就能得到每个区间了，可以取区间角平分线作为特解。

代码比较冗长。

#include <bits/stdc++.h>
#define pt(x) cout << x << endl;
#define Mid ((l + r) / 2)
#define lson (rt << 1)
#define rson (rt << 1 | 1)
using namespace std;
const double Pi = acos(-1.0);
const double eps = 1e-11;
struct Point {
	double x, y;
	Point() : x(0), y(0) {}
	Point(double x, double y) : x(x), y(y) {}
	Point operator+(const Point &a)const{ return Point(x + a.x, y + a.y); }
	Point operator-(const Point &a)const{ return Point(x - a.x, y - a.y); }
} a[109];
typedef Point Vector;
struct line{
	Point a, b;
	line() {}
	line(const Point &a,const Point &b) : a(a), b(b) {}
};
 
struct circle{
	Point c;
	double r;
	circle() : c(Point(0,0)), r(0) {}
	circle(const Point &c, double r) : c(c), r(r) {}
};
Point operator*(double c, const Point &a){ return Point(c * a.x, c * a.y); }
double abs(const Point &a){ return sqrt(a.x * a.x + a.y * a.y); }
Point rot(const Point &a, double theta){ return Point(a.x * cos(theta) - a.y * sin(theta), a.x * sin(theta) + a.y * cos(theta)); }
double arg(const Point &a){
	double t = atan2(a.y, a.x);
	return t < 0 ? t + 2 * Pi:t;
}
int tangent(const circle &C1,const circle &C2,vector<line> &res){
	double d = abs(C1.c - C2.c);
	if(d < eps) return 0;
 
	int c=0;
	if(C1.r + C2.r < d - eps){
		double t = acos((C1.r + C2.r) / d);
		res.push_back(line(C1.c + rot(C1.r / d * (C2.c - C1.c), t), C2.c + rot(C2.r / d * (C1.c - C2.c), t)));
		res.push_back(line(C1.c + rot(C1.r / d * (C2.c - C1.c),-t), C2.c + rot(C2.r / d * (C1.c - C2.c),-t)));
		c += 2;
	} else if(C1.r + C2.r < d + eps){
		Point p = C1.c + C1.r / d * (C2.c - C1.c);
		res.push_back(line(p, p + rot(C2.c - C1.c, Pi / 2)));
		c++;
	}
 
	if(abs(C1.r - C2.r) < d - eps){
		double t1 = acos((C1.r - C2.r) / d), t2 = Pi - t1;
		res.push_back(line(C1.c + rot(C1.r / d * (C2.c - C1.c), t1), C2.c + rot(C2.r / d * (C1.c - C2.c),-t2)));
		res.push_back(line(C1.c + rot(C1.r / d * (C2.c - C1.c),-t1), C2.c + rot(C2.r / d * (C1.c - C2.c), t2)));
		c+=2;
	} else if(abs(C1.r - C2.r) < d + eps){
		Point p = C1.c + C1.r / d * (C2.c - C1.c);
		res.push_back(line(p, p + rot(C2.c - C1.c, Pi / 2)));
		c++;
	}
 
	return c;
}
 
int n;
double maxn, minn, sima, simi;
double calc(double theta) {
	double tt[109] = {0};
	for(int k = 1; k <= n; k++) {
		tt[k] = rot(a[k], theta).y;
	}
	sort(tt + 1, tt + 1 + n);
	double mr = 0, sum = 0;
	sum = 2.0;
	mr = tt[1] + 1.0;
	for(int k = 1; k <= n; k++) {
		if(tt[k] + 1.0 < mr) continue;
		if(tt[k] - 1.0 < mr) {
			sum += tt[k] + 1.0 - mr;
		} else {
			sum += 2.0;
		}
		mr = tt[k] + 1.0;
	}
	return sum;
}
double calc_peak(double phi) {
	pair<double, int> p[109];
	for(int i = 1; i <= n; i++) {
		p[i] = make_pair(rot(a[i], phi).y, i);
	}
	sort(p + 1, p + 1 + n);
	int pre = 1;
	double x = 0, y = 0;
	for(int i = 2; i <= n + 1; i++) {
		if(i == n + 1 || p[i - 1].first + 1 < p[i].first - 1) {
			x += (a[p[i - 1].second].x - a[p[pre].second].x);
			y += (a[p[i - 1].second].y - a[p[pre].second].y);
			pre = i;
		}
	}
	double ans = arg({y, x});
	if(ans > Pi) ans -= Pi;
	return ans;
}
 
void work() {
	if(n == 0) exit(0);
	for(int i = 1; i <= n; i++) 
		scanf("%lf%lf", &a[i].x, &a[i].y);
	maxn = -1e18;
	minn = 1e18;
	vector<line> v;
	for(int i = 1; i <= n; i++) {
		for(int j = i + 1; j <= n; j++) {
			tangent(circle(a[i], 1.0), circle(a[j], 1.0), v);
		}
	}
	vector<double> th;
	for(int i = 0; i < v.size(); i++) {
		line &l = v[i];
		double tt = arg(l.b - l.a);
		if(tt < Pi) tt = Pi - tt;
		else tt = 2 * Pi - tt;
		th.push_back(tt);
	}
	th.push_back(0);
	th.push_back(Pi);
	sort(th.begin(), th.end());
	int tot = 0;
	for(int i = 0; i < th.size(); i++) 
		if(i == 0 || th[i] > th[i - 1] + eps) 
			th[tot++] = th[i];
	th.erase(th.begin() + tot, th.end());
	for(int i = 0; i < tot; i++) {
		double sum = calc(th[i]);
		if(sum > maxn + eps) {
			maxn = sum;
			sima = th[i];
		}
		if(sum < minn - eps) {
			minn = sum;
			simi = th[i];
		}
		if(i - 1 >= 0) {
			double phi0 = (th[i - 1] + th[i]) / 2;
			double phi1 = calc_peak(phi0);
			double sum = calc(phi1);
			if(sum > maxn + eps) {
				maxn = sum;
				sima = phi1;
			}
			if(sum < minn - eps) {
				minn = sum;
				simi = phi1;
			}
		}
	}
	printf("%.15f\n%.15f\n", simi, sima);
}
signed main()
{
	ios :: sync_with_stdio(0);
	cin.tie(0);
	while(~scanf("%d", &n)) work(); 
	return 0;
}

标签：Shadows,1146,sum,th,theta,C2,C1,rot,AIZU
来源： https://www.cnblogs.com/onglublog/p/15519283.html

本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享；
2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关；
3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关；
4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除；
5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。