点散射方法进行混响仿真，信号为单频信号，散射体高斯分布，无多普勒频域

clc;close all;clear all;

%点散射方法进行混响仿真，信号为单频信号，散射体高斯分布，无多普勒频域

M=100;   %阵元数

N=1000;

deltar=1;   %每次扫描的角度

deltasita=1;      

p=(deltar)*(deltasita)/(2*pi);

f0=25000;

fs=150000;

c=1500;

T=0.02;

t=0:1/fs:T;  %信号持续时间  %1*3001

x=exp(j*2*pi*f0.*t);  %单频信号

r=zeros(M,N);

sita=zeros(M,N);

area=zeros(M,N);

A=zeros(M,N);      %振幅 %大概是产生信号的代码，不太懂   

  for kk=1:N   %N=1000           %背景噪声处理之类的程序？        

r(1,kk)=p*exprnd(1); %exprnd函数，生成服从指数分布的随机数%只有第一行有数据，其余地方没有数据        

sita(1,kk)=2*pi*rand(1);        

%R(1,,ii)=[r(1,ii)*cos(sita(1,ii)),r(1,ii)*sin(sita(1,ii)),1]        

area(1,kk)=2*pi*p*r(1,kk);        

A(1,kk)=sqrt(area(1,kk))*abs(randn(1));    

end    

%处理2-100行数据    

for k=2:M     %100       

  for kk=1:N          %1000            

r(k,kk)=r(k-1,kk)+p*exprnd(1);  %把上一行的加过来 100*1000        

     sita(k,kk)=2*pi*rand(1);            

area(k,kk)=2*pi*p*r(k,kk);            

A(k,kk)=sqrt(area(k,kk))*abs(randn(1));        

end    

end    

figure;plot(A);   %A是振幅，100*1000

%===============================================  

   tao=2*r./c;      %过去再回来的时间 100%1000    

m=1;    

%n=1,nn=1000,t1:1*21;n=2,nn=1000,t1:1*35....n=100,nn=1000,t1:1*3884     %tao一直在变    

for n=1:M        %100       

  for nn=1:N       %1000           

  t1=0:1/fs:tao(n,nn);    %tao是一百行，一千列，行代表阵元数，列代表采样一千个点 1*3436           

  s1=0*t1;           

  t2=tao(n,nn)+1/fs:1/fs:T+tao(n,nn);            

s2=A(n,nn)*exp(j*2*pi*f0.*(t2-tao(n,nn)));         

    s=[s1,s2];            

L(m)=length(s);            

m=m+1;                   

end    

end       

  L_max=max(L);  

   r=zeros(1,L_max);   

  for n=1:M     

    for nn=1:N         

    t1=0:1/fs:tao(n,nn);         

    g1=0*t1;            

t2=tao(n,nn)+1/fs:1/fs:T+tao(n,nn);       

      g2=A(n,nn)*exp(j*2*pi*f0.*(t2-tao(n,nn)));      

       g3=zeros(1,L_max-length(t1)-length(t2));         

    g=[g1,g2,g3];    

         r=r+g;                       

        end    

end

    figure;  

   t3=0:1/fs:(L_max-1)/fs;

    plot(t3,r);