;+
; $Id: fov2radec.pro,v 1.3 2008/01/04 21:34:10 colaninn Exp $
;
; Project: STEREO-SECCHI
;
; Name: fov2radec
;
; Purpose: To convert HI pixel positions to RA-Dec pairs
;
; Category: STEREO, SECCHI, HI, Calibration
;
; Explanation: HI pixel positions are converted into a general AZP
;              form, which is converted to cartesian coordintes. These
;              are then rotated from HI pointing, to S/C pointing then
;              finally aligned with the RA-Dec frame. The Cartesian
;              coordinates are finally converted to RA-Dec.
;
; Syntax: radec = fov2radec(x, y, index)
;
; Example: radec = fov2radec(512,512,index)
;
; Inputs: x     - an n-element array of x pixel positions
;         y     - an n-element array of y pixel positions
;         index - the fits header for the HI image in question. This is 
;               used to deduce the pointings and instrument - etc.
;
; Opt. Inputs: None
;
; Outputs: radec - a 2xn array of RA-Dec pairs
;
; Opt. Outputs: None
;
; Keywords: system - which coordinate system to work in 'hpc' or 'gei'
;           note that in a change from the norm, gei is default
;
; Calls: azp2cart, hi2sc, sc2cart, get_hi_params
;
; Common: None
;
; Restrictions: None 
;
; Side effects: None
;
; Prev. Hist.: None
;
; History: 27-Jun-2007, Daniel Brown
;
; Contact: Daniel Brown (dob@aber.ac.uk) 
;
; $Log: fov2radec.pro,v $
; Revision 1.3  2008/01/04 21:34:10  colaninn
; changed fltarr to dblarr,  to 180D
;
; Revision 1.2  2007/12/14 13:21:57  bewsher
; Added _EXTRA keyword
;
; Revision 1.1  2007/11/28 12:05:53  bewsher
; First releases of hi_calib_point and associated code
;
;-

function fov2radec,xv,yv,index,yaw=yaw,pitch=pitch,roll=roll,hi_offsetx=offset_hi,hi_offsety=pitch_hi,hi_roll=roll_hi,mu=mu,fov=d,ccd_offsetx=ccdosx,ccd_offsety=ccdosy,pixmult=pmult,system=sys,_extra=extra

if (n_elements(sys) eq 0) then sys='gei'

if (sys eq 'gei') then begin
   if (n_elements(yaw) eq 0) then yaw = index.sc_yawa
   if (n_elements(pitch) eq 0) then pitch = index.sc_pita
   if (n_elements(roll) eq 0) then roll = index.sc_rolla
endif else begin
   if (n_elements(yaw) eq 0) then yaw = index.sc_yaw
   if (n_elements(pitch) eq 0) then pitch = index.sc_pitch
   if (n_elements(roll) eq 0) then roll = -index.sc_roll
endelse

if (n_elements(ccdosx) eq 0) then ccdosx=0.
if (n_elements(ccdosy) eq 0) then ccdosy=0.
if (n_elements(pmult) eq 0) then pmult=float(index.naxis1)/2.0


get_hi_params,index,pitch_hi,offset_hi,roll_hi,mu,d,_extra=extra


ang=(90. - 0.5*d)*!dpi/180.
rng=(1.0+mu)*cos(ang)/(sin(ang)+mu)

vfov=transpose([[reform(xv)],[reform(yv)]])

nst=n_elements(vfov(0,*))

vv4=dblarr(4,nst)
vv4(0,*) = ((vfov(0,*)-ccdosx)/pmult - 1.0)*rng
vv4(1,*) = ((vfov(1,*)-ccdosy)/pmult - 1.0)*rng
vv4(3,*) = 1.0

vv3=azp2cart(vv4,mu)
vv2=hi2sc(vv3,roll_hi,pitch_hi,offset_hi)
vv=sc2cart(vv2,roll,pitch,yaw)

radec = dblarr(2,nst)

rd = 180D/!dpi ;degrees to radians
for i=0,nst-1 do begin
    th=acos(vv(2,i))
    radec(1,i) = 90. - th*rd

    cphi=vv(0,i)/sin(th)
    sphi=vv(1,i)/sin(th)

    th1=acos(cphi)*rd
    th2=asin(sphi)*rd
    if (th2 gt 0) then begin
        ra = th1
    endif else begin
        ra = -th1
    endelse
    radec(0,i) = ra
endfor

return,radec
end