?? ray tracing.cpp
字號:
//ray tracing
#include"math.h"
#include"iostream"
#include"fstream"
#include "vector"
using namespace std;
class CPointd
{
public:
double x;
double y;
double V;
CPointd(){x=0,y=0;}
};
int main()
{
int i,j,xi,yj,NS,NG,n=0,N=500;
char opfname[20];ifstream ifile;
double f=1.0,df=2,e=100000;int a=0,b=0;
double temp,tangle0,tangle1,Pa,v;
CPointd *p,p0,p1,p2;bool judge=true;
vector <CPointd> pp;
cin>>opfname;
ifile.open(opfname,ios::in|ios::binary);
ifile.read((char*)&NS,4);
ifile.read((char*)&NG,4);
p=new CPointd[NS*NG];
for(i=0;i<NS;i++)
for(j=0;j<NG;j++)
{
ifile.read((char*)&((p+i*NG+j)->x),8);
ifile.read((char*)&temp,8);
ifile.read((char*)&temp,8);
ifile.read((char*)&((p+i*NG+j)->y),8);
ifile.read((char*)&((p+i*NG+j)->V),8);
(p+i*NG+j)->y=-(p+i*NG+j)->y;
(p+i*NG+j)->V=3000;//velocity
}
cin>>xi>>yj;
int kx=0,ky=yj;
while(fabs(e)>0.000000001&&n<=N)
{
pp.clear();
judge=true;
kx=0,ky=yj;
p0.x=0;
p0.y=(p+ky)->y;
p2.x=f*(p+xi*NG)->x;
p2.y=0;
tangle1=(p2.x-p0.x)/(p2.y-p0.y);
Pa=(p+kx*NG+ky)->V*sqrt(pow(p2.x,2)+pow(p0.y,2))/p2.x;
p1.x=(p+kx*NG+ky-1)->x,p1.y=(p+kx*NG+ky-1)->y;
tangle0=(p1.x-p0.x)/(p1.y-p0.y);
pp.push_back(p0);
while(judge==true)
{
if(tangle1<=tangle0)
{
if(ky>=1)
{
ky--;
p2.y=(p+ky)->y;
}
else
{
p2.y=0; judge=false;
}
p2.x=(p2.y-p0.y)*tangle1+p0.x;
if(fabs(tangle1-tangle0)<0.0000000000001)
{ kx++;
if(kx>xi)
judge=false;
}
}
else
{
p2.x=(p+kx*NG)->x;
p2.y=(p2.x-p0.x)/tangle1+p0.y;
kx++;
if(kx>xi)
judge=false;
}
if(ky>=0)
{p1.x=(p+kx*NG)->x,p1.y=(p+ky-1)->y;}
else
{
p1.x=(p+kx*NG)->x,p1.y=0;
}
p0.x=p2.x,p0.y=p2.y;
v=(p+kx*NG+ky)->V;
tangle0=(p1.x-p0.x)/(p1.y-p0.y);
if(sqrt(pow(Pa/v,2)-1)!=0)
tangle1=1.0/sqrt(pow(Pa/v,2)-1);
else
{
judge=false;
}
pp.push_back(p0);
}
if(p0.x<(p+xi*NG)->x&&sqrt(pow(Pa/v,2)-1)>0)
{
if(a!=0)
{
df=df/2;
}
f=f+df;
b=1;
}
else if(p0.y<0)
{
if(b==0)
f=f/2.0;
else
{df=df/2;
f=f-df;}
a=1;
}
n++;
e=sqrt(pow(p0.x-(p+xi*NG)->x,2)+pow(p0.y,2));
}
if(fabs(e)<0.000000001)
{
cout<<pp.size()<<endl;
for(i=0;i<pp.size();i++)
cout<<pp[i].x<<" "<<-pp[i].y<<endl;
}
else if(fabs(e)>0.000000001)
cout<<"這種速度分布兩點不可進行透射追蹤!"<<endl;
p=NULL;
return 0;
}
?? 快捷鍵說明
復(fù)制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -