;+
; NAME :    HR_XY2SPH (procedure)
; PURPOSE :
;   This procedure is to calculate the spherical coordinates of a given point on
;   the image from Huairou vector amgnetograph when its Cartesian coordinates
;   (X,Y) are provided, or vice when set keyword /rev.
; CATEGORY :
;   idl/huairou
; CALLING SEQUENCE :
;   hr_xy2sph,fname,x,y,theta,phi,gx=gx,gy=gy,dir=dir,pix=pix,rev=rev,h=h,$
;		silent=silent
; INPUTS :
;   fname : file name containing Huairou observation data (hbyymmdd.hhmm.fit)
;   X & Y: Cartesian coordinates of the given point (X > 0 for west/right, Y > 0
;          for north/up). The origin is at the disk center. When set /rev, X is
;          the latitude ( >0 for north) and Y the longitude ( >0 for west)
; OPTIONAL INPUT PARAMETERS :
;   none
; OUTPUTS :
;   theta: latitude of the given point (theta>0 for North) or X (>0 for west)
;   phi:   longitude of the given point (phi>0 for West) or Y (>0 for north)
; KEYWORD PARAMETERS :
;   gx=gx  : if set, it provides the X coordinate (arcsec) of FOV center
;               default: get it from file head
;   gy=gy  : if set, it provides the Y coordinate (arcsec) of FOV center
;               default: get it from file head
;   dir=dir: directory where file is located
;   pix=pix: pix scale of the image/magnetogram. default: 1.830"/pixel
;   rev=rev: if set, calculate Cartesian coordinates
;	h=h    : if set, take fname as a variable containing the head info
;	silent : if set, silently excute the program
; COMMON BLOCKS :   none
; SIDE EFFECTS :
; RESTRICTIONS :
;   You are requested to conform with the image orientation assumption, i.e.,
;   the North is at the top and west at the right.
; PROCEDURE :
;
; MODIFICATION HISTORY :
;   Hui LI, 2001/05/24 writing
;	Hui LI, 2001/06/04 add /h and /silent
;--------------------------------------- ---------------------------------
;-
pro hr_xy2sph,fname,x,y,theta,phi,gx=gx,gy=gy,dir=dir,pix=pix,rev=rev,$
	h=h,silent=silent

if keyword_set(pix) then pix1=pix else pix1=1.830

; -------------- check and read file ----------------

if not keyword_set(h) then begin
	if keyword_set(dir) then begin
    	if !version.os_family eq 'Windows' then begin
    	    if strmid(dir,strlen(dir)-1,1) eq '\' then dir=dir else dir=dir+'\'
	    endif else begin
        	if strmid(dir,strlen(dir)-1,1) eq '/' then dir=dir else dir=dir+'/'
    	endelse
	endif else begin
    	dir=''
	endelse
	fname=dir+fname
	ff=findfile(fname)
	if ff(0) eq '' then begin
    	print,'File '+fname+' not found ! Procedure stopped.'
    	return
	endif
	img=readfits(fname,hd,date0,time,/sil)
	s=size(img) & nx=s(1) & ny=s(2)
endif else begin
	hd=fname
	nx=sxpar(hd,'NAXIS1')
	ny=sxpar(hd,'NAXIS2')
endelse

; ---------- calculate solar parameters --------

date0=sxpar(hd,'DATE-OBS')
time=sxpar(hd,'TIME-OBS')
if keyword_set(gx) then gx=gx else gx=sxpar(hd,'ORTH-COX')
if keyword_set(gy) then gy=gy else gy=sxpar(hd,'ORTH-COY')

tmp=strupcase(strmid(date0,0,3))
if tmp eq 'JAN' then mon='01'
if tmp eq 'FEB' then mon='02'
if tmp eq 'MAR' then mon='03'
if tmp eq 'APR' then mon='04'
if tmp eq 'MAY' then mon='05'
if tmp eq 'JUN' then mon='06'
if tmp eq 'JUL' then mon='07'
if tmp eq 'AUG' then mon='08'
if tmp eq 'SEP' then mon='09'
if tmp eq 'OCT' then mon='10'
if tmp eq 'NOV' then mon='11'
if tmp eq 'DEC' then mon='12'
date=strmid(date0,7,4)+mon+strmid(date0,4,2)

sunparam,'UT',date,time,rsun,p,b,cl,doy

; ---------- calculate spherical coordinates -------------------

r0=rsun/pix1
xc=gx/pix1 & yc=gy/pix1
if xc le 0 then xc=fix(xc-0.5) else xc=fix(xc+0.5)
if yc le 0 then yc=fix(yc-0.5) else yc=fix(yc+0.5)
if keyword_set(rev) then begin
    hr_sph2xy,x,y,b,p,r0,theta,phi,/noprint
    theta=theta-xc+nx/2 & phi=phi-yc+ny/2
    xx=theta & yy=phi
    if n_elements(size(xx)) eq 3 and not keyword_set(silent) then begin
	    if xx lt 0 or xx gt nx-1 or yy lt 0 or yy gt ny-1 then begin
    	    print,'The given point is out of the FOV of the Telescope'
    	endif else begin
        	print,'(X,Y) = ('+strcompress(string(theta),/rem) $
            	+','+strcompress(string(phi),/rem)+')'
       		print,'The given point locate at ('+strcompress(string(theta),/rem) $
            	+','+strcompress(string(phi),/rem)+') on the image.'
    	endelse
    endif
endif else begin
    xx=x+xc-nx/2 & yy=y+yc-ny/2
    hr_sph2xy,xx,yy,b,p,r0,theta,phi,/rev,/noprint
	if n_elements(size(theta)) eq 3 and not keyword_set(silent) then begin
    	if theta le 0 then thet=fix(theta-0.5) else thet=fix(theta+0.5)
	    if phi le 0 then ph=fix(phi-0.5) else ph=fix(phi+0.5)
	    if theta ge 0 then sn='N' else sn='S'
    	if phi ge 0 then ew='W' else ew='E'
    	print,'Spherical coordinates of the given point is '+sn+string(abs(thet),$
        	format='(i2.2)')+ew+string(abs(ph),format='(i2.2)')
    endif
endelse

end