;+
; FILE:   emdpb.pro
;* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
; PURPOSE
;	Electron mass density (gram/cm2) per 1.0 Bsun brightness (pB)
;	in the solar corona.
;
; SYNTAX
;	gcm2 = emdpb (r)
;
; INPUTS
;	r	radial distance from solar limb in Rsun units.
;
; OUTPUT
;	gcm2	grams per 1.0 Bsun per cm^2 at the specified distance.
;
; HISTORY
;	Andrew L. Stanger   HAO/NCAR   17 June 1998
;	23 April 1999: documentation & cosmetic changes only.
;* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
;-

FUNCTION emdpb, r

NELEM  = 101		; # elements in pbtab array.
irad   = 0		;
dr     = 0.10		; delta r  (radial increment). */
ro     = 1.00		; r offset (radial offset). */
Ne_val = 0.0		;
Me_val = 0.0		;
emass  = 2.0e-24	; mass (grams) per coronal electron. */
zero   = 0.0		;

;* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
;* # electrons per square cm per Bsun brightness in solar corona 
;* as a function of radius (1.00 -> 11.0 Rsun in increments of 0.1 Rsun.
;* Values computed from program jack (obtained from R. Munro).
;* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
;*/
epbtab = fltarr (NELEM)

epbtab [  0 :  49] = [3.884e+25, 1.957e+25, 1.895e+25, 1.984e+25, 2.137e+25, $
                     2.327e+25, 2.547e+25, 2.790e+25, 3.054e+25, 3.337e+25, $
                     3.639e+25, 3.960e+25, 4.297e+25, 4.652e+25, 5.023e+25, $
                     5.412e+25, 5.817e+25, 6.238e+25, 6.676e+25, 7.130e+25, $
                     7.601e+25, 8.088e+25, 8.591e+25, 9.110e+25, 9.646e+25, $
                     1.020e+26, 1.076e+26, 1.135e+26, 1.195e+26, 1.257e+26, $
                     1.320e+26, 1.385e+26, 1.451e+26, 1.519e+26, 1.589e+26, $
                     1.660e+26, 1.733e+26, 1.807e+26, 1.883e+26, 1.961e+26, $
                     2.040e+26, 2.121e+26, 2.204e+26, 2.288e+26, 2.374e+26, $
                     2.461e+26, 2.550e+26, 2.640e+26, 2.733e+26, 2.826e+26]

epbtab [ 50 : 100] = [2.922e+26, 3.019e+26, 3.117e+26, 3.217e+26, 3.319e+26, $
                     3.422e+26, 3.527e+26, 3.634e+26, 3.742e+26, 3.852e+26, $
                     3.963e+26, 4.076e+26, 4.191e+26, 4.307e+26, 4.425e+26, $
                     4.544e+26, 4.665e+26, 4.788e+26, 4.912e+26, 5.038e+26, $
                     5.165e+26, 5.294e+26, 5.425e+26, 5.557e+26, 5.691e+26, $
                     5.826e+26, 5.964e+26, 6.102e+26, 6.242e+26, 6.384e+26, $
                     6.528e+26, 6.673e+26, 6.819e+26, 6.968e+26, 7.117e+26, $
                     7.269e+26, 7.422e+26, 7.577e+26, 7.733e+26, 7.891e+26, $
                     8.050e+26, 8.211e+26, 8.374e+26, 8.538e+26, 8.704e+26, $
                     8.872e+26, 9.041e+26, 9.212e+26, 9.384e+26, 9.558e+26, $
                     9.733e+26] 

   rtab = fltarr (NELEM)
   rtab = indgen (NELEM) * dr + ro

;   print, rtab

   IF (r LT ro)  THEN RETURN, zero

   irad = FIX ((r - ro) / dr)

;   print, 'r: ', r, ' ro: ', ro, ' irad: ', irad

   Ne_val = (epbtab [irad]			$
	  + (epbtab [irad+1] - epbtab [irad])	$
	  * double ((r - irad*dr - ro) / dr)) 

   Me_val = float (Ne_val) * emass

   RETURN, Me_val
END