概率:频率和信度

频率学派:基于频率的解释,P(A)被认为是无限次重复试验事件A发生的频率,例如:当我们说随机抛硬币出现的头像概率为1/2,是指当重复抛硬币足够多次时,出现头像的频率接近于1/2。

对于一些重复性试验困难的问题中,这种解释就存在一定的挑战。比如:

如果要刻画明天下雨的概率或判断火星上曾经存在生命的可能性,频率的解释就无能为力了,主要是因为,明天的天气或者火星的历史都是不可重复的。

也就是说,当面对一些不可重复性的试验,基于频率的解释就行不通了。

贝叶斯学派:此时,基于信度的解释认为P(A)是观察者认为事件A发生的可信程度。

统计推断:在机器学习、数据挖掘等实际应用中,对于只能观测到有限数据样本的情况下,就需要“逆向工程”,推断数据背后的规律,也就是推断数据产生的过程,这个过程就是统计推断。

统计推断:给定观测数据x_{1},x_{2},...,x_{N}\sim F,推断/估计/学习概率分布F或其数字特征(如均值、方差等)

在统计推断中,主要有频率推断和贝叶斯推断两种方法。

频率推断将参数\theta看成未知但固定的,通过优化目标函数找到最优逼近\hat{\theta },这种估计方法也叫点估计。

贝叶斯推断将未知参数看作随机变量,推断其后验概率分布p(\theta\mid D )

概率:数据产生过程→观测数据

统计/学习:观测数据→数据产生过程

在进行统计推断时,需构建一个统计模型,一般分为参数化模型和非参数化模型。

统计模型:一组分布的集合M。

参数化模型:集合M中的分布可用有限哥参数表示:

M=\left \{ p(x;\theta ):\theta \in \Theta \right \}

其中,\theta为未知参数,\Theta是可行参数空间。

非参数化模型:集合M不能用有限个参数进行描述的模型或参数个数为无限多个。

对于参数化模型,统计推断的目标是估计位置参数\theta \in \Theta;对于非参数化模型,统计推断的目标是直接估计F。由于约束更少,所以后者一般更困难。

理解贝叶斯推断:贝叶斯方法将概率看做对事件(如明天会下雨)发生的信度。因此,可以对很多事情进行概率描述,包括模型的未知参数\theta。此外,当观察到新的数据时,对未知变量的信度也会相应发生变化,例如:当听到明天天气预报后,对明天是否会下雨会有更加确信的判断,并选择适当的行程安排。这个过程就可以用贝叶斯推断。

贝叶斯推断的基本流程:贝叶斯推断将未知参数\theta看作随机变量

(1)用p(\theta)描述在看到数据之前对参数可能取值的信度,成为参数\theta的先验分布;

(2)给定数据集D=\left \{ x_{_{i}} \right \}_{i=1}^{N},假设统计模型p(x\mid\theta )描述在给定参数\theta的情况下,生成数据x_{i}的信度,则p(D\mid \theta )称为参数\theta的似然函数;

(3)利用贝叶斯公式,得到给定数据后参数的概率分布p(\theta\mid D ),成为参数\theta的后验分布:

p(\theta\mid D )=\frac{p(D\mid \theta ))}{p(D)}

其中,p(D)称为证据。对上式两边积分,得到p(D )=\int p\left ( D\mid \theta \right )p(\theta )d\theta

与先验p(\theta)相比,后验分布p(\theta\mid D )蕴含了从数据D中观测到的信息,刻画了关于参数\theta更新后的概率分布。

与频率方法相比,统计推断把\theta看成未知参数,其值通过某个估计如MLE确定;但这个估计本身的不确定性是通过考虑数据集D的分布刻画(如估计的方差)。在贝叶斯推断中,模型的不确定性是通过参数\theta的分布刻画,数据集D是给定的。

%% ================ 雷达系统参数设置 ================ % 基本参数 pd = 0.9; % 检测概率 pfa = 1e-6; % 虚警概率 max_range = 20000; % 最大探测距离 (m) tgt_rcs = 1.5; % 目标RCS (m²) int_pulsenum = 10; % 积累脉冲数 B=7.5e6; fs = B; % 采样率 (Hz) prf = 7500; % 脉冲重复频率 (Hz) pulse_width = 10e-6; % 脉冲宽 (s) fc = 10e9; % 载波频率 (Hz) v = physconst('LightSpeed'); % 传播速 (m/s) lambda = v/fc; % 波长 (m) %% ================ 系统组件初始化 ================ % 1. 波形发生器 load BasicMonostaticRadarExampleData; waveform = phased.LinearFMWaveform(... 'SampleRate', fs,... 'PRF', prf,... 'PulseWidth', pulse_width,... 'SweepBandwidth', B,... 'SweepDirection', 'Up',... 'SweepInterval', 'Positive'); % 2. 天线阵列 antenna = phased.IsotropicAntennaElement(... 'FrequencyRange', [0.8*fc 1.2*fc],... 'BackBaffled', true); ura = phased.URA(... 'Element', antenna,... 'Size', [25 25],... 'ElementSpacing', [lambda/2 lambda/2]); % 3. 发射机 array_gain_db = 10*log10(25*25); % 阵列增益 snr_min = albersheim(pd, pfa, int_pulsenum); noise_power = noisepow(1/waveform.PulseWidth); % 噪声系数3dB, 290K参考温 peak_power = ((4*pi)^3 * noise_power * max_range^4 * db2pow(snr_min)) / ... (db2pow(2*(array_gain_db+20)) * tgt_rcs * lambda^2); transmitter = phased.Transmitter(... 'PeakPower', peak_power/625,... 'Gain', 20,... 'LossFactor', 0,... 'InUseOutputPort', true,... 'CoherentOnTransmit', true); % 4. 接收机 receiver = phased.ReceiverPreamp(... 'SampleRate', fs,... 'Gain', 20,... 'NoiseFigure', 0,... 'ReferenceTemperature', 290,... 'EnableInputPort', true,... 'SeedSource', 'Property',... 'Seed', 2023); % 5. 辐射器和收集器 radiator = phased.Radiator(... 'Sensor', ura,... 'OperatingFrequency', fc,... 'PropagationSpeed', v,... 'CombineRadiatedSignals', true,... 'WeightsInputPort', true); collector = phased.Collector(... 'Sensor', ura,... 'OperatingFrequency', fc,... 'PropagationSpeed', v,... 'Wavefront', 'Plane'); %% ================ 目标模型参数设置 ================ tgtpos = [[13500*cosd(2);13500*sind(2); 0],[13530*cosd(2);13530*sind(2); 0],[13470*cosd(2);13470*sind(2); 0],[18000*cosd(2); 18000*sind(2); 0],[18030*cosd(2); 18030*sind(2); 0],[179700*cosd(2); 17970*sind(2); 0]]; tgtvel = [[-450*cosd(2);-450*sind(2); 0],[-450*cosd(2); -450*sind(2); 0],[-450*cosd(2);-450*sind(2); 0],[-450*cosd(6);-450*sind(6); 0],[-450*cosd(6);-450*sind(6); 0],[-450*cosd(6);-450*sind(6); 0]]; tgtmotion = phased.Platform('InitialPosition',tgtpos,'Velocity',tgtvel); tgtrcs = [1.5 1.5 1.5 1.5 1.5 1.5]; target = phased.RadarTarget('MeanRCS',tgtrcs,'OperatingFrequency',fc); % Calculate the range, angle, and speed of the targets [tgtrng,tgtang] = rangeangle(tgtmotion.InitialPosition,... sensormotion.InitialPosition); numtargets = length(target.MeanRCS); %% ================ 目标模型初始化 ================ tgtmotion = phased.Platform('InitialPosition', tgtpos, 'Velocity', tgtvel); target = phased.RadarTarget('MeanRCS', tgtrcs, 'OperatingFrequency', fc); % 7. 传播信道 channel = phased.FreeSpace(... 'SampleRate', fs,... 'TwoWayPropagation', true,... 'OperatingFrequency', fc); % 8. 波束形成组件 steeringvec = phased.SteeringVector(... 'SensorArray', ura,... 'PropagationSpeed', v); beamformer = phased.PhaseShiftBeamformer(... 'SensorArray', ura,... 'OperatingFrequency', fc,... 'PropagationSpeed', v,... 'DirectionSource', 'Input port'); % 9. 雷达平台 sensormotion = phased.Platform(... 'InitialPosition', [0; 0; 0],... 'Velocity', [0; 0; 0]); %% ================ 扫描参数设置 ================ initialAz = 60; % 起始方位角 () endAz = -60; % 结束方位角 () scanstep = -4; % 扫描步长 () 负值表示从右向左扫描 % 确保扫描范围正确 if initialAz < endAz % 如果起始角小于结束角,交换值 [initialAz, endAz] = deal(endAz, initialAz); scanstep = -scanstep; % 反转扫描方向 end % 计算扫描网格 scangrid = initialAz-2 : scanstep : endAz+2; numscans = length(scangrid); disp(['扫描角范围: ', num2str(scangrid(1)), '° 到 ', num2str(scangrid(end)), '°']); disp(['扫描步数: ', num2str(numscans), ' (步长: ', num2str(scanstep), '°)']); disp(['每个波位脉冲数: ', num2str(int_pulsenum)]); pulsenum = int_pulsenum * numscans; % 总脉冲数 revisitTime = pulsenum / prf; % 扫描周期 % 计算每个PRI的采样点数 pri = 1 / prf; % 脉冲重复间隔 samples_per_pri = round(pri * fs); % 3000 (22.5e6 / 7500) fast_time_grid = (0:samples_per_pri-1)/fs; % 快时间网格 % 预分配接收信号矩阵 rxpulses = zeros(samples_per_pri, pulsenum); % 计算目标的距离、角和速 [tgtrng,tgtang] = rangeangle(tgtmotion.InitialPosition,... sensormotion.InitialPosition); numtargets = length(target.MeanRCS); % 创建导向矢量和波束形成器 steeringvec = phased.SteeringVector('SensorArray',ura,... 'PropagationSpeed',v); % 计算不同方向上的相位偏移权重 % 创建接收端的波束形成器 beamformer = phased.PhaseShiftBeamformer('SensorArray',ura,... 'OperatingFrequency',fc,'PropagationSpeed',v,... 'DirectionSource','Input port'); % 输入端口指定方向 % 定义传播信道 channel = phased.FreeSpace(... 'SampleRate',fs,... 'TwoWayPropagation',true,... 'OperatingFrequency',fc); % 创建快时间网格(从0开始,步长1/fs,上限1/prf) fast_time_grid = unigrid(0, 1/fs, 1/prf, '[)'); % 预分配接收信号矩阵 rxpulses = zeros(numel(fast_time_grid),pulsenum); % 预分配内存存储接收到的脉冲信号 % 创建阻塞式干扰机(有效辐射功率6000W,采样帧数=波形脉冲数×采样率/PRF) jammer = barrageJammer('ERP',6000,... 'SamplesPerFrame',waveform.NumPulses*waveform.SampleRate/waveform.PRF); % 创建干扰信号传播信道(单向传播) jammerchannel = phased.FreeSpace('TwoWayPropagation',false,... 'SampleRate',fs,'OperatingFrequency', fc); % 干扰机位置坐标 jammerloc = [15000; 15000; 30]; % 计算干扰机相对于雷达的角 [~,jamang] = rangeangle(jammerloc); % 主处理循环 - 处理每个脉冲 for m = 1:pulsenum % 更新平台和目标位置 [sensorpos,sensorvel] = sensormotion(1/prf); % 调用运动模型函数更新雷达平台位置/速 [tgtpos,tgtvel] = tgtmotion(1/prf); % 更新时间步长 = 1/PRF(脉冲重复间隔) % 计算目标相对参数 [tgtrng,tgtang] = rangeangle(tgtpos,sensorpos); % 计算目标相对雷达的距离角[方位角; 俯仰角] % 计算当前扫描角和波束权重 scanid = floor((m-1)/int_pulsenum) + 1; % 确定当前扫描角索引 sv = steeringvec(fc,scangrid(scanid)); % 获取当前方向的导向矢量 w = conj(sv); % 计算共轭权重(用于发射波束形成) % 生成和发射脉冲 pulse = waveform(); % 生成线性调频脉冲 [txsig,txstatus] = transmitter(pulse); % 通过发射机放大信号 txsig = radiator(txsig,tgtang,w); % 通过辐射器形成波束 txsig = channel(txsig,sensorpos,tgtpos,sensorvel,tgtvel); % 信号传播到目标 % 信号传播和目标反射 tgtsig = target(txsig); % 目标反射信号 % 生成干扰信号 jamsig = jammer(); % 生成阻塞式干扰 % 传播干扰信号到阵列 jamsig = jammerchannel(jamsig,jammerloc,[0;0;0],[0;0;0],[0;0;0]); % 收集干扰信号 jamsig = collector(jamsig,jamang); % 通过阵列收集干扰信号 % 接收和处理回波信号 rxsig = collector(tgtsig,tgtang); % 收集目标回波信号 rxsig = receiver(rxsig,~(txstatus>0)); % 通过接收机放大信号(发射时关闭接收) % 波束形成处理(接收方向与发射方向相同) rxpulses(:,m) = 625*beamformer(rxsig+jamsig,[scangrid(scanid);0]); end %% ================ 匹配滤波 ================ % 匹配滤波处理 matchingcoeff = getMatchedFilter(waveform); % 获取匹配滤波器系数 matchedfilter = phased.MatchedFilter(... 'Coefficients',matchingcoeff,... 'GainOutputPort',true); % 创建匹配滤波器对象(输出处理增益) [mf_pulses, mfgain] = matchedfilter(rxpulses); % 应用匹配滤波 % 重塑脉冲矩阵为[距离门×积累脉冲数×扫描角] mf_pulses = reshape(mf_pulses,[],int_pulsenum,numscans); % 补偿匹配滤波引入的延迟 matchingdelay = size(matchingcoeff,1)-1; % 计算延迟量 sz_mfpulses = size(mf_pulses); % 保存原始尺寸 % 补偿延迟:移除前matchingdelay个样本,末尾补零 mf_pulses = [mf_pulses(matchingdelay+1:end) zeros(1,matchingdelay)]; mf_pulses = reshape(mf_pulses,sz_mfpulses); % 恢复原始维 % 脉冲积累(非相干积累) int_pulses = pulsint(mf_pulses,'noncoherent'); % 沿脉冲维积累 int_pulses = squeeze(int_pulses); % 移除单一维 %% ================ 检测与距离估计 ================ % 定义距离门 range_gates = v*fast_time_grid/2; % 距离门对应实际距离 % 创建时间可变增益(TVG)补偿距离衰减 tvg = phased.TimeVaryingGain(... 'RangeLoss',2*fspl(range_gates,lambda),... % 双程自由空间路径损失 'ReferenceLoss',2*fspl(max(range_gates),lambda)); % 参考距离处的损失 tvg_pulses = tvg(mf_pulses); % 应用TVG % 脉冲积累(非相干) int_pulses = pulsint(tvg_pulses,'noncoherent'); % 沿脉冲维积累 int_pulses = squeeze(int_pulses); % 移除单一维 % 计算检测阈值 % 在设计系统中,采样率是噪声带宽的两倍 noise_bw = receiver.SampleRate/2; % 计算噪声带宽 % 计算噪声功率(考虑接收机噪声系数和参考温) npower = noisepow(noise_bw,... receiver.NoiseFigure,receiver.ReferenceTemperature); % 计算基于虚警概率的检测阈值 threshold = npower * db2pow(npwgnthresh(pfa,int_pulsenum,'noncoherent')); % 增加匹配滤波处理增益和阵列增益 threshold = threshold * db2pow(mfgain)*625; 在以上添加了干扰机的代码基础上,编写反干扰代码,使其能重新检测到目标
最新发布
07-18
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值