function [mu2,R,T] = trb(mu1,er1,er2,p)

% function [mu2,R,T] = trb(mu1,er1,er2,p)
% Brian Hornbuckle, 3-07-98.  Modified 4-13-06.  Modified 4-24-08.  
% trb = transmissivity, reflectivity for brightness temperatures
% Given an direction cosine incident on medium 2 from medium 1, this function
% calculates the transmission direction cosine and the reflectivity and
% transmissivity for brightness temperatures.  These expressions happen to be
% the same as the power transmissivity and reflectivity.
% Inputs:   mu1 = incident direction cosine = cos(theta_incident)
% 	        er1 = relative dielectric constant of medium 1 = complex
%	        er2 = relative dielectric constant of medium 2 = complex
%	        p = polarization, vertical=1 or horizontal=2
% Outputs:  mu2 = transmission direction cosine = cos(theta_transmitted)
%	        R = reflectivity for brightness temperature
%	        T = transmissivity for brightness temperature
% Note:  if mu1 is a vector, er1 and er2 should be scalars, and all outputs are vectors.
%        if mu1 is a scalar, er1 and er2 can be vectors, and all outputs will be vectors.
% Polarization orientation follows FSA convention, i.e. h=zxk, v=hxk, where
% k is the direction of propagation.
% Only the real part of each relative dielectric constant is considered,
% so medium 1 and 2 must be low loss for good results.  

er1 = real(er1);
n1 = sqrt(er1);				% index of refraction in medium 1
er2 = real(er2);			
n2 = sqrt(er2);				% index of refraction in 2

mu2 = sqrt(1-(n1./n2).^2.*sin(acos(mu1)).^2);	% transmission direction cosine
  
if p == 2					                    % horizontal polarization
  gamma = (n1.*mu1-n2.*mu2)./(n1.*mu1+n2.*mu2);	% Fresnel coefficients
  tau = 1+gamma;
    
elseif p == 1					                % vertical polarization
  gamma = (n2.*mu1-n1.*mu2)./(n2.*mu1+n1.*mu2);	% Fresnel coefficients
  tau = mu1./mu2.*(1-gamma);
    
else
  error('v or h not entered for polarization');
end;

R = abs(gamma).^2;
T = 1-R;