单通道跳频信号分选的子函数，包括实现同步正交和异步非正交跳频网台信号分选附matlab代码.zip

%% 单通道信号分选（异步非正交） %%%%%%%%%%%%%%%%%输入%%%%%%%%%%%%%%%%%%%%% %OriData:单通道跳频信号 %dataout:生成阵列信号的跳频序列（主要用于后续准确率的计算） %NF:频点数 %fs:采样率 %Rs:跳速 %f1-f2:载频范围 %Span:STFT的窗长 %Step:STFT的步长 %%%%%%%%%%%%%%%%%%输出%%%%%%%%%%%%%%%%%%%%% %TEST:分选的正确率 %R:估计的跳速 %Number_est:估计的网台数 %Type:正交/非正交 %HT_Number_est:频率集/起跳时间/结束时间 %Bandwidth:跳频信号带宽 %Duty_ratio:驻留时间的占空比 %HT_Result:画分选后的跳频图案的 %HT_Mosaic:跳频拼接结果 %Dataout_sorting:分选后的跳频序列 %Xishu:画三维时频图 function [TEST,R,Number_est,Type,HT_Number_est,Bandwidth,Duty_ratio,HT_Result,HT_Mosaic,Dataout_sorting,Xishu] = Single_channel_sorting(OriData,dataout,NF,fs,Rs,f1,f2,Span,Step) Rd = round((f2-f1)/(NF-1)/3); period = size(dataout,2); %% 滤波 b = fir1(512,abs(f2-f1)/fs*2); OriData1 = OriData.*exp(-1i*2*pi*rem(f1,fs)*(0:length(OriData)-1)/fs); y = filter(b,1,OriData1); y = y.*exp(1i*2*pi*rem(f1,fs)*(0:length(y)-1)/fs); clear OriData; OriData = y; clear y b; %% STFT Xishu=abs((STFT(OriData,Span,Step))); %% 画频谱图 Len=floor((length(OriData)-Span)/Step); Y=(1:Span)/Span-1/2; %归一化频率，相当于fs=1 XX=(1:Len)*Step/fs; %归一化时间 Z=Xishu/max(max(Xishu)); Z = Z(1:Len,:); figure;contour(XX,Y,Z.'); ylabel('归一化频率f/fs'); xlabel('时间t/s'); title('分选前的信号等高线图'); axis([0 0.01 0 0.5]); %% 聚类每个窗内的频率 deta=(f2-f1)/(NF-1); [f,c]=Single_Cluster1(Xishu,Span,fs,deta,f1,f2); %% 信源个数估计 est = sum(c~=0,2); T_est = sortrows(tabulate(nonzeros(est)),2); if T_est(find(T_est(:,1)==size(c,2)),3) <= 5 T_est_2 = T_est(end-2:end,:); % Number_est = mode(est); Number_est = max(T_est_2(:,1)); else Number_est = size(c,2); end c = c(:,1:Number_est); %% 估计跳变时刻 min_interval=deta*Span/fs; PRI = fs/Rs/Step; [HT,TR]=Single_HT_estimation(fs,c,Span,Step,deta,PRI); K = round(min(f1,f2)/fs)-1/2; HT_test = HT; HT(1,:)=round((HT(1,:)+Span*K-f1/fs*Span)/min_interval); %转换成0-255频点 % HT_test(1,:)=round((HT_test(1,:)+Span/2*K-f1/fs*Span)/min_interval)*deta+f1; %转换成频率 HT_test(1,:) = (HT_test(1,:)+Span*K)*fs/Span; Type = []; if abs(HT(2,1)-HT(2,2)) > 2 %正交非正交的判断 Type = 1; else Type = 0; end %% 精估跳变时刻 if Type == 1 [FH_HT_2_3] = Jump_time_estimation_FIR(OriData,HT_test,Rs,Rd,fs,Step); HT(2:3,:) = FH_HT_2_3/Step; HT_test(2:3,:) = round(HT(2:3,:)); if HT(2,1) <= 0 HT(2,1)=1; end for i = 1:length(HT) if HT(2,i) <= 0 HT(:,i) = 0; end end for i = 1:length(HT) if HT(1,i) < 0 HT(:,i) = 0; end end HT=HT'; HT(all(HT==0,2),:)=[]; HT=HT'; %% 分选 toa_data=[]; %将此变量清空的原因是重新赋值时维度不够会报错 toa_data(:,1)=HT(2,:).'; toa_data(:,2)=1:length(toa_data(:,1)); PRI=fs/Rs/Step; % Error_toa=[0.07 0.07 0.07 0.08 0.08 0.095 0.095 0.10]; Error_toa=ones(1,Number_est)*0.05; Error_lin=ones(1,Number_est)*0.06; %参数可以调节 极大了影响了正确率 % Error_lin=[0.07 0.07 0.07 0.08 0.08 0.095 0.095 0.10]; seq_out =[]; seq_out_temp = []; seq_rest = toa_data; PRI_end = []; PRI_end_temp =[]; for i=1:Number_est [seq_out_temp,num_out,seq_rest,num_rest,PRI_end_temp,PRI_jieju] = fun_Sequence_retrieve(seq_rest,PRI,Error_toa(i),Error_lin(i)); seq_out = [seq_out seq_out_temp]; PRI_end = [PRI_end PRI_end_temp]; end % 加判断以免报错 % if size(seq_out,2) ~= size(dataout,1) % TEST=zeros(1,Number_est);R=0; Type=-1; % HT=zeros(3,period*Number_est); Bandwidth=0; Duty_ratio=0; % HT_Result=0; Dataout_sorting=zeros(Number_est,period); % Xishu=0; HT_Mosaic=0; test_Dataout_sorting=zeros(Number_est,period); c=0; % HT_Number_est = zeros(3*Number_est,10); % disp('无法分选'); % return; % end %% 判断是否从第一个周期开始，排序 %去零处理 for i = 1:size(seq_out,2) for j = 1:size(seq_out,1) if seq_out(1,i)>0 break; else seq_out(1:end-1,i)=seq_out(2:end,i); end end end seq_out=[HT(2,seq_out(1,:));seq_out]; for i=1:size(seq_out,2) % if max(seq_out(2,:)) < Number_est % break; % else if seq_out(1,i) >= PRI*0.9 b=(seq_out(1,i)/(PRI)); if b >=0.8 && b < 1 b = 1; end b = floor(b); seq_out(1,i)=seq_out(1,i)-b*PRI; if seq_out(1,i) < 1 seq_out(1,i) = 1; end end % end end seq_out=seq_out'; seq_out=sortrows(seq_out,1); seq_out=seq_out'; seq_out=seq_out(2:end,:); seq_out=seq_out'; %% 估计正确率 seq_out_test = zeros(Number_est,500); seq_out_test(1:size(seq_out,1),1:size(seq_out,2)) = seq_out; TEST = []; test_Dataout_sorting = zeros(Number_est,period);%含错标记的跳频序列 Dataout_sorting = zeros(Number_est,period);%分选后的跳频序列 % if dataout(1,1) >= 0 for i = 1:Number_est if i > size(dataout,1) | i > size(seq_out,1) break; end Dataout_sorting(i,1:length(nonzeros(seq_out_test(i,1:end)))) = HT(1,nonzeros(seq_out_test(i,1:end))); [rate,test_HT] = compute_rate(dataout(i,1:period),HT(:,nonzeros(seq_out_test(i,:))),PRI); TEST = [TEST rate]; test_Dataout_sorting(i,1:length(test_HT)) = test_HT; end % else % TEST = zeros(); % test_Dataout_sorting = ; % end % plot(1:1:length(HT(2,nonzeros(seq_out(1,1:end)))),HT(2,nonzeros(seq_out(1,1:end))),'r') Duty_ratio = mean(abs(HT(2,Number_est:end)-HT(3,Number_est:end)))/(mean(PRI_end)); % Duty_ratio = strcat(num2str(Duty_ratio*100),'%'); R=fs/(mean(PRI_end)*Step); %跳速估计 HT_Result = cell(1,Number_est); %画时频图 第一行为频率(f1-f2),第二、三行为Span上的采样点数 for i = 1:Number_est HT_Result{i} = HT_test(:,nonzeros(seq_out_test(i,:))); end %% 跳频拼接 Input_Duty_ratio = mean(abs(HT(2,Number_est:end)-HT(3,Number_est:end)))/(mean(PRI_end)); c_seq = Dataout_sorting; Nh = round(fs/Rs); f_seq = c_seq*deta+f1; time_delay = round(abs(HT(2,1)-HT(2,2)))*Step;HT_Mosaic =1; % HT_Mosaic = signal_mosaic(fs,Rs,time_delay,f1,f2,Rd,f_seq(:,1:10),Input_Duty_ratio); HT(1,:) = (HT(1,:)*deta+f1)/1e6; %/1e6即表示频率的单位是MHz HT(2:3,:) = HT(2:3,:)/fs*1e3; %*1e3即表示时间的单位是毫秒(ms)、*1e6即表示时间的单位是微秒(us) Bandwidth = max(HT(1,:))-min(HT(1,:)); % 分选后跳频网台的频率及起跳时间 for ni = 1:size(dataout,1) nll(ni) = (length(nonzeros(seq_out_test(ni,:)))); end nl = min(nll); if nl >= 10 nl = 10; end HT_Number_est = zeros(3*size(dataout,1),nl); % for nn = 1:size(dataout,1) % HT_Number_est(1+(nn-1)*3:nn*3,:) = [f_seq(nn,1:nl)/1e6;HT_test(2:3,nonzeros(seq_out(nn,1:nl)))*Step/fs*1e3]; % end Dataout_sorting = test_Dataout_sorting; if Type == 1 Type = '非正交'; elseif Type == 0 Type = '正交'; end elseif Type == 0 && rem(Number_est,2)~=0 HT(2:3,:) = round(HT(2:3,:)); p=1; beg(p) = 1; flag = 0; for i=2:size(HT,2) if abs(HT(2,i)-HT(2,i-1)) >= fs/Rs/Step*0.7 endl(p)=i-1; p=p+1; beg(p)=i; flag = 0; else flag = flag+1; end if flag > Number_est-1 endl(p)=i-1; p=p+1; beg(p)=i; flag = 0; end if i==size(HT,2) endl(p)=size(HT,2); break; end end pos=endl-beg+1; for i=1:length(pos) if pos(i)<round(Number_est/2) pos(i)=0; beg(i)=0; endl(i)=0; end end pos=nonzeros(pos); beg=nonzeros(beg); endl=nonzeros(endl); pos_save = pos; %% 序列分选 HT_seq = HT; PRI=fs/Step/Rs; Seq_h = zeros(Number_est,length(HT)); for i = 1:length(HT_seq) pos = find(nonzeros(HT_seq(2,:))<=PRI*0.9+PRI*(i-1)); HT_seq(2,pos) = 2*max(HT_seq(2,:)); if length(pos) > Number_est Seq_h(1:Number_est,i) = HT_seq(1,pos(1:Number_est)); else Seq_h(1:length(pos),i) = HT_seq(1,pos); end end Seq_h = Seq_h(:,1:size(dataout,2)); Seque