SUBROUTINE VHVFAC
*
*  Module Number: 15.9.2.1
*
*  Module Name: voltfac
*
*  Keyphrase:
*  ----------
*  Calculate high voltage factors at different detector high voltage settings
*
*  Description:
*  ------------
*  Calculate high voltage factors, i. e. count rate relative to that obtained
*  at the nominal (reference) high voltage setting.  Input data must have the
*  same target, use the same aperture (same detector), and are
*  roughly obtained at the same temperature/epoch.
*  Maximum number of input data points is 2000.
*
*  FORTRAN Name: VHVFAC.FOR
*
*  Keywords of Accessed Files :
*  --------------------------
*  Name                         I/O     Description / Comments
*
*  Modules Called:
*  ---------------
*  CDBS:
*       VHVFGT, VHVFPT, CDRFHV
*  SDAS:
*       UMSPUT
*  OTHERS:
*       none
*
*  History:
*  --------
*  Version      Date        Author          Description
*     1       12-20-85     J.-C. Hsu        Design and coding
*     2       09-15-85     J.-C. Hsu        F77 standard
*     3       09-12-89     J.-C. Hsu        error propagation
*-------------------------------------------------------------------------------
*
*== local:
*                                       --size of input data arrays
        INTEGER                 SIZE
        PARAMETER               (SIZE = 2000)
*                                       --observed quantity (e. g. count rate)
*                                       --and its standard deviation
        REAL                    COUNT(SIZE), DCOUNT(SIZE),
*                                       --high voltage factor and its standard
*                                       --deviation
     :                          HVFAC(SIZE), DHVFAC(SIZE),
*                                       --temperature
     :                          TEMP(SIZE),
*                                       --reference voltage
     :                          REFHV(5),
*                                       --the nominal high voltage setting, at
*                                       --which the high voltage factor is ONE
     :                          HV0,
*                                       --high voltage setting
     :                          HV(SIZE),
*                                       --summation of data taken at reference
*                                       --voltage
     :                          SUM, SUMWT, SUMSQ, WT, REFCNT, ERR
*                                       --epoch of the observation
        DOUBLE PRECISION        EPOCH(SIZE)
*                                       --number of input data points
        INTEGER                 NPTS,
*                                       --number of output points
     :                          NP, NREF,
*                                       --indices corresponding to hvfac
     :                          INDEX(SIZE),
*                                       --detector ID
     :                          DETID,
*                                       --error status
     :                          STATUS, STATOK,
*                                       --loop index
     :                          I, J
        CHARACTER*6             CHAR6
*                                       --data type (analog or digital)
        CHARACTER*7             TYPE
*                                       --aperture name
        CHARACTER*10            APERT
*                                       --target names
        CHARACTER*20            TARGET
*                                       --column name of the temperature
        CHARACTER*16            TNAME
*                                       --output table name
        CHARACTER*128           OFILE
*                                       --error message context
        CHARACTER*130           CONTXT, MESS
*=========================begin hsp.inc=========================================
*                                       --status return code
        INTEGER         OK, ERRNUM(20)
*                                       --message destination and priority
        INTEGER         DEST, PRIO
        DATA OK /0/
        DATA ERRNUM /701, 702, 703, 704, 705, 706, 707, 708, 709, 710,
     :               711, 712, 713, 714, 715, 716, 717, 718, 719, 720/
        DATA DEST, PRIO /1, 0/
*=========================end hsp.inc===========================================
*-------------------------------------------------------------------------------
*
*  input parameters from parameter file and data from input table
*
        CALL VHVFGT (    COUNT, DCOUNT, TEMP, HV, EPOCH, DETID, TYPE, 
     :               APERT, TARGET, NPTS, OFILE, TNAME, STATUS)
        IF (STATUS .NE. OK .OR. NPTS .GT. SIZE) THEN
            WRITE (CHAR6, '(I6)') SIZE
            CONTXT = 'cannot get input parameters/data or more than '
     :               // CHAR6 // ' input points'
            GO TO 999
        END IF
*
*  find the nominal high voltage of the detector
*
        CALL CDRFHV (REFHV,     STATUS)
        IF (STATUS .NE. OK) THEN
            CONTXT = 'error getting nominal voltages'
            GO TO 999
        END IF
*
*  average count rates taken at reference voltage and calculate high
*  voltage factors
*
        HV0 = REFHV (DETID)
*
*  initialization
*
        SUM = 0.
        SUMWT = 0.
        SUMSQ = 0.
        NREF = 0
*
*  pick out the measurements made at the reference voltage
*
        DO 10 I = 1, NPTS
            IF (HV(I) .EQ. HV0) THEN
                WT = 1. / DCOUNT(I)**2
                SUMWT = SUMWT + WT
                SUM = SUM + WT * COUNT(I)
                SUMSQ = SUMSQ + WT * COUNT(I)**2
                NREF = NREF + 1
                J = I
            END IF
   10   CONTINUE
*
*  check if there is any data taken at the nominal voltage
*
        IF (NREF .LE. 0) THEN
            STATUS = ERRNUM(1)
            CONTXT = 'no data at the nominal voltage setting'
            GO TO 999
        END IF
*
*  calculate mean and mean error of count rates at reference voltage
*
        REFCNT = SUM / SUMWT
        IF (NREF .EQ. 1) THEN
            ERR = DCOUNT(J)
        ELSE
            ERR = SQRT((SUMSQ/SUMWT - REFCNT**2) / REAL(NREF-1))
        END IF
*
*  calculate high voltage factors
*
        NP = 0
        DO 20 I = 1, NPTS
            IF (HV(I) .NE. HV0) THEN
                NP = NP + 1
                INDEX(NP) = I
                HVFAC(NP) = COUNT(I) / REFCNT
                DHVFAC(NP) = HVFAC(NP) * SQRT((ERR / REFCNT)**2 +
     :                          (DCOUNT(I) / COUNT(I))**2)
            END IF
   20   CONTINUE
*
*  check if there is any data points besides reference voltage
*
        IF (NP .EQ. 0) THEN
            STATUS = ERRNUM(2)
            CONTXT = 'no data point besides reference voltage'
            GO TO 999
        END IF
*
*  put the result in output table
*
        CALL VHVFPT (OFILE, HV, HVFAC, DHVFAC, TEMP, EPOCH, INDEX,
     :               NP, TYPE, TNAME, DETID, APERT, TARGET,    STATUS)
        IF (STATUS .NE. OK) THEN
            GO TO 999
        END IF
*
        STATUS = OK
        GO TO 1000
*
*  write error message
*
  999   MESS = 'VHVFAC: ' // CONTXT
        CALL UMSPUT (MESS, DEST, PRIO,     STATOK)
*
 1000   RETURN
        END