电子情报侦察分析综合实验

一、实验目的

领会专业的特色，掌握电子情报侦察分析的基本流程，理解雷达原理与雷达截获分析课程的联系。

二、实验要求

1、工具采用MATLAB

2、对实验结果进行理论分析，并给出自己的见解

三、实验计划

1、拟将18名学员分成两组，进行电子情报模拟对抗

2、每名学员编程产生雷达信号，汇集给教员，由教员随机配对，对产生的信号进行分析。

四、实验步骤和内容

4.1产生不同体制的雷达信号

1、实验原理

加强对简单脉冲信号、频率分集信号、频率捷变信号、重频参差信号、PRI抖动和PRI跳变信号、PRI滑变信号、PRI排定信号、脉组PRI变化信号、双脉冲信号、脉冲压缩信号等信号的直观认识；使学员从多角度描述雷达信号。

常见的雷达信号形式内容：

（1）简单脉冲信号

指载频、脉冲重复频率和脉冲宽度等三个参数均固定不变的信号

（2）频率分集信号

指一部雷达同时或在相隔很短的时间内发射多个载频以完成同一任务。频率分集信号中，不同载频的脉冲可以同时发射，也可以顺序发射。频率分集信号的频谱是各个不同频率脉冲频谱的合成。

（3）脉冲重复间隔（PRI）变化信号

PRI变化信号脉冲指脉冲的PRI在相邻脉冲之间或脉组之间变化的信号。

①重频参差信号。重频参差信号指PRI在相邻脉冲之间或脉组之间变化的信号

②PRI抖动信号和PRI跳变信号

PRI抖动指相邻脉冲的PRI在一定范围内抖动，即相邻的脉冲时间间隔不相等。一般PRI抖动的范围小于中心值的5%。

PRI跳变信号是指PRI在其平均值附近的有限几个点上按照某一分布规律跳动，其跳动范围可高达平均值得30%。

PRI抖动信号与PRI跳变信号的差别是，前者可以在无限个点上变化，后者只能在有限个点上变化；前者变化范围小，后者变化范围大。

③PRI滑变信号

这种信号的特点是PRI值周期性单调增加或下降，在两个极值之间周期性变化，PRI滑变信号可以消除盲距，还可以提供恒定高度范围内的仰角扫描，达到最佳性能。在这种应用中，最小PRI的值与最大的值之比大致等于仰角扫描系统的最小作用距离与最大作用距离之比。通常这个值为3~6。

④PRI排定信号

这种信号的特点是PRI值按排定的多个数值周期性的变化，其PRI数值一把为十几至二十几点。PRI排定信号一把用于多功能雷达。

⑤脉组PRI变化

脉组PRI变化信号是指信号在PRI为T1值上发射几个脉冲，然后又在PRI为T2值上发射几个脉冲，并且具有周期性。每组的脉冲数一般不少于四个。

（4） 双脉冲信号

这是指在每个PRI内发射两个脉冲，这两个脉冲的宽度可以相同，也可以不同，两脉冲间的间隔通常为脉冲宽度的几倍。

（5） 脉冲压缩信号

是指单脉冲内的频率呈线性变化。

2、实验要求

（1）、绘制出各种形式的雷达信号；

（2）、频率变化的信号要频率分析各雷达信号；

（3）、总结该雷达信号特征，说明选用该波形的目的以及这样波形的优缺点；

（4）、给出源代码。

4.2雷达信号截获与定位实验

1、实验原理：

不同雷达发射不同体制雷达信号，经电磁空间传播被雷达信号截获接收机接收，本实验，根据侦察方程模拟该过程。

仿真3~5部雷达信号，体制分别为简单体制、重频参差、脉压信号（线性调频）、频率捷变等，脉冲个数至少有100个；雷达信号统一为1；雷达天线增益函数为，设每部雷达波束20度，扫描速度10转/分；不考虑损耗与大气传播影响；

雷达信号截获接收机有1~3个，且同时截获雷达信号；

截获接收机天线采用全向天线；

2、实验要求：

（1）、仿真出雷达信号截获接收机截获的雷达信号；

（2）、给出雷达信号参数（载频、脉宽、到达时刻、脉幅）测量方法，并给出测量值；

（3）、根据截获的雷达信号，给出相应的定位方法，并给出位置信息。

4.3雷达信号分选实验

1、实验原理

对截获接收机的PDW输出流进行去交错处理，即雷达信号脉冲分选。按照分选采用的参数的不同分为多参数分选与PRI分选等。

多参数分选即按照到达方向、载频、脉宽等参数进行分选，也可以看作是多参数聚类分析。常用的聚类算法有K均值、最大最小方法以及神经网络方法等。

PRI分选即在多参数分选后，按照雷达脉冲到达时刻信息，获取脉冲间隔信息，然后采用统计、变换域等方法进行PRI分选。

2、实验内容

（1）选取任一种聚类方法，完成多参数雷达信号多参数分选；

（2）采用CDIF等PRI分选方法，完成PRI分选，重点观察重频参差信号的分选效果。

3、实验要求

（1）以实验二的结果进行该实验。

五、实验结果及分析

5.1产生不同体制的雷达信号

1.简单脉冲信号

产生的简单脉冲如图5.1所示。

图5.1

参数设置如下：

2.线性调频脉冲信号

产生的脉冲如图5.2所示。

图5.2

参数设置如下：

3.重频参差脉冲信号

产生的脉冲如图5.3所示。

图5.3

参数设置如下：

4.重频参差脉冲信号

产生的脉冲如图5.4所示。

图5.4

参数设置如下：

1.仿真出雷达信号截获接收机截获的雷达信号；

（1）仿真两个雷达接收机，其参数如下：

一号接收机：（50km ,62km） 天线有效面积：1

二号接收机：（1.5km ,6km） 天线有效面积：1.2

在信号接收过程中，由于四部雷达的距离不同，会产生不同的时延，且在传输过程中会夹杂着噪声，因此，在仿真接收信号中加入了时延和高斯噪声（信噪比为20，信噪比不能太大，否则有可能噪声将信号淹没，无法进行检测）。

加时延的代码如下：

c=3e8;

shiyan1=rand(1,4)/1e3;

for i=1:4

shiyan2(i)=shiyan1(i)+(R2(i)-R1(i))/c;

end

加噪声的代码如下：

SNR=20;

snr=1/(10^(SNR/10));

Y11=[zeros(1,fs*shiyan1(1)) Y1];len11=length(Y11);

noise11=sqrt(snr/2)*(randn(1,len11));

Y_11=Y11+noise11;

在接收信号时，由于四种信号的点数不同，需要以最长的信号为准，在其他信号后补零，以保证正确的叠加。

（2）两个接收机收到的雷达信号分别如图5.5和图5.6所示。

图5.5

图5.6

2． 给出雷达信号参数（载频、脉宽、到达时刻、脉幅）测量方法，并给出测量值。

（1）从接收到的信号中提取脉冲

思路步骤：先设定一个门限值，然后对信号进行判别，如果大于门限，则判为起始位置，若小于门限，则判为结束位置。其中需要注意的是：

A．门限的设置必须合理；

B．在判断起始位置和结束位置时，需要考虑脉内小于或大于门限的点，因此，通过将后边200个点求和，设定一个门限，判断是否为起始和结束位置。

为了判断门限更简便，首先将信号二值化处理，然后比较容易的设定门限。

信号二值化代码如下：

as=zeros(1,length(s));

for i=1:length(s)

if s(i)＞=T %%如果信号幅度大于等于门限，则判为1

as(i)=1;

else

as(i)=0;%否则判为0

end

end

求脉冲的起始位置和结束位置的代码如下：

for i=2:length(s)-300

summ=sum(as(i:i+300));

if as(i)==1&&count==0

if summ＞=220

b(m,1)=i; % b矩阵的第一列存放脉冲起始位置

m=m+1;

count=1;

end

end

if count==1&&summ＜40

b(n,2)=i; % b矩阵的第二列存放脉冲结束位置

n=n+1;

count=0;

end

end

实验结果如图5.7所示：测得脉冲个数为212个

图5.7

（2）求脉宽

思路步骤：上一步将脉冲的起始位置和结束位置存放在一个矩阵b中，脉宽的测量只需用结束位置减去起始位置，求出点数，然后除以采样频率fs，即可得到脉宽。

求脉宽的代码如下：

a=b(:,2)-b(:,1); %a矩阵存放脉宽

figure;

stem(a);

axis([0,200,0,4000]);

实验结果如图：

图5.8

通过观察可看出，测得的脉宽大量集中于6us，10us，20us，30us，还有部分零散分布在30—40us，与了雷达发射信号的脉宽近似符合，测量结果比较准确。

（3）求载频

思路步骤：求载频时，将每个脉冲做傅里叶变换，从中提取出最大幅度的即为载频位置，然后根据公式fc=k/N*fs得到载频。

需要注意的是，在傅里叶变换后，会有较大的直流分量的影响，因此需要将直流分量去除，然后再求最大值。

求载频的代码如下：

N=10000;

for i=1:length(b)

maichong=[s(b(i,1):b(i,2))];

maichong_FFT=abs(fft(maichong,N));

maichong_FFT(1:60)=0;

maichong_FFT(9940:10000)=0;

[Fcmax,k]=max(maichong_FFT);

Fc(i)=k/N*fs; %Fc矩阵存放载频

End

实验结果如图5.9所示：

图5.9

通过观察可看出，测得的载频大量集中于5MHz，10MHz，15MHz，20MHz，还有部分零散分布在20—25MHz，与了雷达发射信号的载频近似符合，测量结果比较准确。

5.3雷达信号分选实验

思路步骤：根据测得的脉宽和载频，利用K均值聚类法，对测得的脉冲进行聚类分选，其难点在于K均值聚类。

K均值聚类的思想是，先人为选定四个与真实值较为接近的聚类中心，然后设定四个聚类种类，求出其余各点与聚类中心的距离，将这些距离求均值，如果此均值与原聚类中心的距离小于一个门限值，则说明此聚类中心比较合适，聚类结束；若此均值与原聚类中心的距离大于门限，则将此均值当做新的聚类中心，从头开始，在进行聚类，直到聚类结束，最终得到四类。

实验中较大的难点在于门限值的设置，需要根据每次得到的数据进行分析，设定一个合适的门限。

实验结果：

通过上表可以看出，四种体制的脉冲信号的脉宽和载频的测量值与实际值在误差允许的范围内较为精确，基本能够分选出四种体制信号的参数；但脉冲个数的聚类结果与实际值相差较大。

脉冲个数误差较大的原因主要有：a. 脉冲的重叠，产生很大脉宽的脉冲，在筛选脉冲时将其舍弃;

b. 聚类时，由于误差的存在，部分脉冲与聚类中心距离较大，在门限值之外，因此将其舍弃。

六、实验心得

本次实验是一个比较综合的实验，是对以前学习的五门课程知识的综合总结。虽然最终没有完成全部的实验内容，但通过完成的实验内容，我感觉收获已经很多。

实验中产生四种不同体制的雷达信号，这在电子战实验中已经做过了，但需要产生100个脉冲，因此需要理解脉冲重复间隔和时间的关系，根据采样定理，还需要使采样率大于2倍的载频。在这部分中，我对各种体制雷达信号的原理有了一定掌握，之前对于这些信号调制原理根本不懂。此外，在测量脉宽和载频时，加深了对信号点数与脉宽、载频之间的转换关系。

通过实验，我收获最大的是使用matlab编程的能力有了明显提高，这个实验的程序应该是我大学四年编程最长最综合的一次，在编程过程中不断理解原理，思考编程思路，如何将原理转变成程序来实现，不断调试程序、修改程序，最终得到较为满意的实验结果，这为我下学期的毕业设计打下了较好的基础。

实验过程中，我得到了教员以及一些能力水平比较高的同志的帮助和指导，可以看出自己与某些同志之间的差距，这对我的帮助也是很大的，在感谢教员和同志们的同时，我也会从整个实验过程中总结经验，在今后的学习中弥补不足。

七、源代码附录

%% %%%%%%%%简单脉冲雷达信号%%%%%%%%

clc

clear all

close all

PW=6e-6; %脉宽

fs=10e7; %采样频率

fc=10e6; %载频

prt=30e-6;%脉冲重复间隔

pri=1/prt; %脉冲重复频率

t=0:1/fs:100*prt;

y=cos(2*pi*fc*t-0.5);

duty=100*PW/prt; %占空比

s=(square(2*pi*pri*t,duty)+1)/2; %产生方波

Y1=y.*s;

plot(t*1e6,Y1);

xlabel('t/us');ylabel('幅度');title('简单脉冲信号');

%% %%%%%%%%%%%线性调频雷达信号%%%%%%%%%%

pw=20e-6; %发射脉宽

B=10e6; %调频带宽

% fs=10e7; %采样频率

fc=20e6; %载频

prt=70e-6;

pri=1/prt; %脉冲重频

t=0:1/fs:prt;

f0=rectpuls(t-pw/2,pw);

Fc=f0.*exp(j*2*pi*(fc*t+B/pw*t.*t));

% figure

% plot(abs(fft(Fc)));

% figure;

% plot(t*1e6,real(Fc));title('线性调频信号');xlabel('t/us');

Y2=repmat(Fc,1,100);

figure;

plot(real(Y2));title('线性调频脉冲发射信号');

%% %%%%%%%%%PRI参差%%%%%%%%%

fc=5e6; %载频

tao=10e-6; %脉宽

prt1=100e-6;

prt2=200e-6;

prt3=300e-6; %脉冲间隔

% fs=3*fc; %采样频率

prt=prt1+prt2+prt3; %参差脉冲间隔

pri=1/prt; %参差脉冲重频

t=0:1/fs:33*prt;

zhankong=tao/prt*100; %占空比

x=sin(2*pi*fc*t);

y =(square(2*pi*pri*t,zhankong)+1)/2;

s1=(square(2*pi* pri *(t-prt1), zhankong)+1)/2;

s2=(square(2*pi* pri *(t-prt1-prt2), zhankong)+1)/2;

s3=(square(2*pi* pri *(t-prt1-prt2-prt3), zhankong)+1)/2;

Y3=x.*s1+x.*s2+x.*s3;

figure;

plot(t*1e6,abs(Y3));title('PRI参差信号');

xlabel('t/us');ylabel('幅度');

%% % %%%%%%%重频抖动信号%%%%%%

fc=15e6; %载频

pw=30e-6;%脉宽

prix=(rand(1,4)*0.5*50+50)*1e-6;

pri1=prix(1);pri2=prix(2);pri3=prix(3);pri4=prix(4);

pri=pri1+pri2+pri3+pri4;

prt=1/pri;

zhankong=pw/pri*100; %占空比

t=0:1/fs:25*pri;

x=sin(2*pi*fc*t);

y1=(square(2*pi*prt*t,zhankong)+1)/2;

y2=(square(2*pi*prt*(t-pri1),zhankong)+1)/2;

y3=(square(2*pi*prt*(t-pri1-pri2),zhankong)+1)/2;

y4=(square(2*pi*prt*(t-pri1-pri2-pri3),zhankong)+1)/2;

y5=(square(2*pi*prt*(t-pri1-pri2-pri3-pri4),zhankong)+1)/2;

Y4=x.*y2+x.*y3+x.*y4+x.*y5;

figure;

plot(t*1e6,abs(Y4));

title('PRI抖动雷达脉冲信号');

xlabel('t/us');ylabel('幅度');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% 雷达信号截获

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%4

%%%设定雷达和接收机的坐标

% clc;close all;clear all;

leidax=[1e4 2.3e4 3.5e4 7e4];

leiday=[4e4 1.8e4 4.5e4 3e4];

jieshoux=[5e4 6.2e4];

jieshouy=[1.5e3 6e3];

% stem(leidax,leiday);hold on;

% stem(jieshoux,jieshouy,'r');

% axis([0,1e5,0,5e4]);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% 计算接收机与雷达的距离

%%%%%%%%与接收机1的距离%%%%%%

for i=1:4

R1(i)=sqrt((leidax(i)-jieshoux(1))^2+(leiday(i)-jieshouy(1))^2);

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%与接收机2的距离%%%%%%

for i=1:4

R2(i)=sqrt((leidax(i)-jieshoux(2))^2+(leiday(i)-jieshouy(2))^2);

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%

Pt=[1e7 2e7 3e7 4e7];%各雷达的发射功率

Ar=[1 1.2];%两个接收机的有效接收面积

Gt=[30 31 32 33];

Gtt=10.^(Gt./10); %各雷达天线增益

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% 接收机的接收功率

for i=1:4

Pr1(i)=(Pt(i).*Gtt(i)*Ar(1))./(4*pi*(R1(i)^2));%求接收机1的接收功率

end

for i=1:4

Pr2(i)=(Pt(i).*Gtt(i)*Ar(2))./(4*pi*(R1(i)^2));%求接收机2的接收功率

end

%% 加时沿

c=3e8;

shiyan1=rand(1,4)/1e3; %产生接收机1的随机时延

for i=1:4

shiyan2(i)=shiyan1(i)+(R2(i)-R1(i))/c; %接收机2的时延与时延1有关

end

%% 加噪声

SNR=20; %信噪比

snr=1/(10^(SNR/10));

% noise=-30;

%%%%接收机1接收的信号

Y11=[zeros(1,fs*shiyan1(1)) Y1];len11=length(Y11); %%加时延

Y12=[zeros(1,fs*shiyan1(2)) Y2];len12=length(Y12);

Y13=[zeros(1,fs*shiyan1(3)) Y3];len13=length(Y13);

Y14=[zeros(1,fs*shiyan1(4)) Y4];len14=length(Y14);

% noise11=awgn(Y11,noise);

% noise12=awgn(Y12,noise);

% noise13=awgn(Y13,noise);

% noise14=awgn(Y14,noise);

noise11=sqrt(snr/2)*(randn(1,len11));

noise12=sqrt(snr/2)*(randn(1,len12));

noise13=sqrt(snr/2)*(randn(1,len13));

noise14=sqrt(snr/2)*(randn(1,len14));

Y_11=Y11+noise11; %%加噪声

Y_12=Y12+noise12;

Y_13=Y13+noise13;

Y_14=Y14+noise14;

YY1=[Y_11 zeros(1,length(Y_13)-length(Y_11))];

YY2=[Y_12 zeros(1,length(Y_13)-length(Y_12))];

YY4=[Y_14 zeros(1,length(Y_13)-length(Y_14))];

YY3=Y_13;

receive1=YY1+YY2+YY3+YY4; %%接收机1接收的信号

figure;

plot(abs(receive1));

title('接收机1收到的信号');

%%%%接收机2接收的信号

Y21=[zeros(1,fs*shiyan2(1)) Y1];len21=length(Y21); %加时延

Y22=[zeros(1,fs*shiyan2(2)) Y2];len22=length(Y22);

Y23=[zeros(1,fs*shiyan2(3)) Y3];len23=length(Y23);

Y24=[zeros(1,fs*shiyan2(4)) Y4];len24=length(Y24);

% noise21=awgn(Y21,noise);

% noise22=awgn(Y22,noise);

% noise23=awgn(Y23,noise);

% noise24=awgn(Y24,noise);

noise21=sqrt(snr/2)*(randn(1,len21));

noise22=sqrt(snr/2)*(randn(1,len22));

noise23=sqrt(snr/2)*(randn(1,len23));

noise24=sqrt(snr/2)*(randn(1,len24));

Y_21=Y21+noise21;

Y_22=Y22+noise22; %加噪声

Y_23=Y23+noise23;

Y_24=Y24+noise24;

YYY1=[Y_21 zeros(1,length(Y_23)-length(Y_21))];

YYY2=[Y_22 zeros(1,length(Y_23)-length(Y_22))];

YYY4=[Y_24 zeros(1,length(Y_23)-length(Y_24))];

YYY3=Y_23;

receive2=YYY1+YYY2+YYY3+YYY4; %%接收机2接收的信号

figure;

plot(abs(receive2));

title('接收机2收到的信号');

fid = fopen('C:\Documents and Settings\Administrator\桌面\信号分析\receive1.dat','w');

fwrite(fid,receive1,'double');

fid = fopen('C:\Documents and Settings\Administrator\桌面\信号分析\receive2.dat','w');

fwrite(fid,receive2,'double');

%% 信号分析

clc;close all;clear all;

fid=fopen('C:\Documents and Settings\Administrator\桌面\信号分析\receive1.dat');

s=fread(fid,'double');

s=abs(s');

plot(s);

% noise=s(1:10000);

% noisemean=mean(noise);

% noisevar=var(noise);

T=0.3; %%%设定门限

fs=10e7;

% b=zeros(length(s),2);

%%脉宽

m=1;n=1;q=0;summ=0;count=0;%定义参数初值

%% 将信号二值化

as=zeros(1,length(s));

for i=1:length(s)

if s(i)＞=T %%如果信号幅度大于等于门限，则判为1

as(i)=1;

else

as(i)=0;%否则判为0

end

end

%% 求脉冲的起始位置和结束位置

for i=2:length(s)-300

summ=sum(as(i:i+300));

if as(i)==1&&count==0

if summ＞=220

b(m,1)=i; % b矩阵的第一列存放脉冲起始位置

m=m+1;

count=1;

end

end

if count==1&&summ＜40

b(n,2)=i; % b矩阵的第二列存放脉冲结束位置

n=n+1;

count=0;

end

end

%% 求脉宽

a=b(:,2)-b(:,1); %a矩阵存放脉宽

figure;

stem(a);

axis([0,200,0,4000]);

title('脉宽');

%% 载频

N=10000;

for i=1:length(b)

maichong=[s(b(i,1):b(i,2))];

maichong_FFT=abs(fft(maichong,N));

maichong_FFT(1:60)=0;

maichong_FFT(9940:10000)=0;

[Fcmax,k]=max(maichong_FFT);

Fc(i)=k/(2*N)*fs; %Fc矩阵存放载频

end

figure;

stem(Fc);

title('载频');

%% K均值聚类法对脉冲进行分类

% %选择聚类中心

maizu=[a';Fc/10000];

jlzx=[maizu(:,1),maizu(:,2),maizu(:,102),maizu(:,128)]; %聚类中心

jlzl(1,1)=1;

jlzl(2,1)=2;

jlzl(3,1)=102;

jlzl(4,1)=128; %聚类种类

%%

%设置判决条件

panjue=20;

panjue_js=25;

while panjue_js＞panjue

%对每个脉冲进行分类

jljl=zeros(1,4); %聚类距离

m1=2;m2=2;m3=2;m4=2;

M=length(a);

jljl_min=zeros(1,M);

for n=1:M

jljl(:)=sqrt((jlzx(1,:)-maizu(1,n)).^2+(jlzx(2,:)-maizu(2,n)).^2);

jljl_min(n)=min(jljl);

if jljl_min(n)＞500

continue;

end

if jljl(1)==jljl_min(n)

jlzl(1,m1)=n;

m1=m1+1;

elseif jljl(2)==jljl_min(n)

jlzl(2,m2)=n;

m2=m2+1;

elseif jljl(3)==jljl_min(n)

jlzl(3,m3)=n;

m3=m3+1;

else

jlzl(4,m4)=n;

m4=m4+1;

end

end

Mjl(1)=m1-1;Mjl(2)=m2-1;Mjl(3)=m3-1;Mjl(4)=m4-1; %聚类种类每类个数

%对每类求出新的聚类中心

jlzx_n=zeros(2,4);

for n=1:4

sum1=0;sum2=0;

for nn=1:Mjl(n)

sum1=sum1+maizu(1,jlzl(n,nn));

sum2=sum2+maizu(2,jlzl(n,nn));

jlzx_n(1,n)=sum1/Mjl(n);

jlzx_n(2,n)=sum2/Mjl(n);

end

end

%判决条件（新旧两个聚类中心）距离较小

panjue_js=sqrt(sum(sum((jlzx-jlzx_n).^2)));

%设置新的判决中心

jlzx=jlzx_n;

end