clear;
f0=1500;
fs=80*f0;
Ns=2048;
N1=4;
r=1;
c=1500;
fai=0;
w=2*pi*f0;
ph1=-180;                   
ph2=180;                   
phn=1;                  
ax=(ph2-ph1)/phn+1;

for n=1:N1;
    for k=1:Ns;
        p(n,k)=exp(j*w*(k-1)/fs-j*w*r*cos(fai*pi/180-(n-1)*2*pi/N1)/c);
        vx(n,k)=p(n,k)*cos(fai*pi/180);
        vy(n,k)=p(n,k)*sin(fai*pi/180);  
    end
%     plot((1:Ns),s)
%     hold on
end

for i=1:ax;     
    num=ph1+(i-1)*phn;
    fain=num*pi/180;
  Xp1=0;
    Xv1=0;
    for n=1:N1;
        Xp(n,1,i)=p(n,1)*exp(j*w*r*cos(fain-(n-1)*2*pi/N1)/c); 
        Xp1=Xp1+Xp(n,1,i);
        Xv(n,1,i)=Xp(n,1,i)*(cos(fai*pi/180)*cos(fain)+sin(fai*pi/180)*sin(fain));
        Xv1=Xv1+Xv(n,1,i);
    end
    Rp(i)=abs((Xp1+Xv1)*Xv1);
end  
i=1:ax;
maxRp=max(Rp);
Rpp=20*log10(Rp(i)/maxRp);
figure
plot(ph1+(i-1)*phn,Rpp);
grid on
hold on
axis([-180 180 -100 0])
xlabel('角度')
ylabel('dB')