%Code from Leathers et al., 2004, NRL report
% Irradiances Out The Top Of Scattering Medium
% Monte Carlo Matlab script for computing the intensities returned from 
% an isotoprically scattering semi-infinite medium
c = 1; % attenuation coefficient (1/m)
w0 = 0.8; % single-scattering albedo
N = 1e5; % number of photons to trace
ns = 100; % max number of scatters per photon
E = 0; % initialize detector
for i=1:N
z = 0; % initial position depth
muz = 0.1; % incident light direction (collimated)
w = 1; % initial photon weight
for j=1:ns
s = - log(rand)/c; % geometric path length
z = z + muz*s; % move photon
if (z<0) E=E+w; break, end % count photons leaving out top
w=w*w0; % absorb fraction of photon packet
muz= 1 - 2*rand; % isotopic scattering
end
end
E = E/N % normalize result