题目地址:LA3616
直接O(nm)暴力枚举.. 简单题
代码:
#include<iostream>
#include<cmath>
#include<cstdio>
#include<iomanip>
#include<ctime>
#include<cstdlib>
#include<vector>
using namespace std;
const double eps=1e-8;
// 1e-8?
const double PI=acos(-1.0);
int dcmp(double x)
{
return (x>eps)-(x<-eps);
}
struct Point3
{
double x,y,z;
Point3(double x=0,double y=0,double z=0):x(x),y(y),z(z) {}
};
typedef Point3 Vector3;
Vector3 operator+(Point3 A,Point3 B)
{
return Point3(A.x+B.x,A.y+B.y,A.z+B.z);
}
Vector3 operator-(Point3 A,Point3 B)
{
return Point3(A.x-B.x,A.y-B.y,A.z-B.z);
}
Vector3 operator*(Point3 A,double p)
{
return Point3(A.x*p,A.y*p,A.z*p);
}
Vector3 operator/(Point3 A,double p)
{
return Point3(A.x/p,A.y/p,A.z/p);
}
bool operator<(Point3 &A,Point3 &B)
{
return A.x<B.x||(A.x==B.x&&A.y<B.y)||(A.x==B.x&&A.y==B.y&&A.z<B.z);
}
bool operator==(Point3 &A,Point3 &B)
{
return dcmp(A.x-B.x)==0&&dcmp(A.y-B.y)==0&&dcmp(A.z-B.z)==0;
}
ostream & operator<<(ostream & out,Vector3 A)
{
out<<A.x<<' '<<A.y<<' '<<A.z<<endl;
return out;
}
Point3 read_point3()
{
Point3 p;
scanf("%lf%lf%lf",&p.x,&p.y,&p.z);
return p;
}
double Dot(Vector3 A,Vector3 B)
{
return A.x*B.x+A.y*B.y+A.z*B.z;
}
double Length(Vector3 A)
{
return sqrt(Dot(A,A));
}
double Angle(Vector3 A,Vector3 B)
{
return acos(Dot(A,B)/Length(A)/Length(B));
}
Vector3 Cross(Vector3 A,Vector3 B) // 真正的叉积
{
return Vector3(A.y*B.z-A.z*B.y,A.z*B.x-A.x*B.z,A.x*B.y-A.y*B.x);
}
// 面积的两倍
double Area2(Point3 A,Point3 B,Point3 C)
{
return Length(Cross(B-A,C-A));
}
double DistanceToPlane(const Point3 &p,const Point3 & p0,const Vector3 &n)
{
return fabs(Dot(p-p0,n))/Length(n);
}
// n不用管方向
Point3 GetPlaneProjection(Point3 p,Point3 p0,Vector3 n)
{
return p-n*Dot(p-p0,n)/Dot(n,n);
}
// 三点式
Point3 GetPlaneProjection(Point3 p,Point3 p0,Point3 p1,Point3 p2)
{
Vector3 n=Cross(p1-p0,p2-p0);
return p-n*Dot(p-p0,n)/Dot(n,n);
}
Point3 LinePlaneIntersection(Point3 p1,Point3 p2,Point3 p0,Vector3 n)
{
double t=Dot(n,p0-p1)/Dot(n,p2-p1);
return p1+(p2-p1)*t;
}
// 点在三角形内部
bool PointInTri(Point3 p,Point3 p0,Point3 p1,Point3 p2)
{
double area1=Area2(p, p0, p1);
double area2=Area2(p, p1, p2);
double area3=Area2(p, p2, p0);
return dcmp(area1+area2+area3-Area2(p0,p1,p2))==0;
}
// 线段与三角形相交 不考虑共面情形
bool TriSegIntersection(Point3 p0,Point3 p1,Point3 p2,Point3 A,Point3 B,Point3 & p)
{
Vector3 n=Cross(p1-p0, p2-p0);
if(dcmp(Dot(n,B-A))==0) return false;
else{
double t=Dot(n,p0-A)/Dot(n,B-A);
if(dcmp(t)<0||dcmp(t-1)>0) return false;
else
{
p=A+(B-A)*t;
return PointInTri(p,p0,p1,p2);
}
}
}
double DistanceToLine(Point3 P,Point3 A,Point3 B)
{
Vector3 v1=B-A;
Vector3 v2=P-A;
return Length(Cross(v1, v2))/Length(v1);
}
double DistanceToSegment(Point3 P,Point3 A,Point3 B)
{
if(A==B) return Length(P-A);
Vector3 v1=B-A,v2=P-A,v3=P-B;
if(dcmp(Dot(v1,v2))<0) return Length(v2);
else if(dcmp(Dot(v1,v3))>0) return Length(v3);
else return Length(Cross(v1,v2))/Length(v1);
}
double volume6(Point3 A,Point3 B,Point3 C,Point3 D)
{
return Dot(D-A,Cross(B-A,C-A));
}
// 0-1随机数
double rand01() { return rand()/(double)RAND_MAX;}
double randeps() {
//srand((unsigned) time(NULL));
return (rand01()-0.5)*eps;}
Point3 add_noise(Point3 P)
{
return Point3(P.x+randeps(),P.y+randeps(),P.z+randeps());
}
// 三维凸包 增量法
struct Face {
int v[3];
Face(int a, int b, int c) { v[0] = a; v[1] = b; v[2] = c; }
Vector3 Normal(const vector<Point3>& P) const {
return Cross(P[v[1]]-P[v[0]], P[v[2]]-P[v[0]]);
}
// f是否能看见P[i]
int CanSee(const vector<Point3>& P, int i) const {
return Dot(P[i]-P[v[0]], Normal(P)) > 0;
}
};
// 增量法求三维凸包
// 注意:没有考虑各种特殊情况(如四点共面)。实践中,请在调用前对输入点进行微小扰动
vector<Face> CH3D(const vector<Point3>& P) {
int n = P.size();
vector<vector<int> > vis(n);
for(int i = 0; i < n; i++) vis[i].resize(n);
vector<Face> cur;
cur.push_back(Face(0, 1, 2)); // 由于已经进行扰动,前三个点不共线
cur.push_back(Face(2, 1, 0));
for(int i = 3; i < n; i++) {
vector<Face> next;
// 计算每条边的“左面”的可见性
for(int j = 0; j < cur.size(); j++) {
Face& f = cur[j];
int res = f.CanSee(P, i);
if(!res) next.push_back(f);
for(int k = 0; k < 3; k++) vis[f.v[k]][f.v[(k+1)%3]] = res;
}
for(int j = 0; j < cur.size(); j++)
for(int k = 0; k < 3; k++) {
int a = cur[j].v[k], b = cur[j].v[(k+1)%3];
if(vis[a][b] != vis[b][a] && vis[a][b]) // (a,b)是分界线,左边对P[i]可见
next.push_back(Face(a, b, i));
}
cur = next;
}
return cur;
}
double DistanceToPlane( const Point3 &p,const Point3 & p0,const Point3 & p1,const Point3 & p2)
{
return fabs(volume6(p,p0,p1,p2)/Area2(p0, p1, p2));
}
// 四面体质心 就是四个点的平均
Point3 Centroid(const Point3 & A,const Point3 & B,const Point3 & C,const Point3 & D)
{
return (A+B+C+D)/4.0;
}
struct ConvexPolyhedron
{
int n; // 总顶点数
vector<Point3> P,P2; // P是真实顶点, P1是扰乱后的顶点
vector<Face> faces;
void init(int n)
{
this->n=n;
P.resize(n);
P2.resize(n);
for(int i=0;i<n;i++)
{
P[i]=read_point3();
P2[i]=add_noise(P[i]);
}
faces=CH3D(P2);
}
// 本题需要,写在类里面
double mindist(Point3 C)
{
double ans=1e30;
for(int i=0;i<faces.size();i++)
{
ans=min(ans,DistanceToPlane(C,P[faces[i].v[0]], P[faces[i].v[1]], P[faces[i].v[2]])); // 都是P这个数组
}
return ans;
}
Point3 centroid()
{
Point3 C=Point3(0,0,0);
Point3 tot=Point3(0,0,0);
double totv=0;
for(int i=0;i<faces.size();i++)
{
Point3 p0=P[faces[i].v[0]];
Point3 p1=P[faces[i].v[1]];
Point3 p2=P[faces[i].v[2]];
double v=volume6(p0, p1, p2, C);
tot=tot+Centroid(p0,p1,p2,C)*v;
totv+=v;
}
return tot/totv;
}
};
// self
// p0,p1在面A,B,C 的同侧
bool SameSide(Point3 p0,Point3 p1,Point3 A,Point3 B,Point3 C)
{
Vector3 n=Cross(B-A, C-A);
Vector3 v1=p0-A;
Vector3 v2=p1-A;
return Dot(v1, n)*Dot(v2,n)>=0;
}
// 改写读点函数
bool read_point3(Point3 & P)
{
if(scanf("%lf%lf%lf",&P.x,&P.y,&P.z)==3)
{
return 1;
}
else return 0;
}
// 求质心
Point3 Centroid(Point3 &A,Point3 & B,Point3 &C,Point3 &D)
{
return (A+B+C+D)/4.0;
}
Point3 Centroid(Point3 &A,Point3 & B,Point3 &C,Point3 &D,Point3 &E) // 注意A,B,C 为四面体公共平面
{
double v1=volume6(A, B, C, D);
double v2=volume6(A, B, C, E);
v1=fabs(v1);
v2=fabs(v2);
Point3 C1=Centroid(A, B, C, D);
Point3 C2=Centroid(A, B, C, E);
return (C1*v1+C2*v2)/(v1+v2);
}
bool OnSamePlane(Point3 &A,Point3 & B,Point3 &C,Point3 &D)
{
return dcmp(volume6(A, B, C, D))==0;
}
int n,m;
Point3 P[505];
Point3 T[55];
double T_ang[55];
bool cansee(int x)
{
for(int i=0;i<m;i++)
{
if(dcmp(Angle(P[x], T[i])-T_ang[i])<0)
{
return 1;
}
}
return 0;
}
int main()
{
while(cin>>n)
{
if(!n) break;
for(int i=0;i<n;i++)
{
P[i]=read_point3();
}
cin>>m;
for(int i=0;i<m;i++)
{
T[i]=read_point3();
scanf("%lf",&T_ang[i]);
}
int cnt=0;
for(int i=0;i<n;i++)
{
if(cansee(i))
{
cnt++;
}
}
cout<<cnt<<endl;
}
}