;-----------------------------------------------------------------------------
; figure_14b.pro:
; Reproduces the results shown in Figure 14b of Cranmer & Saar (2011).
;
; Calls:   IDL procedure boreas.pro
; Inputs:  ASCII file titled 'evol_track_1Msun.dat' (in current directory)
; Outputs: Plot to current device
;-----------------------------------------------------------------------------

; First read in theoretical 2D grid of ZAMS models.

    FILE1 = 'evol_track_1Msun.dat'

    get_lun,unit & openr,unit,FILE1 & readf,unit,npts
    aage    = dblarr(npts)
    Mstar   = dblarr(npts)
    Rstar   = dblarr(npts)
    Lstar   = dblarr(npts)
    Protday = dblarr(npts)
    for i=0,npts-1 do begin
      readf,unit,x1,x2,x3,x4,x5
      aage(i)    = x1
      Mstar(i)   = x2
      Rstar(i)   = x3
      Lstar(i)   = x4
      Protday(i) = x5
    endfor
    close,unit & free_lun,unit

; For each star, call BOREAS to compute the theoretical mass loss rate.

    Mdot_theory = dblarr(npts)
    FeHnow      = 0.0

    for i=0,npts-1 do begin
      BOREAS,Mstar(i),Rstar(i),Lstar(i),Protday(i),FeHnow, $
        Rossby,fstar,Mdot,Mdot_hot,Mdot_cold,MATR
      Mdot_theory(i) = Mdot
    endfor

; Make a simple plot similar to Figure 14b of Cranmer & Saar (2011).

    xr = [5.8, 10.4]
    yr = [-15., -9.]
    plot,xr,yr, /nodata, $
      xstyle=1,ystyle=1,xrange=xr,yrange=yr

    oplot,aage,alog10(Mdot_theory)

    oplot,[9.66,9.66],yr, linestyle=1

    end