标签:des blog os io for div amp log
直线关于球的多次反射,求最后一次反射点
#include <iostream> #include <cstdio> #include <cstring> #include <algorithm> #include <cmath> using namespace std; const double inf=1e10; const double eps=1e-8; struct point { double x,y,z; // point (double _x,double _y,double _z){ x=_x; y=_y; z=_z; }; }; struct sphe { point cent; double r; }; struct vect { point st,des; }; sphe cir[110]; vect livc; int n; point operator -(const point &u,const point &v){ point ret; ret.x=u.x-v.x; ret.y=u.y-v.y; ret.z=u.z-v.z; return ret; } double dot(point x,point y){ return x.x*y.x+x.y*y.y+x.z*y.z; } point xmulti(point u,point v){ point ret; ret.x=(u.y*v.z-v.y*u.z); ret.y=(u.z*v.x-u.x*v.z); ret.z=(u.x*v.y-u.y*v.x); return ret; } double dis(point x,point y){ return sqrt((x.x-y.x)*(x.x-y.x)+(x.y-y.y)*(x.y-y.y)+(x.z-y.z)*(x.z-y.z)); } double vlen(point x){ return sqrt(x.x*x.x+x.y*x.y+x.z*x.z); } point construct(){ point crop; crop.x=crop.y=crop.z=0; double stoc=inf; point tmpcrop; point foot,tmpfoot; bool flag; point tmp; int k; while(true){ flag=false; stoc=inf; for(int i=0;i<n;i++){ if(dot(livc.des-livc.st,cir[i].cent-livc.st)>=-eps){//判断圆是否与直线同向,通过点积判方向 double D=vlen(xmulti(livc.des-livc.st,cir[i].cent-livc.st))/dis(livc.st,livc.des); // cout<<D<<‘ ‘<<i<<endl; if(D-cir[i].r<=eps){ //半径小于D,相交 flag=true; // cout<<"YES"<<endl; double u=dot(cir[i].cent-livc.st,livc.des-livc.st)/(dis(livc.st,livc.des)*dis(livc.st,livc.des)); //计算垂足。可通过向量的比例所得方程,联合垂直点积为0的方程解得 tmpfoot=livc.st; tmpfoot.x+=u*(livc.des.x-livc.st.x); tmpfoot.y+=u*(livc.des.y-livc.st.y); tmpfoot.z+=u*(livc.des.z-livc.st.z); // cout<<tmpfoot.x<<‘ ‘<<tmpfoot.y<<‘ ‘<<tmpfoot.z<<‘ ‘<<endl; u=sqrt((cir[i].r*cir[i].r-D*D))/dis(livc.st,livc.des); //计算交点。垂足到圆上交点方向与直线反方向相同 //通过两者距离比计算出向量的转化 tmpcrop=tmpfoot; tmp=livc.st-livc.des; tmpcrop.x+=tmp.x*u; tmpcrop.y+=tmp.y*u; tmpcrop.z+=tmp.z*u; D=dis(tmpcrop,livc.st); // cout<<D<<endl; if(D<stoc){ //若与多个圆相交,选取较近的一个 stoc=D; crop=tmpcrop; k=i; } } } } if(!flag) return crop; double tu=dot(livc.st-cir[k].cent,crop-cir[k].cent)/(dis(crop,cir[k].cent)*dis(crop,cir[k].cent)); tmpfoot=cir[k].cent; //计算反射线。直线st点关于交点与球心的直线 对称点作为反射线的des点 tmpfoot.x+=tu*(crop.x-cir[k].cent.x); tmpfoot.y+=tu*(crop.y-cir[k].cent.y); tmpfoot.z+=tu*(crop.z-cir[k].cent.z); //知直线st点到反射线des点的方向与st点到关于对称线垂足方向相同且为两倍 livc.des.x=((tmpfoot.x-livc.st.x)*2+livc.st.x); //通过这样可以求对称点 livc.des.y=((tmpfoot.y-livc.st.y)*2+livc.st.y); livc.des.z=((tmpfoot.z-livc.st.z)*2+livc.st.z); livc.st=crop; // cout<<livc.des.x<<‘ ‘<<livc.des.x<<‘ ‘<<livc.des.x<<endl; } } int main(){ point tmp; double r; while(scanf("%d",&n),n){ livc.st.x=livc.st.y=livc.st.z=0; scanf("%lf%lf%lf",&tmp.x,&tmp.y,&tmp.z); livc.des=tmp; for(int i=0;i<n;i++){ scanf("%lf%lf%lf%lf",&cir[i].cent.x,&cir[i].cent.y,&cir[i].cent.z,&cir[i].r); } tmp=construct(); printf("%.4lf %.4lf %.4lf\n",tmp.x,tmp.y,tmp.z); } }
标签:des blog os io for div amp log
原文地址:http://www.cnblogs.com/jie-dcai/p/3902389.html