| |
| subroutine thom (u,v,w,mySeed) |
| c*********************************************************************** |
| c selects a scattering angle for a thomson collision and then |
| c generates new direction cosines. |
| c*********************************************************************** |
| implicit double precision (a-h,o-z) |
| integer mySeed(*) |
| real*8 ranf |
| external ranf |
| dimension csn (128) |
| |
| data(csn(i),i=1,67)/ |
| | 1. d0,.98952871d0,.9789466 d0,.96825135d0,.95744053d0, |
| x .94651168d0, |
| | .93546224d0,.92428964d0,.91299119d0,.90156417d0,.89000577d0, |
| x .87831314d0, |
| | .86648335d0,.8545134d0, .84240023d0,.83014072d0,.81773168d0, |
| x .80516986d0, |
| | .79245197d0,.77957465d0,.76653448d0,.75332801d0,.73995175d0, |
| x .72640218d0, |
| | .71267574d0,.69876886d0,.68467798d0,.67039954d0,.65592999d0, |
| x .64126583d0, |
| | .62640361d0,.61133997d0,.59607164d0,.58059548d0,.5649085 d0, |
| x .5490079d0, |
| | .53289112d0,.51655582d0,.5d0, .48322199d0,.46622053d0, |
| x .44899477d0, |
| | .43154441d0,.41386968d0,.39597143d0,.3778512 d0,.35951125d0, |
| x .34095466d0, |
| | .32218536d0,.30320817d0,.2840289 d0,.26465437d0,.24509241d0, |
| x .22535195d0, |
| | .20544297d0,.18537655d0,.1651648 d0,.14482089d0,.12435893d0, |
| x .10379394d0, |
| | .083141761d0,.062418936d0,.041642596d0,.020830322d0,.0d0, |
| | -.020830322d0,-.041642596d0/ |
| |
| data(csn(i),i=68,128)/ |
| | -.062418936d0,-.083141761d0,-.10379394d0,-.12435893d0, |
| x -.14482089d0, |
| | -.1651648d0,-.18537655d0,-.20544297d0,-.22535195d0, |
| x -.24509241d0, |
| | -.26465437d0,-.2840289d0,-.30320817d0,-.32218536d0, |
| x -.34095466d0,-.35951125d0, |
| | -.3778512d0,-.39597143d0,-.41386968d0,-.43154441d0, |
| x -.44899477d0,-.46622053d0, |
| | -.48322199d0,-.5d0,-.51655582d0,-.53289112d0,-.5490079d0, |
| x -.5649085d0, |
| | -.58059548d0,-.59607164d0,-.61133997d0,-.62640361d0, |
| x -.64126583d0,-.65592999d0, |
| | -.67039954d0,-.68467798d0,-.69876886d0,-.71267574d0, |
| x -.72640218d0,-.73995175d0, |
| | -.75332801d0,-.76653448d0,-.77957465d0,-.79245197d0, |
| x -.80516986d0,-.81773168d0, |
| | -.83014072d0,-.84240023d0,-.8545134d0,-.86648335d0, |
| x -.87831314d0, |
| | -.89000577d0,-.90156417d0,-.91299119d0, |
| | -.92428964d0,-.93546224d0, |
| | -.94651168d0,-.95744053d0,-.96825135d0,-.9789466d0, |
| x -.98952871d0/ |
| |
| index = ranf(mySeed)*128 + 1 |
| sctang = csn (index) |
| t3 = 6.2831853d0 * ranf(mySeed) |
| t4 = sqrt(1.0d0 - sctang*sctang) |
| t5 = sqrt(1.0d0 - w*w) |
| |
| c--------------------update scattered photon direction |
| c--------------------(from tid-26607) |
| temp1 = t4/t5 |
| temp2 = sin (t3) |
| temp3 = cos (t3) |
| temp4 = w*temp3 |
| t1 = u |
| |
| c--------------------new u, v, w |
| |
| u = u*sctang + temp1*(u*temp4 - v*temp2) |
| v = v*sctang + temp1*(v*temp4 + t1*temp2) |
| w = w*sctang - (t4*t5)*temp3 |
| |
| return |
| end |
