function Radiation_Decay_Coefficient_Modified(d,lamda,Bg,phi,fx,n)
dd = 10e-6; % depth under the WG
D = d; % the depth in which delta_eps is 1/e of the maximum
n1e = 2.3278;
n2e = 2.2028;
n1o = 2.2460;
n2o = 2.2860;
n3 = 1.00;
Vs = 3500; % Velocity of the SAW
c = 3e8;
mu0 = 4*pi*1e-7;
eps0 = 8.854e-12;
k = 2*pi/lamda;
omega = k*c;
% Before SAW
h1 = sqrt(n1e^2*k^2-Bg^2);
a2 = sqrt(Bg^2-n2e^2*k^2);
a3 = sqrt(Bg^2-n3^2*k^2);
dx1 = d/n;
dx2 = (dd-d)/n;
dE = d+1/a2+1/a3;
CTE = sqrt(omega*mu0/(Bg*dE));
En1 = [];
En2 = [];
for x1 = linspace(0,-d,n);
E= CTE*cos(h1*x1-phi);
En1 = [En1;E];
end;
for x2 = linspace(-d,-dd,n);
E = CTE*cos(phi+h1*d)*exp(a2*(x2+d));
En2 = [En2;E];
end;
alpha = [];
fs0 = Vs/(2*pi)*(Bg-n1o*k)+1;
for fs = linspace(fs0,900e6,100);
Enu1 = [];
Enu2 = [];
% After SAW
X = (2*pi/Vs)*fs; % SAW period
Bv = Bg - X;
sigma = sqrt(n1o^2*k^2-Bv^2);
rho = sqrt(n2o^2*k^2-Bv^2);
Delta = sqrt(n3^2*k^2-Bv^2); % Should be imaginary
for x1 = linspace(0,-d,n);
E = 2/(pi*1)*sin(1*pi/2)*cos(sigma*1*x1-phi);
Enu1 = [Enu1;E];
end;
for x2 = linspace(-d,-dd,n);
E = 2/(pi*1)*sin(pi*1/2)*cos(sigma*d)*cos(rho*1*(x2+d)-phi);
Enu2 = [Enu2;E];
end;
delta_eps1 = [];
delta_eps2 = [];
for x1 = linspace(0,-d,n);
E = fs*exp(x1/D);
% E = 1e-6*(exp(x1/D));
delta_eps1 = [delta_eps1;E];
end;
for x2 = linspace(-d,-dd,n);
E = fs*exp(x2/D);
% E = 1e-6*(exp(x2/D));
% E = 1e-6*exp(x2-1.5e-6/D);
delta_eps2 = [delta_eps2;E];
end;
A = (En1.*delta_eps1.*Enu1)*dx1+(En2.*delta_eps2.*Enu2)*dx2;
B = sum(A);
couple = omega*eps0/4*B;
E = pi.*Bv./rho.*couple^2*1e-21/3; %%%% Normalization
alpha = [alpha;E];
end;
fs = linspace(fs0/1e6,fx,100);
plot(fs,alpha,'color',rand(1,3))
xlabel('Acoustic Frequency (MHz)');
ylabel('\alpha_\nu (arb. units)');
title('Radiation Decay for TE0 mode for d = 1.0 \mum');
end