subroutine aptscaz (au, av, aw, cosph, bu, bv, bw, nerr)

ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c                             SUBROUTINE APTSCAZ
c
c     call aptscaz (au, av, aw, cosph, bu, bv, bw, nerr)
c
c     Version:  aptscaz  Updated    1991 August 13 17:30.
c               aptscaz  Originated 1991 July 25 16:30.
c
c     Author:   Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c     Purpose:  To randomly sample a unit direction vector b = (bu, bv, bw),
c               at an angle phi from a unit direction vector a = (au, av, aw)
c               in 3-D space, where cosph = cos (phi).
c
c     Input:    au, av, aw, cosph.
c
c     Output:   bu, bv, bw, nerr.
c
c     Glossary:
c
c     au,av,aw  Input    The u, v and w components of the unit vector "a".
c                          Coordinates u, v and w may be any three
c                          orthogonal coordinates, such as x, y and z.
c
c     bu,bv,bw  Output   The u, v and w components of a unit vector randomly
c                          sampled from a cone at angle phi from unit vector
c                          "a", where cosph = cos (phi).
c                          To get vector "b" isotropic in 3-D space, randomly
c                          sample cosph between -1.0 and 1.0.
c
c     cosph     Input    The cosine of the angle phi.
c
c     nerr      Output   Indicates an input error, if not 0.
c                          1 if cosph is not between -1.0 and 1.0.
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.

c.... Local variables.

c---- Absolute value of aw.
      common /laptscaz/ awabs
c---- Component u of a random vector.
      common /laptscaz/ cu
c---- Component v of a random vector.
      common /laptscaz/ cv
c---- Has magnitude 1.0, same sign as aw.
      common /laptscaz/ faxis
c---- Is sin (phi) / normalization factor.
      common /laptscaz/ fnorm
c---- Rotation factor.
      common /laptscaz/ frota
c---- Rotation factor.
      common /laptscaz/ frotb
c---- Normalization factor for cu, cv.
      common /laptscaz/ sumsq
cbugc***DEBUG begins.
cbug      common /laptscaz/ sumsqs
cbug 9901 format (/ 'aptscaz scattering by cosph around vector A:' /
cbug     &  '  cosph=   ',1pe22.14 /
cbug     &  '  au,av,aw=',1p3e22.14)
cbug      write (3, 9901) cosph, au, av, aw
cbugc***DEBUG ends.

c.... Test input for errors.

      nerr = 0

      if ((cosph .lt. -1.0) .or. (cosph .gt. 1.0)) then
        nerr = 1
        go to 210
      endif

c.... Randomly sample a unit vector at angle phi from the (+/-)w axis,
c....   and do rotation that moves (+/-)w axis onto vector "a".

  100 cu    = 1.0 - 2.0 * ranf( )
      cv    = 1.0 - 2.0 * ranf( )
      sumsq = cu**2 + cv**2
      if (sumsq .gt. 1.0) go to 100

      fnorm = sqrt ((1.0 - cosph**2) / sumsq)
      cu    = cu * fnorm
      cv    = cv * fnorm

      faxis = sign (1.0, aw)
      awabs = faxis * aw
      frota = cu * au + cv * av
      frotb = cosph - frota / (1.0 + awabs)

      bu    = cu + au * frotb
      bv    = cv + av * frotb
      bw    = cosph * aw - faxis * frota
cbugc***DEBUG begins.
cbug 9902 format (/ 'aptscaz results:')
cbug 9903 format ('  bu,bv,bw=',1p3e22.14 /
cbug     &  '  sumsq=   ',1pe22.14)
cbug      write ( 3, 9902)
cbug      sumsqs = bu**2 + bv**2 + bw**2
cbug      write ( 3, 9903) bu, bv, bw, sumsqs
cbugc***DEBUG ends.

  210 return

c.... End of subroutine aptscaz.      (+1 line.)
      end

UCRL-WEB-209832