请选择 进入手机版 | 继续访问电脑版

Hello Mat

 找回密码
 立即注册

QQ登录

只需一步,快速开始

查看: 1926|回复: 23

51基于布谷鸟算法的导航波束形成问题

[复制链接]

662

主题

781

帖子

4245

积分

管理员

Rank: 9Rank: 9Rank: 9

积分
4245
发表于 2017-7-9 15:08:16 | 显示全部楼层 |阅读模式
基于布谷鸟算法的导航波束形成问题:
百度网盘链接:
链接:http://pan.baidu.com/s/1hsmpdJ2

具体链接在halcom.cn论坛,联系人QQ:3283892722
该论坛是一个学习交流平台,我会逐一的和大家分享学习。
欢迎大家录制视频,你可在论坛进行打赏分享。
视频专用播放器:http://halcom.cn/forum.php?mod=viewthread&tid=258&extra=page%3D1\

具体的代码如下:
  1. % 布谷鸟函数优化分析
  2. % cuckoo_search -- CS算法
  3. %% 清空环境
  4. clc          % 清屏
  5. clear all;   % 删除workplace变量
  6. close all;   % 关掉显示图形窗口
  7. warning off; % 消除警告
  8. tic
  9. %% 参数初始化
  10. maxiter = 20;    % 迭代次数
  11. sizepop = 20;    % 种群规模
  12. pa=0.25;         % 发现新巢的概率Discovery rate of alien eggs/solutions
  13. %        M   K  phi
  14. popmin = [2, 0, 331, 31, 91, 151, 211, 271];   % x 最小值
  15. popmax = [6, 1, 390, 90, 150, 210, 270, 330];  % x 最大值
  16. nvar = 8;  % 未知量个数

  17. % 导航-波束形成问题
  18. theta_signal=50;    % 俯仰角
  19. phi_signal=200;     % 方位角

  20. % 产生初始粒子
  21. for i=1:sizepop
  22.     nest(i,:)= popmin+(popmax-popmin).*rand(size(popmin)) ;
  23.     fitness(i)= fun( [nest(i,1), nest(i,2)*50000, nest(i,3:nvar)], theta_signal,phi_signal );  % 适应度函数值--目标函数值--求解全局的极大值
  24. end
  25. % 找当前最好的个体
  26. [fmax,bestindex]=max(fitness);
  27. bestnest=nest(bestindex,:);   % 全局最佳
  28. %% 寻优
  29. for i=1:maxiter
  30.     disp(['Iteration ' num2str(i)]);
  31.     % 粒子位置更新
  32.     new_nest=get_cuckoos(nest,bestnest,popmin,popmax); % Levy flights 粒子位置更新
  33.     % 比较更新
  34.     [fnew,best,nest,fitness]=get_best_nest(nest,new_nest,fitness, theta_signal,phi_signal);
  35.     % 发现新巢、随机粒子的位置
  36.     new_nest=empty_nests(nest,popmin,popmax,pa) ;
  37.     % 继续评估新巢的适应度值
  38.     [fnew,best,nest,fitness]=get_best_nest(nest,new_nest,fitness, theta_signal,phi_signal);
  39.     % 寻找全局最优解  
  40.     if fnew>fmax,
  41.         fmax=fnew;
  42.         bestnest=best;
  43.     end
  44.    
  45.     CS_yy(i)=fmax;   % 保存极大值
  46. end  % 结束迭代循环
  47. CS_fitness_iter = CS_yy;
  48. time = toc
  49. save CS_fitness_iter.mat
  50. %% 结果分析
  51. plot(CS_yy,'bo-','Linewidth',2)
  52. title(['适应度曲线  ' '终止代数=' num2str(maxiter)]);
  53. grid on
  54. xlabel('进化代数');ylabel('适应度');
  55. % 结果输出
  56. disp('最优解')
  57. zbest_result = fix( [ bestnest(1),bestnest(2)*50000, bestnest(3:nvar)] );
  58. disp(zbest_result)   % 最佳个体值
  59. fprintf('\n')

  60. %% 结果
  61. tic
  62. [fitness_zbest, F]= fun(zbest_result, theta_signal, phi_signal);
  63. time1 = toc
  64. figure,
  65. plot(0:359, F(:,50)); grid on   %截面图:
  66. xlabel('方位角');
  67. ylabel('增益(dB)');
  68. title('俯仰角50度时截面')
复制代码
寻找新的解:
  1. function new_nest=empty_nests(nest,Lb,Ub,pa)
  2. % A fraction of worse nests are discovered with a probability pa
  3. n=size(nest,1);
  4. % Discovered or not -- a status vector
  5. K=rand(size(nest))>pa;
  6. %% New solution by biased/selective random walks
  7. stepsize=rand*(nest(randperm(n),:)-nest(randperm(n),:));
  8. new_nest= nest+stepsize.*K;
  9. for j=1:size(new_nest,1)

  10.     for k=1:length(nest(j,:))
  11.        if new_nest(j,k)>Ub(k)
  12.            new_nest(j,k)=Ub(k);
  13.        end
  14.        if new_nest(j,k)<Lb(k)
  15.            new_nest(j,k)=Lb(k);
  16.        end
  17.     end
  18.   
  19. end
复制代码
寻找当前最优的解:
  1. function [fmax,best,nest,fitness]=get_best_nest(nest,newnest,fitness, theta_signal,phi_signal)
  2. % Evaluating all new solutions
  3. for j=1:size(nest,1),
  4. %     fnew=Ackley_Fun(newnest(j,:));  % 适应度值
  5.     fnew= fun( [newnest(j,1), newnest(j,2)*50000,newnest(j, 3:end) ], theta_signal,phi_signal );  % 适应度函数值--目标函数值--求解全局的极大值
  6.    
  7.     if fnew>=fitness(j),
  8.        fitness(j)=fnew;
  9.        nest(j,:)=newnest(j,:);
  10.     end
  11. end
  12. % Find the current best
  13. [fmax,K]=max(fitness) ;
  14. best=nest(K,:);
复制代码
位置更新:
  1. % 位置更新
  2. function nest=get_cuckoos(nest,best,Lb,Ub)
  3. % Levy flights
  4. n=size(nest,1);
  5. % Levy exponent and coefficient
  6. beta=3/2;
  7. sigma=(gamma(1+beta)*sin(pi*beta/2)/(gamma((1+beta)/2)*beta*2^((beta-1)/2)))^(1/beta);

  8. for j=1:n,
  9.     s=nest(j,:);
  10.     % implementing Levy flights
  11.     % For standard random walks, use step=1;
  12.     %% Levy flights by Mantegna's algorithm
  13.     u=randn(size(s))*sigma;
  14.     v=randn(size(s));
  15.     step=u./abs(v).^(1/beta);  
  16.     stepsize=0.01*step.*(s-best);
  17.     nest(j,:) = s+stepsize.*randn(size(s));
  18.    % Apply simple bounds/limits
  19.    for k=1:length(nest(j,:))
  20.        if nest(j,k)>Ub(k)
  21.            nest(j,k)=Ub(k);
  22.        end
  23.        if nest(j,k)<Lb(k)
  24.            nest(j,k)=Lb(k);
  25.        end
  26.    end
  27.    
  28. end
复制代码


适应度函数:
  1. function [fitness, F]= fun(x, theta_signal, phi_signal)
  2. % 输入:
  3. %      theta_signal=50;    % 俯仰角
  4. %      phi_signal=200;     % 方位角
  5. %      x:
  6. %         N=7;      % 阵元数目
  7. %         M=6;      % 延迟数目
  8. %         K=5000;   % 快拍数
  9. % 输出:
  10. %   fitness:适应度值,主瓣-旁瓣的差值最大
  11. %      F:输出的beampattern

  12. % x=fix(x); % 取整
  13. N = 7;                   % 阵元数目
  14. M = fix(x(1));                   % 延迟数目
  15. K = fix(x(2));                   % 快拍数
  16. m=[0:M-1]';               

  17. c=3*10^8;                   % 光速
  18. f0=1.26852*10^9;            % 信号中心频率
  19. fs=5*10^9;                  % 采样频率
  20. Ts=1/fs;                    % 单位时延
  21. t = 0:Ts:(K*Ts-Ts);   
  22. wl=c/f0;                    % 信号波长
  23. r=wl/2;                     % 圆阵半径
  24. % phi1=[0, 60, 120, 180, 240, 300];  % 圆周上阵元角度位置 6围1圆阵
  25. % phi1 = 0:360/(N-1):360-360/(N-1);
  26. phi1 = fix( x(3:end) );
  27. P=[zeros(2,1), r*[cos(phi1*pi/180);sin(phi1*pi/180)]]; % 位置

  28. % theta_signal=50;    % 俯仰角
  29. % phi_signal=200;     % 方位角

  30. theta_inter=40;     % 干扰1
  31. phi_inter=50;

  32. theta_inter2=60;    % 干扰2
  33. phi_inter2=250;

  34. theta_inter3=50;    % 干扰3
  35. phi_inter3=60;

  36. theta_inter4=30;    % 干扰4
  37. phi_inter4=120;

  38. theta_inter5=30;    % 干扰5
  39. phi_inter5=80;

  40. waven_signal=[sin(theta_signal*pi/180)*cos(phi_signal*pi/180), sin(theta_signal*pi/180)*sin(phi_signal*pi/180) ];  % 信号波数(wavenumber)

  41. waven_inter=[sin(theta_inter*pi/180)*cos(phi_inter*pi/180), sin(theta_inter*pi/180)*sin(phi_inter*pi/180)];
  42. waven_inter2=[sin(theta_inter2*pi/180)*cos(phi_inter2*pi/180), sin(theta_inter2*pi/180)*sin(phi_inter2*pi/180)];
  43. waven_inter3=[sin(theta_inter3*pi/180)*cos(phi_inter3*pi/180), sin(theta_inter3*pi/180)*sin(phi_inter3*pi/180)];
  44. waven_inter4=[sin(theta_inter4*pi/180)*cos(phi_inter4*pi/180), sin(theta_inter4*pi/180)*sin(phi_inter4*pi/180)];
  45. waven_inter5=[sin(theta_inter5*pi/180)*cos(phi_inter5*pi/180), sin(theta_inter5*pi/180)*sin(phi_inter5*pi/180)];

  46. snr1=-20;                   % 信噪比

  47. inr1=30;                    % 干噪比
  48. inr2=30;                    
  49. inr3=30;
  50. inr4=30;
  51. inr5=30;

  52. s_a0=exp(j*2*pi/wl*(waven_signal*P)');    % 期望信号1
  53. s_a1=exp(j*2*pi/wl*(waven_inter*P)');     % 干扰1
  54. s_a2=exp(j*2*pi/wl*(waven_inter2*P)');    % 干扰2
  55. s_a4=exp(j*2*pi/wl*(waven_inter3*P)');
  56. s_a5=exp(j*2*pi/wl*(waven_inter4*P)');
  57. s_a6=exp(j*2*pi/wl*(waven_inter5*P)');

  58. t_a=exp(j*2*pi*f0*Ts*m);                          

  59. steer_signal=kron(t_a,s_a0);           % Kronecker积,期望及干扰导向矢量

  60. steer_jamming1=kron(t_a,s_a1);
  61. steer_jamming2=kron(t_a,s_a2);
  62. steer_jamming3=kron(t_a,s_a4);
  63. steer_jamming4=kron(t_a,s_a5);
  64. steer_jamming5=kron(t_a,s_a6);

  65. C=[steer_signal];      % 构造线性约束  C'W=G
  66. g=[1]; % 两个目标的无失真约束

  67. amp_signal=sqrt(2*10^(snr1/10))*cos(2*pi*f0*t+randn);     % 信号加上个随机相位

  68. %白噪声干扰,同时也属于宽带干扰
  69. amp_jamming1=sqrt(10^(inr1/10))*randn(1,K);            
  70. amp_jamming2=sqrt(10^(inr2/10))*randn(1,K);
  71. amp_jamming3=sqrt(10^(inr3/10))*randn(1,K);
  72. amp_jamming4=sqrt(10^(inr4/10))*randn(1,K);           
  73. amp_jamming5=sqrt(10^(inr5/10))*randn(1,K);
  74. %amp_jamming1=sqrt(10^(inr1/10))*wgn(1,K,0,'complex');
  75. %amp_jamming2=sqrt(10^(inr2/10))*wgn(1,K,0,'complex');
  76. %amp_jamming3=sqrt(10^(inr3/10))*wgn(1,K,0,'complex');
  77. %amp_jamming4=sqrt(10^(inr4/10))*wgn(1,K,0,'complex');
  78. %amp_jamming5=sqrt(10^(inr5/10))*wgn(1,K,0,'complex');

  79. % 单音干扰
  80. %amp_jamming1= sqrt(2*10^(inr1/10))*sin(2*pi*f0*t+randn);
  81. %amp_jamming2= sqrt(2*10^(inr2/10))*sin(2*pi*f0*t+randn);
  82. %amp_jamming3= sqrt(2*10^(inr3/10))*sin(2*pi*f0*t+randn);
  83. %amp_jamming4= sqrt(2*10^(inr4/10))*sin(2*pi*f0*t+randn);
  84. %amp_jamming5= sqrt(2*10^(inr5/10))*sin(2*pi*f0*t+randn);

  85. noise=(randn(N*M,K));

  86. jam=steer_jamming1*amp_jamming1+steer_jamming2*amp_jamming2+steer_jamming3*amp_jamming3...
  87.    +steer_jamming4*amp_jamming4+steer_jamming5*amp_jamming5;
  88. %jam=steer_jamming2*amp_jamming2;
  89. SS=steer_signal*amp_signal;
  90. X=SS+jam+noise;                             % 接收信号(信号+干扰+噪声)

  91. OO=jam+noise;
  92. Rs = 1/K*SS*SS';
  93. Rni=1/K*OO*OO';
  94. Rn=1/K*noise*noise';
  95. Ri=1/K*jam*jam';
  96. Rr=Rn+Ri;
  97. Rx = 1/K*X*X';
  98. invR=inv(Rx);

  99. W_LCMV=g'*inv(C'*invR*C)*C'*invR;

  100. theta=(0:1:89)*pi/180;         
  101. phi =(0:1:359)*pi/180;
  102. beam = zeros(360, 90);
  103. for ii=1:length(phi)
  104.     tempphi=phi(ii);
  105.     waven_sweep=[(sin(theta).*cos(tempphi))', (sin(theta)*sin(tempphi))'];
  106.     svsweep=exp(j*2*pi/wl*(waven_sweep*P)');
  107.     steer = kron(t_a, svsweep);
  108.     beam(ii,:)=W_LCMV*steer;
  109. end
  110.    
  111. F=20*log10(abs(beam)/max((max(abs(beam)))));    % beampattern

  112. % ref_point = 200;  % 该坐标位置值越大越好
  113. % fitness1 = fun_F(F, 70, ref_point);
  114. % fitness2 = fun_F(F, 60, ref_point);
  115. % fitness3 = fun_F(F, 50, ref_point);
  116. % fitness4 = fun_F(F, 40, ref_point);
  117. % fitness5 = fun_F(F, 30, ref_point);
  118. % fitness = fitness1 + fitness2 + fitness3 + fitness4 + fitness5;

  119. ref_point = phi_signal;  % 该坐标位置值越大越好
  120. fitness = fun_F(F, theta_signal, ref_point);

  121. clear invR jam C c beam amp_jamming1 amp_jamming2 amp_jamming3 amp_jamming4 amp_jamming5 SS
  122. clear steer_jamming1 steer_jamming2 steer_jamming3 steer_jamming4 steer_jamming5
复制代码
funF极值求解函数:
游客,如果您要查看本帖隐藏内容请回复






参考文献:【1】Symmetric Linear Antenna Array Geometry Synthe
【2】THINNING OF PLANAR CIRCULAR ARRAY ANTENNAS USING FIREFLY ALGORITHM

本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有帐号?立即注册

x
程序:算法QQ  3283892722
群智能算法视频售价:38元,链接http://halcom.cn/forum.php?mod=forumdisplay&fid=73
回复

使用道具 举报

0

主题

45

帖子

241

积分

中级会员

Rank: 3Rank: 3

积分
241
发表于 2017-7-22 20:48:29 | 显示全部楼层
太牛了,会这么多算法
回复 支持 反对

使用道具 举报

0

主题

21

帖子

72

积分

注册会员

Rank: 2

积分
72
发表于 2017-9-18 16:39:11 | 显示全部楼层
非常好的代码,谢谢
回复 支持 反对

使用道具 举报

0

主题

50

帖子

112

积分

注册会员

Rank: 2

积分
112
发表于 2017-10-29 23:13:49 | 显示全部楼层
基于布谷鸟算法的导航波束形成问题
回复 支持 反对

使用道具 举报

1

主题

7

帖子

25

积分

新手上路

Rank: 1

积分
25
发表于 2017-11-28 16:58:22 | 显示全部楼层
太棒了,我正找这个
回复 支持 反对

使用道具 举报

0

主题

41

帖子

91

积分

注册会员

Rank: 2

积分
91
发表于 2017-12-1 11:36:24 | 显示全部楼层
学习了,谢谢楼主分享~
回复 支持 反对

使用道具 举报

0

主题

15

帖子

33

积分

新手上路

Rank: 1

积分
33
发表于 2017-12-7 15:36:06 | 显示全部楼层
基于布谷鸟算法的导航波束形成问题,万分感谢
回复 支持 反对

使用道具 举报

0

主题

42

帖子

114

积分

注册会员

Rank: 2

积分
114
发表于 2017-12-21 10:23:03 | 显示全部楼层

学习了,谢谢楼主分享~
回复 支持 反对

使用道具 举报

0

主题

4

帖子

12

积分

新手上路

Rank: 1

积分
12
发表于 2017-12-25 14:47:56 | 显示全部楼层
感谢 研究研究 !!
回复 支持 反对

使用道具 举报

0

主题

7

帖子

55

积分

注册会员

Rank: 2

积分
55
发表于 2018-1-5 10:43:22 | 显示全部楼层
谢谢分享,学习一下!
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则


Python|Opencv|MATLAB|Halcom.cn  

GMT+8, 2019-2-22 09:45 , Processed in 0.137218 second(s), 25 queries .

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.

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