# How to vectorize this piece of code and why doesn't e come out to be zero though it must come out to be zero because u and b are equal?

2 views (last 30 days)
Sadiq Akbar on 10 Dec 2022
Commented: Bruno Luong on 11 Dec 2022
u=[30 50 75 -30 -50 -75];
b=u;
Noise=5;
M = 6;%Constant1
N = 6;%Constant2
d = 0.5;%Constant3
K = length(u)/2; %Constant4
alpha=ones(1,K);%[1 1 1 1 1];
a=zeros(M,K);%aTx
f=zeros(N,K);%fRx
c=zeros(M*N, length(u)-K);% Extended Matrix
%%%%%%%%%%%%%%%%%%%%
% Swapping vector b
%%%%%%%%%%%%%%%%%%%%
[~, ix] = sort(u); % u is my desired vector
[~, ix1(ix)] = sort(b);
b = b(ix1);
%%%%%%%%%%%%%%%%%%%%
% Observed response
%%%%%%%%%%%%%%%%%%%
for i=1:K
for h=1:M
a(h,i)=exp(j*2*pi*(h-1)*d*sind(u(i))); %%%%% a
end
for p=1:N
f(p,i)=exp(j*2*pi*(p-1)*d*sind(u(K+i))); %%%%% f
end
end
for g= 1:K
c(:,g)=kron(a(:,g),f(:,g));% Extended vector
end
yo=c*alpha'; % Observed Response (Noise Free)
yo=awgn(yo,Noise);% Uncomment for Noise consideration
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Estimated response
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
ae=zeros(M,K);
fe=zeros(N,K);
ce=zeros(M*N, length(u)-K);
for i=1:K
for h=1:M
ae(h,i)=exp(j*2*pi*(h-1)*d*sind(b(i))); %%%%% ae
end
for p=1:N
fe(p,i)=exp(j*2*pi*(p-1)*d*sind(b(K+i))); %%%%% fe
end
end
for g= 1:K
ce(:,g)=kron(ae(:,g),fe(:,g));
end
ye=ce*alpha';% Estimated Response
%%%%%%%%%%%%%%%%%
% MSE
%%%%%%%%%%%%%%%%%
e=0.0;
for s=1:M*N
e=e+(abs(yo(s,1)-ye(s,1))).^2;
end
e=e/M*N;

Bruno Luong on 10 Dec 2022
Not tested but this
for i=1:K
for h=1:M
ae(h,i)=exp(j*2*pi*(h-1)*d*sind(b(i))); %%%%% ae
end
for p=1:N
fe(p,i)=exp(j*2*pi*(p-1)*d*sind(b(K+i))); %%%%% fe
end
end
can be replaced by
h=(1:M)';
p=(1:N)';
i=1:K;
br = b(:).';
ae=exp(j*2*pi*(h-1)*d.*sind(br(i)));
fe=exp(j*2*pi*(p-1)*d.*sind(br(K+i)));
And
for g= 1:K
ce(:,g)=kron(ae(:,g),fe(:,g));
end
By
ce = reshape(reshape(ae,K,1,[]).*fe,K,[]);
Bruno Luong on 11 Dec 2022
You are welcome.

Walter Roberson on 10 Dec 2022
Edited: Walter Roberson on 10 Dec 2022
format long g
u=[30 50 75 -30 -50 -75];
b=u;
Noise=5;
M = 6;%Constant1
N = 6;%Constant2
d = 0.5;%Constant3
K = length(u)/2; %Constant4
alpha=ones(1,K);%[1 1 1 1 1];
a=zeros(M,K);%aTx
f=zeros(N,K);%fRx
c=zeros(M*N, length(u)-K);% Extended Matrix
%%%%%%%%%%%%%%%%%%%%
% Swapping vector b
%%%%%%%%%%%%%%%%%%%%
[~, ix] = sort(u); % u is my desired vector
[~, ix1(ix)] = sort(b);
b = b(ix1);
%%%%%%%%%%%%%%%%%%%%
% Observed response
%%%%%%%%%%%%%%%%%%%
for i=1:K
for h=1:M
a(h,i)=exp(j*2*pi*(h-1)*d*sind(u(i))); %%%%% a
end
for p=1:N
f(p,i)=exp(j*2*pi*(p-1)*d*sind(u(K+i))); %%%%% f
end
end
for g= 1:K
c(:,g)=kron(a(:,g),f(:,g));% Extended vector
end
yo=c*alpha'; % Observed Response (Noise Free)
yo=awgn(yo,Noise);% Uncomment for Noise consideration
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Estimated response
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
ae=zeros(M,K);
fe=zeros(N,K);
ce=zeros(M*N, length(u)-K);
for i=1:K
for h=1:M
ae(h,i)=exp(j*2*pi*(h-1)*d*sind(b(i))); %%%%% ae
end
for p=1:N
fe(p,i)=exp(j*2*pi*(p-1)*d*sind(b(K+i))); %%%%% fe
end
end
for g= 1:K
ce(:,g)=kron(ae(:,g),fe(:,g));
end
ye=ce*alpha';% Estimated Response
%%%%%%%%%%%%%%%%%
% MSE
%%%%%%%%%%%%%%%%%
e=0.0;
for s=1:M*N
e=e+(abs(yo(s,1)-ye(s,1))).^2;
end
e=e/(M*N);
e2 = sum( abs(yo - ye).^2 )
e2 =
11.7479934461356
yo - ye
ans =
0.334761462628025 - 0.251604299511244i 0.178519020660132 - 0.11790756943532i 0.284315581572755 - 0.408074305354306i 0.436075962283236 + 0.0206136540004817i 0.148113323815692 - 0.1129397164173i 0.148601229842305 - 0.684227855074489i -0.045229189659987 + 1.06816976506834i 0.0449597505407739 + 0.0474860403089703i -0.117067171711324 - 0.335272836828924i -0.170838058879666 - 0.485416275775231i -1.01036972409615 + 0.0987920670194429i 0.17724058327333 + 0.358430065054272i -0.296057798011593 + 0.268159218550655i 0.219874574050923 - 0.303912184656327i -0.196397184200888 + 0.579342311435272i 0.362391793015109 - 0.00637536507210101i -0.198930152531317 + 0.17081587160563i 0.00126195142908353 + 0.196629724720236i 0.451323400274059 - 0.19338404893942i 0.0163549272291462 + 0.211263099020562i 0.153134628152902 - 0.360072343716162i 0.125618160624793 - 0.447663575852295i 0.18682261208064 - 0.440681989383234i 0.56658020901609 - 0.0413352507512756i 0.216864956628124 + 0.402128673330828i -0.650894787871801 + 0.195828199364969i 0.526482743315456 + 0.161261991699198i -0.0222655470394779 - 0.303132876336713i 0.249479128033312 - 0.125089115081018i 0.23011373904232 + 0.168827365022589i 0.37643644331261 + 0.00726229738016904i -0.933652368281796 - 0.958904637492813i -0.552407757174997 - 0.128690424103339i 0.0864282827509827 + 0.565038735206412i 0.617124996371567 - 1.11650348178558i 0.390001758508875 + 0.358139661751764i
Sadiq Akbar on 10 Dec 2022
Thanks dear Walter Roberson for your kind response. Yes, you are right that awgn() has been added to yo and not to ye. So that's the reason. But what about vectorization of the code? Can you replace all the for loops by vectorization?
Regards,