标签存档：计算几何

【算法分析】
以(x0,y0)为起点向x轴正方向作一水平射线，如果中间穿过奇数个多边形上的点，应该是在多边形内，否则应该是在多边形外。
但是有不少特殊情况。
WA了好多次，找了下资料，原来遵循下面几个原则即可：
1对于多边形的水平边不作考虑；
2对于多边形的顶点和L相交的情况，如果该顶点是其所属的边上纵坐标较大的顶点，则计数，否则忽略；
3对于P在多边形边上的情形，直接可判断P属于多边形
【CODE】http://xiudaima.appspot.com/code/detail/3434005

【题目链接】http://acm.sgu.ru/problem.php?contest=0&problem=120

【题目大意】

给定一个N边形，这个N边形上的点都被按顺时针顺序标号为1,2,3…。然后给出标号为m1和m2的点的坐标，最后按标号顺序输出这个多边形所有的点。

【算法分析】

大概思路是先算出中心，然后用向量一直转，就能得到所有的点。

我们先连接m1,m2两个点得到线段L，然后作这条线段的中垂线L’。然后可以很轻松算出图中的角α的值（就是pi/n*(abs(m1-m2))）。然后搞到0~pi之间就可以了。然后我们知道中心点Mid到L的距离d。

然后对向量m1->m2旋转个90°/270°，再把长度调成d。就可以得出中心点的位置了。

当然，如果abs(m1-m2)*2==n要特判一下。因为取不到tan。。。= =，应该有更好的方法。

然后可以得出中心Mid的坐标。然后构造Mid->a[m1]的向量，然后把这个向量每次旋转2*pi/n，就可以了得出所有坐标了噢。

^_^

【关于向量旋转】

学过极坐标的话。。。这个太简单了。

下面推导一下~~

对于向量(x,y)，如果要顺时针旋转的角度β的话，那么这样弄。

转化为极坐标:

x=L*cos(α)

y=L*sin(α)

然后旋转以后变成

x=L*cos(α-β)=L*cos(α)*cos(β)+L*sin(α)*sin(β)

y=L*sin(α-β)=L*sin(α)*cos(β)-L*cos(α)*sin(β)

然后变回直角坐标系

x=x*cos(β)+y*sin(β)

y=y*cos(β)-x*sin(β)

然后就完成了旋转

【CODE】

C++ CODE   :SGU120 1 2 3 4 5 6 7 8 91011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586 #include

【题目大意】
给定n个球，与一条射线，然后让你模拟这个射线在这些球之间弹，输出弹到得球的标号。如果弹到<=10个就有多少输出多少。否则，输出前10个+ 'etc.'

【算法分析】
我们用三维向量来解决。
1、首先根据点到球心距离等于半径 + 我们的运动向量列出一个方程，解这个一元二次方程可以得到我们什么时候碰到这个球。
2、然后碰到以后呢。就把射线起点，球上的碰触点，球心这三个点，搞成一个平面，然后再在上面以向量：碰触点->球心作为法线，考虑镜面反射即可。
然后这里我再次用了一个解方程+向量伸缩的方法得到了他反射以后的向量。

关于第二步的反射部分，我们所知道的条件有：入射向量，法线向量。我的做法具体就是先构造出向量tmp=发现向量-入射向量。并使得|tmp|=|入射向量|。然后就可以通过一些几何知识来搞了。

【其它】WA了几遍，因为我解一元二次方程的时候判断无解写得是delta<1e-6，然后如果1e-6<=delta<0的话，我就先将它置为0再进行后面的计算。
我改成delta<0就为无解就AC了。T_T。。。精度什么的，最讨厌了。

【CODE】
#include #include #include #include #include #define dis(A,B) sqrt( (A.x-B.x)*(A.x-B.x)+(A.y-B.y)*(A.y-B.y)+(A.z-B.z)*(A.z-B.z) )
#define Sqr(x) ((x)*(x))
#define max(x,y) ((x)>(y)?(x):(y))
#define min(x,y) ((x)<(y)?(x):(y))
using namespace std;
const int N=1005;
const double eps=1e-6;
int n;
double res1,res2;
struct Circle{double x,y,z,r;}a[N];
struct Point{double x,y,z;}p,fx;

void init(){
     scanf("%d",&n);
     int i,j;
     for (i=1;i<=n;i++)
       cin >> a[i].x >> a[i].y >> a[i].z >> a[i].r;
     cin >> p.x >> p.y >> p.z >> fx.x >> fx.y >> fx.z;
     fx.x-=p.x;
     fx.y-=p.y;
     fx.z-=p.z;
}

Point Get_Point(double rate){
      Point res;
      res.x=fx.x*rate+p.x;
      res.y=fx.y*rate+p.y;
      res.z=fx.z*rate+p.z;
      return res;
}

double Get_Meet_Rate(Circle Cir){
       double a,b,c,delta,x1,x2;
       Point tmp;
       tmp.x=Cir.x-p.x;
       tmp.y=Cir.y-p.y;
       tmp.z=Cir.z-p.z;
       a=Sqr(fx.x)+Sqr(fx.y)+Sqr(fx.z);
       b=-2*(fx.x*tmp.x + fx.y*tmp.y + fx.z*tmp.z);
       c=-Sqr(Cir.r)+Sqr(tmp.x)+Sqr(tmp.y)+Sqr(tmp.z);
       delta=b*b-4*a*c;
       if (delta<0) return 1e22;
       delta=sqrt(delta);
       x1=(-b+delta)/(2*a);
       x2=(-b-delta)/(2*a);
       if (x1<=eps) x1=1e22;
       if (x2<=eps) x2=1e22;
       res1=x1; res2=x2;
       return min(x1,x2);
}

Point opposite(Point tmp){
      Point res;
      res.x=-tmp.x;
      res.y=-tmp.y;
      res.z=-tmp.z;
      return res;
}

void Reflect(int i,double rate){
     Point cross=Get_Point(rate);
     Get_Meet_Rate(a[i]);
     if ( fabs(res1-res2) < eps ){
       p=cross;
       return;
     }
     Point normal_line,tmp;
     normal_line.x=(a[i].x-cross.x)/2;
     normal_line.y=(a[i].y-cross.y)/2;
     normal_line.z=(a[i].z-cross.z)/2;
     double t=Sqr(normal_line.x)+Sqr(normal_line.y)+Sqr(normal_line.z);
     t/=(normal_line.x*fx.x + normal_line.y*fx.y + normal_line.z*fx.z);
     fx.x*=t; fx.y*=t; fx.z*=t;
     normal_line.x*=2;
     normal_line.y*=2;
     normal_line.z*=2;
    
     tmp.x=normal_line.x-fx.x;
     tmp.y=normal_line.y-fx.y;
     tmp.z=normal_line.z-fx.z;
    
     fx=opposite(tmp);
     p=cross;
    
     bool flag=false;
     if (fx.x>10) flag=true;
     if (fx.y>10) flag=true;
     if (fx.z>10) flag=true;
     if (!flag){
       fx.x*=10;
       fx.y*=10;
       fx.z*=10;
     }
}

void work(){
     int i,best,cnt=0;
     double rate,tmp;
     while (1){
           rate=1e20;
           best=0;
           for (i=1;i<=n;i++){
               tmp=Get_Meet_Rate(a[i]);
               if (tmp                 rate=tmp;
                 best=i;
               }
           }
           if (best==0) break;
           Reflect(best,rate);
           cnt++;
           if (cnt>10) break;
           if (cnt!=1) cout << " ";
           cout << best;
     }
     if (cnt>10) cout << " etc." << endl;
            else cout << endl;
}

int main(){
    init();
    work();
    return 0;
}

【题目大意】
求平面距离最远的点对。
【算法分析】
凸包。
旋转卡壳求出对踵点。
就是弄两个向量，一个贴着凸包上某条边，另一个于该向量平行并与凸包相交，并使得这两个向量把整个图都夹在中间。
【其它】
以前n^2水过。。。于是重写。
另外我这种转法不是按最小角度转的，于是凸包上要注意去掉共线3点的情况。
去这种情况有点小技巧。。。。就是先判他们凸包连续的3点是否共线。然后在构造两个向量。
P:L[i-1]->L[i]
Q:L[i]->L[i+1]
然后用点积判他们方向是否相同，如果相同那就删，否则是不能删的。。。
（YY一下只有两个点的凸包）
【CODE】
#include #include #include #include #define dis(A,B) ((A.x-B.x)*(A.x-B.x)+(A.y-B.y)*(A.y-B.y))
using namespace std;
const int N=100005;
struct Point{int x,y;}a[N],L[N],Std;
int n,tot,v[N];

int Del_cmp(const void *A,const void *B){
Point P=*((Point*)A),Q=*((Point*)B);
if (P.x!=Q.x) return P.x-Q.x;
return P.y-Q.y;
}

void Delsame(){
qsort(a+1,n,sizeof(Point),Del_cmp);
int i,t=1;
for (i=2;i<=n;i++)
if (a[i].x!=a[i-1].x || a[i].y!=a[i-1].y)
a[++t]=a[i];
n=t;
}

int CJ(Point Z,Point A,Point B){
return (A.x-Z.x)*(B.y-Z.y)-(A.y-Z.y)*(B.x-Z.x);
}

int cmp(const void *A,const void *B){
int res=CJ(Std,*((Point*)A),*((Point*)B));
if (res) return res;
return (((Point*)A)->x+((Point*)A)->y-((Point*)B)->x-((Point*)B)->y);
}

void TB(){
Point P,Q,tmp,M;
int i,j;
M.x=M.y=0x7FFFFFFF;
for (i=1;i<=n;i++)
if (a[i].x M=a[i];
j=i;
}
tmp=a[1]; a[1]=a[j]; a[j]=tmp;
Std=a[1];
qsort(a+2,n-1,sizeof(Point),cmp);
L[1]=a[1];
a[n+1]=a[1];
for (tot=1,i=2;i<=n+1;i++){
while (tot>1 && CJ(L[tot-1],L[tot],a[i])>0) tot–;
L[++tot]=a[i];
}
tot–;
}

void Del_Useless_Point(){
if (tot<3) return;
Point P,Q;
int i,j;
for (i=1;i<=tot;i++) v[i]=0;
L[0]=L[tot]; L[tot+1]=L[1];
for (i=1;i<=tot;i++)
if (CJ(L[i-1],L[i],L[i+1])==0){
P.x=L[i].x-L[i-1].x;
P.y=L[i].y-L[i-1].y;
Q.x=L[i+1].x-L[i].x;
Q.y=L[i+1].y-L[i].y;
if (P.x*Q.x+P.y*Q.y>0) v[i]=1;
}
n=0;
for (i=1;i<=tot;i++)
if (!v[i])
a[++n]=L[i];
}

int max(int x,int y){return x>y?x:y;}
int min(int x,int y){return x int i,j,ans=0;
Point P,Q;
bool check;
Std.x=Std.y=0;
a[0]=a[n]; a[n+1]=a[1];
for (i=j=1;i<=n;i++){
for (;;j=j==n?1:j+1){
check=true;
P.x=a[i+1].x-a[i].x;
P.y=a[i+1].y-a[i].y;

Q.x=a[j+1].x-a[j].x;
Q.y=a[j+1].y-a[j].y;
if (CJ(Std,P,Q)<0) check=false; Q.x=a[j-1].x-a[j].x;
Q.y=a[j-1].y-a[j].y;
if (CJ(Std,P,Q)<0) check=false;
if (check) break;
}
ans=max(ans,dis(a[i],a[j-1]));
ans=max(ans,dis(a[i],a[j]));
ans=max(ans,dis(a[i],a[j+1]));
}
printf("%dn",ans);
}

int main(){
int i,j;
while (scanf("%d",&n)!=EOF){
for (i=1;i<=n;i++) scanf("%d%d",&a[i].x,&a[i].y);
Delsame();
TB();
Del_Useless_Point();
switch (n){
case 0:puts("0"); break;
case 1:puts("0"); break;
default:RC(); break;
}
}
return 0;
}

【题目大意】
给定一个线段，和n个圆。
如果线段上某点被某一个圆覆盖，那么就算被覆盖。
最终问被覆盖点占线段总长的百分比。
【算法分析】
圆与线段交点有很多很快的方法，但是我介绍一种异常沙茶，巨慢，却不大容易错的方法。（实际上就是不用分类什么的而已。。。）
首先，对线段两点建立一个向量，然后就可以利用向量伸缩来很方便的表示该线段。
对于一个线段和一个圆。我们先用点到圆心距离为函数。由于他的单峰性质，二分得那个离圆最近的点。
（好像人们都比较喜欢三分，在此表示不解）
现在我们得到了圆覆盖掉的线段上的一个点。
于是我们现在就可以在这个点的基础上，向线段的两个端点爬山。
（就是设一个步长，如果爬出去就停止，并不停将步长除以2）
然后爬完了。得到那两个边缘点就是圆覆盖了线段的两个边缘点。
整个算法复杂度O(lg K) K为坐标距离除以eps。
【一些小技巧】
设向量为(x1,y1) -> (x2,y2)
那么我们就可以用一个数rate表示线段上的点。（其中0<=rate<=1）
具体就是( x1+rate*(x2-x1) , y1+rate*(y2-y1) )。
然后就比较方便了。
【对于本题】
本题的话，就是可以按每个圆覆盖点对的前端rate值排序。
然后就很容易贪心得到每部分极长被覆盖距离。
【其它】
程序有点货不对板。。。求圆与线段的一个交点时用的直接暴力枚举，不是二分。。。
【CODE】
#include #include #include #include #include #define sqr(x) ((x)*(x))
using namespace std;
const int N=125;
const int Times=3000;
const double Start_Step=100;
const double eps=1e-6;
struct Point{int x,y;}p[N];
int n,cnt;
int R[N];
struct Line{double St,Ed;}L[N];

void op(){
for (int i=1;i<=n+2;i++) swap(p[i].x,p[i].y);
}

inline bool Cover(double x,int i){
if (x double y=p[1].y+(x-p[1].x)/(p[2].x-p[1].x)*(p[2].y-p[1].y);
return sqr(R[i])>=sqr(x-p[i].x)+sqr(y-p[i].y);
}

void work(){
n+=2;
double sx=(double)(p[2].x-p[1].x)/Times,sy=(double)(p[2].y-p[1].y)/Times;
int i,j;
cnt=0;
for (i=3;i<=n;i++){
double x=p[1].x,y=p[1].y;
for (j=0;j<=Times;j++){
if (Cover(x,i)){
cnt++;
double tx=x,ty=y,step;
for (step=Start_Step;step>eps;step/=2)
while (Cover(x-step,i)) x-=step;
L[cnt].St=x;
x=tx;
for (step=Start_Step;step>eps;step/=2)
while (Cover(x+step,i)) x+=step;
L[cnt].Ed=x;
break;
}
x+=sx; y+=sy;
}
}
}

int cmp(const void *A,const void *B){
Line AA=*((Line*)A),BB=*((Line*)B);
if (AA.St if (AA.St>BB.St) return 1;
return 0;
}

void solve(){
qsort(L+1,cnt,sizeof(Line),cmp);
double now=0,Sum=0,last;
for (int i=1;i<=cnt;i++){
if (i==1 || last Sum+=now;
now=L[i].Ed-L[i].St;
last=L[i].Ed;
}
else{
if (L[i].Ed>last){
now+=L[i].Ed-last;
last=L[i].Ed;
}
}
}
Sum+=now;
printf("%.2lfn",Sum/(p[2].x-p[1].x)*100);
}

int main(){
while (1){
scanf("%d",&n); if (!n) return 0;
scanf("%d%d%d%d",&p[1].x,&p[1].y,&p[2].x,&p[2].y);
for (int i=3;i<=n+2;i++)
scanf("%d%d%d",&p[i].x,&p[i].y,&R[i]);
if (abs(p[1].x-p[2].x) if (p[1].x>p[2].x) swap(p[1],p[2]);
work();
solve();
}
}