pro euvi_stereopair,file_a,file_b,index_a,index_b,image_pair,para
;+
; Project     : STEREO
;
; Name        : EUVI_STEREOPAIR 
;
; Category    : Data analysis 
;
; Explanation : This routine rotates STEREO image A into the AB plane,
;		and adjusts STEREO image B to the same Sun center position,
;		Sun distance, roll angle, and orientation in AB plane,
;		for proper stereoscopic viewing (using routines by JP Wuelser).
;
; Syntax      : IDL> euvi_stereopair,file_a,file_b,index_a,index_b,image_pair,para
;
; Inputs      : FILE_A(nfiles_a) =image array A 
;               FILE_B(nfiles_b) =image array B 
;		INDEX_A(nfiles_a)=structure array of FITS header of image A
;		INDEX_B(nfiles_b)=structure array of FITS header of image A
;
; Outputs     : IMAGE_PAIR	=fltarr(2,nx,ny) in DN/s, flatfield-subtracted
;		PARA		=structure containing image parameters
;
; History     : 31-Aug-2007, Version 1 written by Markus J. Aschwanden
;		             using routines from J.-P. Wuelser
;		May 1, 2008: add stacking of images from all files of A or B
;
; Contact     : aschwanden@lmsal.com
;-

;NUMBER OF FILES____________________________________________
nfiles_a =n_elements(file_a)
nfiles_b =n_elements(file_b)
nfiles	=nfiles_a<nfiles_b
x0	=index_a(0).crpix1
y0	=index_a(0).crpix2
for ifile=0,nfiles-1 do begin
 print,'X,Y-OFFSET=',index_a(ifile).crpix1-x0,index_a(ifile).crpix2-y0
endfor
print,'X,Y-OFFSET is not corrected (chose file selection without offsets)'

;IMAGE PARAMETERS___________________________________________
nx	=index_a(0).naxis1		;STEREO A
ny	=index_a(0).naxis2
crpix1	=index_a(0).crpix1
crpix2	=index_a(0).crpix2
cdelt1	=index_a(0).cdelt1
cdelt2	=index_a(0).cdelt2
rsun	=index_a(0).rsun
l_a     =index_a(0).hgln_obs 
b_a     =index_a(0).hglt_obs 
da      =index_a(0).dsun_obs ;distance from Sun [m]

l_b     =index_b(0).hgln_obs 	;STEREO_B
b_b     =index_b(0).hglt_obs 
db      =index_b(0).dsun_obs ;distance from Sun [m]

;IMAGE A_____________________________________________________
delt	=rsun/cdelt1
naxis	=[nx,ny]
rpix	=[crpix1,crpix2]
disp	={naxis:fix(naxis),delt_1:float(delt),rpix:rpix,rval:[0.,0.,0.],carr:0}
dmy	=sccreadfits(file_b(0),iio,/nodata)
swcs	=fitshead2wcs(iio)
print,file_a(0)
image_aa=sb_getimage(file_a(0),disp,/sc,swcs=swcs)	;,/despike) ;--> DN/s
print,minmax(image_aa)
if (nfiles_a ge 2) then begin
 for ifile=1,nfiles_a-1 do begin
  image   =sb_getimage(file_a(ifile),disp,/sc,swcs=swcs) ;,/despike) ;--> DN/s
  image_aa=image_aa+image
  print,minmax(image_aa)
 endfor
 image_aa=image_aa/float(nfiles_a)
 print,minmax(image_aa)
 print,'Stacking ',nfiles_a,' images of STEREO/A'
endif

;IMAGE B_____________________________________________________
dmy	=sccreadfits(file_a(0),iio,/nodata)
swcs	=fitshead2wcs(iio)
image_bb=sb_getimage(file_b(0),disp,/sc,swcs=swcs)	;,/despike) ;--> DN/s
print,minmax(image_bb)
if (nfiles_b ge 2) then begin
 for ifile=1,nfiles_b-1 do begin
  image   =sb_getimage(file_b(ifile),disp,/sc,swcs=swcs) ;,/despike) ;--> DN/s
  image_bb=image_bb+image
  print,minmax(image_bb)
 endfor
 image_bb=image_bb/float(nfiles_b)
 print,minmax(image_bb)
 print,'Stacking ',nfiles_b,' images of STEREO/B'
endif

;STEREO IMAGE PAIR___________________________________________
image_pair=intarr(2,nx,ny)
image_pair(0,*,*)=fix(image_bb)
image_pair(1,*,*)=fix(image_aa)

;STEREO ANGLE________________________________________________
sepa    =sqrt((l_a-l_b)^2+(b_a-b_b)^2)  ;approximation for small angles
dl_rad  =(l_a-l_b)*!pi/180.
db_rad  =(b_a-b_b)*!pi/180.
sep_rad =acos(cos(dl_rad)*cos(db_rad))
ab_rad  =atan(tan(db_rad)/sin(dl_rad))
sep_deg =(180./!pi)*sep_rad
ab_deg  =(180./!pi)*ab_rad
print,'separation angle        A-B =',sep_deg
print,'angle of AB plane       AB  =',ab_deg

;PARAMETER STRUCTURE__________________________________________
para	={nx:nx,$
	  ny:ny,$
	  cdelt1:cdelt1,$
	  cdelt2:cdelt2,$
	  crpix1:crpix1,$
	  crpix2:crpix2,$
	  da:da,$
	  db:db,$
          rsun:rsun,$
	  sep_deg:sep_deg,$
	  ab_deg:ab_deg}
help,para,/str
end
