cod matlab

l1=4
l2=2


for n = 1:20
 %valores tomados de ads para la primera matriz, todas las frecuencias
p11= [9.943E-1-2.440E-2i -1.335E-2-9.768E-1i -9.664E-1-1.607E-2i -4.688E-2+9.572E-1i 9.479E-1+7.871E-2i 1.107E-1-9.374E-1i -9.258E-1-1.420E-1i -1.741E-1+9.142E-1i 9.017E-1+2.059E-1i 2.380E-1-8.865E-1i -8.719E-1-2.685E-1i -2.989E-1+8.552E-1i 8.388E-1+3.291E-1i 3.576E-1-8.201E-1i -8.055E-1-3.870E-1i -4.176E-1+7.808E-1i 7.606E-1+4.411E-1i 4.684E-1-7.377E-1i -7.190E-1-4.924E-1i -5.210E-1+6.935E-1i]
p12= [1.177E-3+1.399E-3i 1.698E-2-6.154E-3i 1.325E-2-1.417E-2i 7.072E-3+1.124E-2i 3.773E-3+3.942E-2i -3.779E-3+1.933E-2i -1.440E-3-1.722E-2i 1.820E-2+1.728E-3i 2.502E-2+4.499E-2i -4.648E-3+2.765E-2i -3.648E-2-1.671E-2i -2.644E-2+7.443E-3i -2.145E-3+4.947E-2i -8.167E-3-4.354E-3i -2.844E-2-9.388E-2i -1.154E-2-5.945E-2i 2.519E-2+4.548E-2i 1.989E-2+2.915E-2i -1.089E-2-8.090E-2i -5.444E-3-7.246E-2i]
p21= [-6.173E-4-4.506E-3i -9.470E-3+2.845E-3i 1.288E-2+1.351E-2i 2.239E-2-8.180E-3i 3.692E-3-3.678E-2i -1.323E-2-1.993E-2i 2.401E-3+1.572E-2i 2.791E-2-1.034E-3i 2.560E-2-4.177E-2i -1.014E-2-2.567E-2i -3.428E-2+1.619E-2i -2.368E-2-8.102E-3i -5.696E-3-5.126E-2i -4.624E-3+5.300E-3i -2.011E-2+9.279E-2i -1.751E-2+5.668E-2i 7.542E-3-4.712E-2i 2.514E-2-2.783E-2i 1.120E-2+8.183E-2i -1.291E-2+7.168E-2i]
p22= [1.004E0+3.148E-2i -1.468E-2+1.023E0i -1.035E0+1.653E-2i -4.967E-2-1.043E0i 1.049E0-8.529E-2i 1.234E-1+1.053E0i -1.055E0+1.604E-1i -1.995E-1-1.056E0i 1.057E0-2.399E-1i 2.815E-1+1.053E0i -1.049E0+3.212E-1i -3.628E-1-1.043E0i 1.036E0-4.045E-1i 4.453E-1+1.025E0i -1.018E0+4.880E-1i -5.322E-1-1.000E0i 9.866E-1-5.707E-1i 6.118E-1+9.683E-1i -9.527E-1+6.522E-1i -6.939E-1-9.283E-1i]

%se señala que valor usar de los vectores anteriores
a11=p11(1,n)
a12=p12(1,n)
a21=p21(1,n)
a22=p22(1,n)

%valores tomados de ads para la segunda matriz, todas las frecuencias
p33= [9.990E-1-9.636E-3i 6.932E-1-7.046E-1i 5.405E-3-9.826E-1i -6.774E-1-7.060E-1i -9.742E-1-3.528E-2i -7.216E-1+6.497E-1i -6.590E-2+9.649E-1i 6.216E-1+7.372E-1i 9.561E-1+9.722E-2i 7.529E-1-5.917E-1i 1.297E-1-9.450E-1i -5.617E-1-7.668E-1i -9.344E-1-1.613E-1i -7.824E-1+5.305E-1i -1.927E-1+9.218E-1i 5.015E-1+7.942E-1i 9.108E-1+2.238E-1i 8.097E-1-4.675E-1i 2.570E-1-8.942E-1i -4.370E-1-8.183E-1i]
p43= [-6.812E-4-3.937E-3i -1.014E-2-1.345E-3i -1.197E-2+8.551E-3i -1.041E-2+2.096E-2i -8.317E-3+2.769E-2i -6.542E-3+2.157E-2i -4.789E-3+1.228E-3i -2.117E-3-2.494E-2i 3.191E-3-4.201E-2i 9.261E-3-3.806E-2i 1.453E-2-1.164E-2i 1.528E-2+2.388E-2i 1.045E-2+4.724E-2i -1.136E-3+4.064E-2i -1.395E-2+5.801E-3i -2.110E-2-4.369E-2i -2.140E-2-7.747E-2i -1.261E-2-7.410E-2i -1.785E-3-3.127E-2i 1.611E-3+2.972E-2i]
p34=[-9.302E-4-4.222E-4i 9.669E-3-9.194E-4i 8.531E-3-1.097E-2i 1.553E-3-2.304E-2i -7.628E-3-2.860E-2i -1.426E-2-2.107E-2i -1.550E-2-4.224E-4i -8.528E-3+2.504E-2i 3.560E-3+4.138E-2i 1.499E-2+3.623E-2i 2.009E-2+9.030E-3i 1.590E-2-2.585E-2i 6.106E-3-4.674E-2i -7.301E-3-3.932E-2i -1.625E-2-1.847E-3i -1.652E-2+4.573E-2i -1.061E-2+7.783E-2i -5.574E-3+7.285E-2i -4.828E-3+2.923E-2i -1.145E-2-3.174E-2i]
p44=[1.002E0+1.759E-2i 7.093E-1+7.221E-1i 4.686E-3+1.018E0i -7.091E-1+7.364E-1i -1.026E0+3.551E-2i -7.653E-1-6.903E-1i -6.911E-2-1.032E0i 6.704E-1-7.925E-1i 1.038E0-1.047E-1i 8.214E-1+6.471E-1i 1.412E-1+1.039E0i -6.232E-1+8.488E-1i -1.042E0+1.785E-1i -8.767E-1-5.959E-1i -2.159E-1-1.040E0i 5.710E-1-9.017E-1i 1.042E0-2.551E-1i 9.308E-1+5.389E-1i 2.954E-1+1.034E0i -5.096E-1+9.506E-1i]

%se señala que valor usar de los vectores anteriores
b11=p33(1,n)
b12=p34(1,n)
b21=p43(1,n)
b22=p44(1,n)


M1= [a11 a12; a21 a22];%matriz Rd
M2= [b11 b12; b21 b22]%matriz thru (Rt)
inv_M2= inv(M2)% inverso de la matriz thru
T = M1*inv_M2 % matriz T
t11 = T(1,1)
t12 = T(1,2)
t21 = T(2,1)
t22 = T(2,2)

%Solución de la ecuacion cuadratica 18 del articulo
a=t21;
b=t22-t11;
c=-t12;
d=b^2-4*a*c;
if d>0
    x1=(-b+sqrt(d))/(2*a);
    x2=(-b-sqrt(d))/(2*a);
    disp('valor positivo de la raiz')
elseif d==0
    x1=-b/(2*a);
    x2=-b/(2*a);
    disp('la raiz vale 0')
else
    x1=(-b+i*sqrt(-d))/(2*a);
    x2=(-b-i*sqrt(-d))/(2*a);
    disp('la raiz es negativa')
end
disp('valores de la ecuación de segundo grado:')
x1 , x2 % soluciones una es bm y la otra dm=(am/cm)
mag_sol1= abs(x1)
mag_sol2= abs(x2)

%comprobando y asignando que valor le corresponde a bm y a (am/cm)
%discriminante
if mag_sol1>mag_sol2
    bm= x2;
    dm= x1;
else
    bm= x1;
    dm= x2
end

%solución ecuacion 17 articulo
lamda= (dm-bm)/((dm*t22)+(bm*dm*t21)-(bm*t11)-(t12))
cte_prop=(log(lamda)/(l2-l1))
alfa= real(cte_prop)
beta=imag(cte_prop)
  disp(alfa)
end
%n=3%posición del valor de cada parametro respecto a la frecuencia (frecuencia 1)...