; file: plotsolrot.idl
; init: Dec 7 2013 Rob Rutten Sac Peak 
; note: coefficients from http://en.wikipedia.org/wiki/Solar_rotation

a=14.713 
b=-2.936
c=-1.787
phi=indgen(9)*10.
omega=a+b*sin(phi*!pi/180.)^2+c*sin(phi*!pi/180.)^4
omega=omega/24.
print,' phi, degrees/hour '
print,ntostr(phi,format='(F14.2)')
print,ntostr(omega,format='(F16.2)')
cgplot,phi,omega,/window,/ynozero,charsize=2,$
  xtitle='latitude [degrees]',ytitle='solar rotation  [degrees/hour]'

end
