您现在的位置:首页 >> ★免费资源 >> 源码下载 >> 内容

MATLAB代做|FPGA代做--点散射体团簇体积散射蒙特卡罗模拟主程序

时间:2018-12-12 0:29:33 点击:

  核心提示:点散射体团簇体积散射蒙特卡罗模拟主程序...
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                            ptscext.m                               %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Main program for the Monte Carlo simulations of volumetric 
% scattering by clusters of point scatterers
% -- Part of the Electromagnetic Wave MATLAB Library (EWML)--
%    <http://www.emwave.com/>

% Original: by K.H. Ding, November 1998

% input parameters
nc=1%input('Enter number of clusters : ');
npc=2790%input('Enter number of point scatterers per cluster : ');
wavel=0.0319%input('Enter incident wavelenght (meter) : ');
L=1%input('Enter cubic box length (# of wavelengths) : ');
lc=1%input('Enter cluster length (# of wavelengths) : ');
fp=6.2887e-006%input('Enter real part of scattering amplitude : ');
tai=pi/2%input('Enter polar incidence angle (degree) : ');
phi=0%input('Enter azimuthal incidence angle (degree) : ');
nrlz=20%input('Enter number of realizations : ');
seed=123245%input('Enter seed for random numbers : ');

ntas=10;
gmatas=(1:ntas-1)./sqrt((2*(1:ntas-1)).^2-1);
[evtas,xatas]=eig(diag(gmatas,-1)+diag(gmatas,1));
[xatas,katas]=sort(diag(xatas));
watas=2*evtas(1,katas)'.^2;

nphs=20;
gmaphs=(1:nphs-1)./sqrt((2*(1:nphs-1)).^2-1);
[evphs,xaphs]=eig(diag(gmaphs,-1)+diag(gmaphs,1));
[xaphs,kaphs]=sort(diag(xaphs));
waphs=2*evphs(1,kaphs)'.^2;

dtr=pi/180;
twopi=2*pi;

N=nc*npc;
ncps=L/lc;
k=twopi/wavel
fpp=twopi*fp*fp;
f=(fp+i*fpp)*wavel
hsl=L*wavel/2;
csl=lc*wavel;
hcsl=csl/2;
V=(2*hsl)^3;
n0=N/V
Kei=n0*4*pi*imag(f)/k

% scattering directions
nsang=ntas*nphs;
kxs=ones(nsang,1);
kys=ones(nsang,1);
kzs=ones(nsang,1);
wf=ones(nsang,1);
iang=0;
for ita=1:ntas
  wt=watas(ita);
  tas=acos(xatas(ita));
  for iph=1:nphs
    iang=iang+1;
    wp=waphs(iph);
    wf(iang)=wt*wp;
    phs=pi+pi*xaphs(iph);
    cps=cos(phs);
    sps=sin(phs);
    cts=cos(tas);
    sts=sin(tas);
    kxs(iang)=-i*k*sts*cps;
    kys(iang)=-i*k*sts*sps;
    kzs(iang)=-i*k*cts;
    tdeg(iang)=tas/dtr;
    pdeg(iang)=phs/dtr;
  end
end

% incident direction
tai=tai*dtr;
phi=phi*dtr;
cti=cos(tai);
sti=sin(tai);
cpi=cos(phi);
spi=sin(phi);
kxi=i*k*sti*cpi;
kyi=i*k*sti*spi;
kzi=-i*k*cti;

fext=fopen('ptscext.dat','w+');

rand('seed',seed);

xrow=zeros(1,N);
yrow=zeros(1,N);
zrow=zeros(1,N);
xc=zeros(npc,nc);
yc=zeros(npc,nc);
zc=zeros(npc,nc);
xccol=zeros(N,1);
yccol=zeros(N,1);
zccol=zeros(N,1);
rx=zeros(N,N);
ry=zeros(N,N);
rz=zeros(N,N);
rr=ones(N,N);
rr1=ones(N,N);
b=ones(N,1);
Z=ones(N,N);
Psi=ones(N,1);
erow=ones(1,N);
ecol=ones(N,1);
emat=ones(npc,nc);
xcmat=eye(nc);
ycmat=eye(nc);
zcmat=eye(nc);

FF=ones(nsang,1);
FFc=zeros(nsang,1);
Pt=zeros(nsang,1);
Pc=zeros(nsang,1);
P=zeros(nsang,1);

% Monte Carlo Simulation
for ir=1:nrlz
  ir1=ir-1;

% positions of point scatterers
  if nc == 1
    xrow=hsl*(erow-2.0*rand(1,N));
    yrow=hsl*(erow-2.0*rand(1,N));
    zrow=hsl*(erow-2.0*rand(1,N));

  else

    %QQ1224848052

    xcmat=diag(hcsl+csl*floor(ncps*rand(1,nc)));
    ycmat=diag(hcsl+csl*floor(ncps*rand(1,nc)));
    zcmat=diag(hcsl+csl*floor(ncps*rand(1,nc)));
    xc=emat*xcmat;
    yc=emat*ycmat;
    zc=emat*zcmat;
    xccol=xc(:);
    yccol=yc(:);
    zccol=zc(:);
    xrow=xccol'+hcsl*(erow-2.0*rand(1,N));
    yrow=yccol'+hcsl*(erow-2.0*rand(1,N));
    zrow=zccol'+hcsl*(erow-2.0*rand(1,N));
  end

% incident
  b=exp(kxi*xrow'+kyi*yrow'+kzi*zrow');

% calculate separation between pairs of point scatterers
  rx=xrow'*erow-ecol*xrow;
  ry=yrow'*erow-ecol*yrow;
  rz=zrow'*erow-ecol*zrow;
  rr=sqrt(rx.^2+ry.^2+rz.^2);
  rr1=rr;
  for j=1:N
    rr1(j,j)=1.0;
  end

% impedance matrix
  Z=-f*exp(i*k*rr)./rr1;
  for j=1:N
    Z(j,j)=1.0;
  end

% exciting field
  Psi=Z\b;

% total scattering amplitude
  FF=f*exp(kxs*xrow+kys*yrow+kzs*zrow)*Psi;
  Pt=(ir1*Pt+abs(FF).*abs(FF))/ir;

% cohenent scattering amplitude
  FFc=(ir1*FFc+FF)/ir;
  Pc=abs(FFc).*abs(FFc);

% incoherent scattered power
  P=Pt-Pc;

% extinction
  Kec=pi*wf'*P/V;
  fprintf(fext,'%5u %9.6e %9.6e %9.6e \n',ir,Kec*wavel,Kei*wavel,Kec/Kei);
  fprintf('\n     realization = %6u \n',ir);
end
% N realizations done
fclose(fext);

% plot
clear;
load ptscext.dat -ascii;
figure;
plot(ptscext(:,1),ptscext(:,4));
xlabel('number of realizations ');
ylabel('ext. coeff. (cluster)/(uniform)');
title('      Convergence Test');
clear;

联系:highspeedlogic

QQ :1224848052

微信:HuangL1121

邮箱:1224848052@qq.com

网站:http://www.mat7lab.com/

网站:http://www.hslogic.com/

微信扫一扫:

作者:不详 来源:网络
  • 您是如何找到本站的?
  • 百度搜索
  • Google搜索
  • 查阅资料过程中
  • 论坛发现
  • 百度贴吧发现
  • 朋友介绍
本站最新成功开发工程项目案例
相关评论
发表我的评论
  • 大名:
  • 内容:
  • matlab代做|matlab专业代做|matlab淘宝代做|matlab代写(www.hslogic.com) © 2019 版权所有 All Rights Reserved.
  • Email:highspeed_logic@163.com 站长QQ: 1224848052

    专业代做/代写/承接、MATLAB、SIMULINK、FPGA项目、博士/硕士/本科毕业设计、课题设计、论文,毕业论文,Coursework、Eassy、Assignment