;+
; FILE:   emdsun.pro
;* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
; PURPOSE
;	Electron mass density (gram/cm2) per 1.0 Bsun brightness
;	in the solar corona.

; SYNTAX
;	gcm2 = emdsun (r)
;
; INPUTS
;	r	distance from solar limb in Rsun units.
;
; OUTPUTS
;	gcm2	grams per Bsun per cm^2 at the specified distance.
;
; HISTORY
;	Andrew L. Stanger   HAO/NCAR   4 June 1998
;	 2 October 1998: increase etab from 100 to 311 elements,
;		         to increase the range from 11.0 to 30.0 Rsun.
;	23 April 1999: ro was incorrectly left set to 1.0 when additional
;		       elements were added on 2 Oct 1998; changed it to 0.0
;* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
;-
FUNCTION emdsun, r

NELEM  = 311		; # elements in etab array.
irad   = 0		;
dr     = 0.10		; delta r  (radial increment). */
ro     = 0.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).
;* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
;*/
etab = fltarr (NELEM)

;etab [0:49] = [3.602e+24, 5.954e+24, 7.910e+24, 9.960e+24, 1.215e+25, $
;                 1.448e+25, 1.697e+25, 1.962e+25, 2.242e+25, 2.539e+25, $
;                 2.851e+25, 3.180e+25, 3.524e+25, 3.885e+25, 4.261e+25, $
;                 4.654e+25, 5.062e+25, 5.487e+25, 5.928e+25, 6.384e+25, $
;                 6.857e+25, 7.346e+25, 7.851e+25, 8.371e+25, 8.908e+25, $
;                 9.461e+25, 1.003e+26, 1.062e+26, 1.122e+26, 1.183e+26, $
;                 1.247e+26, 1.312e+26, 1.378e+26, 1.446e+26, 1.516e+26, $
;                 1.587e+26, 1.660e+26, 1.735e+26, 1.811e+26, 1.889e+26, $
;                 1.968e+26, 2.049e+26, 2.131e+26, 2.216e+26, 2.301e+26, $
;                 2.389e+26, 2.478e+26, 2.568e+26, 2.660e+26, 2.754e+26]

;etab [50:100] = [2.849e+26, 2.946e+26, 3.045e+26, 3.145e+26, 3.247e+26, $
;                 3.350e+26, 3.455e+26, 3.562e+26, 3.670e+26, 3.780e+26, $
;                 3.891e+26, 4.004e+26, 4.119e+26, 4.235e+26, 4.353e+26, $
;                 4.472e+26, 4.593e+26, 4.716e+26, 4.840e+26, 4.966e+26, $
;                 5.093e+26, 5.222e+26, 5.353e+26, 5.485e+26, 5.619e+26, $
;                 5.755e+26, 5.892e+26, 6.030e+26, 6.170e+26, 6.312e+26, $
;                 6.456e+26, 6.601e+26, 6.747e+26, 6.896e+26, 7.046e+26, $	
;                 7.197e+26, 7.350e+26, 7.505e+26, 7.661e+26, 7.819e+26, $
;                 7.978e+26, 8.140e+26, 8.302e+26, 8.466e+26, 8.632e+26, $
;                 8.800e+26, 8.969e+26, 9.140e+26, 9.312e+26, 9.486e+26, $
;                 9.662e+26]

; - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
;--- etab: Ne/Bsun-cm2
; - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

; 0.0 - 0.9 Rsun
etab [  0:  9] = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]

; 1.0 - 1.9 Rsun
etab [ 10: 19] = [3.602e+24, 5.954e+24, 7.910e+24, 9.960e+24, 1.215e+25, $
                  1.448e+25, 1.697e+25, 1.962e+25, 2.242e+25, 2.539e+25]

; 2.0 - 2.9 Rsun
etab [ 20: 29] = [2.851e+25, 3.180e+25, 3.524e+25, 3.885e+25, 4.261e+25, $
		  4.654e+25, 5.062e+25, 5.487e+25, 5.928e+25, 6.384e+25]

; 3.0 - 3.9 Rsun
etab [ 30: 39] = [6.857e+25, 7.346e+25, 7.851e+25, 8.371e+25, 8.908e+25, $
                  9.461e+25, 1.003e+26, 1.062e+26, 1.122e+26, 1.183e+26]

; 4.0 - 4.9 Rsun
etab [ 40: 49] = [1.247e+26, 1.312e+26, 1.378e+26, 1.446e+26, 1.516e+26, $
                  1.587e+26, 1.660e+26, 1.735e+26, 1.811e+26, 1.889e+26]

; 5.0 - 5.9 Rsun
etab [ 50: 59] = [1.968e+26, 2.049e+26, 2.131e+26, 2.216e+26, 2.301e+26, $
                  2.389e+26, 2.478e+26, 2.568e+26, 2.660e+26, 2.754e+26]

; 6.0 - 6.9 Rsun
etab [ 60: 69] = [2.849e+26, 2.946e+26, 3.045e+26, 3.145e+26, 3.247e+26, $
                  3.350e+26, 3.455e+26, 3.562e+26, 3.670e+26, 3.780e+26]

; 7.0 - 7.9 Rsun
etab [ 70: 79] = [3.891e+26, 4.004e+26, 4.119e+26, 4.235e+26, 4.353e+26, $
                  4.472e+26, 4.593e+26, 4.716e+26, 4.840e+26, 4.966e+26]

; 8.0 - 8.9 Rsun
etab [ 80: 89] = [5.093e+26, 5.222e+26, 5.353e+26, 5.485e+26, 5.619e+26, $
                  5.755e+26, 5.892e+26, 6.030e+26, 6.170e+26, 6.312e+26]

; 9.0 - 9.9 Rsun
etab [ 90: 99] = [6.456e+26, 6.601e+26, 6.747e+26, 6.896e+26, 7.046e+26, $
                  7.197e+26, 7.350e+26, 7.505e+26, 7.661e+26, 7.819e+26]

; 10.0 - 10.9 Rsun
etab [100:109] = [7.978e+26, 8.140e+26, 8.302e+26, 8.466e+26, 8.632e+26, $
                  8.800e+26, 8.969e+26, 9.140e+26, 9.312e+26, 9.486e+26]

; 11.0 - 11.9 Rsun
etab [110:119] = [9.662e+26, 9.838e+26, 1.002e+27, 1.020e+27, 1.038e+27, $
                  1.056e+27, 1.075e+27, 1.093e+27, 1.112e+27, 1.131e+27]

; 12.0 - 12.9 Rsun
etab [120:129] = [1.150e+27, 1.170e+27, 1.189e+27, 1.209e+27, 1.229e+27, $
                  1.249e+27, 1.269e+27, 1.289e+27, 1.309e+27, 1.330e+27]

; 13.0 - 13.9 Rsun
etab [130:139] = [1.351e+27, 1.372e+27, 1.393e+27, 1.414e+27, 1.435e+27, $
                  1.457e+27, 1.479e+27, 1.501e+27, 1.523e+27, 1.545e+27]

; 14.0 - 14.9 Rsun
etab [140:149] = [1.567e+27, 1.590e+27, 1.612e+27, 1.635e+27, 1.658e+27, $
                  1.681e+27, 1.705e+27, 1.728e+27, 1.752e+27, 1.776e+27]

; 15.0 - 15.9 Rsun
etab [150:159] = [1.800e+27, 1.824e+27, 1.848e+27, 1.873e+27, 1.897e+27, $
                  1.922e+27, 1.947e+27, 1.972e+27, 1.997e+27, 2.023e+27]

; 16.0 - 16.9 Rsun
etab [160:169] = [2.048e+27, 2.074e+27, 2.100e+27, 2.126e+27, 2.152e+27, $
                  2.178e+27, 2.205e+27, 2.231e+27, 2.258e+27, 2.285e+27]

; 17.0 - 17.9 Rsun
etab [170:179] = [2.312e+27, 2.340e+27, 2.367e+27, 2.395e+27, 2.423e+27, $
                  2.450e+27, 2.478e+27, 2.507e+27, 2.535e+27, 2.564e+27]

; 18.0 - 18.9 Rsun
etab [180:189] = [2.593e+27, 2.622e+27, 2.651e+27, 2.680e+27, 2.709e+27, $
                  2.739e+27, 2.769e+27, 2.799e+27, 2.829e+27, 2.859e+27]

; 19.0 - 19.9 Rsun
etab [190:199] = [2.890e+27, 2.920e+27, 2.951e+27, 2.981e+27, 3.012e+27, $
                  3.044e+27, 3.075e+27, 3.107e+27, 3.139e+27, 3.170e+27]

; 20.0 - 20.9 Rsun
etab [200:209] = [3.203e+27, 3.234e+27, 3.267e+27, 3.298e+27, 3.331e+27, $
                  3.365e+27, 3.397e+27, 3.430e+27, 3.463e+27, 3.497e+27]

; 21.0 - 21.9 Rsun
etab [210:219] = [3.530e+27, 3.564e+27, 3.599e+27, 3.633e+27, 3.666e+27, $
                  3.701e+27, 3.735e+27, 3.769e+27, 3.804e+27, 3.841e+27]

; 22.0 - 22.9 Rsun
etab [220:229] = [3.874e+27, 3.911e+27, 3.948e+27, 3.982e+27, 4.016e+27, $
                  4.053e+27, 4.091e+27, 4.126e+27, 4.164e+27, 4.200e+27]

; 23.0 - 23.9 Rsun
etab [230:239] = [4.236e+27, 4.274e+27, 4.309e+27, 4.348e+27, 4.385e+27, $
                  4.420e+27, 4.460e+27, 4.498e+27, 4.537e+27, 4.574e+27]

; 24.0 - 24.9 Rsun
etab [240:249] = [4.613e+27, 4.651e+27, 4.691e+27, 4.728e+27, 4.769e+27, $
                  4.807e+27, 4.845e+27, 4.885e+27, 4.928e+27, 4.965e+27]

; 25.0 - 25.9 Rsun
etab [250:259] = [5.005e+27, 5.044e+27, 5.087e+27, 5.129e+27, 5.170e+27, $
                  5.205e+27, 5.250e+27, 5.288e+27, 5.332e+27, 5.374e+27]

; 26.0 - 26.9 Rsun
etab [260:269] = [5.416e+27, 5.458e+27, 5.499e+27, 5.540e+27, 5.581e+27, $
                  5.621e+27, 5.669e+27, 5.709e+27, 5.757e+27, 5.798e+27]

; 27.0 - 27.9 Rsun
etab [270:279] = [5.840e+27, 5.882e+27, 5.924e+27, 5.968e+27, 6.012e+27, $
                  6.057e+27, 6.104e+27, 6.144e+27, 6.193e+27, 6.235e+27]

; 28.0 - 28.9 Rsun
etab [280:289] = [6.279e+27, 6.324e+27, 6.372e+27, 6.412e+27, 6.464e+27, $
                  6.508e+27, 6.555e+27, 6.596e+27, 6.648e+27, 6.693e+27]

; 29.0 - 29.9 Rsun
etab [290:299] = [6.732e+27, 6.783e+27, 6.827e+27, 6.875e+27, 6.927e+27, $
                  6.971e+27, 7.019e+27, 7.071e+27, 7.116e+27, 7.165e+27]

; 30.0 - 30.9 Rsun
etab [300:309] = [7.207e+27, 7.265e+27, 7.304e+27, 7.348e+27, 7.397e+27, $
                  7.452e+27, 7.500e+27, 7.553e+27, 7.599e+27, 7.651e+27]

; 31.0 Rsun
etab [310] = 7.694e+27

; - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

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

;   print, rtab

; - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

   IF (r LT ro)  THEN RETURN, zero	;

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

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

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

   Me_val = float (Ne_val) * emass	;

   RETURN, Me_val			;
END