      function BRDFbotm(rmup,sphip,umu,vphi)
c
c     user routine on file ..\maincode\BRDFminn.f
c
c     This routine is an example of a non-Lambertian BRDF, which
c     could be used in place of the default Lambertian BRDF on file
c     ..\maincode\BRDFLamb.f
c
c     This routine defines the Minnaert BRDF (originally used for modeling 
c     the lunar surface).  The BRDF is given by
c
c     BRDF(mu',phi',mu,phi) = (rho/pi) * (mu' * mu)**k
c
c     where rho is the "bottom reflectivity" rflbot (= R = Eu/Ed for a Lambertian
c     surface, which is k=0).  The parameter k is a measure of how
c     much this BRDF deviates from a Lambertian BRDF.

c*****Note:  For simplicity, this example of a non-Lambertian BRDF uses 
c     rho = R = rbottom as is used for the Lambertian bottom on file BRDFLamb.f.
c     However, the reflectivity rho is not equal to the irradiance reflectance R
c     for non-Lambertian bottoms.

c*****NOTE:  For non-Lambertian bottoms, increase the default values of
c     nsubmu = 1 and nsubphi = 1 in routine qaBRRF.  The default values are 
c     intended for use only with Lambertian BRDFs, which do not depend on angles.
c
      COMMON /CMISC/ imisc(30),fmisc(30) 

      data kall/0/
c
      save
c
c     initialization on first call:
c
      if(kall.eq.0) then
         pi = fmisc(1)
         ibotm = imisc(12)
         wavel = fmisc(13)
         rho = rbottom(ibotm,wavel)  ! initialize rbottom

c        use k = 0 to get a Lambertian surface:
         k = 2
         write(10,100) k
         kall = 1
      end if
c
c     subsequent calls start here:
c
      BRDFbotm = (rbottom(ibotm,wavel)/pi) * (rmup * umu)**k
c
      return
c
  100 format(/' The bottom BRDF is given by the Minnaert BRDF:'//
     1"   BRDF(mu',phi',mu,phi) = (rho/pi) * (mu' * mu)**k"//
     2'   with k =',i3)
      end     