pro aia_loop_manu,fileset,wave_,iwave,fov,nsig,vers
;+
; Project     : AIA/SDO
;
; Name        : AIA_LOOP_MANU 
;
; Category    : Data analysis   
;
; Explanation : Manual loop tracingand storing coordinates of 
;               traced loop in LOOPFILE
;
; Syntax      : IDL>aia_loop_manu,fileset,wave_,iwave,fov,nsig
;
; Inputs      : fileset  = beginning of filename 
;               wave_(6) = wavelengths
;               iwave    = number of selected wavelength for tracing
;	        fov[4]   = [x1,y1,x2,y2] field-of-view in solar radii
;		nsig	 = contrast (number of sigmas)
;               vers     = label in filename (e.g., vers='a')
;
; Outputs     : loopfile   = x,y coordinates of traced loops
;
; History     : 09-Mar-2011, Version 1 written by Markus J. Aschwanden
;               18-May-2011, change units of FOV from pixels to solar radii
;		11-Jun-2011, add VERS to specify filename
;
; Contact     : aschwanden@lmsal.com
;-

;____________________READ DATA_______________________________________
loopfile=fileset+vers+'_'+wave_(iwave)+'_manu.dat'
if (max(abs(fov)) gt 1.3) then STOP,'Specify FOV in units of solar radii'
searchfile=fileset+'*'+wave_(iwave)+'.fits'
datafiles=file_search(searchfile)
datafile =datafiles(0)
read_sdo,datafile,index0,data0
aia_prep,index0,data0,index,data
cdelt1  =index.cdelt1
crpix1  =index.crpix1
crpix2  =index.crpix2
rsun    =index.rsun_obs
rpix    =rsun/cdelt1
x1	=fov[0]
x2	=fov[2]
y1	=fov[1]
y2	=fov[3]
i1      =long(x1*rpix+crpix1+0.5)
i2      =long(x2*rpix+crpix1+0.5)
j1      =long(y1*rpix+crpix2+0.5)
j2      =long(y2*rpix+crpix2+0.5)
print,'FOV in pixels [i1,j1,i2,j1]=',i1,j1,i2,j2
image1   =data(i1:i2,j1:j2)
nx	 =i2-i1+1
ny	 =j2-j1+1

;_____________________MANUAL LOOP TRACING______________________________
window,0,xsize=nx,ysize=ny
!p.position=[0,0,1,1]
!x.range=[x1,x2]
!y.range=[y1,y2]
!x.style=1
!y.style=1
loadct,3
plot,[0,0] 
statistic,image1,zav,zsig
tv,bytscl(image1,min=zav-nsig*zsig,max=zav+nsig*zsig)
read,'Number of loops         = ',nl
read,'Number of spline points = ',np
openw,2,loopfile
for il=0,nl-1 do begin
 x_	=fltarr(np)
 y_	=fltarr(np)
 f_	=fltarr(np)
 s_	=fltarr(np)
 for ip=0,np-1 do begin
  print,'click point with cursor ip=',ip
  cursor,x,y
  print,'x,y=',x,y
  wait,0.2
  oplot,[x,x],[y,y],psym=1,symsize=3
  i	=long(x*rpix+crpix1-i1+0.5)
  j	=long(y*rpix+crpix2-j1+0.5)
  x_(ip)=x
  y_(ip)=y
  f_(ip)=image1(i,j)
  if (ip ge 1) then s_(ip)=s_(ip-1)+$
       sqrt((x_(ip)-x_(ip-1))^2+(y_(ip)-y_(ip-1))^2)
  printf,2,il,x_(ip),y_(ip),f_(ip),s_(ip) 
 endfor
endfor
close,2
print,'Loop coordinates written into file = ',loopfile
end
