% Stochastic Process
% LOS/NLOS Identification
clc
clear all
close all
load channel1;
h1=h;
load channel2;
h2=h;
dt=t(2);
h1_thr = (h1>0.1).*h1;
h2_thr = (h2>0.1).*h2;
h1_thr_e = exp(-t)*ones(1,length(h1(1,:))).*h1_thr;
h2_thr_e = exp(-t)*ones(1,length(h1(1,:))).*h2_thr;
h1_thr_t = (t/t(end))*ones(1,length(h1(1,:))).*h1_thr;
h2_thr_t = (t/t(end))*ones(1,length(h1(1,:))).*h2_thr;
%---------------------Mean Excess Delay---------------------------
t1_m = t'*(h1_thr_e.^2)./sum(h1_thr_e.^2);
t2_m = t'*(h2_thr_e.^2)./sum(h2_thr_e.^2);
% t1_m = t'*(h1_thr.^2)./sum(h1_thr.^2);
% t2_m = t'*(h2_thr.^2)./sum(h2_thr.^2);
[p1,x1]=ksdensity(t1_m);
[p2,x2]=ksdensity(t2_m);
[indx indx11] = min(abs(x1'*ones(1,length(t1_m))-ones(length(x1),1)*t1_m));
[indx indx21] = min(abs(x2'*ones(1,length(t1_m))-ones(length(x1),1)*t1_m));
pM_LOS = sum(p1(indx11)>p2(indx21))/length(indx)
[indx indx12] = min(abs(x1'*ones(1,length(t2_m))-ones(length(x2),1)*t2_m));
[indx indx22] = min(abs(x2'*ones(1,length(t2_m))-ones(length(x2),1)*t2_m));
pM_NLOS = sum(p2(indx22)>p1(indx12))/length(indx)
%---------------------RMS Delay Spread---------------------------
t1_rms = zeros(1,length(h1(1,:)));
for i = 1:length(h1(1,:))
t1_rms(i) = ((t'-t1_m(i)).^2)*(h1_thr_t(:,i).^2)./sum(h1_thr_t(:,i).^2);
end
t2_rms = zeros(1,length(h1(1,:)));
for i = 1:length(h1(1,:))
t2_rms(i) = ((t'-t2_m(i)).^2)*(h2_thr_t(:,i).^2)./sum(h2_thr_t(:,i).^2);
end
% t1_rms = zeros(1,length(h1(1,:)));
% for i = 1:length(h1(1,:))
% t1_rms(i) = ((t'-t1_m(i)).^2)*(h1_thr(:,i).^2)./sum(h1_thr(:,i).^2);
% end
%
% t2_rms = zeros(1,length(h1(1,:)));
% for i = 1:length(h1(1,:))
% t2_rms(i) = ((t'-t2_m(i)).^2)*(h2_thr(:,i).^2)./sum(h2_thr(:,i).^2);
% end
[p3,x3]=ksdensity(t1_rms);
[p4,x4]=ksdensity(t2_rms);
%find the indexes of t1_rms to compute the probability
[indx indx31] = min(abs(x3'*ones(1,length(t1_m))-ones(length(x1),1)*t1_rms));
[indx indx41] = min(abs(x4'*ones(1,length(t1_m))-ones(length(x1),1)*t1_rms));
pR_LOS = sum(p3(indx31)>p4(indx41))/length(indx)
[indx indx32] = min(abs(x3'*ones(1,length(t2_m))-ones(length(x2),1)*t2_rms));
[indx indx42] = min(abs(x4'*ones(1,length(t2_m))-ones(length(x2),1)*t2_rms));
pR_NLOS = sum(p4(indx42)>p3(indx32))/length(indx)
%---------------------Peak-to-Lead Time---------------------------
% [X1,I1] = max(h1_thr);
% [X2,I2] = max(h2_thr);
[X1,I1] = max(h1_thr_e);
[X2,I2] = max(h2_thr_e);
P2L_1 = I1*dt;
P2L_2 = I2*dt;
[p5,x5]=ksdensity(P2L_1);
[p6,x6]=ksdensity(P2L_2);
[indx indx51] = min(abs(x5'*ones(1,length(t1_m))-ones(length(x1),1)*P2L_1));
[indx indx61] = min(abs(x6'*ones(1,length(t1_m))-ones(length(x1),1)*P2L_1));
pP_LOS = sum(p5(indx51)>p6(indx61))/length(indx)
[indx indx52] = min(abs(x5'*ones(1,length(t2_m))-ones(length(x2),1)*P2L_2));
[indx indx62] = min(abs(x6'*ones(1,length(t2_m))-ones(length(x2),1)*P2L_2));
pP_NLOS = sum(p6(indx62)>p5(indx52))/length(indx)
%------------------------Kurtosis---------------------------
kurtosis_1=kurtosis(h1);
kurtosis_2=kurtosis(h2);
[p7,x7]=ksdensity(kurtosis_1);
[p8,x8]=ksdensity(kurtosis_2);
[indx indx71] = min(abs(x7'*ones(1,length(t1_m))-ones(length(x1),1)*kurtosis_1));
[indx indx81] = min(abs(x8'*ones(1,length(t1_m))-ones(length(x1),1)*kurtosis_1));
pK_LOS = sum(p7(indx71)>p8(indx81))/length(indx)
[indx indx72] = min(abs(x7'*ones(1,length(t2_m))-ones(length(x2),1)*kurtosis_2));
[indx indx82] = min(abs(x8'*ones(1,length(t2_m))-ones(length(x2),1)*kurtosis_2));
pK_NLOS = sum(p8(indx82)>p7(indx72))/length(indx)
%----------------------Number of peaks-------------------------------------
for i=1:1000
Z11(i) = length(findpeaks(h1_thr(:,i)));
end
for i=1:1000
Z12(i) = length(findpeaks(h2_thr(:,i)));
end
[p9,x9]=ksdensity(Z11);
[p10,x10]=ksdensity(Z12);
[indx indx91] = min(abs(x9'*ones(1,length(t1_m))-ones(length(x1),1)*Z11));
[indx indx101] = min(abs(x10'*ones(1,length(t1_m))-ones(length(x1),1)*Z11));
pZ_LOS = sum(p9(indx91)>p10(indx101))/length(indx)
[indx indx92] = min(abs(x9'*ones(1,length(t2_m))-ones(length(x2),1)*Z12));
[indx indx102] = min(abs(x10'*ones(1,length(t2_m))-ones(length(x2),1)*Z12));
pZ_NLOS = sum(p10(indx102)>p9(indx92))/length(indx)
%----------------------Number of nonzero samples-------------------------------------
Z21 = sum(sign(h1_thr)+1)/2;
Z22 = sum(sign(h2_thr)+1)/2;
[p11,x11]=ksdensity(Z21);
[p12,x12]=ksdensity(Z22);
[indx indx111] = min(abs(x11'*ones(1,length(t1_m))-ones(length(x1),1)*Z21));
[indx indx121] = min(abs(x12'*ones(1,length(t1_m))-ones(length(x1),1)*Z21));
pZ2_LOS = sum(p11(indx111)>p12(indx121))/length(indx)
[indx indx112] = min(abs(x11'*ones(1,length(t2_m))-ones(length(x2),1)*Z22));
[indx indx122] = min(abs(x12'*ones(1,length(t2_m))-ones(length(x2),1)*Z22));
pZ2_NLOS = sum(p12(indx122)>p11(indx112))/length(indx)
%----------------------Maximum Peak-------------------------------------
% Z31 = max(h1_thr);
% Z32 = max(h2_thr);
Z31 = max(h1_thr_e);
Z32 = max(h2_thr_e);
[p13,x13]=ksdensity(Z31);
[p14,x14]=ksdensity(Z32);
[indx indx131] = min(abs(x13'*ones(1,length(t1_m))-ones(length(x1),1)*Z31));
[indx indx141] = min(abs(x14'*ones(1,length(t1_m))-ones(length(x1),1)*Z31));
pZ3_LOS = sum(p13(indx131)>=p14(indx141))/length(indx)
[indx indx132] = min(abs(x13'*ones(1,length(t2_m))-ones(length(x2),1)*Z32));
[indx indx142] = min(abs(x14'*ones(1,length(t2_m))-ones(length(x2),1)*Z32));
pZ3_NLOS = sum(p14(indx142)>=p13(indx132))/length(indx)
%----------------------Energy of h(t) -------------------------------------
% Z41 = sum(h1_thr.^2)*dt;
% Z42 = sum(h2_thr.^2)*dt;
Z41 = sum(h1_thr_e.^2)*dt;
Z42 = sum(h2_thr_e.^2)*dt;
[p15,x15]=ksdensity(Z41);
[p16,x16]=ksdensity(Z42);
[indx indx151] = min(abs(x15'*ones(1,length(t1_m))-ones(length(x1),1)*Z41));
[indx indx161] = min(abs(x16'*ones(1,length(t1_m))-ones(length(x1),1)*Z41));
pZ4_LOS = sum(p15(indx151)>p16(indx161))/length(indx)
[indx indx152] = min(abs(x15'*ones(1,length(t2_m))-ones(length(x2),1)*Z42));
[indx indx162] = min(abs(x16'*ones(1,length(t2_m))-ones(length(x2),1)*Z42));
pZ4_NLOS = sum(p16(indx162)>p15(indx152))/length(indx)
%-----------------------JOINT--------------------------------
PJ1_LOS = sum(p1(indx11).*p5(indx51)>p2(indx21).*p6(indx61))*100/length(indx)
PJ1_NLOS = sum(p1(indx12).*p5(indx52)<p2(indx22).*p6(indx62))*100/length(indx)
PJ2_LOS = sum(p1(indx11).*p13(indx131)>p2(indx21).*p14(indx141))*100/length(indx)
PJ2_NLOS = sum(p1(indx12).*p13(indx132)<p2(indx22).*p14(indx142))*100/length(indx)
PJ3_LOS = sum(p13(indx131).*p5(indx51)>p14(indx141).*p6(indx61))*100/length(indx)
PJ3_NLOS = sum(p13(indx132).*p5(indx52)<p14(indx142).*p6(indx62))*100/length(indx)
PJ4_LOS = sum(p1(indx11).*p5(indx51).*p13(indx131)>p2(indx21).*p6(indx61).*p14(indx141))*100/length(indx)
PJ4_NLOS = sum(p1(indx12).*p5(indx52).*p13(indx132)<p2(indx22).*p6(indx62).*p14(indx142))*100/length(indx)
PJ5_LOS = sum(p1(indx11).*p5(indx51).*p13(indx131).*p15(indx151)>p2(indx21).*p6(indx61).*p14(indx141).*p16(indx161))*100/length(indx)
PJ5_NLOS = sum(p1(indx12).*p5(indx52).*p13(indx132).*p15(indx152)<p2(indx22).*p6(indx62).*p14(indx142).*p16(indx162))*100/length(indx)
%-----------------------PLOT-----------------------------
figure (1)
plot(x1,p1,x2,p2)
legend('LOS','NLOS')
title('PDF of Mean Excess Delay')
xlabel('Mean Excess Delay')
ylabel('')
figure (2)
plot(x3,p3,x4,p4)
legend('LOS','NLOS')
title('PDF of RMS Delay Spread')
xlabel('RMS Delay Spread')
ylabel('')
figure (3)
plot(x5,p5,x6,p6)
legend('LOS','NLOS')
title('PDF of Peak-to-Lead Time')
xlabel('Peak-to-Lead Time')
ylabel('')
figure (4)
plot(x7,p7,x8,p8)
legend('LOS','NLOS')
title('PDF of Kurtosis')
xlabel('Kurtosis')
ylabel('')
figure (5)
plot(x9,p9,x10,p10)
legend('LOS','NLOS')
title('PDF of Number o
评论2