1.软件版本
matlab2015a
2.系统程序
clc;
clear;
close all;
warning off;
addpath 'func\'
addpath 'func\mfiles\'
load data.mat
fs = 100e6;
Nfft = 4096;
fc = 5.8e9;
% [x1,t1,ssf1,yy1] = plotspec(x,1/fs);
% figure;
% plot(ssf1,yy1) % plot magnitude spectrum
% xlabel('frequency');
% ylabel('magnitude') % label the axes
frameNumber = floor(length(x)/Nfft);
txBlockFFT = zeros(frameNumber,Nfft);
%对信号进行短时傅里叶变换,得到跳频图案
for i = 0:frameNumber-1
start = i*Nfft;
txBlockFFT(i+1,:) = fftshift(fft(x(start+1:start+Nfft)));
%找出跳频图案中横轴上每个时间上对应的最大峰值,从时间开始到结束写个循环就行,
[maxValue(i+1),maxIndex(i+1)] = max(abs(txBlockFFT(i+1,:)));
end
%计算峰值均值
Level = 10*mean(maxValue);
delta_f = fs/Nfft;
f = delta_f: delta_f: fs;
f = f - fs/2;
%detect vaalid signal
validIndexCount = 1;
validIndex = zeros(1,1);
validIndex(validIndexCount) = 0;
validFrameCount = 0;
figure;
for i = 1:frameNumber-1
i
%峰值搜索
if(max(abs(txBlockFFT(i+1,:))> Level))
validFrameCount = validFrameCount + 1;
validFrame(validFrameCount) = i+1;
if(abs(maxIndex(i+1) - validIndex(validIndexCount)) > 40)
validIndexCount = validIndexCount +1;
validIndex(validIndexCount) = maxIndex(i+1) ;
detectFHFc = validIndex(2:end)*fs/Nfft- fs/2;
end
detectFHResult(i+1) = maxIndex(i+1)*fs/Nfft- fs/2+fc;
t = (0:length(detectFHResult)-1)*Nfft/fs;
plot(t,detectFHResult,'c*');
hold on
end
axis([0,0.1983,5.74e9,5.86e9]);
pause(0.001);
end
hold on;
y_label = fc*ones(1,length(detectFHResult));
plot(t,y_label,'r');
ylim([fc-50e6 fc+50e6]);
xlabel('时间s')
ylabel('频率Hz');
grid on;
FT = zeros(1000,length(t));
%归一化频率
indx=find(detectFHResult==0);
detectFHResult1=detectFHResult;
detectFHResult1(indx)=[];
detectFHResult2 = detectFHResult-min(detectFHResult1);
indx=find(detectFHResult2<0);
detectFHResult2(indx)=0;
detectFHResult3 = floor(1000*detectFHResult2/max(detectFHResult2));
for i = 1:length(t)
if detectFHResult3(i)>0
FT(detectFHResult3(i),i) = 1;
end
end
SE = strel('disk',24);
FT1 = imdilate(FT,SE);
SE1 = strel('disk',4);
FT1a = imerode(FT1,SE1);
SE1 = strel('disk',8);
FT2a = imerode(FT1,SE1);
SE1 = strel('disk',12);
FT3a = imerode(FT1,SE1);
SE1 = strel('disk',16);
FT4a = imerode(FT1,SE1);
SE1 = strel('disk',20);
FT5a = imerode(FT1,SE1);
SE1 = strel('disk',23);
FT6a = imerode(FT1,SE1);
FTf = FT1a+FT2a+FT3a+FT4a+FT5a+FT6a;
[X,Y]=meshgrid(t,0.001:0.001:1); %横坐标为1953~2006,纵坐标1~54
figure;
contour(X,Y,FTf);
xlabel('time/s');
ylabel('归一化频率f/s');
3.仿真效果
A01-160