function index_nat,wave_um

;-----------------------------------------------------------------
; Purpose:
;   Returns the refractive index for crystalline "beta" NAT at 196 K, 
;   from Toon et al., JGR, 99, pp25,631, 1994.  The data are stored in 
;   the routine for ease and simplicity.
;
; Input:
;   wave_um.......input wavelength, microns
;
; Output:
;   result........complex refractive index
;
; Source: Mark Hervig 
;-----------------------------------------------------------------

     numw = 109          ; # wavelengths in Toon et al. 
     wave = fltarr(numw) ; 1/cm
     dxr  = fltarr(numw) ; real part of index
     dxi  = fltarr(numw) ; imaginary part of index

     waves = $
    '482.000 494.000 505.000 517.000 525.000 528.000 555.000 571.000 582.000 602.000 '+$
    '617.000 629.000 636.000 648.000 671.000 687.000 702.000 717.000 721.000 729.000 '+$
    '733.000 741.000 748.000 779.000 806.000 814.000 822.000 829.000 837.000 845.000 '+$
    '872.000 891.000 906.000 941.000 976.000 1034.00 1049.00 1057.00 1065.00 1111.00 '+$
    '1169.00 1184.00 1196.00 1219.00 1277.00 1300.00 1323.00 1354.00 1362.00 1373.00 '+$
    '1400.00 1415.00 1423.00 1435.00 1454.00 1469.00 1504.00 1539.00 1589.00 1608.00 '+$
    '1631.00 1662.00 1682.00 1701.00 1712.00 1743.00 1786.00 1828.00 1878.00 1909.00 '+$
    '1940.00 1979.00 2087.00 2148.00 2349.00 2480.00 2519.00 2600.00 2700.00 2750.00 '+$
    '2796.00 2920.00 2989.00 3047.00 3089.00 3128.00 3170.00 3201.00 3228.00 3267.00 '+$
    '3313.00 3363.00 3379.00 3448.00 3494.00 3533.00 3568.00 3595.00 3618.00 3626.00 '+$
    '3629.00 3656.00 3737.00 3853.00 3996.00 4108.00 4420.00 4891.00 5805.00'
     reads,waves,wave

     dxrs = $
    '1.95000 1.83200 1.79200 1.82300 1.85400 1.86900 1.93600 1.94500 1.93000 1.85000 '+$
    '1.79600 1.82500 1.85400 1.89800 1.95000 1.90600 1.77200 1.65200 1.64600 1.57200 '+$
    '1.50800 1.44300 1.46500 1.52900 1.61800 1.64800 1.60100 1.53700 1.52900 1.53800 '+$
    '1.58200 1.60100 1.60300 1.64600 1.68700 1.74700 1.75800 1.76100 1.77200 1.83000 '+$
    '1.89100 1.89700 1.88800 1.86800 1.87400 1.85100 1.65800 1.47600 1.45800 1.39000 '+$
    '1.18200 1.09500 1.12200 1.19600 1.30000 1.36000 1.43300 1.46200 1.48600 1.49500 '+$
    '1.48900 1.48500 1.48200 1.49900 1.51100 1.52900 1.51300 1.47500 1.42800 1.40800 '+$
    '1.40400 1.41300 1.42400 1.42400 1.42600 1.43900 1.44800 1.45200 1.45000 1.44000 '+$
    '1.42900 1.41600 1.42400 1.43800 1.45100 1.46400 1.46900 1.45500 1.42900 1.41100 '+$
    '1.40100 1.34600 1.31100 1.22400 1.23400 1.26000 1.28500 1.30300 1.31600 1.31900 '+$
    '1.32100 1.33300 1.35700 1.37900 1.39600 1.40500 1.42400 1.44000 1.45400'
     reads,dxrs,dxr

     dxis = $
    '0.222000 0.254000 0.199000 0.157000 0.150000 0.153000 0.212000 0.273000 0.319000 0.371000 '+$
    '0.306000 0.267000 0.263000 0.282000 0.403000 0.525000 0.616000 0.552000 0.559000 0.600000 '+$
    '0.581000 0.456000 0.374000 0.241000 0.218000 0.265000 0.327000 0.293000 0.242000 0.214000 '+$
    '0.166000 0.165000 0.150000 0.120000 0.119000 0.130000 0.144000 0.141000 0.141000 0.178000 '+$
    '0.273000 0.318000 0.351000 0.387000 0.529000 0.655000 0.798000 0.692000 0.700000 0.710000 '+$
    '0.625000 0.424000 0.328000 0.242000 0.190000 0.183000 0.196000 0.212000 0.215000 0.225000 '+$
    '0.235000 0.232000 0.221000 0.210000 0.213000 0.241000 0.279000 0.296000 0.280000 0.257000 '+$
    '0.231000 0.210000 0.188000 0.184000 0.154000 0.136000 0.137000 0.147000 0.159000 0.166000 '+$
    '0.162000 0.141000 0.129000 0.126000 0.129000 0.144000 0.174000 0.202000 0.212000 0.208000 '+$
    '0.221000 0.262000 0.259000 0.161000 8.49000E-02 4.74000E-02 2.80000E-02 2.06000E-02 1.57000E-02 '+$
    '1.52000E-02 1.35000E-02 1.25000E-02 1.01000E-02 9.04000E-03 9.44000E-03 0. 0. 0. 0.'
     reads,dxis,dxi

;---- interpolate index to wave_cm 

      wave_cm = 10000.0/wave_um   ; turn microns into 1/cm

      if (wave_cm ge wave(0) and wave_cm le wave(numw-1) ) then begin

        for i = 0,(numw-2) do begin

          if (wave(i) le wave_cm and wave(i+1) gt wave_cm) then begin

            term   = (wave_cm-wave(i)) /(wave(i+1)-wave(i))
            dr     = term * (dxr(i+1)-dxr(i)) + dxr(i)
            di     = term * (dxi(i+1)-dxi(i)) + dxi(i)
            result = complex(dr,di)
            goto,jump1

          endif
        endfor

      endif

      print  ,wave_um,' not in range in index_NAT'
      result = complex(0.0,0.0)
      jump1:
      return,result
      end