声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 3876|回复: 12

[声学基础] 刚性球对平面波脉冲的散射过程演示

[复制链接]
发表于 2009-2-4 15:27 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
声学二维谱能很好的模拟物体的脉冲声散射过程
能仿真声散射试验声光显示技术方法的效果
这里给出一个简单的例子,供参考。
下面给出主程序,说明思路:
结果见附件,可以看到散射波清晰可见。
clear;
clc;
ka=3;
a=0.5;
k=ka/a;
x=-4:0.1:4;              %采用直角坐标系
y=-4:0.1:4;
n=length(x);
m=length(y);
p=zeros(m,n);
df=300000/1024;
fo=5;
b=5;
omnume=128;
for iso=1:omnume
    sou(1,iso)=exp(-(iso*df/1000-fo)^2/b^2);
end   
f=1:omnume;
f=f*df;
t=1:omnume;
ka=2*pi*f/1500*a;
for mm=1:m
    for nn=1:n
        r=sqrt(y(mm)^2+x(nn)^2);
             ps=0;
         if (r>a)
           theta=acos(x(nn)/r);
           for l=0:80
                Pl=legendre(l,cos(theta));
                Pl=Pl(1);
                ps=ps-(2*l+1)*i^l*dsphericalbessel(l,ka)./dsphericalbesselh(l,1,ka)*sphericalbesselh(l,1,k*r)*Pl;
           end
           frek=sou.*(exp(i*2*pi*f/1500*(x(nn)+1.0))+ps);
           fre=[frek,0,conj(fliplr(frek))];
           frek=fre(1,1:2*omnume);
           temp=fft(frek);
           p(mm,nn)=temp(80);
         end
%*********************************
    end
end   
pcolor(x,y,abs(p));
colormap(gray);
shading interp;
xlabel('x');
ylabel('y');
感兴趣的朋友可参见张海澜的文章
望大家多多交流,共同进步!
sphere_scattering.jpg

评分

2

查看全部评分

回复
分享到:

使用道具 举报

发表于 2009-2-12 10:53 | 显示全部楼层

请教几个函数

dsphericalbessel(l,ka)是不是l阶球贝塞尔函数,是不是等于sqrt(pi/2/ka)*besselj(l+0.5,ka);dsphericalbesselh(l,1,ka)是不是第一类l阶球汉克尔函数,是不是等于besselj(l,ka)+i*bessely(l+0.5,ka);sphericalbesselh(l,1,k*r)是不是也是第一类l阶球汉克尔函数。

[ 本帖最后由 evans_xu 于 2009-2-12 10:56 编辑 ]
发表于 2009-2-12 16:50 | 显示全部楼层

回复 沙发 evans_xu 的帖子

呵呵,那两个函数分别是球贝塞尔函数的一阶导数、第一类球汉克尔函数的一阶导数.
 楼主| 发表于 2009-2-12 16:53 | 显示全部楼层
对的,下面是具体函数:
function y=sphericalbessel(n,z)
y=sqrt(pi/2./z).*bessel(n+0.5,z);
function y=dsphericalbessel(n,z)
y=(n*sphericalbessel(n-1,z)-(n+1)*sphericalbessel(n+1,z))/(2*n+1);
function y=sphericalbesselh(n,k,z)
y=sqrt(pi/2./z).*besselh(n+0.5,k,z);
function y=dsphericalbesselh(n,k,z)
y=(n*sphericalbesselh(n-1,k,z)-(n+1)*sphericalbesselh(n+1,k,z))/(2*n+1);
发表于 2009-2-13 13:08 | 显示全部楼层

回复 地板 Sommerfeld 的帖子

我用你给的程序计算了,结果跟你的不一样啊,是不是我有哪些地方没注意到啊?
发表于 2009-2-17 09:23 | 显示全部楼层
楼主指点一下啊,小弟很想学习一下,先谢了:handshake
发表于 2009-2-22 15:36 | 显示全部楼层
有没有人用楼主的程序计算过啊,请高手出来指点一下,谢谢了!
发表于 2009-3-2 23:26 | 显示全部楼层
确实,我也用了你给的程序算了,可是算的不对。你能留给邮箱,以便我向你请教下吗?非常感谢!
发表于 2009-3-3 09:36 | 显示全部楼层
早上试跑下, 需要跑很久, 结果与附图好像的确差异很大!
具体不懂, 所以帮不忙!
发表于 2010-11-17 02:42 | 显示全部楼层
回复 1 # Sommerfeld 的帖子

F:\Temporary Files\MATLABuntitled.jpg
这是我安装楼主的程序跑出的结果,不懂~~~
发表于 2012-2-27 15:47 | 显示全部楼层
语句中先关的注释不多,直接读程序语句理解程序不太方便。
发表于 2012-4-22 12:10 | 显示全部楼层
我也运行了,但没有得到理想的结果,楼主的程序还是有些不懂
发表于 2017-10-31 15:36 | 显示全部楼层
将近十年前的帖子了,不知楼主还在不在。
参考文献能否给点看看。
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-6-10 02:04 , Processed in 0.067917 second(s), 24 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表