
      SUBROUTINE CALBFK(phabkg, varbkg, breasc, region, brgion)

      IMPLICIT NONE

      REAL phabkg(*), varbkg(*), breasc(*), region(*), brgion(*)

C---
C Fakes a new background based on Poisson or Gaussian deviates from the 
c current one.
C PHABKG    I/R  Input with the expected background contribution to
C                the count rate, returned is the new rate
C VARBKG    I/R  The variance on the background - only set if Gaussian
C                randomization used
C BREASC    I    The background AREASCAL values
C REGION    I    The source BACKSCAL values
C BRGION    I    The background BACKSCAL values
C---

      REAL grasp, expos, gaunum
      INTEGER ifile, ipha

      REAL dgbtme
      INTEGER dgndst, dgnbnb, dgnbne
      LOGICAL dgqfst, dgdtps
      EXTERNAL dgndst, dgnbnb, dgnbne, dgqfst, dgbtme, dgdtps

C loop over files

      DO ifile = 1, dgndst()

         expos = dgbtme(ifile)

         IF ( expos .LE. 0 ) THEN
            DO ipha = dgnbnb(ifile), dgnbne(ifile)
               phabkg(ipha) = 0.0
            ENDDO
            GOTO 100
         ENDIF

c Calculate the background counts

         DO ipha = dgnbnb(ifile), dgnbne(ifile)
            grasp = expos * breasc(ipha)*brgion(ipha)/region(ipha)
            phabkg(ipha) = grasp * phabkg(ipha) 
         ENDDO

c If the randomization is required then get Poisson or Gaussian numbers based
c on these means

         IF ( dgqfst() ) THEN

            IF ( dgdtps(ifile, 'b') ) THEN

               CALL gtpoir( phabkg(dgnbnb(ifile)),
     &                      (dgnbne(ifile)-dgnbnb(ifile)+1) )

            ELSE

               DO ipha = dgnbnb(ifile), dgnbne(ifile)
                  CALL gtgaus(gaunum, 1)
                  grasp = expos * breasc(ipha)*brgion(ipha)/region(ipha)
                  varbkg(ipha) = grasp*grasp*varbkg(ipha)
                  phabkg(ipha) = phabkg(ipha) 
     &                + gaunum*sqrt(varbkg(ipha))
               ENDDO

            ENDIF

         ENDIF

c Convert back to count/s/cm^2
 
         DO ipha = dgnbnb(ifile), dgnbne(ifile)
            IF ( breasc(ipha) .NE. 0. ) THEN
               grasp = expos * breasc(ipha)
               phabkg(ipha) = phabkg(ipha) / grasp
            ENDIF
         ENDDO

         IF ( dgqfst() ) THEN
            DO ipha = dgnbnb(ifile), dgnbne(ifile)
               IF ( breasc(ipha) .NE. 0. ) THEN
                  grasp = expos * breasc(ipha)
                  varbkg(ipha) = varbkg(ipha) / grasp / grasp
               ENDIF
            ENDDO
         ENDIF


C end loop over files

 100     CONTINUE

      ENDDO

      RETURN
      END




