C     Bolometric corrections (and in turn, color indices) using alpha_std 
C     models are computed for a sample of stars. The sample is read from the 
C     file "input.sample", where each row identifies a star of known 
C     logg, [Fe/H], Teff.
C 
C     Up to five filters can be chosen in selectbc.data (1 = ialf, for standard 
C     alpha variation).
C     After setting the desired parameters in selectbc.data, "getinputbcs" and 
C     "getbctable" must be compiled and run as usual (see README file). The 
C     desired value of E(B-V) must also be specified at this stage, and it 
C     will be used subsequently for all stars in "input.sample". 
C     bcstars.for can be compiled and linked as follows (e.g., ifort compiler):
C
C     ifort bcstars.for bcutil.for -o bcstars.exe
C
C     For each star the program returns BCs in the filters selected above 
C     (unit 6, i.e. on the monitor, from left to right, in the same order 
C     specified in selectbc.data, one star per row in the same order as 
C     "input.sample"). 
C     The user can then compute any color index (or combination of those) 
C     by simply doing the difference between the BCs in different bands.
C     NOTICE: since bolometric corrections -BCs- are returned, then a colour 
C     index i-j = BC_j - BC_i
C
C     If [alpha/Fe] of the input stars is also known, the interpolation using 
C     MARCS models with m04, p00 and p04 can also be done, instead of assuming 
C     standard alpha variation. In this case bcm04(5),bcp00(5),bcp04(5) must 
C     be uncommented below, and getbc_m04, getbc_p00.for, getbc_p04.for called. 
C     See colstar.for for an example of interpolation in [alpha/Fe].
C     In selectbc.data is also necessary to cycle through ialf, from 1 to 4. 

      real*4 a(3),x(3),bcstd(5)!,bcm04(5),bcp00(5),bcp04(5)
      character*6 fil(5)
      integer status

 1001 format(F4.2,X,F5.2,X,F6.0)
 1002 format(3x,'BC_1',4x,'BC_2',4x,'BC_3',4x,'BC_4',4x,'BC_5')
 1003 format(5(F7.3,X))

      write(6,1002) 

c *** read input file in loop till the end
      open(3,file='input.sample')
      readloop: do
      read(3,1001,IOSTAT=status) gv,fe,temp
      if (status /= 0) exit

      teff=alog10(temp)

c *** interpolate the bc_std tables to [Fe/H] 
      call getbc_std(0,fe,gv,teff,ebv,bcstd,fil,nbc)
c *** now perform the gv and Teff interpolation
      call getbc_std(1,fe,gv,teff,ebv,bcstd,fil,nbc)

      bc1 = bcstd(1)
      bc2 = bcstd(2)
      bc3 = bcstd(3)
      bc4 = bcstd(4)
      bc5 = bcstd(5)

c *** output is written on the monitor. User can send the output to a file 
c *** by simply modifying the code below
      write(6,1003) bc1,bc2,bc3,bc4,bc5

      end do readloop
      stop
      end