caiyu/PKU-Millimeter-Wave-Radar-Tutorial
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Delete 雷达点云处理流程仿真 directory

Browse Source
This commit is contained in:
St. Zhang 2024-03-11 21:03:08 +08:00 committed by GitHub
parent f5bca4cb61
commit c6cb4edf45
No known key found for this signature in database
GPG Key ID: B5690EEEBB952194
10 changed files with 0 additions and 452 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

19
雷达点云处理流程仿真/RAfftMatrix.m
View File

@ -1,19 +0,0 @@
%% 2DFFT得到距离-角度热力图
function fft2dData = RAfftMatrix(rawData)
rawData = squeeze(rawData);
[angleBin,rangeBin] = size(rawData);
angleWin = hanning(angleBin);
angleWin2D = repmat(angleWin,1,rangeBin);
rangeWin = hanning(rangeBin)';
rangeWin2D = repmat(rangeWin,angleBin,1);
rawDataWin = rawData .* angleWin2D;
fft1dData = fftshift(fft(rawDataWin,angleBin,1));
fft1dDataWin = fft1dData .* rangeWin2D;
fft2dData = fft(fft1dDataWin,rangeBin,2);
end

14
雷达点云处理流程仿真/RDfftMatrix.m
View File

@ -1,14 +0,0 @@
%% 2DFFT得到距离-速度热力图
function fft2dData = RDfftMatrix(rawData)
rawData = squeeze(rawData);
[rangeBin,dopplerBin] = size(rawData);
rangeWin = hanning(rangeBin);
rangeWin2D = repmat(rangeWin,1,dopplerBin);
dopplerWin = hanning(dopplerBin)';
dopplerWin2D = repmat(dopplerWin,rangeBin,1);
rawDataWin = rawData .* rangeWin2D;
fft1dData = fft(rawDataWin,rangeBin,1);
fft1dDataWin = fft1dData .* dopplerWin2D;
fft2dData = fftshift(fft(fft1dDataWin,dopplerBin,2),2);
end

67
雷达点云处理流程仿真/cfar.m
View File

@ -1,67 +0,0 @@
%% 请了解：什么是CFAR，以及什么是2D-CFAR
function [pointList,cfarRD] = cfar(parameter,accumulateRD)
[rangeLen,dopplerLen] = size(accumulateRD);
%% doppler维度搜索
dopplerMethod = parameter.dopplerMethod;
dopplerSNR = parameter.dopplerSNR;
dopplerWinGuardLen = parameter.dopplerWinGuardLen;
dopplerWinTrainLen = parameter.dopplerWinTrainLen;
dopplerLeft = accumulateRD(:,dopplerLen - dopplerWinGuardLen - dopplerWinTrainLen + 1:dopplerLen);
dopplerRight = accumulateRD(:,1:dopplerWinGuardLen+dopplerWinTrainLen);
dopplercfar = [dopplerLeft accumulateRD dopplerRight];
dopplerCfarList = [];
for rangeIdx = 1:rangeLen
for dopplerIdx = 1:dopplerLen
dopplerCfarIdx = dopplerIdx + dopplerWinGuardLen + dopplerWinTrainLen;
leftCell = dopplercfar(rangeIdx,dopplerIdx:dopplerIdx+dopplerWinTrainLen-1);
rightCell = dopplercfar(rangeIdx,dopplerCfarIdx+dopplerWinGuardLen:dopplerCfarIdx+dopplerWinGuardLen+dopplerWinTrainLen-1);
leftNoise = mean(leftCell);
rightNoise = mean(rightCell);
noise = (leftNoise + rightNoise) / 2;
indexdb = dopplercfar(rangeIdx,dopplerCfarIdx);
targetSnr = indexdb - noise;
if targetSnr > dopplerSNR
dopplerCfarList = [dopplerCfarList dopplerIdx];
cfarRDdoppler(rangeIdx,dopplerIdx) = indexdb;
end
end
end
dopplerCfarList = unique(dopplerCfarList);
%% range维度搜索
rangeMethod = parameter.rangeMethod;
rangeSNR = parameter.rangeSNR;
rangeWinGuardLen = parameter.rangeWinGuardLen;
rangeWinTrainLen = parameter.rangeWinTrainLen;
rangeUp = accumulateRD(rangeLen - rangeWinGuardLen - rangeWinTrainLen + 1:rangeLen,:);
rangeDown = accumulateRD(1:rangeWinGuardLen+rangeWinTrainLen,:);
rangecfar = [rangeUp;accumulateRD;rangeDown];
rangeCfarList = [];
cfarRD = zeros(rangeLen,dopplerLen);
for dopplerIdx = dopplerCfarList
for rangeIdx = 1:rangeLen
rangeCfarIdx = rangeIdx + rangeWinGuardLen + rangeWinTrainLen;
upCell = rangecfar(rangeIdx:rangeIdx+rangeWinTrainLen-1,dopplerIdx);
downCell = rangecfar(rangeCfarIdx+rangeWinGuardLen:rangeCfarIdx+rangeWinGuardLen+rangeWinTrainLen-1,dopplerIdx);
upNoise = mean(upCell);
downNoise = mean(downCell);
if rangeMethod == 1
noise = (upNoise + downNoise) / 2;
indexdb = rangecfar(rangeCfarIdx,dopplerIdx);
targetSnr = indexdb - noise;
if targetSnr > rangeSNR
rangeCfarList = [rangeCfarList [rangeIdx;dopplerIdx]];
cfarRD(rangeIdx,dopplerIdx) = indexdb;
end
elseif dopplerMethod == 2
elseif dopplerMethod == 3
else
end
end
end
pointList = rangeCfarList;
end

10
雷达点云处理流程仿真/chan_Accumulate.m
View File

@ -1,10 +0,0 @@
%% 通道间的非相干累计
function accumulateRD = chan_Accumulate(fft2dDataDB)
[channelNum,rangeBin,dopplerBin] = size(fft2dDataDB);
accumulateRD = zeros(rangeBin,dopplerBin);
for channelId = 1:channelNum
accumulateRD = accumulateRD + (abs(squeeze(fft2dDataDB(channelId,:,:))));
end
end

83
雷达点云处理流程仿真/doa.m
View File

@ -1,83 +0,0 @@
%% 三种角度估计的方法
function [angle,doa_abs] = doa(parameter,antVec)
doaMethod = parameter.doaMethod;
if doaMethod == 1
[angle,doa_abs] = dbfMethod(parameter,antVec);
elseif doaMethod == 2
[angle,doa_abs] = fftMethod(parameter,antVec);
elseif doaMethod == 3
[angle,doa_abs] = caponMethod(parameter,antVec);
else
end
end
function [angle,doa_abs] = dbfMethod(parameter,antVec)
txAntenna = parameter.txAntenna;
rxAntenna = parameter.rxAntenna;
virtualAntenna = parameter.virtualAntenna;
lambda = parameter.lambda;
txNum = length(txAntenna);
rxNum = length(rxAntenna);
dx = parameter.dx;
deg = -90:0.1:90;
weightVec = zeros(virtualAntenna,1);
doa_dbf = zeros(length(deg),1);
kk = 1;
for degscan = deg
for txId = 1:txNum
for rxId = 1:rxNum
dphi = ((txId-1) * rxNum + rxId - 1) * 2 * pi / lambda * dx * sind(degscan);
weightVec((txId-1) * rxNum + rxId) = exp(-1i * dphi);
end
end
doa_dbf(kk) = antVec'*weightVec;
kk = kk + 1;
end
doa_abs = abs(doa_dbf);
[pk,loc]=max(doa_abs);
angle = deg(loc);
end
function [angle,doa_abs] = fftMethod(parameter,antVec)
angleIndex = asin((-512:1:512-1)/512) * 180 / pi;
doa_fft=fftshift(fft(antVec,1024));
doa_abs=abs(doa_fft);
[pk,loc]=max(doa_abs);
angle = angleIndex(loc);
end
function [angle,doa_abs] = caponMethod(parameter,antVec)
txAntenna = parameter.txAntenna;
rxAntenna = parameter.rxAntenna;
txNum = length(txAntenna);
rxNum = length(rxAntenna);
lambda = parameter.lambda;
dx = parameter.dx;
virtualAntenna = parameter.virtualAntenna;
Rx = antVec * antVec' / virtualAntenna;
deg = -90:0.1:90;
a = zeros(virtualAntenna,1);
% a = zeros(1,virtualAntenna);
kk = 1;
for degscan = deg
for txId = 1:txNum
for rxId = 1:rxNum
virtualAntennaId = (txId-1) * rxNum + rxId - 1;
dphi = 2 * pi / lambda * dx * virtualAntennaId * sind(degscan);
a((txId-1) * rxNum + rxId) = exp(-1i * dphi);
end
end
RxInv = inv(Rx);
P_out(kk) = 1/(a'*RxInv*a);
kk = kk + 1;
end
doa_abs = abs(P_out);
[pk,loc]=max(doa_abs);
angle = deg(loc);
figure;
plot(deg,20*log10(abs(P_out)));
end

14
雷达点云处理流程仿真/generateCfarParameter.m
View File

@ -1,14 +0,0 @@
%% 生成CFAR所需的参数
function parameter = generateCfarParameter()
parameter.dopplerMethod = 1; %1：ca-cfar 2：so-cfar 3：go-cfar
parameter.dopplerSNR = 10;
parameter.dopplerWinGuardLen = 2;
parameter.dopplerWinTrainLen = 8;
parameter.rangeMethod = 1; %1：ca-cfar 2：so-cfar 3：go-cfar
parameter.rangeSNR = 10;
parameter.rangeWinGuardLen = 2;
parameter.rangeWinTrainLen = 4;
end

36
雷达点云处理流程仿真/generateParameter.m
View File

@ -1,36 +0,0 @@
% 参数设置
function parameter = generateParameter()
parameter.c = 3e8; %光速
parameter.stratFreq = 76.5e9; %起始频率
parameter.Tr = 10e-6; %扫频时间 也就是周期
parameter.Samples = 256; %采样点
parameter.Fs = 25.6e6; %采样率
parameter.rangeBin = parameter.Samples ; %rangebin
parameter.Chirps = 512; %chirp数
parameter.dopplerBin = parameter.Chirps; %dopplerbin
parameter.Slope = 30e12; %chirp斜率
parameter.Bandwidth = parameter.Slope * parameter.Tr ; %发射信号有效带宽
parameter.BandwidthValid = parameter.Samples/parameter.Fs*parameter.Slope; %发射信号带宽
parameter.centerFreq = parameter.stratFreq + parameter.Bandwidth / 2; %中心频率
parameter.lambda = parameter.c / parameter.centerFreq; %波长
parameter.txAntenna = ones(1,3); %发射天线个数
parameter.rxAntenna = ones(1,4); %接收天线个数
parameter.txNum = length(parameter.txAntenna);
parameter.rxNum = length(parameter.rxAntenna);
parameter.virtualAntenna = length(parameter.txAntenna) * length(parameter.rxAntenna);
parameter.dz = parameter.lambda / 2; %接收天线俯仰间距
parameter.dx = parameter.lambda / 2; %接收天线水平间距
parameter.doaMethod = 2; %测角方法选择 1-dbf 2-fft 3-capon
parameter.target = [
100 -20 0; %target1 range speed angle
0 10 -30; %target2 range speed angle
0 20 30; %target2 range speed angle
];
end

67
雷达点云处理流程仿真/generateSignal.m
View File

@ -1,67 +0,0 @@
%% 生成混频信号
function rawData = generateSignal(Parameter)
c = Parameter.c; %光速
stratFreq = Parameter.stratFreq; %起始频率
Tr = Parameter.Tr; %扫频时间
samples = Parameter.Samples; %采样点
fs = Parameter.Fs; %采样率
rangeBin = Parameter.rangeBin; %rangeBin
chirps = Parameter.Chirps; %chirp数
dopplerBin = Parameter.dopplerBin; %dopplerBin
slope = Parameter.Slope; %chirp斜率
bandwidth = Parameter.Bandwidth; %发射信号带宽
centerFreq = Parameter.centerFreq; %中心频率
lambda = Parameter.lambda;
txAntenna = Parameter.txAntenna; %发射天线
txNum = length(txAntenna); %发射天线数
rxAntenna = Parameter.rxAntenna; %接收天线
rxNum = length(rxAntenna); %接收天线数
dz = Parameter.dz; %俯仰间距
dx = Parameter.dx; %水平间距
target = Parameter.target; %目标
targetNum = size(target,1); %目标数
rawData = zeros(txNum*rxNum,rangeBin,dopplerBin);
t = 0:1/fs:Tr-(1/fs); %chirp采样的时间序列
for chirpId = 1:chirps
for txId = 1:txNum
St = exp((1i*2*pi)*(centerFreq*(t+(chirpId-1)*Tr)+slope/2*t.^2)); %发射信号
for rxId = 1:rxNum
Sif = zeros(1,rangeBin);
for targetId = 1:targetNum
%%连续帧 目标设置，如果不需要连续帧，令Parameter.frame=0，即可。
if targetId==1
targetRange = target(targetId,1)-Parameter.frame;
targetSpeed = target(targetId,2);
targetAngle = target(targetId,3);
elseif targetId==2
targetRange = target(targetId,1)+0.5*Parameter.frame;
targetSpeed = target(targetId,2);
targetAngle = target(targetId,3);
elseif targetId==3
targetRange = target(targetId,1)+Parameter.frame;
targetSpeed = target(targetId,2);
targetAngle = target(targetId,3);
end
tau = 2 * (targetRange + targetSpeed * (txId - 1) * Tr) / c;
fd = 2 * targetSpeed / lambda;
wx = ((txId-1) * rxNum + rxId) / lambda * dx * sind(targetAngle);
Sr = 10*exp((1i*2*pi)*((centerFreq-fd)*(t-tau+(chirpId-1) * Tr)+slope/2*(t-tau).^2 -wx)); %回波信号
Sif = Sif + St .* conj(Sr);
%叠加20dB高斯白噪声
Sif = awgn(Sif,20);
end
rawData((txId-1) * rxNum + rxId,:,chirpId) = Sif;
end
end
end
end

112
雷达点云处理流程仿真/main.m
View File

@ -1,112 +0,0 @@
%% 功能：毫米波雷达点云生成流程MATLAB仿真
%% 思考：如何将IWR1843EVM+DCA1000采集得到的数据移植到这份代码上
clc;clear;close all;
Frame =1; %帧数设置；
for frame =1:Frame %帧数设置
%% 雷达参数设置
parameter = generateParameter();
parameter.frame = frame;
%% 雷达回波信号建模
rawData = generateSignal(parameter);
firstChirp = rawData(1,:,1);
figure(1);
plot(real(firstChirp));
hold on;
plot(imag(firstChirp));
xlabel('采样点数'); ylabel('幅值');title('原始数据实虚部');%第1个chirp。
%% 雷达信号处理
rangeRes = parameter.c / (2 * parameter.BandwidthValid); %距离分辨率 有效带宽
rangeIndex = (0:parameter.rangeBin-1) * rangeRes;
speedRes = parameter.lambda / (2 * parameter.dopplerBin * parameter.Tr);
dopplerIndex = (-parameter.dopplerBin/2:1:parameter.dopplerBin/2 - 1) * speedRes;
angleRes = parameter.lambda / (parameter.virtualAntenna * parameter.dx) * 180 / pi;
angleIndex = (-parameter.virtualAntenna/2:1:parameter.virtualAntenna/2 - 1) * angleRes;
%% 1D FFT
fft1dData = fft(firstChirp);
figure(2);
plot(rangeIndex,db(abs(fft1dData)./max(abs(fft1dData))));
xlabel('距离(m)'); ylabel('幅值(dB)');title('距离维FFT');
%% 2D FFT
%% 距离-多普勒谱
channelNum = size(rawData,1);
rangebinNum = size(rawData,2);
dopplerbinNum = size(rawData,3);
fft2dDataPower= zeros(size(rawData));
fft2dDataDB = zeros(size(rawData));
fftRADataPower= zeros(size(rawData));
for chanId = 1:1:channelNum
fft2dDataPower(chanId,:,:) = RDfftMatrix(rawData(chanId,:,:));
end
figure(3);
mesh(dopplerIndex',rangeIndex,db(abs(squeeze(fft2dDataPower(chanId,:,:)))));
view(2);
xlabel('速度(m/s)'); ylabel('距离(m)'); zlabel('幅值');
title('距离-多普勒谱');
mesh(abs(squeeze(fft2dDataPower(chanId,:,:))));
display_static=[30,100,400,300];
set(gcf,'Position',display_static); % [左下角x,左下角y,宽度,高度]
%% 距离-角度谱
for dopplerId = 1:1:dopplerbinNum
fftRADataPower(:,:,dopplerId) = RAfftMatrix(rawData(:,:,dopplerId));
end
figure(4);
imagesc(rangeIndex,angleIndex,(abs(squeeze(fftRADataPower(:,:,dopplerId)))));
view(2);
xlabel('距离(m)'); ylabel('角度'); zlabel('幅值');
title('距离-角度谱');
%% 多通道非相干积累
accumulateRD = chan_Accumulate((fft2dDataPower));
figure(5);
imagesc(dopplerIndex',rangeIndex,db(accumulateRD));
view(2);
xlabel('速度(m/s)'); ylabel('距离(m)'); zlabel('幅值');
title(['通道积累 第',num2str(frame),'帧']);
pause(0.01);
%% CFAR检测
cfarParameter = generateCfarParameter(); %生成cfar数据
[pointList,cfarRD] = cfar(cfarParameter,db(accumulateRD));
figure(6);
mesh(dopplerIndex',rangeIndex,cfarRD);
xlabel('速度(m/s)'); ylabel('距离(m)'); zlabel('幅值');
title('cfar');
display_static=[30,100,400,300];
set(gcf,'Position',display_static); % [左下角x,左下角y,宽度,高度]
%% peakSearch
[RD_pearkSearch,peakSearchList] = peakSearch(cfarRD,pointList);
detectPointNum = size(peakSearchList,2);
%% DOA估计
for targetIdx = 1:detectPointNum
rangeBin = peakSearchList(1,targetIdx);
speedBin = peakSearchList(2,targetIdx);
range = (rangeBin - 1) * rangeRes;
speed = (speedBin - parameter.dopplerBin/2 - 1) * speedRes;
ant = squeeze(fft2dDataPower(:,rangeBin,speedBin));
[angle,doa_abs] = doa(parameter,ant);
figure(7);
angleIndex = asin((-512:1:512-1)/512) * 180 / pi;
hold on;
plot(angleIndex,doa_abs);grid on
title('测角结果');
xlabel('角度');ylabel('幅值');
fprintf('目标%d的距离为%f,速度为%f,角度为%f\n',targetIdx,range,speed,angle);
end
end

30
雷达点云处理流程仿真/peakSearch.m
View File

@ -1,30 +0,0 @@
%% 寻找峰值
function [RD_pearkSearch,peakSearchList] = peakSearch(RD_cfar,cfarTargetList)
peakSearchList = [];
[rangeLen, dopplerLen] = size(RD_cfar);
RD_pearkSearch = zeros(rangeLen, dopplerLen);
length = size(cfarTargetList,2);
for targetIdx = 1:length
rangeIdx = cfarTargetList(1,targetIdx);
dopplerIdx = cfarTargetList(2,targetIdx); %坐标
if rangeIdx > 1 && rangeIdx < rangeLen && dopplerIdx > 1 && dopplerIdx < dopplerLen %边界点不考虑
if RD_cfar(rangeIdx,dopplerIdx) > RD_cfar(rangeIdx - 1,dopplerIdx) && ...
RD_cfar(rangeIdx,dopplerIdx) > RD_cfar(rangeIdx + 1,dopplerIdx) && ...
RD_cfar(rangeIdx,dopplerIdx) > RD_cfar(rangeIdx,dopplerIdx - 1) && ...
RD_cfar(rangeIdx,dopplerIdx) > RD_cfar(rangeIdx,dopplerIdx + 1)
RD_pearkSearch(rangeIdx,dopplerIdx) = RD_cfar(rangeIdx,dopplerIdx);
cfarTarget = [rangeIdx ; dopplerIdx];
peakSearchList = [peakSearchList cfarTarget];
end
end
end
end