文章目录
- 概述
- 初始化环境参数
概述
无线电层析成像是一种通过获取一定区域内多对相对固定的无线通信节点间的某种测量数据后,按照一定的数学处理方法,对区域内的障碍物目标以图像的形式
展现出来的成像技术。
开山之作:
J. Wilson and N. Patwari, “Radio tomographic imaging with wireless
networks,” IEEE Trans. Mobile Comput., vol. 9, no. 5, pp. 621–632,
May 2010.
Radio tomographic imaging with wireless networks
优势:
- 不要求物体佩戴任何设备;
- 无线电波可以穿透物理结构。
初始化环境参数
代码是初始化环境参数,并生成用于深度学习的训练数据和测试数据的代码。主要包括以下部分:
输入参数:包括移动射频节点数量、每个节点收集数据的航点数、信噪比范围、路径损耗指数范围、全球平均发射功率水平范围、信号角度范围、Kappa参数(用于生成空间相关的高斯随机场)、节点遍历区域大小、图像原点坐标、图像大小、像素大小、白色像素的衰减、垂直、水平和方形墙的最小/最大值等。
数据预分配:预先分配一些变量的空间,包括像素数、UWB链路数、每个链路的UWB数据数等。
节点位置计算:根据输入参数计算每个节点的所有坐标。
数据存储:预先分配存储数据的空间,包括节点位置、UWB节点间距离、SLF权重、不匹配实验的SLF权重和偏置的权重。
生成权重和距离:根据节点位置和输入参数生成SLF权重和UWB节点间距离。
生成空间相关的高斯随机场:根据图像大小、像素大小和Kappa参数生成空间相关的高斯随机场,并进行Cholesky分解得到上三角矩阵。
%M:移动RF节点的数量
%Tw:每个节点收集数据的航点数
%sig_epsilon_range:噪声信噪比的范围
%alpha_range:路径损耗指数的范围
%b_range:全局平均发射功率水平的范围
%sig_theta_range:SLF图像的标准差范围
%kappa:空间相关性高斯随机场的参数
%areabounds:节点遍历区域的大小
%imgorg:图像的左下角坐标
%imgdims:图像的大小(像素)
%pixelsize:每个像素的长度/宽度 %dBwhite:白色像素的衰减值(分贝)
%wallcount_vec:垂直、水平和方形墙壁的最小/最大数量
clear; clc; close all; % 清空命令窗口、清空工作区、关闭所有图窗
M = 4; %number of mobile RF nodes 移动射频节点的数量
Tw = 6; %number of waypoints for each node to gather data from 每个节点收集数据的航点数,一条边上有6个航点
sig_epsilon_range = [0.3,1,3]; % SNR is about 30,40,50 信噪比约为30、40、50
%对于信噪比(SNR),通常使用以分贝(dB)为单位的对数比值来表示。在这里,sig_epsilon_range = [0.3,1,3] 是指信噪比的范围,其中的值是对应于SNR的对数比值。
%SNR = 10 * log10(sig_epsilon_range^2)
%当 sig_epsilon_range = 0.3 时,SNR ≈ 30 dB 当 sig_epsilon_range = 1 时,SNR ≈ 40 dB 当 sig_epsilon_range = 3 时,SNR ≈ 50 dB
%低噪、中噪、高噪
alpha_range = [0.9 1]; % path loss exponent 路径损耗指数
b_range = [90 100]; %global avg TX "level" (db) 全球平均发射功率水平(分贝)
sig_theta_range = [0.01,0.03,0.09]; % 信号角度范围
%sig_theta_range中的三个值分别为0.01、0.03和0.09,表示了三种不同的信号角度范围。较小的值表示信号角度范围较窄,即信号的传输方向相对集中。较大的值表示信号角度范围较宽,即信号的传输方向相对分散。
%通过在生成数据时使用不同的sig_theta值,可以模拟不同的信号角度环境。这样可以更好地评估和比较算法在不同信号角度条件下的性能和鲁棒性。
kappa = 0.21; % kappa参数,用于生成空间相关的高斯随机场
%kappa参数用于生成空间相关的高斯随机场(GRF)。GRF是一种随机过程,用于建模空间中的随机变量,如图像或地理数据。
%在这里,kappa参数被用于生成空间相关的高斯随机场的协方差矩阵。协方差矩阵描述了随机变量之间的相关性和方差。kappa的值决定了GRF的空间相关性的程度。
%较小的kappa值会导致GRF的相关性较弱,即相邻像素之间的相关性较低。这意味着生成的随机场在空间上具有较大的变化和不规则性。
%较大的kappa值会导致GRF的相关性较强,即相邻像素之间的相关性较高。这意味着生成的随机场在空间上具有较为平滑和连续的变化。
%通过调整kappa值,可以控制生成的随机场的空间相关性,从而适应不同的应用需求和模拟场景。
areabounds = [0 5 0 5]; %size of the node traverse area (m) [xmin xmax ymin ymax] 节点遍历区域的大小(米)
imgorg = [0.5 0.5]; %lower left corner of the image [x y] (m) 图像的左下角坐标(米)
imgdims = [40 40]; %image size in pixels 图像的大小(像素)
pixelsize = 0.1; %l/w of each square pixel (m) 每个正方形像素的长度/宽度(米)
dBwhite = 1; %dB atten of a white pixel 白色像素的衰减(分贝)
wallcount_vec = [0 2 0 2 1 3]; %min/max for vertical , horizontal and square walls 垂直、水平和方形墙的最小/最大值
%表示垂直、水平和方形墙壁的最小和最大数量。这个向量中的值用来限制在区域内生成的墙壁的数量。
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%preallocate data 预分配数据
K = imgdims(1)*imgdims(2); %pixel count 像素数
N_Link = M*(M-1); %Number of UWB Links UWB链路数
N_Tw = Tw^2; %number of UWB data per link 每个链路的UWB数据数
%每条链路连接连接两个节点,每个节点可以在6个航点上移动
x_1=0; y_1=0;
x_2 = imgdims(2)*pixelsize+2*imgorg(1); %UWB node range UWB节点范围
y_2 = imgdims(1)*pixelsize+2*imgorg(2);
x_distance = (imgdims(2)-1)*pixelsize/(Tw-1); %distance of two position of per node 每个节点两个位置之间的距离
y_distance = (imgdims(1)-1)*pixelsize/(Tw-1);
x_org = imgorg(1)+0.5*pixelsize; %coordinate of first UWB node 第一个UWB节点的坐标
y_org = imgorg(2)+0.5*pixelsize;
% all coordinates of each node 每个节点的所有坐标
c = NaN(M,Tw);
for Twi = 1:Tw %Tw point per node
pos_node(1,Twi) = x_org + (Twi-1)*x_distance + 1j * y_1;
pos_node(2,Twi) = x_2 + 1j * (y_org + (Twi-1)*y_distance);
pos_node(3,Twi) = x_org + (Twi-1)*x_distance + 1j * y_2;
pos_node(4,Twi) = x_1 + 1j * (y_org + (Twi-1)*y_distance);
end
%通过 for 循环,对于每个节点,根据给定的参数计算节点的坐标,并将其存储在 pos_node 矩阵中。
%具体来说,对于每个节点,根据节点的索引 Twi 和给定的参数(x_org、y_org、x_distance、y_distance),计算出节点在不同路径点上的坐标。
%对于节点 1,其 x 坐标为 x_org + (Twi-1)*x_distance,y 坐标为 y_1。在这里,路径点是沿着 x 轴移动的。
%(0.55,0) (1.33,0) (2.11,0) (2.89,0) (3.67,0) (4.45,0)
%对于节点 2,其 x 坐标为 x_2,y 坐标为 y_org + (Twi-1)*y_distance。在这里,路径点是沿着 y 轴移动的。
%(5,0.55) (5,1.33) (5,2.11) (5,2.89) (5,3.67) (5,4.45)
%对于节点 3,其 x 坐标为 x_org + (Twi-1)*x_distance,y 坐标为 y_2。在这里,路径点是沿着 x 轴移动的。
%(0.55,5) (1.33,5) (2.11,5) (2.89,5) (3.67,5) (4.45,5)
%对于节点 4,其 x 坐标为 x_1,y 坐标为 y_org + (Twi-1)*y_distance。在这里,路径点是沿着 y 轴移动的。
%(0,0.55) (0,1.33) (0,2.11) (0,2.89) (0,3.67) (0,4.45)
%store the data 存储数据
nodepos = NaN(2,N_Tw,N_Link);
%nodepos 是一个大小为 2 x N_Tw x N_Link 的矩阵,用于存储每个测量的UWB节点的位置。N_Tw代表每个节点的路径点数,N_Link代表链路的数量。
d = NaN(N_Link*N_Tw,1); %distance of UWB node for each measurements 每个测量的UWB节点间距离
W = NaN(N_Link*N_Tw,K); %weight of slf slf的权重
W_ENR = NaN(N_Link*N_Tw,K); %weight of slf for mismatched experiments 用于不匹配实验的slf权重
Z = zeros(N_Link*N_Tw,N_Link); %weight of bias 偏置的权重
%node position calculation 节点位置计算
N_Linki = 0;
for i1 = 1:M
for i2 = 1:M
if i1 ~= i2
N_Linki = N_Linki + 1;
for Twi = 1:Tw
nodepos(1,(Twi-1)*Tw+1:Twi*Tw,N_Linki)=pos_node(i1,Twi)*ones(1,Tw);
nodepos(2,(Twi-1)*Tw+1:Twi*Tw,N_Linki)=pos_node(i2,:);
%分别记录下对于每条链路,节点的位置
end
end
end
end
%(Twi-1)Tw+1:TwiTw表示节点i1在nodepos矩阵中的位置索引范围。例如,当Twi=1时,索引范围为1:Tw;当Twi=2时,索引范围为Tw+1:2*Tw,依此类推。
%通过pos_node(i1,Twi)*ones(1,Tw)的赋值,将节点i1在路径点Twi上的坐标赋给了nodepos矩阵的第一个维度中指定的索引范围。这样,nodepos矩阵中的相应位置就存储了节点i1在不同路径点上的坐标信息。
%generate weight and distance 生成权重和距离
for N_Linki = 1:N_Link
W((N_Linki-1)*N_Tw+1:N_Linki*N_Tw,:) = ...
f_gen_Wi_ellipse(nodepos(:,:,N_Linki),pixelsize,imgorg,imgdims,1,2);
% Inverse Area Elliptical Model:
%B. R. Hamilton, X. Ma, R. J. Baxley, and S. M. Matechik, ‘‘Propagation
%modeling for radio frequency tomography in wireless networks, IEEE J.
%Sel. Topics Signal Process., vol. 8, no. 1, pp. 55–65, Feb. 2014.
% W((N_Linki-1)*N_Tw+1:N_Linki*N_Tw,:) = ...
% f_gen_Wi(nodepos(:,:,N_Linki),pixelsize,imgorg,imgdims,1,2); % Normalized Ellipse Model
di = abs( nodepos(1,:,N_Linki) - nodepos(2,:,N_Linki) )'; %link lengths
% 使用abs函数计算向量差的绝对值,得到节点之间的距离(欧式距离)。
d((N_Linki-1)*N_Tw+1:N_Linki*N_Tw,1) = 20 * log10( di );%将距离转换为分贝单位
Z((N_Linki-1)*N_Tw+1:N_Linki*N_Tw,N_Linki) = ones(N_Tw,1);%将Z矩阵中对应链路的列设置为1,表示偏置的权重。
end
C_theta = f_gen_C(imgdims,pixelsize,kappa); % spatial correlated Gaussian random field
R_theta = chol(C_theta); % 生成C_theta的上三角矩阵
%根据图像大小、像素大小和Kappa参数生成空间相关的高斯随机场,并进行Cholesky分解得到上三角矩阵。
Inverse Area Elliptical Model:
Propagation Modeling for Radio Frequency Tomography in Wireless Networks
function [ Wi ] = f_gen_Wi_ellipse( nodepos, pixelsize, imgorg, imgdims, i1, i2 )
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% f_gen_Wi : GENERATE A WEIGHT MATRIX FOR THE ITH LINK (BIN) IN THE SYSTEM
% This function calls f_pixel_weights to compile all the relevant weights for all measurements taken
% between a given pair of nodes
% nodepos - array of node positions with time
% pixelsize - l/w of each square pixel (m)
% imgorg - [x y] position of lower left corner of image(m)
% imgdims - [rows cols] of the image
% i1 - node index of the TX, i2 - node index of the RX (link i)
% Wi ( ni x K ) - weights for link i
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
ni = size(nodepos,2); %# of time steps
K = imgdims(1)*imgdims(2); %pixel count
Wi = NaN(ni,K);
if i1 ~= i2
for j = 1:ni %loop over individual link's measurements
p1 = [real(nodepos(i1,j)) imag(nodepos(i1,j))];
p2 = [real(nodepos(i2,j)) imag(nodepos(i2,j))];
Wi(j,:) = f_pixel_weights_3(imgdims,pixelsize,imgorg,p1,p2);
end
else
Wi = zeros(ni,K); %invalid link
end
end
function Xij = f_pixel_weights_3(imgdims, pixelsize, imgorg, p1, p2)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% COMPUTES THE PIXEL WEIGHTS FOR SINGLE MEASUREMENT J OF LINK I (tomographic projection) %
% Xij (1 x numpixels) - row vector of weights for the given points. inverse area ellipse %
% imgdims (1x2) - number of [rows cols] in the image. the image is assumed to reside in %
% the first quadrant, with the lower left corner at the origin %
% pixelsize - the PHYSICAL length of a pixel's sides, assumed to be square pixels (m) %
% imgorg (1x2) - physical [x y] location of the lower left corner of the image (m) %
% p1, p2 (1x2) - [x y] locations of the two nodes of interest (m) %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% inverse area ellipse model
freq_RF = 2.4e9; % RF 2.4 GHz
lanpt = 3e8/freq_RF; %wave length
p1 = p1' - imgorg';
p2 = p2' - imgorg'; %scale points to "single units"
fai = (p1-p2)/norm(p1-p2); %projection vector
c = (p1+p2)/2; % center of the ellipse
d = norm(p1-p2);
K = imgdims(1)*imgdims(2); %pixel count
bb = zeros(1,K);
u = NaN(1,K);
u_max = 0.5*sqrt(1*lanpt*d);
NI = 0;
for j=1:imgdims(2) %col
for i=1:imgdims(1) %row
NI = NI + 1; % number of pixels
sx = (j-0.5)*pixelsize;
sy = (imgdims(1)-i+0.5)*pixelsize;
px = fai(1)*(sx-c(1))+fai(2)*(sy-c(2));
py = -fai(2)*(sx-c(1))+fai(1)*(sy-c(2));
a0 = 1; b0 = -(px^2+py^2-d^2); c0 = -py^2*d^2;
u_2 = (-b0 + sqrt(b0^2-4*a0*c0))/(2*a0);
u(1,NI) = sqrt(u_2);
if u(1,NI)<u_max
if u(1,NI)<0.05
u(1,NI) = 0.05;
end
bb(1,NI) = 1/(pi*u(1,NI)*sqrt(u(1,NI)^2+d^2));
end
end
end
Xij = bb;
% imshow(reshape(Xij,[40,40]),[0 max(Xij)])
% %DEBUG plot code
end