        REAL*8 FUNCTION MODEL(MASS, IFILE, DPARAM, NPAR, FLAG)
C--
C  Returns the probability, normalised to the total number of 
C  events expected for the model.
C
C       MASS    Mass value of datum
C       IFILE   Channel number (1 or 2)
C       DPARAM  REAL*8 list of fit parameters
C       NPAR    Number of fit parameters
C       FLAG    = 0 to calculate probability (no other function at
C                   the moment)
C--
        IMPLICIT NONE
C
        REAL*8          MASS
        REAL*8          DPARAM(*)
        INTEGER         IFILE, NPAR, FLAG
C--
C  Switch on FLAG
C--
        MODEL = 0.0D0
C
        IF (FLAG.EQ.0) THEN
C--
C  Calculate probability.  For this test, it is just a gaussian
C  signal on a 4rth-order polynomial background:
C
C  P(x) = p4 + p5*m + p6*m**2 + p7*m**3 + p8*m**4
C
C  This is very crude but shows the basic philosophy.  Note 
C  normalisations should be defined carefully...
C--
          MODEL = DPARAM(4) +
     >            DPARAM(5) * MASS +
     >            DPARAM(6) * MASS**2 +
     >            DPARAM(7) * MASS**3 +
     >            DPARAM(8) * MASS**4
C
C  Force this sucker to be 
          IF (MODEL.LT.0.0D0) THEN   
C  +ve definite....
            MODEL = 0.0D0            
          ENDIF
C
C  Add signal
          IF (IFILE.EQ.1) THEN  
            MODEL = MODEL + 
     >         DPARAM(1)*DEXP(-0.5D0*((MASS-DPARAM(2))/DPARAM(3))**2)/
     >                  (DPARAM(3)*DSQRT(2.0D0*3.14159D0))
          ENDIF
C
        ELSE
          WRITE(*,*) 'MODEL:  FLAG value is incorrect - HALT!'
          STOP
        ENDIF
C
        RETURN
        END
        REAL*8 FUNCTION INTMODEL(IFILE, DPARAM, NPAR)
C--
C  Returns the integral of the above model function, performed
C  crudely by doing a brute-force integration of the MODEL function.
C
C       IFILE   Channel number (1 or 2)
C       DPARAM  REAL*8 list of fit parameters
C       NPAR    Number of fit parameters
C--
        IMPLICIT NONE
C
        REAL*8          DPARAM(*)
        INTEGER         IFILE, NPAR
C--
C  Region of integration.
C--
        REAL*8          LOWEDG 
        PARAMETER       (LOWEDG=0.6D0)
        REAL*8          HIGEDG
        PARAMETER       (HIGEDG=1.2D0)
        REAL*8           DELTAX
        PARAMETER       (DELTAX = HIGEDG-LOWEDG)
        INTEGER          NSTEPS
        PARAMETER       (NSTEPS=1000)
C--
C  Local variables
C--
        REAL*8          MASS
        REAL*8          BINSIZ
        INTEGER         I
C--
C  External functions
C--
        REAL*8          MODEL
C--
C  Initialise
C--
        INTMODEL = 0.0D0
C--
C  Step over the mass interval.
C--
        BINSIZ = DELTAX/FLOAT(NSTEPS)
        DO I=1,NSTEPS
          MASS = LOWEDG + (FLOAT(I)-0.5D0)*BINSIZ
          INTMODEL = INTMODEL + 
     >          MODEL(MASS, IFILE, DPARAM, NPAR, 0)*BINSIZ
        ENDDO
C
        RETURN
        END
