1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121
|
clear;
NbEch=100;
MargeEch=50
k=0:NbEch+MargeEch;
Fsignal=1000;
Fsignal2=1500;
Fech = 96000;
Tech = 1/Fech;
Angle = -40;
Angle2 = -35;
ECART_MIC_M = 0.1;
VitSon=1500;
Lambda = VitSon / Fsignal;
Lambda2 = VitSon / Fsignal2;
retard = sin(3.141592*Angle/180) * ECART_MIC_M/VitSon;
retard2 = sin(3.141592*Angle2/180) * ECART_MIC_M/VitSon;
Mic1Cplx = exp( - 2 * %i * 3.141592 * Fsignal *( k(:)*Tech-0*retard) );
Mic1 = exp( - 2 * %i * 3.141592 * Fsignal *( k(:)*Tech-0*retard) ) + exp( - 2 * %i * 3.141592 * Fsignal2 *( k(:)*Tech-0*retard2) ) ;
Mic2 = exp( - 2 * %i * 3.141592 * Fsignal *( k(:)*Tech-1*retard) ) + exp( - 2 * %i * 3.141592 * Fsignal2 *( k(:)*Tech-1*retard2) ) ;
Mic3 = exp( - 2 * %i * 3.141592 * Fsignal *( k(:)*Tech-2*retard) ) + exp( - 2 * %i * 3.141592 * Fsignal2 *( k(:)*Tech-2*retard2) ) ;
Mic4 = exp( - 2 * %i * 3.141592 * Fsignal *( k(:)*Tech-3*retard) ) + exp( - 2 * %i * 3.141592 * Fsignal2 *( k(:)*Tech-3*retard2) ) ;
Mic1 = real ( Mic1 );
Mic2 = real ( Mic2 );
Mic3 = real ( Mic3 );
Mic4 = real ( Mic4 );
QuartPeriodeK = (Fech/Fsignal) / 4;
for i=1 : (NbEch + MargeEch - QuartPeriodeK)
Mic1(i) = Mic1(i) + %i * Mic1(i+QuartPeriodeK);
Mic2(i) = Mic2(i) + %i * Mic2(i+QuartPeriodeK);
Mic3(i) = Mic3(i) + %i * Mic3(i+QuartPeriodeK);
Mic4(i) = Mic4(i) + %i * Mic4(i+QuartPeriodeK);
end
QuartPeriodeK = (Fech/Fsignal2) / 4;
for i=1 : (NbEch + MargeEch - QuartPeriodeK)
Mic1(i) = Mic1(i) + %i * Mic1(i+QuartPeriodeK);
Mic2(i) = Mic2(i) + %i * Mic2(i+QuartPeriodeK);
Mic3(i) = Mic3(i) + %i * Mic3(i+QuartPeriodeK);
Mic4(i) = Mic4(i) + %i * Mic4(i+QuartPeriodeK);
end
Mic1 = Mic1(1:NbEch);
Mic2 = Mic2(1:NbEch);
Mic3 = Mic3(1:NbEch);
Mic4 = Mic4(1:NbEch);
X=[Mic1, Mic2, Mic3, Mic4]';
R=X*conj(X)';
[vect, val]=spec (R);
vp1=[vect(:,1)];
vp2=[vect(:,2)];
vp3=[vect(:,3)];
vp4=[vect(:,4)];
ValX = 0;
ValY1 = 0;
Index = 1;
maxDOA1 = -100000;
for angle_degre=-90:90
SinAngle = sin ( angle_degre *3.141592/180 );
arg = - 2.0 * 3.141592654 * ECART_MIC_M * SinAngle / Lambda;
A = [ exp(0); exp(%i*arg); exp(%i*2*arg); exp(%i*3*arg)];
//PowMUSIC = conj(A)'*vp2*conj(vp2)'*A +conj(A)'*vp3*conj(vp3)'*A +conj(A)'*vp4*conj(vp4)'*A ;
PowMUSIC = conj(A)'*vp3*conj(vp3)'*A +conj(A)'*vp4*conj(vp4)'*A ;
ValY1(Index) = 1 / abs (PowMUSIC);
if ( ValY1(Index) > maxDOA1 )
maxDOA1 = ValY1(Index);
end
ValX(Index) = angle_degre;
Index = Index + 1;
end
Index = 1;
for angle_degre=-90:90
ValY1(Index) = 65536 * ValY1(Index) / maxDOA1;
Index = Index + 1;
end
//plot (ValX, ValY1);
ValX = 0;
ValY2 = 0;
Index = 1;
maxDOA2 = -100000;
for angle_degre=-90:90
SinAngle = sin ( angle_degre *3.141592/180 );
arg = - 2.0 * 3.141592654 * ECART_MIC_M * SinAngle / Lambda2;
A = [ exp(0); exp(%i*arg); exp(%i*2*arg); exp(%i*3*arg)];
//PowMUSIC = conj(A)'*vp2*conj(vp2)'*A +conj(A)'*vp3*conj(vp3)'*A +conj(A)'*vp4*conj(vp4)'*A ;
PowMUSIC = conj(A)'*vp3*conj(vp3)'*A +conj(A)'*vp4*conj(vp4)'*A ;
ValY2(Index) = 1 / abs (PowMUSIC);
if ( ValY2(Index) > maxDOA2 )
maxDOA2 = ValY2(Index);
end
ValX(Index) = angle_degre;
Index = Index + 1;
end
Index = 1;
for angle_degre=-90:90
ValY2(Index) = 65536 * ValY2(Index) / maxDOA2;
Index = Index + 1;
end
//plot (ValX, ValY2);
Index = 1;
ValY_tot = 0;
for angle_degre=-90:90
ValY_tot(Index) = ValY1(Index) + ValY2(Index);
Index = Index + 1;
end
plot (ValX, ValY_tot);
abs(val) |
Partager