c ======================================================================
c     usage: partdemo1 <np>
c
c     Generates some particle data in sdf form.
c
c     Specifically, generates <np> particles randomly distributed on 
c     surface of sphere, and outputs <np> x 5 data where 2 fcns 
c     associated with (x,y,z) data are
c
c     1) max(abs(xj),abs(yj),abs(zj))
c     2) min(abs(xj),abs(yj),abs(zj))
c ======================================================================
      program      partdemo1

      implicit     none

      character*9  cdnm
      parameter  ( cdnm = 'partdemo1' )

		real*8       r8arg
      integer      iargc,  i4arg

      integer      np
      real*8       a  

      integer      npmax
      parameter  ( npmax = 10 000 000 )

      real*8       r(npmax,6)
      real*8       t

      if( iargc() .lt. 1  .or.  iargc() .gt. 2  ) go to 900

      np = i4arg(1,-1)
      if( np .eq. -1 ) go to 900
      if( np .lt. 1  .or.  np .gt. npmax ) then
         write(0,*) 'cdnm: np=', np
         write(0,*) 'cdnm: Exceeds range 1 to ', npmax
         go to 900
      end if
		a = r8arg(2,0.1d0)

		write(0,*) cdnm//': a=', a
      call genpart(r,np,a)
      t = 0.0d0 

      call gft_out_part_xyzf('partdemo1',t,r,np,3)

      stop

 900  continue
         write(0,*) 'usage: '//cdnm//' <np> [<a>]'
      stop
      
      end

      subroutine genpart(r,np,a)
         implicit     none

         real*8       drand48

         integer      np
         real*8       r(np,6)
			real*8       a

         integer      j
         real*8       pi, th1, th2, xj, yj, zj, nrm

         call srand(np)

			pi = 4.0d0 * atan(1.0d0)
         do j = 1 , np 
				th1 = pi * (2.0d0 * drand48() - 1.0d0)
				th2 = pi * (2.0d0 * drand48() - 1.0d0)
            xj = cos(th1) * (1.0d0 + a * cos(th2))
            yj = sin(th1) * (1.0d0 + a * cos(th2))
            zj = a * sin(th2)
            r(j,1) = xj
            r(j,2) = yj
            r(j,3) = zj
            r(j,4) = sin(th1)
            r(j,5) = sin(th2)
         end do

         return

      end