;+
;
; computes doppler profile given 
;
; ny0 - line frequency at rest
; ny - frequency grid
; te - temperature of formation of line
; u(x,y) - velocity
; chi: non-thermal line broadening in km/s
; awgt: atomic weight of atom
;
; output: phiny(ny,x,y)
; 
; Oivind Wikstol     12.02.04  ver 1.0
;
;-
pro phiny,ny0,ny,u,te, awgt, chi, phi
; define constants:
cc=2.997929e10
bk=1.38046e-16

chi = chi*1.e5  ; km/s to cm/s
u = u*1.e5      ; km/s to cm/s

dopp_width=double(ny0/cc*sqrt(2.*bk*te/awgt + chi*chi))
dny = double(ny-ny0)
sz = size(u)
nny = n_elements(ny)
nx = sz[1]
if sz[0] gt 1 then xny = sz[2] else xny = 1
phi = fltarr(nny, nx, xny)

; non 
for i = 0, nny-1 do begin
  lshift = (dny(i)-ny0*u/cc)/dopp_width
  lshift = lshift*lshift
  phi(i, *, *) = 1./sqrt(!pi)/dopp_width*exp(-lshift)
endfor

return
end