pro ramaty_gysy,delta,bb,theta,nv,freq,fi,rc,omega,tau $
    ,j1=j1,j2=j2,kf=kf,etr=etr,razin=razin $
    ,e1d=e1d,e2d=e2d,a1d=a1d,a2d=a2d,rccd=rccd $
    ,errflag=errflag

;+
; NAME:
;       ramaty_gysy
;
; PURPOSE:
;       Derive gyrosyncrotron radiation from power-law electrons.
;        (c.f. Ramaty 1969, ApJ, 158, 753)
;
; CALLING SEQUENCE:
;       ramaty_gysy,delta,bb,theta,nv,freq,fi,rc
;       ramaty_gysy,delta,bb,theta,nv,freq,fi,rc,omega,tau
;
; INPUTS:
;       delta   : power-law index of non-thermal e- distribution
;       bb      : magnetic field strength [G]
;       theta   : angle [degree] between Line-of-sight and magnetic vector
;       nv      : [particles] total number (E>10keV) of non-thermal e-
;       omega   : source size in surface angle [steradian]
;
; OUTPUTS:
;       freq    ; [GHz] frequency
;       fi      ; [SFU] flux density  (need 'nv' input)
;       rc      : polarization degree
;       tau     ; optical thickness
;
; OPTIONAL OUTPUTS:
;       errflag ; 
;           the number of omitted Bessel functions which should
;           be added to the results. The larger, the greater the error is.
;
; HISTORY:
;       Writtern 2000-05-03 T. Yokoyama
;       Revised  2000-10-18 T. Yokoyama nv,omega,tau
;       Revised  2001-08-28 T. Yokoyama 
;           Option 'errflag' is added since new versions of IDL
;           does NOT support the Bessel function which has
;           an order higher than 20. Array 'errflag' returns
;           the number of omitted Bessel functions which should
;           be added to the results.
;
;-

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; Universal constant

c = 3e10         ; [cm/s] speed of light
kb = 1.38*10.^(-16)  ; [cgs] Bolzmann constant
cm2_per_omega = (7.e7/(!dtor/3600.))^2
cgs2sfu = 1.e19 ; SFU/[erg/s/cm^2/Hz]

;
; assume 
;

if not keyword_set(omega) then omega= 1.e9^2/cm2_per_omega

;
; change scaling from 10keV (Dulk) to 1MeV (Ramaty)
;
nv_over_e_norm=(delta-1)*nv
anor=nv_over_e_norm*(100)^(1-delta)  ; used 1MeV/10keV = 100

cost=cos(theta * !dtor)

  src_size=sqrt(omega) * !radeg / 3600.

ramaty_gysy_core,delta,anor,bb,cost,src_size $
    ,ffbd,phi1,phi2,tau,rcd=rcd $
    ,j1=j1,j2=j2,kf=kf,etr=etr,alpha=razin $
    ,e1d=e1d,e2d=e2d,a1d=a1d,a2d=a2d,rccd=rccd $
    ,errflag=errflag

fgyro = bb * 2.8e6
freq=fgyro*ffbd*1.e-9
fi = (phi1+phi2)*cgs2sfu
rc=rcd


return
end