!+GCORRMF ! ----------------- subroutine gcorrf ! ----------------- ! --- DESCRIPTION ------------------------------------------------------ ! This task shifts the channels of an OGIP standard RMF file, ! by a specified amount. ! ---------------------------------------------------------------------- ! --- VARIABLES -------------------------------------------------------- ! IMPLICIT NONE character(180) infile, outfile, ebdfile, shftfile integer errflg, chatter logical killit ! ! --- VARIABLE DIRECTORY ----------------------------------------------- ! ! infile char : input RMF filename ! ebdfile char : input EBD filename - usually same as RMF filename ! outfile char : output filename ! shftfile char : input shift filename ! ! --- AUTHORS/MODIFICATION HISTORY ------------------------------------- ! ! Rehana Yusaf (1.0.0: August 1995) ! ! Banashree M Seifert (1.1.0: July 1996) ! . bug fixed -- DMA allocation for N_chan was not done ! but was declared as integer p_N_chan ! ! Banashree M Seifert (1.2.0: Aug, 1996) ! . fractional shift introduced ! . replaced by screen display routines ! Banashree M Seifert (1.3.0: Oct 1996) ! . replaced by rdebd3,rdrmf3 in place of rdebd1,rdrmf1 ! . filename character(180) ! ! Banashree M Seifert (1.4.0: Nov21 1996) ! . replaced by call to grprmf in place of call to grp_rmf ! change in grprmf contains passing parameter first channel ! no. (0 or 1) ! ! Banashree M Seifert (1.4.1: Mar 25, 1997) ! . clobber parameter was not active ! ! Peter D Wilson (1.4.2: Jul 01, 1998) ! . Updated for new FCPARS behavior ! Ning Gan (1.4.3: Nov 01, 1999) ! . Fixed bugs in gcorr_rmf. Now the output file will ! have the same group structure as the input file. ! MFC (1.4.4 Apr 16 2020) - f90 conversion ! ------------------------------------------------------------------------ character(5) version parameter (version = '1.4.4') character(8) taskname taskname = 'gcorrmf' !- ! ---------------------------------------------------------------------- ! ! --- GET PARAMETERS --- ! errflg = 0 call gcorr_gp(infile, ebdfile, outfile, shftfile, & killit, chatter, errflg) IF (errflg.NE.0) THEN call wterrm(taskname, version, 'returning from gcorr_gp') goto 100 ENDIF call wtbegm(taskname, version, chatter) ! --- COMPRESS RMF DATA --- call gcorr_shft(infile, ebdfile, outfile, shftfile, & version, killit, errflg, chatter) IF (errflg.NE.0) THEN call wterrm(taskname, version, 'returning from gcorr_shft') ENDIF 100 call wtendm(taskname, version, errflg, chatter) return end ! ---------------------------------------------------------------- ! END OF GCORRMF ! ---------------------------------------------------------------- !+GCORR_GP ! --------------------------------------------------------- subroutine gcorr_gp(infile, ebdfile, outfile, shftfile, & killit, chatter, errflg) ! --------------------------------------------------------- ! --- DESCRIPTION ------------------------------------------------------ ! Gets parameters. ! --- VARIABLES -------------------------------------------------------- ! IMPLICIT NONE character*(*) infile, outfile, shftfile, ebdfile character(80) ebdfilename, rmffilename character(70) desc, defval integer errflg, chatter integer status, rmfextno, ebdextno integer flen, fcstln logical killit, ext ! ! --- VARIABLE DIRECTORY ----------------------------------------------- ! ! Arguments ... ! ! infile char : input file name ! ebdfile char : ebdfile input file, usually same as infile ! outfile char : Output filename ! nchan int : final number of channels, entered directly/from pha ! chatter int : Chattiness flag, >20 verbose ! errflg int : Error flag ! ! --- CALLED ROUTINES ------------------------------------------------- ! ! UCLGST : (HOST) Get string input ! UCLGSI : (HOST) Get integer input ! FCECHO : (FTOOLS) Screen write ! ! --- COMPILATION/LINKING --------------------------------------------- ! ! CALTOOLS, FTOOLS ! ! --- AUTHORS/MODIFICATION HISTORY ------------------------------------ ! ! Rehana Yusaf (1995 August 9) 1.0.0; ! ! Banashree M Seifert (1996 Nov20) 1.1.0: ! . substituted screen display routine, e.g., wterrm etc ! Banashree M Seifert (1997 Mar25) 1.2.0: ! . clobber parameter was not checking if clobber=yes is given ! on command line. ! Peter D Wilson (1998 Jul 01) 1.2.1: ! . Drop INQUIRE tests for input files. Let FTOPEN check ! ------------------------------------------------------------- character(9) subname parameter (subname = 'gcorr_gp') character(5) version parameter (version = '1.2.1') ! ------------------------------------------------------------- ! GET CLOBBER status = 0 call uclgsb('clobber', killit, status) IF (status.NE.0) THEN desc = 'getting clobber parameter !' call wterrm(subname, version, desc) errflg = 1 return ENDIF ! GET CHATTER status = 0 call uclgsi('chatter', chatter, status) IF (status.NE.0) THEN desc = 'getting chatter parameter !' call wterrm(subname, version, desc) errflg = 1 return ENDIF ! GET INFILE status = 0 call uclgst('infile', infile, status) IF (status.NE.0) THEN desc = 'getting infile parameter !' call wterrm(subname, version, desc) errflg = 1 return ENDIF call crmvlbk(infile) IF (infile(1:2) .EQ.' ') THEN desc = 'input file has to be entered !' call wterrm(subname, version, desc) errflg = 1 return ENDIF ! PDW 7/1/98: Don't bother! Let FTOPEN determine if file exists ! Strip off extension with ftrtnm, for later status = 0 call ftrtnm(infile, rmffilename, status) ! call fcpars(infile,rmffilename,rmfextno,status) ! ext = .true. ! flen = fcstln(rmffilename) ! INQUIRE(FILE=rmffilename(:flen),EXIST=ext) ! IF (.NOT.ext) THEN ! errflg = 1 ! desc = 'File does not exist :'//rmffilename ! call wterrm(subname,version,desc) ! return ! ENDIF ! GET EBDFILE status = 0 call uclgst('ebdfile', ebdfile, status) IF (status.NE.0) THEN desc = 'getting ebdfile parameter !' call wterrm(subname, version, desc) errflg = 1 return ENDIF call crmvlbk(ebdfile) IF (ebdfile.EQ.'%') THEN ebdfile = rmffilename ! PDW 7/1/98: Don't bother! Let FTOPEN determine if file exists ! ELSE ! call fcpars(ebdfile,ebdfilename,ebdextno,status) ! IF (ebdfilename.NE.rmffilename) THEN ! ext = .true. ! flen = fcstln(ebdfilename) ! INQUIRE(FILE=ebdfilename(:flen),EXIST=ext) ! IF (.NOT.ext) THEN ! errflg = 1 ! desc = 'File does not exist :'//ebdfilename ! call wterrm(subname,version,desc) ! return ! ENDIF ! ENDIF ENDIF ! PUT DEFAULT EBDFILE defval = '%' status = 0 call uclpst('ebdfile', defval, status) IF (status.NE.0) THEN desc = 'putting default ebdfile parameter !' call wterrm(subname, version, desc) errflg = 1 return ENDIF ! GET SHFTFILE status = 0 call uclgst('shftfile', shftfile, status) IF (status.NE.0) THEN desc = 'getting shftfile parameter !' call wterrm(subname, version, desc) errflg = 1 return ENDIF ext = .true. call crmvlbk(shftfile) INQUIRE(FILE = shftfile, EXIST = ext) IF (.NOT.ext) THEN errflg = 1 desc = 'File does not exist :' // shftfile call wterrm(subname, version, desc) return ENDIF ! GET OUTFILE status = 0 call uclgst('outfile', outfile, status) IF (status.NE.0) THEN desc = 'getting outfile parameter !' call wterrm(subname, version, desc) errflg = 1 return ENDIF call crmvlbk(outfile) IF (outfile.EQ.' ') THEN desc = 'outfile must be entered !' call wterrm(subname, version, desc) errflg = 1 return ENDIF ext = .false. INQUIRE(FILE = outfile, EXIST = ext) IF (ext) THEN if(.not. killit) then errflg = 1 desc = 'File exists :' // outfile call wterrm(subname, version, desc) desc = 'try with !' // outfile call wtinfo(chatter, 0, 2, desc) return endif ENDIF return end ! --------------------------------------------------------------------- ! END OF GCORR_GP ! --------------------------------------------------------------------- !+GCORR_SHFT ! ----------------------------------------------------------------- subroutine gcorr_shft(infile, ebdfile, outfile, shftfile, & mnver, killit, errflg, chatter) ! ----------------------------------------------------------------- ! --- DESCRIPTION -------------------------------------------------- ! ! This routine shifts the channels, in the EBOUNDS, and SPECRESP ! extensions of a user defined RMF file, and writes the rebinned ! data to an output file. ! ! --- VARIABLES ---------------------------------------------------- ! IMPLICIT NONE character*(*) infile, outfile, mnver, shftfile, ebdfile integer errflg, chatter logical killit ! ! --- VARIABLE DIRECTORY ------------------------------------------- ! ! Arguments ... ! ! infile char : input filename ! outfile char : output filename ! phatlscop char : Telescope name, read from Phafile ! phainstrum char : Instrument name, read from phafile ! phafile char : phafile name ! phachan int : Array of pha channels ! mnver char : Main task version number ! errflg int : error flag, 0 is ok ! chatter int : Chattines flag ! ! --- COMPILATION/LINKING ------------------------------------------ ! ! FTOOLS,CALLIB,FITSIO ! ! --- AUTHORS/MODIFICATION HISTORY --------------------------------- ! ! Rehana Yusaf (1995 August 9) 1.0.0; ! ! Banashree Mitra Seifert (Aug, 1996) 1.1.0: ! . replaced screen display routines ! ! Banashree Mitra Seifert (Oct, 1996) 1.2.0: ! . replaced rdebd1 by rdebd3 ! . replaced rcrmf1 by rdrmf3 ! . replaced wtrmf1 by wtrmf3 ! . replaced wtebd1 by wtebd3 ! . filename character(180) ! . initialised status=0 (LINUX problem) ! . formatting i8 instead of i for print out ! ! Peter D Wilsom (Jul 01, 1998) 1.2.1: ! . Update filename handling for new FCPARS/FTOPEN behavior ! ! Alex M. ( Feb 16, 1999) 1.2.2: ! . Replaced RDRMF3/WTRMF3 by RDRMF4/WTRMF4. Added 'call RMFSZ' ! ------------------------------------------------------------------- character(11) subname parameter (subname = 'gcorr_shft') character(5) version parameter (version = '1.2.2') ! ------------------------------------------------------------------ ! **** DYNAMIC MEMORY ALLOCATION **** ! the following MEM common block definition is in the system iraf77.inc ! file LOGICAL MEMB(100) ! INTEGER*2 MEMS(100) ! INTEGER*4 MEMI(100) ! INTEGER*4 MEML(100) INTEGER (kind=2) :: MEMS(100) INTEGER (kind=4) :: MEMI(100) INTEGER (kind=4) :: MEML(100) REAL MEMR(100) DOUBLE PRECISION MEMD(100) COMPLEX MEMX(100) EQUIVALENCE (MEMB, MEMS, MEMI, MEML, MEMR, MEMD, MEMX) COMMON /MEM/ MEMD ! ************************************ ! MAX ARRAY SIZES ( Alex -- added maxne, maxelt, and numelt ) ! increased size of keys array (MJT 17Feb2010) integer maxne, maxchan, maxgrp, maxelt, numelt ! INTERNALS ! ... parameters & static arrays integer maxshft parameter (maxshft = 256) integer nshfts real sfacts(maxshft) integer stbins(maxshft), endbins(maxshft) character(16) ebdtlscop, ebdinstrum, ebdfilt, ebddet character(16) rmftlscop, rmfinstrum, rmffilt, rmfdet real ebdarea integer iunit, ounit, n_ebd, ebderr, cerr, lchan integer nkeys, nmore, nk_comm, nk_hist integer imove, next(50), ninstr, nfound, nsearch character(20) instr(50), outhdu(9, 50), extname(50), outver(9, 50) integer imaxgrp, rmfchan, ienerg, detchans, fchan INTEGER numelt_cor, imaxgrp_cor INTEGER imaxchan, imaxelt, imaxne integer keyspres, j, i, rmfextno, ebdextno, iunit2 real lo_thresh, rmfarea integer block, rw, status, htype character(4) cstart, cend character(6) compchar character(70) bincomm, exhist, errinfo, desc, message character(80) rec, hist(5), comm(5) logical man_key, qorder, isorder character(5) rsp_rmfversn, ebd_rmfversn character(20) hduclas3, taskname character(180) ebdfilename, rmffilename character(8) task, dmode, keys(100), matext character(80) rmf_extname, rmf_extcomm character(80) ebd_extname, ebd_extcomm character(8) ebdchantype, rmfchantype real febdchan, frmfchan ! ... pointers to "arrays" to be dynamically allocated ! integer*8 p_ngrp, p_fmatrix, p_F_chan, p_N_chan, p_order ! integer*8 p_energ_lo, p_energ_hi, p_e_min, p_e_max ! integer*8 p_ebdchan, p_en_slice, p_buffer, p_fmatrix_cor ! integer*8 p_ngrp_cor, p_F_ch, p_N_ch integer (kind=8) :: p_ngrp, p_fmatrix, p_F_chan, p_N_chan, p_order integer (kind=8) :: p_energ_lo, p_energ_hi, p_e_min, p_e_max integer (kind=8) :: p_ebdchan, p_en_slice, p_buffer, p_fmatrix_cor integer (kind=8) :: p_ngrp_cor, p_F_ch, p_N_ch ! ... "arrays" to be dynamically allocated ! integer ngrp(maxne) real fmatrix(maxelt) ! integer F_chan(maxgrp) integer N_chan(maxgrp) ! real energ_lo(maxne) real energ_hi(maxne) ! real e_min(maxchan) real e_max(maxchan) ! integer ebdchan(maxchan) real en_slice(maxchan) ! real buffer(maxchan) integer N_ch(maxgrp) ! real p_fmatrix_cor(2*maxelt) integer F_ch(maxgrp) ! integer p_ngrp_cor(maxne) ! ! ! --- USER INFO --- ! desc = 'using ' // subname // ' Ver ' // version call wtinfo(chatter, 10, 2, desc) ! --- Alex -- initialize parameters isorder = .false. qorder = .false. maxne = 0 maxelt = 0 maxchan = 0 maxgrp = 0 ! --- FIND EBOUNDS EXTENSION --- status = 0 call fcpars(ebdfile, ebdfilename, ebdextno, status) ! PDW 7/1/98: Call ftrtnm to really strip off extension # call ftrtnm(ebdfile, ebdfilename, status) call crmvlbk(ebdfilename) status = 0 block = 2880 rw = 0 call cgetlun(iunit) ! PDW 7/1/98: Call ftopen with original filename ! call ftopen(iunit,ebdfilename,rw,block,status) call ftopen(iunit, ebdfile, rw, block, status) IF (status.NE.0) THEN errinfo = 'opening EBOUNDS file !' call wterrm(subname, version, errinfo) errflg = 1 goto 482 ENDIF ! ! --- FIND EXTNUM if EXTNUM IS NOT SPECIFIED AT THE END OF EBDFILE --- ! IF (ebdextno.LE.0) THEN ninstr = 2 instr(1) = 'RESPONSE' instr(2) = 'EBOUNDS' nsearch = 50 call fndhdu(chatter, iunit, ninstr, instr, nsearch, nfound, next, & outhdu, outver, extname, status) ! --- CHECK FOR (old-style) EXTNAME if HDUCLASS not found --- IF (nfound.LE.0) THEN errinfo = 'Ext with allowed HDUCLASn keywords not found' call wtwarm(subname, version, chatter, 20, errinfo) errinfo = 'searching for EXTNAME= EBOUNDS' call wtwarm(subname, version, chatter, 20, errinfo) call fndext(chatter, iunit, 'EBOUNDS', & nsearch, nfound, next, outhdu, outver, extname, status) ENDIF IF (nfound.LE.0) THEN errinfo = 'EBOUNDS extension found ' call wterrm(subname, version, errinfo) errflg = 1 goto 482 ENDIF IF (nfound.GT.1) THEN errinfo = 'Input file contains >1 EBOUNDS datasets' call wterrm(subname, version, errinfo) write(errinfo, '(a,i12,a)')'... ', nfound, ' extensions found' call wterrm(subname, version, errinfo) do i = 1, nfound write(errinfo, '(a,i12,a)')'Ext ', next(i), ':' call wterrm(subname, version, errinfo) write(errinfo, '(4X,a,a)')'EXTNAME = ', extname(i) call wterrm(subname, version, errinfo) do j = 1, 4 write(errinfo, '(4X,a,i2,2a)')& 'HDUCLAS', j, ' = ', outhdu(j, i) call wterrm(subname, version, errinfo) enddo enddo errinfo = & 'Extension number must be specified via ebdfile parameter' call wterrm(subname, version, errinfo) errinfo = 'for example INPUT.EBD[1]' call wterrm(subname, version, errinfo) errinfo = 'ebdfile is a hidden parameter' call wterrm(subname, version, errinfo) errflg = 1 goto 482 ENDIF ENDIF ! --- MOVE TO APPROPRIATE PLACE IN FILE --- IF (ebdextno.LE.0) THEN IF (next(1).GT.0) THEN imove = next(1) status = 0 call ftmrhd(iunit, imove, htype, status) errinfo = ' Problem moving to EBOUNDS xtens' IF (status.NE.0) THEN call wtferr(subname, version, status, errinfo) errflg = 1 goto 482 ENDIF ENDIF ebdextno = next(1) ELSE status = 0 call ftmahd(iunit, ebdextno + 1, htype, status) errinfo = 'Problem moving to EBOUNDS extension' IF (status.NE.0) THEN errflg = 1 call wterrm(subname, version, errinfo) goto 482 ENDIF ENDIF ! Alex --- Get the maxchan and maxne array sizes status = 0 call ftgkyj(iunit, 'DETCHANS', maxchan, errinfo, status) errinfo = 'reading RMF DETCHANS value' call wtferr(subname, version, status, errinfo) maxchan = maxchan + 2 write(errinfo, '(a,i8)') ' Maxchan array :', maxchan call wtinfo(chatter, 30, 2, errinfo) IF (status.NE.0) THEN errflg = 1 goto 482 ENDIF call ftgkyj(iunit, 'NAXIS2', maxne, errinfo, status) errinfo = 'reading RMF NAXIS2 value' call wtferr(subname, version, status, errinfo) write(errinfo, '(a,i8)') ' Maxne array :', maxne call wtinfo(chatter, 30, 2, errinfo) IF (status.NE.0) THEN errflg = 1 goto 482 ENDIF ! ALLOCATE MEMORY p_e_min = 0 p_e_max = 0 p_ebdchan = 0 status = 0 call udmget(maxne, 6, p_e_min, status) IF (status.NE.0) THEN goto 50 ENDIF call udmget(maxne, 6, p_e_max, status) IF (status.NE.0) THEN goto 50 ENDIF call udmget(maxne, 4, p_ebdchan, status) 50 if(status.NE.0) then message = & 'Failed to allocate Dynamic Memory for EBOUNDS data' call wterrm(subname, version, message) errflg = -1 goto 482 endif ! ! --- READ EBOUNDS EXTENSION --- ! ebderr = 0 detchans = 0 febdchan = 0. ebd_extname = ' ' ebd_extcomm = ' ' call rdebd4(iunit, chatter, ebd_extname, ebd_extcomm, & maxchan, ebdtlscop, ebdinstrum, & ebddet, ebdfilt, ebdarea, ebdchantype, febdchan, & n_ebd, MEMI(p_ebdchan), & MEMR(p_e_min), MEMR(p_e_max), & ebd_rmfversn, ebderr) IF (ebderr.NE.0) THEN errflg = 1 errinfo = 'reading EBOUNDS extension' call wterrm(subname, version, errinfo) goto 482 ENDIF ! ! --- GET SHIFT FACTOR AND CHANNELS --- ! call gt_fchan(MEMI(p_ebdchan), fchan, & lchan, n_ebd, errflg, chatter) task = 'GCORRMF' dmode = ' ' IF (detchans.EQ.0) THEN detchans = n_ebd ENDIF call shft_fact(shftfile, fchan, lchan, nshfts, & stbins, endbins, sfacts, maxshft, errflg, chatter) IF (errflg.NE.0) THEN errinfo = '..determining shift factor' call wterrm(subname, version, errinfo) goto 482 ENDIF ! ! --- SHIFT EBOUNDS DATA --- ! ! ebdcerr = 0 ! call gcorr_ebd(n_ebd,MEMI(p_ebdchan),MEMR(p_e_min),MEMR(p_e_max), ! & nshfts,stbins,endbins,sfacts,ebdcerr,chatter) ! ! --- READ RMF DATA --- ! ! --- FIND RMF EXTENSION --- call fcpars(infile, rmffilename, rmfextno, status) ! PDW 7/1/98: Call ftrtnm to really strip off extension # call ftrtnm(infile, rmffilename, status) call crmvlbk(rmffilename) IF (rmffilename.EQ.ebdfilename) THEN iunit2 = iunit call ftmahd(iunit2, 1, htype, status) ELSE status = 0 call ftclos(iunit, status) status = 0 block = 2880 rw = 0 call cgetlun(iunit2) ! PDW 7/1/98: Call ftopen with original filename call ftopen(iunit2, infile, rw, block, status) IF (status.NE.0) THEN errinfo = 'opening RMF file !' call wterrm(subname, version, errinfo) errflg = 1 goto 482 ENDIF ENDIF ! ! --- FIND EXTNUM if EXTNUM IS NOT SPECIFIED AT THE END OF INFILE --- ! IF (rmfextno.LE.0) THEN ninstr = 2 instr(1) = 'RESPONSE' instr(2) = 'RSP_MATRIX' nsearch = 50 call fndhdu(chatter, iunit2, ninstr, instr, nsearch, nfound, next, & outhdu, outver, extname, status) ! --- CHECK FOR (old-style) EXTNAME if HDUCLASS not found --- IF (nfound.LE.0) THEN errinfo = & 'Ext with allowed HDUCLASn keywords not found' call wtinfo(chatter, 20, 1, errinfo) errinfo = 'searching for EXTNAME= MATRIX' call wtinfo(chatter, 20, 1, errinfo) call fndext(chatter, iunit2, 'MATRIX', & nsearch, nfound, next, outhdu, outver, extname, status) ENDIF IF (nfound.LE.0) THEN errinfo = ' ... searching for EXTNAME=SPECRESP' call wtinfo(chatter, 20, 1, errinfo) call fndext(chatter, iunit2, 'SPECRESP MATRIX', & nsearch, nfound, next, outhdu, outver, extname, status) ENDIF IF (nfound.LE.0) THEN errinfo = 'RSP_MATRIX extension not found ' call wtinfo(chatter, 20, 1, errinfo) errflg = 1 goto 482 ENDIF IF (nfound.GT.1) THEN errinfo = 'Input file contains >1 MATRIX datasets' call wterrm(subname, version, errinfo) write(errinfo, '(a,i12,a)')' ', nfound, ' extensions found' call wterrm(subname, version, errinfo) do i = 1, nfound write(errinfo, '(a,i12,a)')'Ext ', next(i), ':' call wterrm(subname, version, errinfo) write(errinfo, '(4X,a,a)')'EXTNAME = ', extname(i) call wterrm(subname, version, errinfo) do j = 1, 4 write(errinfo, '(4X,a,i2,2a)')& 'HDUCLAS', j, ' = ', outhdu(j, i) call wterrm(subname, version, errinfo) enddo enddo errinfo = & '... Extension number must be specified via infile parameter' call wterrm(subname, version, errinfo) errinfo = 'for example INPUT.RMF[1]' call wterrm(subname, version, errinfo) errflg = 1 goto 482 ENDIF ENDIF ! --- MOVE TO APPROPRIATE PLACE IN FILE --- IF (rmfextno.LE.0) THEN IF (next(1).GT.0) THEN imove = next(1) status = 0 call ftmrhd(iunit2, imove, htype, status) errinfo = 'Problem moving to EBOUNDS xtens' IF (status.NE.0) THEN call wtferr(subname, version, status, errinfo) errflg = 1 goto 482 ENDIF ENDIF rmfextno = next(1) ELSE status = 0 call ftmahd(iunit2, rmfextno + 1, htype, status) errinfo = 'Problem moving to MATRIX extension' IF (status.NE.0) THEN errflg = 1 call wtferr(subname, version, status, errinfo) goto 482 ENDIF ENDIF ! Alex --- Get the array sizes maxne, maxgrp, maxelt status = 0 call rmfsz(iunit2, chatter, maxne, maxgrp, maxelt, status) errinfo = 'reading RMF array sizes maxne, maxgrp, and maxelt' call wtferr(subname, version, status, errinfo) IF (status.NE.0) THEN errflg = 1 goto 482 ENDIF ! -------------------------------------------------------------------- ! ALLOCATE DMA p_energ_lo = 0 p_energ_hi = 0 p_ngrp = 0 p_fmatrix = 0 p_F_chan = 0 p_N_chan = 0 p_ngrp_cor = 0 p_fmatrix_cor = 0 p_F_ch = 0 p_N_ch = 0 p_order = 0 p_en_slice = 0 p_buffer = 0 status = 0 call udmget(maxne, 6, p_energ_lo, status) IF (status.NE.0) THEN goto 70 ENDIF call udmget(maxne, 6, p_energ_hi, status) IF (status.NE.0) THEN goto 70 ENDIF status = 0 call udmget(maxne, 4, p_ngrp, status) IF (status.NE.0) THEN goto 70 ENDIF call udmget(maxelt, 6, p_fmatrix, status) IF (status.NE.0) THEN goto 70 ENDIF call udmget(maxgrp, 4, p_F_chan, status) IF (status.NE.0) THEN goto 70 ENDIF call udmget(maxgrp, 4, p_N_chan, status) IF (status.NE.0) THEN goto 70 ENDIF call udmget(maxne, 4, p_ngrp_cor, status) IF (status.NE.0) THEN goto 70 ENDIF call udmget(2 * maxelt, 6, p_fmatrix_cor, status) IF (status.NE.0) THEN goto 70 ENDIF call udmget(2 * maxgrp, 4, p_F_ch, status) IF (status.NE.0) THEN goto 70 ENDIF call udmget(2 * maxgrp, 4, p_N_ch, status) IF (status.NE.0) THEN goto 70 ENDIF call udmget(2 * maxgrp, 4, p_order, status) IF (status.NE.0) THEN goto 70 ENDIF IF (status.NE.0) THEN goto 70 ENDIF call udmget(maxchan, 6, p_en_slice, status) IF (status.NE.0) THEN goto 70 ENDIF call udmget(maxchan, 6, p_buffer, status) IF (status.NE.0) THEN goto 70 ENDIF 70 if(status.NE.0) then message = 'Failed to allocate dynamic Memory for RMF data' call wterrm(subname, version, message) errflg = -1 goto 482 endif ! --- Alex --- Replaced RDRMF3 by RDRMF4 ! --- Ziqin --- Replaced RDRMF4 by RDRMF5 rmfchan = 0 ienerg = 0 imaxgrp = 0 numelt = 0 rmf_extname = ' ' rmf_extcomm = ' ' call rdrmf5(iunit2, chatter, rmf_extname, rmf_extcomm, & qorder, maxne, maxgrp, maxelt, & rsp_rmfversn, hduclas3, rmftlscop, rmfinstrum, rmfdet, rmffilt, & rmfarea, rmfchantype, frmfchan, rmfchan, ienerg, imaxgrp, numelt, & MEMR(p_energ_lo), MEMR(p_energ_hi), MEMI(p_ngrp), MEMI(p_F_chan), & MEMI(p_N_chan), isorder, MEMI(p_order), MEMR(p_fmatrix), & lo_thresh, errflg) IF (errflg.NE.0) THEN errinfo = 'reading RMF extension ' call wterrm(subname, version, errinfo) goto 482 ENDIF status = 0 ! ! --- SHIFT RMF DATA --- ! ! --- Alex: add MEMI(p_ngrp), MEMI(p_F_chan), and MEMI(p_N_chan) arrays ! --- GCORR_RMF returns the matrix with the shifted channels, ! --- FMATRIX_COR, and new arrays NGRP_COR(), F_CH(), and N_CH() numelt_cor = 0 imaxgrp_cor = 0 call gcorr_rmf(maxchan, maxne, maxelt, maxgrp, rmfchan, ienerg, & MEMR(p_fmatrix), MEMI(p_ngrp), MEMI(p_F_chan), MEMI(p_N_chan), & numelt_cor, imaxgrp_cor, MEMR(p_fmatrix_cor), & MEMI(p_ngrp_cor), MEMI(p_F_ch), MEMI(p_N_ch), MEMR(p_en_slice), & MEMR(p_buffer), lo_thresh, nshfts, stbins, endbins, sfacts, & fchan, lchan, chatter, errflg) ! maxchan i ! maxne i ! maxelt i ! maxgrp i ! rmfchan i ! ienerg i ! MEMR(p_fmatrix) i ! MEMI(p_ngrp) i ! MEMI(p_F_chan) i ! MEMI(p_N_chan) i ! numelt_cor o <=== number of output matrix elements ! imaxgrp_cor o <=== total number of groups in the output matrix ! MEMR(p_fmatrix_cor) o <=== output RMF matrix ! MEMR(p_ngrp_cor) o <=== number of channel groups output array ! MEMI(p_F_ch) o <=== first channels output array ! MEMI(p_N_ch) o <=== number of channels output array ! MEMR(P_en_slice) i ! MEMR(p_buffer) i ! lo_thresh i ! nshfts i ! stbins i ! endbins i ! sfacts i ! fchan i ! lchan i ! chatter i ! errflg o IF (errflg.NE.0) THEN errinfo = 'rebinning RMF data' call wterrm(subname, version, errinfo) goto 482 ENDIF call ftclos(iunit2, status) ! ! --- WRITE DATA --- ! call copyphd(rmffilename, outfile, killit, errflg, chatter) ! Ziqin --- copy the extension before rmfextno & eboundextno IF (ebdfilename.EQ.rmffilename) THEN call gcorr_copyrem(rmffilename, outfile, ebdextno, rmfextno, & -1, errflg, chatter) ENDIF call cgetlun(iunit) block = 2280 ! PDW 7/1/98: Use original filename instead of rmffilename: ! call ftopen(iunit,rmffilename,0,block,status) call ftopen(iunit, infile, 0, block, status) status = 0 block = 2280 call cgetlun(ounit) call ftopen(ounit, outfile, 1, block, status) nk_hist = 0 nk_comm = 0 ebderr = 0 ! ! --- WRITE RMF EXTENSION ( Alex --- Replaced WTRMF3 by WTRMF4 ) ! ! --- WRITE RMF EXTENSION ( Ziqin --- Replaced WTRMF4 by WTRMF5 ) IF (ebdfilename.EQ.rmffilename) THEN IF(rmfextno.gt. ebdextno) THEN goto 504 ENDIF ENDIF 502 CONTINUE errflg = 0 call wtrmf5(ounit, chatter, rmf_extname, rmf_extcomm, & nk_hist, hist, nk_comm, & comm, rsp_rmfversn, hduclas3, rmftlscop, rmfinstrum, rmfdet, & rmffilt, rmfarea, rmfchantype, frmfchan, numelt_cor, rmfchan, & ienerg, imaxgrp_cor, MEMR(p_energ_lo), MEMR(p_energ_hi), & MEMI(p_ngrp_cor), MEMI(p_F_ch), MEMI(p_N_ch), qorder, & MEMI(p_order), MEMR(p_fmatrix_cor), lo_thresh, errflg) IF (errflg.NE.0) THEN errinfo = 'writing RMF extension !' call wterrm(subname, version, errinfo) goto 482 ENDIF ! ! --- WRITE REMAINING KEYWORDS IN RMF EXTENSION --- ! call ftmahd(iunit, rmfextno + 1, htype, status) call getkeys(ounit, keyspres, keys, chatter, status) call ftghsp(iunit, nkeys, nmore, status) do i = 1, nkeys status = 0 call ftgrec(iunit, i, rec, status) errinfo = 'reading RMF record' call wtferr(subname, version, status, errinfo) man_key = .false. do j = 1, keyspres IF (rec(1:8).EQ.keys(j)) THEN man_key = .true. ENDIF enddo IF (.NOT.man_key) THEN call ftprec(ounit, rec, status) ENDIF errinfo = 'writing non-mandatory record ' call wtferr(subname, version, status, errinfo) enddo ! ! --- WRITE EXTRA HISTORY --- ! status = 0 call ftpdat(ounit, status) taskname = task // mnver status = 0 call ftpkys(ounit, 'CREATOR', taskname, & ' s/w task which wrote this dataset', status) exhist = ' The original file was ' // infile status = 0 call ftphis(ounit, exhist, status) ! ! --- WRITE REBIN COMMENT --- ! do i = 1, nshfts write(cstart, 100, IOSTAT = cerr) stbins(i) write(cend, 100, IOSTAT = cerr) endbins(i) write(compchar, 110, IOSTAT = cerr) sfacts(i) bincomm = ' data shifted from' // cstart // ' TO '& // cend // ' with a shift factor : ' // compchar status = 0 call rmvexsp(bincomm, exhist) call ftpcom(ounit, exhist, status) enddo IF (ebdfilename.EQ.rmffilename) THEN IF(rmfextno.gt. ebdextno) THEN goto 505 ENDIF ENDIF 503 CONTINUE IF (ebdfilename.EQ.rmffilename) THEN call gcorr_copyrem(rmffilename, outfile, ebdextno, rmfextno, & 0, errflg, chatter) ENDIF IF (ebdfilename.EQ.rmffilename) THEN IF(rmfextno.gt. ebdextno) THEN goto 502 ENDIF ENDIF 504 CONTINUE ! ! --- WRITE EBOUNDS EXTENSION --- ! call wtebd4(ounit, chatter, ebd_extname, ebd_extcomm, & nk_hist, hist, nk_comm, & comm, ebd_rmfversn, ebdtlscop, ebdinstrum, ebddet, & ebdfilt, ebdarea, ebdchantype, febdchan, n_ebd, & MEMR(p_e_min), MEMR(p_e_max), ebderr) IF (ebderr.NE.0) THEN errflg = 2 errinfo = 'writing EBOUNDS extension' call wterrm(subname, version, errinfo) goto 482 ENDIf ! ! --- WRITE REMAINING KEYWORDS IN EBOUNDS EXTENSION --- ! IF (ebdfilename.EQ.rmffilename) THEN call ftmahd(iunit, ebdextno + 1, htype, status) iunit2 = iunit ELSE status = 0 call ftclos(iunit, status) status = 0 call cgetlun(iunit2) ! PDW 7/1/98: Use original filename instead of ebdfilename: ! call ftopen(iunit2,ebdfilename,0,block,status) call ftopen(iunit2, ebdfile, 0, block, status) status = 0 call ftmahd(iunit2, ebdextno + 1, htype, status) ENDIF status = 0 call getkeys(ounit, keyspres, keys, chatter, status) call ftghsp(iunit2, nkeys, nmore, status) do i = 1, nkeys status = 0 call ftgrec(iunit2, i, rec, status) errinfo = 'readind record from EBOUNDS ext' call wtferr(subname, version, status, errinfo) man_key = .false. do j = 1, keyspres IF (rec(1:8).EQ.keys(j)) THEN man_key = .true. ENDIF enddo IF (.NOT.man_key) THEN call ftprec(ounit, rec, status) ENDIF errinfo = 'writing non-mandatory record ' call wtferr(subname, version, status, errinfo) enddo ! ! --- WRITE EXTRA HISTORY --- ! call ftpdat(ounit, status) taskname = task // mnver status = 0 call ftpkys(ounit, 'CREATOR', taskname, & ' s/w task which wrote this dataset', status) exhist = ' The original file was ' // ebdfile call ftphis(ounit, exhist, status) ! ! --- WRITE shift COMMENT --- ! ! do i=1,nshfts ! write(cstart,100,IOSTAT = cerr) stbins(i) ! write(cend,100,IOSTAT = cerr) endbins(i) ! write(compchar,100,IOSTAT = cerr) sfacts(i) ! bincomm = ' data shifted from'//cstart//' TO ' ! &//cend//' with a shift factor: '//compchar ! call rmvexsp(bincomm,exhist) ! status = 0 ! call ftpcom(ounit,exhist,status) ! enddo status = 0 call ftclos(iunit2, status) status = 0 call ftclos(ounit, status) IF (ebdfilename.EQ.rmffilename) THEN IF(rmfextno.gt. ebdextno) THEN goto 503 ENDIF ENDIF 505 CONTINUE ! ! --- MAKE COMPLETE COPIES OF ANY OTHER EXTENSIONS IN INFILE --- ! IF (ebdfilename.EQ.rmffilename) THEN call gcorr_copyrem(rmffilename, outfile, ebdextno, rmfextno, & 1, errflg, chatter) ENDIF 100 FORMAT(I4) 110 format (f6.2) ! Check for errors 482 if(errflg.ne.0) then errflg = 1 message = ' Fatal' call wterrm(subname, version, message) endif ! Delete the Allocated Memory call udmfre(p_e_min, 6, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_e_max, 6, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_ebdchan, 4, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_energ_lo, 6, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_energ_hi, 6, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_ngrp, 4, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_fmatrix, 6, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_F_chan, 4, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_N_chan, 4, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_ngrp_cor, 4, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_fmatrix_cor, 6, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_F_ch, 4, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_N_ch, 4, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_order, 4, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_en_slice, 6, status) IF (status.NE.0) THEN goto 485 ENDIF call udmfre(p_buffer, 6, status) IF (status.NE.0) THEN goto 485 ENDIF 485 if(status.NE.0) then message = 'Failed to deallocate Dynamic Memory' call wterrm(subname, version, message) errflg = 99 endif ! ---------------------------------------------------------------------- return end ! ---------------------------------------------------------------------- ! END OF GCORR_SHFT ! ---------------------------------------------------------------------- !+GCORR_RMF ! ------------------------------------------------------- subroutine gcorr_rmf(maxchan, maxne, maxelt, maxgrp, & ichan, ienerg, fmatrix, NGRP, F_CHAN, N_CHAN, & numelt_corr, imaxgrp_corr, fmatrix_corr, & NGRP_CORR, F_CH, N_CH, en_slice, buffer, thresh, nshfts, & stbins, endbins, sfacts, fchan, lchan, chatter, err) ! ------------------------------------------------------- ! --- DESCRIPTION -------------------------------------------------- ! ! This routine shifts RMF channels, and the energy arrays are ! changed accordingly. ! ! --- VARIABLES ---------------------------------------------------- ! IMPLICIT NONE integer ichan, chatter, err, nshfts, maxelt, maxgrp integer fchan, lchan, ienerg, numelt, imaxgrp_corr integer stbins(*), endbins(*), maxchan, maxne, numelt_corr real sfacts(*), en_slice(maxchan) REAL BUFFER(maxchan), THRESH REAL fmatrix(maxelt), fmatrix_corr(2 * maxelt) INTEGER NGRP(maxne), NGRP_CORR(ienerg), F_CHAN(maxgrp) INTEGER N_CHAN(maxgrp), F_CH(2 * maxgrp), N_CH(2 * maxgrp) ! ! --- LOCALS --- ! integer newch, ch, stbin, endbin, ie, i, ch_mv real s_shft, frac1 integer shft1 character(100) desc logical curshift INTEGER Igr, Irmf, Imtrx, INDEX, J, Jmax, K, Kmax LOGICAL COUNT ! ! --- MAIN VARIABLES ----- ----------------------------------------- ! ! MAXCHAN i maximum number of channels for each energy row ! MAXNE i maximum number of energy rows ! MAXELT i maximum number of matrix elements ! MAXGRP i maximum number of channel groups ! ICHAN i actual number of channels for each energy row in ! rmf matrix ! IENERG i actual number of energy rows in rmf matrix ! FMATRIX i input rmf matrix array ! NGRP i input array of number of groups ! F_CHAN i input array of first channels ! N_CHAN i input array of number of channels ! EN_SLICE i input array of size = maximum number of channels ! BUFFER i input array of size = maximum number of channels ! THRESH i input value of threshold ! ------------------------------------------------------------------ ! NUMELT_CORR o size of the output corrected matrix array ! IMAXGRP_CORR o total number of groups in the output matrix ! FMATRIX_CORR o output corrected matrix array ! NGRP_CORR o number of groups as a function of energy in ! the output matrix ! F_CH o output array of first channels ! N_CH o ouput array of number of channels !------------------------------------------------------------------- ! ! --- AUTHORS/MODIFICATION HISTORY --------------------------------- ! ! Rehana Yusaf (1995 Aug 11) ! ! Banashree Mitra Seifert (Aug, 1996) 1.1.0: ! . replaced with screen display routines ! . fractional shift introduced instead of integer ! ! AM (March 1999) 1.1.1: ! . modified to take into account the fact that the arrays ! are now 1-dimensional. Added a block that generates ! the output corrected matrix array, calculates the number ! ! Ning Gan (1.1.2: Nov 01, 1999) ! . Fixed the bugs in gcorr_rmf. Now the output file will ! have the same group structure as the input file. ! ------------------------------------------------------------------ character(10) subname parameter (subname = 'gcorr_rmf') character(5) version parameter (version = '1.1.2') integer lch, fch, igr2 logical found_fch, found_lch ! ------------------------------------------------------------------ desc = 'using ' // subname // ' Ver ' // version call wtinfo(chatter, 10, 1, desc) IF (fchan.EQ.0) THEN ch_mv = 1 ELSE ch_mv = 0 ENDIF INDEX = 0 Irmf = 0 Imtrx = 0 Igr = 0 igr2 = 0 numelt_corr = 0 imaxgrp_corr = 0 ! Alex ---- Initialize the output F_CH() and N_CH() arrays DO J = 1, 2 * MAXGRP F_CH(J) = 0 N_CH(J) = 0 ENDDO ! Alex ---- Initialize the output FMATRIX_CORR() array DO J = 1, MAXELT FMATRIX_CORR(J) = 0.E0 ENDDO ! -------------------------------------------------------- DO IE = 1, IENERG COUNT = .FALSE. ! Alex ---- Initialize the BUFFER and EN_SLICE arrays DO CH = 1, ICHAN EN_SLICE(CH) = 0.E0 BUFFER(CH) = 0.E0 ENDDO ! Alex ---- "Uncompress" matrix array and copy its elements for each energy ! row to the EN_SLICE() array Jmax = NGRP(IE) DO J = 1, Jmax Igr = Igr + 1 Kmax = N_CHAN(Igr) DO K = 1, Kmax Imtrx = Imtrx + 1 CH = F_CHAN(Igr) + K - 1 + ch_mv EN_SLICE(CH) = EN_SLICE(CH) + FMATRIX(Imtrx) ENDDO ENDDO ! Alex --- BEGIN LOOPING THROUGH THE CHANNEL SHIFTS do i = 1, nshfts ! Alex --- SHIFT THE CHANNELS stbin = stbins(i) endbin = endbins(i) s_shft = abs(sfacts(i)) shft1 = int(s_shft + 1) frac1 = s_shft - int(s_shft) curshift = .true. ch = stbin ! ----------- case when shift is +ve --------------------------------- IF(sfacts(i).ge.0) then do WHILE(curshift) newch = ch + shft1 IF ((newch.GE.fchan).AND.(newch.LE.lchan)) THEN BUFFER(newch + ch_mv) = BUFFER(newch + ch_mv) + & frac1 * en_slice(ch + ch_mv) ENDIF IF ((newch - 1 .GE.fchan).AND.(newch - 1 .LE.lchan)) THEN BUFFER(newch + ch_mv - 1) = BUFFER(newch + ch_mv - 1) + & (1 - frac1) * en_slice(ch + ch_mv) ENDIF ch = ch + 1 IF (ch.GT.endbin) THEN curshift = .false. ENDIF enddo ! -------------- case when shift < 0 ------------------------------------ else shft1 = -shft1 DO WHILE(curshift) newch = ch + shft1 IF ((newch.GE.fchan).AND.(newch.LE.lchan)) THEN BUFFER(newch + ch_mv) = BUFFER(newch + ch_mv) + & frac1 * en_slice(ch + ch_mv) ENDIF IF ((newch + 1 .GE.fchan).AND.(newch + 1 .LE.lchan)) THEN BUFFER(newch + ch_mv + 1) = BUFFER(newch + ch_mv + 1) + & (1 - frac1) * en_slice(ch + ch_mv) ENDIF ch = ch + 1 IF (ch.GT.endbin) THEN curshift = .false. ENDIF ENDDO ENDIF ! --------------------------------------------------------------------------- ENDDO ! ------------------- THE CHANNELS ARE SHIFTED ----------------------------- ! Alex ***** CREATE CORRECTED MATRIX ARRAY, FMATRIX_CORR(), ***** ! **** FIND NGRP_CORR(), N_CH(), and F_CH() ARRAYS **** ! Alex ---- Find the channel groups, number of channels, and the first channel ! in each group of the output RMF matrix NGRP_CORR(IE) = 0 COUNT = .FALSE. ! Alex ---- COMPRESS THE BUFFER ARRAY AND COPY IT TO THE CORRESPONDING ! ELEMENTS OF THE MATRIX FMATRIX_CORR() ! Ning Gan ---- Find the corrected first channel and last channel in each ! group. DO J = 1, Jmax igr2 = igr2 + 1 Kmax = N_CHAN(igr2) fch = f_chan(igr2) lch = f_chan(igr2) + Kmax - 1 found_fch = .false. found_lch = .false. do i = 1, nshfts shft1 = nint(sfacts(i)) if(.not.found_fch .and.& fch .ge. stbins(i) & .and. & fch .le. endbins(i)) then fch = fch + shft1 ! check the previous channel IF (fch - 1.ge.fchan) then if(BUFFER(fch - 1 + ch_mv) .GT.THRESH) fch = fch - 1 endif if(fch .lt. fchan) fch = fchan found_fch = .true. endif if(.not.found_lch .and. & lch .ge. stbins(i) & .and. & lch .le. endbins(i)) then lch = lch + shft1 ! check the next channel IF (lch + 1.le.lchan) then if(BUFFER(lch + 1 + ch_mv) .GT. THRESH) lch = lch + 1 endif if(lch .gt. lchan) lch = Lchan found_lch = .true. endif if(found_lch.and.found_fch) goto 213 enddo 213 if(fch .gt. lchan) fch = lchan if(lch .lt. fch) lch = fch ! Alex ---- Save the channel group and first channel for each group ! in the corresponding arrays INDEX = INDEX + 1 IF (INDEX.GT.(2 * MAXGRP)) THEN DESC = 'GCORR_RMF: INSUFFICIENT MEMORY FOR F_CH, N_CH' CALL FCECHO(DESC) ERR = 1 RETURN ENDIF NGRP_CORR(IE) = NGRP_CORR(IE) + 1 ! THE BUFFER IS COMPRESSED AND COPIED ! The first channels array: F_CH(INDEX) = fch ! -------------------------------------------- ! The number of channels array: N_CH(INDEX) = lch - fch + 1 ! -------------------------------------------------------- do i = fch, lch Irmf = Irmf + 1 ! Alex --- Copy the elements (with the values exceeding the threshold) ! of the BUFFER array to the "compressed" FMATRIX array FMATRIX_CORR(Irmf) = BUFFER(i + ch_mv) ENDDO enddo ! The number of channel groups array: IMAXGRP_CORR = IMAXGRP_CORR + NGRP_CORR(IE) ! ----------------------------------------------- ENDDO ! Alex ---- THE CHANNEL SHIFT IS DONE FOR ALL ENERGIES ! Save the number of matrix elements NUMELT_CORR = Irmf ! ----------------------------------- desc = subname // ' : matrix is shifted' call wtinfo(chatter, 20, 1, desc) return end ! ----------------------------------------------------------------------- ! << END OF GCORR_RMF >> ! ----------------------------------------------------------------------- !+GCORR_EBD ! ------------------------------------------------------- subroutine gcorr_ebd(n_ebd, ebdchan, e_min, e_max, nshfts, & stbins, endbins, sfacts, err, chatter) ! ------------------------------------------------------- ! --- DESCRIPTION -------------------------------------------------- ! ! This routine rebins EBOUNDS channels, and the energy arrays are ! changed accordingly. ! ! --- VARIABLES ---------------------------------------------------- ! IMPLICIT NONE integer n_ebd, chatter, err, nshfts integer stbins(*), endbins(*) real sfacts(*) integer ebdchan(*) real e_min(*), e_max(*) ! ! --- LOCALS --- ! integer new, i, stbin, stnew, endbin, nbin, endnew logical curbin character(70) desc ! ! --- VARIABLE DIRECTORY ------------------------------------------- ! ! ! --- AUTHORS/MODIFICATION HISTORY --------------------------------- ! ! Rehana Yusaf (1993 July 16) ! character(5) version parameter (version = '1.0.0') !- ! ------------------------------------------------------------------ ! ! --- USER INFO --- ! IF (chatter.GE.10) THEN desc = ' ... using GCORR_EBD Ver ' // version call fcecho(desc) ENDIF ! ! --- COMPRESS --- ! new = 1 do i = 1, nshfts stbin = stbins(i) endbin = endbins(i) stnew = stbin nbin = sfacts(i) curbin = .true. do WHILE(curbin) endnew = stnew + nbin - 1 ebdchan(new) = new e_min(new) = e_min(stnew) e_max(new) = e_max(endnew) new = new + 1 stnew = endnew + 1 IF (stnew.GE.endbin) THEN curbin = .false. ENDIF enddo enddo n_ebd = new - 1 return end ! ----------------------------------------------------------------------- ! END OF GCORR_EBD ! ----------------------------------------------------------------------- !+GCORR_COPYREM ! ------------------------------------------------------ subroutine gcorr_copyrem(infile, outfile, ebdextno, specextno, & flg, errflg, chatter) ! ------------------------------------------------------ ! --- DESCRIPTION ---------------------------------------------------- ! This routine copies any extensions, other than the EBOUNDS, and ! SPECRESP from infile to outfile. ! -------------------------------------------------------------------- ! --- VARIABLES ------------------------------------------------------ ! IMPLICIT NONE character*(*) infile, outfile integer chatter, errflg, ebdextno, specextno, flg ! ! --- INTERNALS ------------------------------------------------------ ! integer ignore1, ignore2, iunit, block, i logical endfile, pkey, hist character(70) ap_in character(4) cnum integer htype, nhdu, status, tothd, iend, end_num, ierr character(32) errstr, wrnstr character(70) desc, errinfo ! ! --- VARIABLE DIRECTORY --------------------------------------------- ! ! Arguments ... ! ! infile char : input filename ! outfile char : output filename ! errflg int : error flag ! chatter int : Chattiness flag ! ! --- AUTHORS/MODIFICATION ------------------------------------------- ! ! Rehana Yusaf (1993 August 1) ! Rehana Yusaf (1993 Nov 22) 1.0.1; spexextno and ebdextno now passed character(5) version parameter (version = '1.0.1') !- ! -------------------------------------------------------------------- ! ! --- USER INFO --- ! errstr = ' ERROR : GCORR_COPYREM Ver ' // version // ':' wrnstr = ' WARNING : GCORR_COPYREM Ver ' // version // ':' IF (chatter.GE.10) THEN desc = ' ... using GCORR_COPYREM Ver ' // version call fcecho(desc) ENDIF ! ! --- DETERMINE SPECRESP AND EBOUNDS EXT Nos --- ! status = 0 call cgetlun(iunit) call ftopen(iunit, infile, 0, block, status) errinfo = errstr // ' opening infile !' call wt_ferrmsg(status, errinfo) IF (status.NE.0) THEN errflg = 1 return ENDIF endfile = .false. tothd = 0 nhdu = 1 do WHILE (.NOT.endfile) status = 0 call ftmahd(iunit, nhdu, htype, status) IF ((status.EQ.107).OR.(status.EQ.207)) THEN endfile = .true. tothd = nhdu - 2 ENDIF nhdu = nhdu + 1 enddo ! ! --- COPY EXTENSIONS, OTHER THAN EBOUNDS OR SPECRESP --- ! IF (specextno.LT.ebdextno) THEN ignore1 = specextno ignore2 = ebdextno ELSE ignore1 = ebdextno ignore2 = specextno ENDIF ! ... COPY EXTENSIONS BEFORE IGNORE1 IF (flg.eq.-1) THEN iend = index(infile(1:), ' ') - 1 pkey = .false. hist = .false. do i = 1, (ignore1 - 1) write(cnum, 100, IOSTAT = ierr) i call crmvlbk(cnum) end_num = index(cnum(1:), ' ') - 1 IF (end_num.EQ.0) THEN end_num = 4 ENDIF ap_in = infile(1:iend) // '[' // cnum(1:end_num) // ']' status = 0 call cfappend(ap_in, outfile, pkey, hist) enddo ENDIF ! ... COPY EXTENSIONS INBETWEEN IGNORE1 AND IGNORE2 IF(flg.eq. 0) THEN IF ((ignore1 + 1).NE.ignore2) THEN do i = (ignore1 + 1), ignore2 - 1 write(cnum, 100, IOSTAT = ierr) i call crmvlbk(cnum) end_num = index(cnum(1:), ' ') - 1 IF (end_num.EQ.0) THEN end_num = 4 ENDIF ap_in = infile(1:iend) // '[' // cnum(1:end_num) // ']' status = 0 call cfappend(ap_in, outfile, pkey, hist) enddo ENDIF ENDIF ! ... COPY EXTENSIONS AFTER IGNORE2 IF(flg.eq.1) THEN IF (ignore2.NE.tothd) THEN do i = 1, (ignore2 + 1), tothd write(cnum, 100, IOSTAT = ierr) i call crmvlbk(cnum) end_num = index(cnum(1:), ' ') - 1 IF (end_num.EQ.0) THEN end_num = 4 ENDIF ap_in = infile(1:iend) // '[' // cnum(1:end_num) // ']' status = 0 call cfappend(ap_in, outfile, pkey, hist) enddo ENDIF ENDIF 100 FORMAT(I4) return end ! --------------------------------------------------------------- ! END OF GCORR_COPYREM ! --------------------------------------------------------------- !+SHFT_FACT ! ------------------------------------------------------ subroutine shft_fact(shftfile, & fchan, lchan, nshfts, stbins, endbins, sfacts, & maxshft, errflg, chatter) ! ------------------------------------------------------ ! --- DESCRIPTION ------------------------------------------------- ! This subroutine reads binning info from an ascii file ! --- VARIABLES --------------------------------------------------- ! IMPLICIT NONE character*(*) shftfile integer fchan, maxshft, nshfts, lchan integer stbins(maxshft), endbins(maxshft) real sfacts(maxshft) integer errflg, chatter ! ! --- VARIABLE DIRECTORY ------------------------------------------ ! ! fchan int : starting channel ! shftfile char : shift filename ! nshfts int : Number of compression factors ! sfacts int : Array of shift factors ! stbins int : Array of shift starting points for each sfact ! endbins int : Array of shift ending points for each sfact ! errflg int : Error flag ! chatter int : chatter flag ! ! --- AUTHORS/MODIFICATION HISTORY --------------------------------- ! Rehana Yusaf (1995 August 10) 1.0.0; ! ! Banashree Mitra Seifert (Aug, 1996)1.1.0: ! . replaced by screen display routines ! . accomodated for shft, so that it can be fractional shift ! ------------------------------------------------------------------ character(10) subname parameter (subname = 'shft_fact') character(5) version parameter (version = '1.0.0') ! !- ! ------------------------------------------------------------------- ! --- INTERNALS --- character(70) subinfo, errinfo character(5) chanchar, chanchar2 character(6) comp real s_fact integer st_bin, end_bin, i, iunit, status integer prev_stbin, prev_endbin ! ! --- USER INFO --- ! subinfo = 'using ' // subname // ' Ver ' // version call wtinfo(chatter, 10, 2, subinfo) ! ! --- OPEN BINFILE --- ! call ftgiou(iunit, errflg) IF (errflg.NE.0) THEN errinfo = 'problem getting free lun' call wterrm(subname, version, errinfo) return ENDIF open(unit = iunit, file = shftfile, status = 'old') ! --- READ COMPRESSION DATA FROM BINFILE --- nshfts = 0 status = 0 prev_stbin = 0 prev_endbin = 0 do i = 1, maxshft + 1 read(iunit, *, IOSTAT = status, end = 100)st_bin, end_bin, s_fact IF (status.NE.0) THEN subinfo = 'invalid number in file' call wterrm(subname, version, subinfo) errflg = 1 return ENDIF IF ((nshfts + 1).GT.maxshft) THEN subinfo = ' The shift array sizes have been exceeded' call wterrm(subname, version, subinfo) write(subinfo, '(a,i12)')' The max array dimension is ', maxshft call wterrm(subname, version, subinfo) errflg = 3 goto 200 ENDIF IF (nshfts.GE.1) THEN IF ((endbins(nshfts) + 1).LT.st_bin) THEN nshfts = nshfts + 1 stbins(nshfts) = endbins(nshfts - 1) + 1 endbins(nshfts) = st_bin - 1 sfacts(nshfts) = 0. ENDIF ENDIF IF ((nshfts.EQ.0).AND.(st_bin.GT.fchan)) THEN nshfts = nshfts + 1 stbins(nshfts) = fchan endbins(nshfts) = st_bin - 1 sfacts(nshfts) = 0. ENDIF nshfts = nshfts + 1 stbins(nshfts) = st_bin endbins(nshfts) = end_bin sfacts(nshfts) = s_fact prev_stbin = stbins(nshfts) prev_endbin = endbins(nshfts) enddo 100 IF (endbins(nshfts).LT.lchan) THEN IF ((nshfts + 1).GT.maxshft) THEN subinfo = ' The shift array sizes have been exceeded' call wtinfo(chatter, 9, 1, subinfo) write(subinfo, '(a,i12)')' The max array dimension is ', maxshft call wtinfo(chatter, 9, 1, subinfo) errflg = 3 goto 200 ENDIF nshfts = nshfts + 1 stbins(nshfts) = endbins(nshfts - 1) + 1 endbins(nshfts) = lchan sfacts(nshfts) = 0. ENDIF ! USER INFO IF (chatter.GE.20) THEN do i = 1, nshfts write(chanchar, 300, IOSTAT = status) stbins(i) write(chanchar2, 300, IOSTAT = status) endbins(i) write(comp, 310, IOSTAT = status) sfacts(i) errinfo = chanchar // ' - ' // chanchar2& // ' are shifted with shift factor ' // comp call rmvexsp(errinfo, subinfo) call wtinfo(chatter, 20, 2, subinfo) enddo ENDIF 200 close(unit = iunit) status = 0 call ftfiou(iunit, status) 300 FORMAT(i4) 310 format (f6.2) 350 FORMAT(' starting channel:', I4, ' ending channel:', I4) return end ! ------------------------------------------------------------- ! END OF SHFT_FACT ! -------------------------------------------------------------