function [Rs,Rp,R] = intrc(n1,n2,k2,theta)
%
% [Rs,Rp,R] = intrc(n1,n2,k2,theta) 
%
% calculates the reflection coefficients for s- and p-polarized components 
% and circular polarized (or unpolarized) reflectivity for light incident 
% from a perfect dielectric medium with refractive index n1 to an absorbing 
% medium with complex refractive index n2+i*k2 with angle of incidence 
% theta (given in radians). 
%
% Input:  n1    - refractive index (medium 1)
%         n2    - refractive index (medium 2)
%         k2    - extinction coefficient (medium 2)
%         theta - angle of incidence
%
% Output: Rs    - reflection coefficient of s-polarized component 
%         Rp    - reflection coefficient of p-polarized component
%         R     - circular polarized (or unpolarized) reflectivity
%
% Last updated: 2009-02-24 (David Bergström)
%

p = sqrt(1/2*(sqrt((n2.^2-k2.^2-n1^2*sin(theta).^2).^2+4*n2.^2.*k2.^2)+(n2.^2-k2.^2-n1^2*sin(theta).^2)));
q = sqrt(1/2*(sqrt((n2.^2-k2.^2-n1^2*sin(theta).^2).^2+4*n2.^2.*k2.^2)-(n2.^2-k2.^2-n1^2*sin(theta).^2)));

Rs = ((n1*cos(theta)-p).^2+q.^2)./((n1*cos(theta)+p).^2+q.^2);
Rp = ((p-n1*sin(theta).*tan(theta)).^2+q.^2)./((p+n1*sin(theta).*tan(theta)).^2+q.^2).*Rs;
R = (Rs+Rp)./2;