fromeigencomposition_lastsolved.m

music高阶谱分析算法研究

        function  [Pmu,search_azimuth,AOA_solved,Pmu_solved,signal_number_juded,big_index]=fromeigencomposition_lastsolved(R,antenna_number,samplepoint_number,antenna_space,lambda_bochang)
                  [eigenvector,D]=eig(R,'nobalance');              %对R进行特征分解。
                   R_eigenvalue=diag(D)';                          %求协方差矩阵R的特征值向量R_eigenvalue
                  [eigenvalue,index]=sort(R_eigenvalue);           %将R_eigenvalue升序排列,Ieig 中前MV个元素含有噪声特征向量信息。
                   R_eigenvalue_decreasing=fliplr(eigenvalue);     %再将R_eigenvalue降序排列。
                   eigenvector=eigenvector(:,index);
                   eigenvector_decreasing=fliplr(eigenvector);
                   
%用MDL准则判断信号源个数signal_number_juded
            for p=antenna_number:-1:1;            
                    a=1.0;
                    b=0.0;
                for i=(antenna_number-p+1):antenna_number;
                       a=a*R_eigenvalue_decreasing(i);
                       b=b+R_eigenvalue_decreasing(i);
                end;
                     a=a^(1.0/p);
                     b=b/p;
                     lr(p)=(a/b)^(samplepoint_number*p/2);
                     MDL(p)=-log(lr(p))+0.25*(antenna_number-p)*(antenna_number+p+1)*log(samplepoint_number);
            end;
                 [eigenvalue_min ,location]=min(MDL);
                  signal_number_juded=antenna_number-location;
                  En=eigenvector_decreasing(:,signal_number_juded+1:antenna_number);         %得到噪声子空间En
                  
%产生空间搜索矩阵search_matrice
                     search_number=501;
                     search_azimuth=linspace(-25,25,search_number);
               for   search_row=1:antenna_number
                     for   search_column=1:search_number
                           search_phase=(search_row-1)*2*pi*antenna_space/lambda_bochang*sin(search_azimuth(search_column)*pi/180);
                           search_matrice(search_row,search_column)=exp(-j*search_phase);
                     end;
               end;

%求空间谱函数值Pmu/dB   
               for   search_column=1:search_number   
                     Pmu(search_column)=search_matrice(:,search_column)'*search_matrice(:,search_column)/(search_matrice(:,search_column)'*En*En'*search_matrice(:,search_column));
               end;
               Pmu=10*log10(abs(Pmu));

%求空间谱函数中的峰值Pmu_solved，以确定信号源入射方位角AOA_solved
               big_index=0;
               for   Pmu_index=2:search_number-1;        
                     if   Pmu(Pmu_index)>Pmu(Pmu_index-1)&Pmu(Pmu_index)>Pmu(Pmu_index+1);
                          big_index=big_index+1;
                          Pmu_big(big_index)=Pmu(Pmu_index);
                          AOA_judged(big_index)=search_azimuth(Pmu_index);
                     end;
               end;
                  [Pmu_big_increasing,big_index_increasing]=sort(Pmu_big);   
                  if   big_index-signal_number_juded+1>=1
                       big_index_signal=big_index_increasing(big_index-signal_number_juded+1:big_index);
                       AOA_solved=AOA_judged(big_index_signal);
                       Pmu_solved=Pmu_big_increasing(big_index-signal_number_juded+1:big_index);
                   else
                       AOA_solved=AOA_judged(big_index_increasing);
                       Pmu_solved=Pmu_big_increasing;
                   end;