C *** This program is used in conjunction with "selectbc.data" to retrieve
C     BC data for the photometric system and filter bandpasses of interest,
C     for all possible values of log g, and for a particular value of E(B-V)
C     (the choices being 0.00, 0.08, 0.14, 0.28, 0.44, and 0.56.)  These 
C     data, which were derived directly from MARCS models, may be compared
C     with the interpolated BC values to check that the interpolation code
C     accurately reproduces the MARCS results.  The user must select the
C     E(B-V) value from the aforementioned choices.
C
C *** The output of this program consists of 13 data files (one for each
C     log g value between -0.5 and +5.5) with names sphm05.data, sph00.data,
C     ..., ppl55.data.  If these files already exist, they will be over-
C     written.  The program chkbctable.for compares the results in these
C     data files with those obtained by interpolating in tables for E(B-V)
C     values from 0.00 to 0.72.     
C
C *** NOTE that "slash" as defined in the data statement (below) **MUST**
C     be defined to be a forward slash (i.e., "slash/'/'/" if this program
C     is executed on a linux or UNIX computer.  (A backslash is used in the
C     case of a WINDOWS operating system.)
C--------1---------2---------3---------4---------5---------6---------7--
      real*4 bc(800,13,5)
      integer*4 itmm(13),ndat(13),ngv(4)
      character*52 fnme,fnme0,fzero
      character*16 bcdata(800,13)
      character*12 filter(41),xfi(5),finme
      character*9 alpha(4),alnme
      character*8 system(7),xps(5),psnme
      character*6 gravity(13),fgrav
      character*5 ftype
      character*4 dtest
      character*3 ebmv(6),rednme
      character*1 slash,fundr
      data alpha/'alpha_std','alpha_p04','alpha_p00','alpha_m04'/,
     1  gravity/'sphm05','sph00 ','sph05 ','sph10 ','sph15 ','sph20 ',
     2  'sph25 ','sph30 ','sph35 ','ppl40 ','ppl45 ','ppl50 ','ppl55 '/,
     3  system/'2mass   ','hst_ab  ','hst_st  ','hst_vega','sdss    ',
     4  'ubvri12 ','ubvri90 '/,ebmv/'r00','r08','r14','r28','r44',
     5  'r56'/,dtest/'test'/,ftype/'.data'/,fundr/'_'/,ngv/13,12,12,12/,
     6  fzero/'                                                    '/,
C *** on linux/UNIX machines, slash must be defined as '/' instead of '\' 
     7  slash/'/'/
      data filter/'2mass_J     ','2mass_H     ','2mass_K     ',
     1  'acs_f435w   ','acs_f475w   ','acs_f555w   ','acs_f606w   ',
     2  'acs_f814w   ','wfc3_f218w  ','wfc3_f225w  ','wfc3_f275w  ',
     3  'wfc3_f336w  ','wfc3_f350lp ','wfc3_f390m  ','wfc3_f390w  ',
     4  'wfc3_f438w  ','wfc3_f475w  ','wfc3_f547m  ','wfc3_f555w  ',
     5  'wfc3_f606w  ','wfc3_f625w  ','wfc3_f775w  ','wfc3_f814w  ',
     6  'wfc3_f850lp ','wfc3ir_f098m','wfc3ir_f110w','wfc3ir_f125w',
     7  'wfc3ir_f140w','wfc3ir_f160w','sdss_u      ','sdss_g      ',
     8  'sdss_r      ','sdss_i      ','sdss_z      ','jc_U        ',
     9  'jc_B        ','jc_V        ','jc_R        ','jc_I        ', 
     *  'jc_UX       ','jc_BX       '/
 1000 format(' enter 1 for E(B-V) = 0.00, 2 for 0.08, 3 for 0.14',/,
     1  7x,'4 for E(B-V) = 0.28, 5 for 0.44, 6 for 0.56')
 1001 format(9x,f5.2)
 1002 format(' E(B-V) =',f5.2)
 1003 format(a16,5f8.4)
 1004 format(i3,' = number of data points')
 1005 format(' files sphm05.data, sph00.data, sph05.data, etc., have'
     1  1x,'been created',/,' for E(B-V) =',f5.2,' to test',5(/,10x,
     2  a8,': ',a12))
c *** the input file "selectbc.data" is assigned to unit 9
      fnme(1:13)='selectbc.data'
      call iofile(9,'in',fnme,13)
      read(9,*) ialf
      ng=ngv(ialf)
      ig=0
      if(ialf.eq.1) go to 5
      ig=1
c *** choose the E(B-V) value to be used in the comparison of MARCS and
c *** interpolated results
    5 write(6,1000)
      read(5,*) ired 
      rednme=ebmv(ired)
c *** the 13 output files (for different log g) are assigned to units 10-22 
      nf1=9
      do 10 i=1,ng
      nf1=nf1+1
      jig=i+ig
      fnme=fzero
      fgrav=gravity(jig)
      fnme(1:6)=fgrav(1:6)
      call str_trim(6,fnme,nch)
      fnme(nch+1:nch+5)=ftype(1:5)
      ntot=nch+5
      call iofile(nf1,'oo',fnme,ntot)
   10 continue
c *** the directory relevant to the desired variation of [alpha/Fe] with [Fe/H]
      alnme=alpha(ialf)
c *** the number of filters for which BC tables are requested is read (unit 9)
      read(9,*) numbc
c *** the main DO-loop is over the number of individual system/filter names
      do 35 k=1,numbc
c *** integers selecting the photometric system and filter are read (unit 9)
      read(9,*) isys,ifil
      psnme=system(isys)
      finme=filter(ifil)
      xps(k)=psnme
      xfi(k)=finme
c *** set up the directory portion of the input file name
      fnme=fzero
      fnme(1:9)=alnme(1:9)
      fnme(10:10)=slash(1:1)
      fnme(11:14)=dtest(1:4)
      fnme(15:15)=slash(1:1)
      fnme(16:18)=rednme(1:3)
      fnme(19:19)=slash(1:1)
      fnme(20:27)=psnme(1:8)
      call str_trim(52,fnme,nch)
      nch=nch+1
      fnme(nch:nch)=slash(1:1)
      nch0=nch
      fnme0=fnme
      nf2=24
      do 30 j=1,ng
      nf2=nf2+1
      jig=j+ig
      fgrav=gravity(jig)
      nch=nch0
      fnme=fnme0 
      fnme(nch+1:nch+6)=fgrav(1:6)
      nch=nch+6
      if(jig.eq.1) go to 15
      nch=nch-1
   15 fnme(nch+1:nch+1)=fundr(1:1)
      nch=nch+1
      fnme(nch+1:nch+12)=finme(1:12)
      call str_trim(52,fnme,nch)
      fnme(nch+1:nch+5)=ftype(1:5)
      nch=nch+5
      call iofile(nf2,'in',fnme,nch)
c *** copy over the header.data information into the unit 10-22 data files
      read(nf2,1001) redden
      do 20 i=1,800
      read(nf2,1003,end=25) bcdata(i,j),bc(i,j,k)
      n=i
   20 continue
      stop 'dimensions of bcdata(i) need to be increased'
   25 ndat(jig)=n
   30 continue
   35 continue
      nf1=9
      do 45 j=1,ng
      jig=j+ig
      nf1=nf1+1
      n=ndat(jig)
      write(nf1,1002) redden
      write(nf1,1004) n
      do 40 i=1,n
      write(nf1,1003) bcdata(i,j),(bc(i,j,k),k=1,numbc)
   40 continue
   45 continue
c *** all of the input and output files are kept on disk
      nf1=9
      nf2=24
      close(unit=9,status='keep')
      do 50 i=1,ng
      nf1=nf1+1
      nf2=nf2+1
      close(unit=nf1,status='keep')
      close(unit=nf2,status='keep')
   50 continue
      write(6,1005) redden,(xps(i),xfi(i),i=1,numbc)
      stop
      end
      subroutine str_trim(nl,str,nch)
      character*(*) str
      nch=nl
      do while (str(nch:nch).eq.' '.and.nch.gt.0)
        nch=nch-1
      end do
      return 
      end
      subroutine iofile(nunit,fstat,fnme,nch)
      integer*4 nunit,nch
      character*52 fnme
      character*2 fstat
 1000 format(' unit',i3,' input file does not exist',/,'name = ',a52)
      if(fstat.eq.'in') go to 10
c *** the named file is apparently an output file
      open(unit=nunit,iostat=ierr,file=fnme(1:nch),status='new')
      if(ierr.eq.0) go to 15
c *** if the named file already exists it is deleted and opened as a
c *** new file; i.e., the existing file will be overwritten
      open(unit=nunit,file=fnme(1:nch),status='old')
      close(unit=nunit,status='delete')
      open(unit=nunit,file=fnme(1:nch),status='new')
      go to 15
c *** the named file is an input file 
   10 open(unit=nunit,iostat=ierr,file=fnme(1:nch),status='old')
      if(ierr.eq.0) go to 15
c *** the named input file apparently does not exist: execution terminated
      write(6,1000) nunit,fnme
      stop ' specified input file does not exist'
   15 return
      end
