NewsEtcto next pagearticle front page Dec. issue NewsEtc. Contents Home | Monte Carlo:Photon spins into new trajectory due to scatteringNewsEtc., Dec., 1998. Steven Jacques, Oregon Medical Laser Center |
The remaining photon weight is now scattered into a new trajectory. If the anisitropy (g) equals 0, then the scattering is isotropic. If g > 0, then the program uses the Henyey-Greenstein scattering function originally introduced by Henyey and Greenstein to mimic the scattering of light from distant stars by galactic dust. The same function approximately mimics the scattering function experimentally observed in biological tissues (see HG function).
The four sections of SPIN are:
/**** SPIN Scatter photon into new trajectory defined by theta and psi. Theta is specified by cos(theta), which is determined based on the Henyey-Greenstein scattering function. Convert theta and psi into cosines ux, uy, uz. *****/ /* Sample for costheta */ rnd = RandomNum; if (g == 0.0) costheta = 2.0*rnd - 1.0; else { double temp = (1.0 - g*g)/(1.0 - g + 2*g*rnd); costheta = (1.0 + g*g - temp*temp)/(2.0*g); } sintheta = sqrt(1.0 - costheta*costheta); /* sqrt() is faster than sin(). */ /* Sample psi. */ psi = 2.0*PI*RandomNum; cospsi = cos(psi); if (psi < PI) sinpsi = sqrt(1.0 - cospsi*cospsi); /* sqrt() is faster than sin(). */ else sinpsi = -sqrt(1.0 - cospsi*cospsi); /* New trajectory. */ if (1 - fabs(uz) <= ONE_MINUS_COSZERO) { /* close to perpendicular. */ uxx = sintheta * cospsi; uyy = sintheta * sinpsi; uzz = costheta * SIGN(uz); /* SIGN() is faster than division. */ } else { /* usually use this option */ temp = sqrt(1.0 - uz * uz); uxx = sintheta * (ux * uz * cospsi - uy * sinpsi) / temp + ux * costheta; uyy = sintheta * (uy * uz * cospsi + ux * sinpsi) / temp + uy * costheta; uzz = -sintheta * cospsi * temp + uz * costheta; } /* Update trajectory */ ux = uxx; uy = uyy; uz = uzz;