SUBROUTINE VGNFAC * * Module Number: 15.9.2.2 * * Module Name: gainfac * * Keyphrase: * ---------- * calculate the analog gain factor * * Description: * ------------ * This routine calculates gain factors of either high gain settings (GAIND = * 0, 2, 3) and low gain settings (VGAIND = 6, 7). For high settings the gain * factor, at the observed temperature and epoch, is directly obtained from * simultaneously obtained analog and digital data. For low settings the gain * factor is calculated by comparing measurements obtained with two adjacent * gain settings. These two observations have to have the same target, same * aperture (same detector), and roughly the same epoch and temperature. The * gain factor of the lower setting is then calculated from the ratio of these * two (analog) measurements and the formerly determined gain factor of the * higher setting . * Maximum number of rows in the input table is 2000. * * FORTRAN Name: VGNFAC.FOR * * Keywords of Accessed Files : * -------------------------- * Name I/O Description / Comments * * Modules Called: * --------------- * CDBS: * VGNGT, VGNPT, CDRFHV * SDAS: * UDMGET, UDMFRE, UMSPUT * OTHERS: * none * * History: * -------- * Version Date Author Description * 1 12-20-85 J.-C. Hsu Design and coding * 2 09-15-87 J.-C. Hsu F77 standard * 3 09-12-89 J.-C. Hsu error propagation *------------------------------------------------------------------------------- * *== local: * --size of the input data arrays INTEGER SIZE PARAMETER (SIZE = 2000) * --real array pointers: * --analog measurement and its standard * --deviation INTEGER CNTA, DCNTA, * --digital count rate and its standard * --deviation : CNTD, DCNTD, * --gain factor and its standard deviation : GNFAC, DGNFAC, * --temperature : TEMP, * --gain setting : GAIN, * --high voltage setting : HV, * --pointers to double precision array: * --epoch of the observation : EPOCH, * --pointers to integer arrays: * --detector ID : DETID, * --indices corresponding to hvfac : INDEX * --reference voltage REAL REFHV(5) * --number of input data points INTEGER NPTS, * --number of output points : N, NROWS, * --indices corresponding to hvfac : HIINDX, LOINDX, * --error status : STATUS, STATOK, STAT(20), * --loop index : I, J, K * --gain setting level keyword DOUBLE PRECISION A, D, G, DG, DA, DD, GHI, DGHI, : ALO, DALO, AHI, DAHI CHARACTER*4 LEVEL CHARACTER*6 CHAR6 * --aperture name CHARACTER*10 APER(SIZE) * --target names CHARACTER*20 TARGET(SIZE) * --column name of the temperature CHARACTER*16 TNAME * --output table name CHARACTER*128 OFILE * --error message context CHARACTER*130 CONTXT, MESS *==========================begin iraf77.inc (without INTEGER*2)================= * Include file for the iraf77 FORTRAN interface to the IRAF VOS * Get IRAF common into main program * LOGICAL MEMB(1) INTEGER MEMI(1) INTEGER MEML(1) REAL MEMR(1) DOUBLE PRECISION MEMD(1) COMPLEX MEMX(1) EQUIVALENCE (MEMB, MEMI, MEML, MEMR, MEMD, MEMX) COMMON /MEM/ MEMD * * File I/O access modes * INTEGER RDONLY PARAMETER (RDONLY = 1) INTEGER RDWRIT PARAMETER (RDWRIT = 2) INTEGER WRONLY PARAMETER (WRONLY = 3) INTEGER APPEND PARAMETER (APPEND = 4) INTEGER NEWFIL PARAMETER (NEWFIL = 5) INTEGER TMPFIL PARAMETER (TMPFIL = 6) INTEGER NEWCPY PARAMETER (NEWCPY = 7) INTEGER NEWIMG PARAMETER (NEWIMG = 5) * * codes for data types * INTEGER TYBOOL PARAMETER (TYBOOL = 1) INTEGER TYCHAR PARAMETER (TYCHAR = 2) INTEGER TYSHOR PARAMETER (TYSHOR = 3) INTEGER TYINT PARAMETER (TYINT = 4) INTEGER TYLONG PARAMETER (TYLONG = 5) INTEGER TYREAL PARAMETER (TYREAL = 6) INTEGER TYDOUB PARAMETER (TYDOUB = 7) INTEGER TYCPLX PARAMETER (TYCPLX = 8) INTEGER TYUSHT PARAMETER (TYUSHT = 11) INTEGER TYUBYT PARAMETER (TYUBYT = 12) * * TYTEXT is a special code for the iraf77 interface; it is not in the VOS * INTEGER TYTEXT PARAMETER (TYTEXT = 13) *========================end iraf77.inc========================================= *=========================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=========================================== *------------------------------------------------------------------------------- * * initialize number of input data points to zero * NROWS = 0 * * allocate dynamic memory for input data arrays * CALL UDMGET (SIZE, TYREAL, CNTA, STAT(1)) CALL UDMGET (SIZE, TYREAL, DCNTA, STAT(2)) CALL UDMGET (SIZE, TYREAL, CNTD, STAT(3)) CALL UDMGET (SIZE, TYREAL, DCNTD, STAT(4)) CALL UDMGET (SIZE, TYREAL, GNFAC, STAT(5)) CALL UDMGET (SIZE, TYREAL, DGNFAC, STAT(6)) CALL UDMGET (SIZE, TYREAL, TEMP, STAT(7)) CALL UDMGET (SIZE, TYREAL, GAIN, STAT(8)) CALL UDMGET (SIZE, TYREAL, HV, STAT(9)) * CALL UDMGET (SIZE, TYDOUB, EPOCH, STAT(10)) CALL UDMGET (SIZE, TYINT, DETID, STAT(11)) CALL UDMGET (SIZE, TYINT, INDEX, STAT(12)) * DO 10 I = 1, 12 IF (STAT(I) .NE. OK) THEN STATUS = ERRNUM(1) CONTXT = 'cannot allocate dynamic memory' GO TO 999 END IF 10 CONTINUE * * input parameters from parameter file and data from input table * CALL VGNGT (MEMR(GAIN), MEMR(CNTA), MEMR(DCNTA), MEMR(CNTD), : MEMR(DCNTD), MEMR(TEMP), MEMD(EPOCH), : APER, MEMR(GNFAC), MEMR(DGNFAC), TARGET, : MEMR(HV), MEMI(DETID), : NPTS, OFILE, LEVEL, TNAME, STATUS) * * number of input data points cannot be over SIZE * 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 * * read nominal high voltage setting table * CALL CDRFHV (REFHV, STATUS) IF (STATUS .NE. OK) THEN CONTXT = 'error getting nominal voltages' GO TO 999 END IF * * process high gain setting case * IF (LEVEL(1:1) .EQ. 'h') THEN DO 20 I = 0, NPTS-1 MESS = ' ' * * measurements must be made at the nominal voltage setting. if not, the data * point will be skipped * IF (MEMR(HV+I) .NE. REFHV(MEMI(DETID+I))) THEN MESS = 'data not obtained at nominal voltage' * * measured quantities must be positive * ELSE IF (MEMR(CNTA+I) .LE. 0. .OR. : MEMR(CNTD+I) .LE. 0.) THEN MESS = 'data are not positive' * * gain setting must be 0, 2, or 3 * ELSE IF (MEMR(GAIN+I) .NE. 0. .AND. MEMR(GAIN+I) .NE. : 2. .AND. MEMR(GAIN+I) .NE. 3.) THEN MESS = 'incorrect gain setting' END IF * * if a row in the input table has error, issue a message and skip the row * IF (MESS .NE. ' ') THEN WRITE (CHAR6, '(I6)') I+1 CONTXT = 'VGNFAC: at input point #' // CHAR6 // : ' : ' // MESS CALL UMSPUT (CONTXT, DEST, PRIO, STATOK) ELSE NROWS = NROWS + 1 N = NROWS - 1 MEMI(INDEX+N) = I + 1 * * calculate gain factor and its error estimate * A = MEMR(CNTA+I) D = MEMR(CNTD+I) DA = MEMR(DCNTA+I) DD = MEMR(DCNTD+I) G = A / D DG = G * SQRT((DA/A)**2 + (DD/D)**2) MEMR(GNFAC+N) = REAL(G) MEMR(DGNFAC+N) = REAL(DG) END IF 20 CONTINUE * * process low gain setting case * the input table for low gain setting case must have pairs of input data, * each pair consists of analog observations made toward the same target at * different gain settings (must be 3 and 6, or 6 and 7) * ELSE IF (LEVEL(1:1) .EQ. 'l') THEN DO 30 I = 1, NPTS / 2 J = I * 2 - 2 K = I * 2 - 1 MESS = ' ' * IF (MEMI(DETID+J) .NE. MEMI(DETID+K)) THEN MESS = 'different detectors' * * measurements must be made at the nominal voltage setting. if not, the data * point will be skipped * ELSE IF (MEMR(HV+J) .NE. REFHV(MEMI(DETID+J)) .OR. : MEMR(HV+K) .NE. REFHV(MEMI(DETID+K))) THEN MESS = 'voltage is different from the ' : // 'nominal high voltage' ELSE IF (TARGET(J+1) .NE. TARGET(K+1)) THEN MESS = 'different targets, target 1 = ' // : TARGET(J+1) // 'target 2 = ' // : TARGET(K+1) ELSE IF (APER(J+1) .NE. APER(K+1)) THEN MESS = 'different apertures, aperture 1 = ' // : APER(J+1) // 'aperture 2 = ' // APER(K+1) ELSE IF (MEMR(CNTA+J) .LE. 0. .OR. MEMR(CNTA+K) .LE. : 0.) THEN MESS = 'non-positive count rate' ELSE IF (MEMR(GAIN+J) .EQ. 3 .AND. MEMR(GAIN+K) .EQ. 6 : .OR. MEMR(GAIN+J) .EQ. 6 .AND. MEMR(GAIN+K) : .EQ. 7) THEN LOINDX = K HIINDX = J ELSE IF (MEMR(GAIN+J) .EQ. 6 .AND. MEMR(GAIN+K) .EQ. 3 : .OR. MEMR(GAIN+J) .EQ. 7 .AND. MEMR(GAIN+K) : .EQ. 6) THEN LOINDX = J HIINDX = K ELSE MESS = 'incorrect gain settings' END IF * IF (MESS .NE. ' ') THEN WRITE (CHAR6, '(I6)') I CONTXT = ' at input data pair #' // CHAR6 // ' : ' : // MESS CALL UMSPUT (CONTXT, DEST, PRIO, STATOK) ELSE NROWS = NROWS + 1 N = NROWS - 1 MEMI(INDEX+N) = LOINDX + 1 * * calculate gain factor and its error estimate * GHI = MEMR(GNFAC+HIINDX) DGHI = MEMR(DGNFAC+HIINDX) ALO = MEMR(CNTA+LOINDX) AHI = MEMR(CNTA+HIINDX) DALO = MEMR(DCNTA+LOINDX) DAHI = MEMR(DCNTA+HIINDX) G = GHI * ALO / AHI DG = G * SQRT((DALO/ALO)**2 + (DAHI/AHI)**2 + : (DGHI/GHI)**2) MEMR(GNFAC+N) = REAL(G) MEMR(DGNFAC+N) = REAL(DG) END IF 30 CONTINUE ELSE STATUS = ERRNUM(1) CONTXT = 'unknown gain setting case ' // LEVEL GO TO 999 END IF * * put the result in output table * CALL VGNPT (OFILE, MEMR(GAIN), MEMR(GNFAC), MEMR(DGNFAC), : MEMR(TEMP), MEMR(HV), MEMD(EPOCH), MEMI(INDEX), : NROWS, TNAME, MEMI(DETID), APER, TARGET, : STATUS) IF (STATUS .NE. OK) THEN GO TO 999 END IF * * free the dynamic memory * CALL UDMFRE (CNTA, TYREAL, STAT(1)) CALL UDMFRE (DCNTA, TYREAL, STAT(2)) CALL UDMFRE (CNTD, TYREAL, STAT(3)) CALL UDMFRE (DCNTD, TYREAL, STAT(4)) CALL UDMFRE (GNFAC, TYREAL, STAT(5)) CALL UDMFRE (DGNFAC, TYREAL, STAT(6)) CALL UDMFRE (TEMP, TYREAL, STAT(7)) CALL UDMFRE (GAIN, TYREAL, STAT(8)) CALL UDMFRE (HV, TYREAL, STAT(9)) * CALL UDMFRE (EPOCH, TYDOUB, STAT(10)) CALL UDMFRE (DETID, TYINT, STAT(11)) CALL UDMFRE (INDEX, TYINT, STAT(12)) * DO 40 I = 1, 12 IF (STAT(I) .NE. OK) THEN STATUS = ERRNUM(1) CONTXT = 'cannot free dynamic memory' GO TO 999 END IF 40 CONTINUE * STATUS = OK GO TO 1000 * * write error message * 999 MESS = 'VGNFAC: ' // CONTXT CALL UMSPUT (MESS, DEST, PRIO, STATOK) * 1000 RETURN END