<code>
 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


</code>


{{:radiation_decay_20151011_nov_11_.jpg|}}

{{:te1.jpg|}}

{{:rgb_plot.jpg|}}