SUBROUTINE rddtg1(unit, ifile, filenm, ftype, ispec, itype, & phabuf, ndata, nungc, qsysdt, llchan, ldchan, & hlchan, hdchan, ierr) INCLUDE '../../inc/xspec.inc' INTEGER nungc REAL phabuf(nungc,5) INTEGER unit, ndata INTEGER llchan, ldchan, hlchan, hdchan INTEGER ierr, itype, ftype, ispec, ifile CHARACTER filenm*(*) LOGICAL qsysdt c kaa 1/24/95 c XSPEC subroutine called by RDDTOG to read from a spectral data file c in OGIP FITS format. The file is assumed to be opened and stationed c at the correct extension. c c unit i i: The I/O unit for the PHA file c ifile i i: The dataset number c filenm c i: The dataset filename c ftype i i: The file type (I or II) c ispec i i: If type II then set to spectrum to read c itype i i: The type of data, c 0 = counts, c 1 = counts/sec, c 2 = counts/sec/area c phabuf r r: The data, the first row of which contains c the data, the second contains the variance, c the third the systematics values, the c fourth the backscal values and the fifth c the areascal values. c ndata i r: The actual number of data points read c nungc i i: Max number of ungrouped channels c qsysdt l r: there was a systematics packages. c llchan i r: Lowest legal data channel c ldchan i r: Lowest actual data channel c hlchan i r: Highest legal data channel c hdchan i r: Highest actual data channel c ierr i r: Error flag (0= no error) REAL syserr, junk, backscal, areascal INTEGER index, filestat, ibin, ijunk, nsize CHARACTER comment*80, context*80, wrtstr*255 LOGICAL qanyf, qpoiss INTEGER lenact EXTERNAL lenact filestat = 0 context = 'RDDTG1 : FITSIO error' c Find the number of stored data points. For the type 2 files we need to c use different methods depending on whether we have variable or constant c length vectors. IF ( ftype .EQ. 1 ) THEN CALL ftgkyj(unit, 'NAXIS2', ndata, comment, filestat) context = 'RDDTG1: Failed to read NAXIS2' IF (filestat .NE. 0) GOTO 100 ELSEIF ( ftype .EQ. 2 ) THEN IF ( itype .EQ. 0 ) THEN CALL ftgcno(unit, .FALSE., 'COUNTS', index, filestat) context = 'RDDTG1: Failed to find COUNTS column' IF (filestat .NE. 0) GOTO 100 ELSEIF ( itype .EQ. 1 ) THEN CALL ftgcno(unit, .FALSE., 'RATE', index, filestat) context = 'RDDTG1: Failed to find RATE column' IF (filestat .NE. 0) GOTO 100 ENDIF CALL ftgdes(unit, index, ispec, ndata, ijunk, filestat) IF ( filestat .NE. 0 ) THEN filestat = 0 CALL ftgtdm(unit, index, 1, ijunk, ndata, filestat) context = & 'RDDTG1: Failed to find size of COUNTS or RATE column' IF (filestat .NE. 0) GOTO 100 ENDIF ENDIF c Check whether the TLMIN1 and TLMAX1 keywords are set. If so then c read the start and end channels in the data. CALL ftgkyj(unit, 'TLMIN1', llchan, comment, filestat) IF (filestat .NE. 0) THEN filestat = 0 llchan = -1 ENDIF CALL ftgkyj(unit, 'TLMAX1', hlchan, comment, filestat) IF (filestat .NE. 0) THEN filestat = 0 hlchan = -1 ENDIF IF (llchan .NE. -1 .AND. hlchan .NE. -1) THEN IF ( ftype .EQ. 1 ) THEN CALL ftgcvj(unit, 1, 1, 1, 1, 0, ldchan, qanyf, filestat) CALL ftgcvj(unit, 1, ndata, 1, 1, 0, hdchan, qanyf, & filestat) ELSEIF ( ftype .EQ. 2 ) THEN CALL ftgcno(unit, .FALSE., 'CHANNEL', index, filestat) IF (filestat .EQ. 0 .AND. index .NE. 0 ) THEN CALL ftgcvj(unit, index, ispec, 1, 1, 0, ldchan, qanyf, & filestat) CALL ftgcvj(unit, index, ispec, ndata, 1, 0, hdchan, & qanyf, filestat) context = 'RDDTG1: Failed to read CHANNEL data' IF (filestat .NE. 0) GOTO 100 ENDIF filestat = 0 ENDIF ENDIF c If counts data then extract into phabuf IF (itype .EQ. 0) THEN CALL ftgcno(unit, .FALSE., 'COUNTS', index, filestat) IF (filestat .NE. 0) THEN context = 'RDDTG1: Failed to find COUNTS column' GOTO 100 ENDIF IF ( ftype .EQ. 1 ) THEN CALL ftgcve(unit, index, 1, 1, ndata, 1.e-32, phabuf, qanyf, & filestat) ELSEIF ( ftype .EQ. 2 ) THEN CALL ftgcve(unit, index, ispec, 1, ndata, 1.e-32, phabuf, & qanyf, filestat) ENDIF context = 'RDDTG1: Failed to read COUNTS data' IF (filestat .NE. 0) GOTO 100 c otherwise read the `pha per sec' data ELSEIF (itype .EQ. 1) THEN CALL ftgcno(unit, .FALSE., 'RATE', index, filestat) IF (filestat .NE. 0) THEN context = 'RDDTG1: Failed to find RATE column' GOTO 100 ENDIF IF ( ftype .EQ. 1 ) THEN CALL ftgcve(unit, index, 1, 1, ndata, 1.e-32, phabuf, & qanyf, filestat) ELSEIF ( ftype .EQ. 2 ) THEN CALL ftgcve(unit, index, ispec, 1, ndata, 1.e-32, phabuf, & qanyf, filestat) ENDIF context = 'RDDTG1: Failed to read RATE data' IF (filestat .NE. 0) GOTO 100 ENDIF c check for the Poisson error keyword CALL ftgkyl(unit, 'POISSERR', qpoiss, comment, filestat) IF (filestat .NE. 0) THEN CALL xwrite('POISSERR keyword not found, assuming FALSE', & 10) qpoiss = .FALSE. filestat = 0 ENDIF c if the Poisson error is set but the data was RATE then write a warning c message IF ( qpoiss .AND. itype .EQ. 1 ) THEN CALL xwrite('*Error*: POISSERR=T used with RATE column',10) CALL xwrite(' Data as RATE requires errors in STAT_ERR', & 10) ENDIF c if the Poisson error keyword is set then errors are square roots c otherwise load them from the stat_err column IF (qpoiss) THEN DO ibin = 1, ndata phabuf(ibin,2) = SQRT(MAX(phabuf(ibin,1),0.)) ENDDO CALL ftgcno(unit, .FALSE., 'STAT_ERR', index, filestat) IF ( filestat .EQ. 0 .AND. index .NE. 0 ) THEN WRITE(wrtstr,'(a,a,a)') 'Warning : ', & filenm(:lenact(filenm)), & ' has both POISSERR=T and STAT_ERR' CALL xwrite(wrtstr, 10) CALL xwrite('Using Poisson errors', 10) ENDIF filestat = 0 ELSE CALL ftgcno(unit, .FALSE., 'STAT_ERR', index, filestat) IF (filestat .NE. 0) THEN filestat = 0 CALL ftgkye(unit, 'STAT_ERR', phabuf(1,2), & comment, filestat) context = 'RDDTG1: Failed to find STAT_ERR' IF ( filestat .NE. 0 ) GOTO 100 DO ibin = 2, ndata phabuf(ibin,2) = phabuf(1,2) ENDDO ELSE IF ( ftype .EQ. 1 ) THEN CALL ftgcve(unit, index, 1, 1, ndata, 1.e-32, & phabuf(1,2), qanyf, filestat) ELSEIF ( ftype .EQ. 2 ) THEN CALL ftgcve(unit, index, ispec, 1, ndata, 1.e-32, & phabuf(1,2), qanyf, filestat) ENDIF context = 'RDDTG1: Failed to read STAT_ERR data' IF (filestat .NE. 0) GOTO 100 CALL ftgkye(unit, 'STAT_ERR', junk, comment, filestat) IF ( filestat .EQ. 0 ) THEN WRITE(wrtstr, '(a,i4,a)') 'Warning : file ', ifile, & ' has both STAT_ERR keyword and column' CALL xwrite(wrtstr, 10) CALL xwrite('Using data from STAT_ERR column', 10) ENDIF filestat = 0 ENDIF ENDIF c check for systematic errors qsysdt = .FALSE. CALL ftgcno(unit, .FALSE., 'SYS_ERR', index, filestat) IF ( filestat .NE. 0 ) THEN filestat = 0 CALL ftgkye(unit, 'SYS_ERR', syserr, comment, filestat) IF (filestat .EQ. 0 .AND. syserr .NE. 0. ) THEN DO ibin = 1, ndata phabuf(ibin,3) = syserr ENDDO qsysdt = .TRUE. ENDIF filestat = 0 ELSEIF ( filestat .EQ. 0 .AND. index .GT. 0 ) THEN IF ( ftype .EQ. 1 ) THEN CALL ftgcve(unit, index, 1, 1, ndata, 1.e-32, & phabuf(1,3), qanyf, filestat) ELSEIF ( ftype .EQ. 2 ) THEN CALL ftgcve(unit, index, ispec, 1, ndata, 1.e-32, & phabuf(1,3), qanyf, filestat) ENDIF context = 'RDDTG1: Failed to read SYS_ERR data' IF ( filestat .NE. 0 ) GOTO 100 qsysdt = .TRUE. CALL ftgkye(unit, 'SYS_ERR', junk, comment, filestat) IF ( filestat .EQ. 0 ) THEN WRITE(wrtstr, '(a,i4,a)') 'Warning : file ', ifile, & ' has both SYS_ERR keyword and column' CALL xwrite(wrtstr, 10) CALL xwrite('Using data from SYS_ERR column', 10) ENDIF filestat = 0 ENDIF c read the BACKSCAL value(s). Note that the type II case is more c complicated now because we need to have the option of a vector c BACKSCAL on each row. CALL ftgcno(unit, .FALSE., 'BACKSCAL', index, filestat) IF ( filestat .NE. 0 ) THEN filestat = 0 CALL ftgkye(unit, 'BACKSCAL', backscal, comment, filestat) IF (filestat .EQ. 0 .AND. backscal .NE. 0. ) THEN DO ibin = 1, ndata phabuf(ibin,4) = backscal ENDDO ELSE WRITE(wrtstr, '(a,i4,a)') 'Warning : file ', ifile, & ' has no BACKSCAL keyword or column' CALL xwrite(wrtstr, 10) CALL xwrite('Assuming BACKSCAL = 1', 10) DO ibin = 1, ndata phabuf(ibin,4) = 1.0 ENDDO ENDIF filestat = 0 ELSEIF ( filestat .EQ. 0 .AND. index .GT. 0 ) THEN IF ( ftype .EQ. 1 ) THEN CALL ftgcve(unit, index, 1, 1, ndata, 1.e-32, & phabuf(1,4), qanyf, filestat) ELSEIF ( ftype .EQ. 2 ) THEN CALL ftgdes(unit, index, ispec, nsize, ijunk, filestat) IF ( filestat .NE. 0 ) THEN filestat = 0 CALL ftgtdm(unit, index, 1, ijunk, nsize, filestat) ENDIF IF ( nsize .GT. 1 ) THEN CALL ftgcve(unit, index, ispec, 1, ndata, 1.e-32, & phabuf(1,4), qanyf, filestat) ELSE CALL ftgcve(unit, index, ispec, 1, 1, 1.e-32, & backscal, qanyf, filestat) DO ibin = 1, ndata phabuf(ibin,4) = backscal ENDDO ENDIF ENDIF context = 'RDDTG1: Failed to read BACKSCAL data' IF ( filestat .NE. 0 ) GOTO 100 CALL ftgkye(unit, 'BACKSCAL', junk, comment, filestat) IF ( filestat .EQ. 0 ) THEN WRITE(wrtstr, '(a,i4,a)') 'Warning : file ', ifile, & ' has both BACKSCAL keyword and column' CALL xwrite(wrtstr, 10) CALL xwrite('Using data from BACKSCAL column', 10) ENDIF filestat = 0 ENDIF c read the AREASCAL value(s). Note that the type II case is more c complicated now because we need to have the option of a vector c AREASCAL on each row. CALL ftgcno(unit, .FALSE., 'AREASCAL', index, filestat) IF ( filestat .NE. 0 ) THEN filestat = 0 CALL ftgkye(unit, 'AREASCAL', areascal, comment, filestat) IF (filestat .EQ. 0 .AND. areascal .NE. 0. ) THEN DO ibin = 1, ndata phabuf(ibin,5) = areascal ENDDO ELSE WRITE(wrtstr, '(a,i4,a)') 'Warning : file ', ifile, & ' has no AREASCAL keyword or column' CALL xwrite(wrtstr, 10) CALL xwrite('Assuming AREASCAL = 1', 10) DO ibin = 1, ndata phabuf(ibin,5) = 1.0 ENDDO ENDIF filestat = 0 ELSEIF ( filestat .EQ. 0 .AND. index .GT. 0 ) THEN IF ( ftype .EQ. 1 ) THEN CALL ftgcve(unit, index, 1, 1, ndata, 1.e-32, & phabuf(1,5), qanyf, filestat) ELSEIF ( ftype .EQ. 2 ) THEN CALL ftgdes(unit, index, ispec, nsize, ijunk, filestat) IF ( filestat .NE. 0 ) THEN filestat = 0 CALL ftgtdm(unit, index, 1, ijunk, nsize, filestat) ENDIF IF ( nsize .GT. 1 ) THEN CALL ftgcve(unit, index, ispec, 1, ndata, 1.e-32, & phabuf(1,5), qanyf, filestat) ELSE CALL ftgcve(unit, index, ispec, 1, 1, 1.e-32, & areascal, qanyf, filestat) DO ibin = 1, ndata phabuf(ibin,5) = areascal ENDDO ENDIF ENDIF context = 'RDDTG1: Failed to read AREASCAL data' IF ( filestat .NE. 0 ) GOTO 100 CALL ftgkye(unit, 'AREASCAL', junk, comment, filestat) IF ( filestat .EQ. 0 ) THEN WRITE(wrtstr, '(a,i4,a)') 'Warning : file ', ifile, & ' has both AREASCAL keyword and column' CALL xwrite(wrtstr, 10) CALL xwrite('Using data from AREASCAL column', 10) ENDIF filestat = 0 ENDIF RETURN c Check for any FITSIO errors 100 CONTINUE IF (filestat .NE. 0) THEN wrtstr = filenm(:lenact(filenm))//' : '//context CALL xwrite(wrtstr, 10) WRITE(wrtstr,'('' FITSIO = '', i4)') filestat CALL xwrite(wrtstr, 10) ENDIF RETURN END