PRO makefinder, ra, dec, ra_pm, dec_pm, $
                ps=ps, title=title, circle=circle, $
                mark=mark, box=box, secscale=secscale, $
                file=file, cscale=cscale, sexi=sexi, $
                DSS=ID

;Creates finder chart from SDSS images.
;Input: ra, dec- Right ascension and declination of target you want to
;                highlight and centre in image (degrees, unless the
;                /sexi switch is set).
;       ra_pm, dec_pm- Size of region to show in plot (degrees, or
;                      minutes of arc if /sexi is set).
;                      (Plus/minus values; i.e. 0.1 will produce a
;                      region 0.2 degreed on its side)
;       OPTIONAL:
;       ps='picture.ps'- Filename of postscript file you want to output image
;                        to.
;       title='Andromeda'- Title to put at the top of the image.
;       circle= rad OR [x, y, rad] (Note, currently x and y must be in degrees) 
;       mark=mark a cross at ra and dec (i.e. the centre of the image)
;       box=xysize|[xsize,ysize]|[xsize,ysize,angle]|[xsize,ysize,xpos,ypos]|
;           [xsize,ysize,xpos,ypos,angle]
;           xysize=box is square of side sysize
;           xsize, ysize=specifies dimensions of a rectangle
;           angle=rotates box by specified angle
;           xpos, ypos = position of box if not at center (MUST be in
;                        degrees, not sexigesimal)
;       secscale=arcsec; places a scale on the plot of 'arcsec' length
;       file=filename of fits file to produce image from (file should
;            be in the current working directory, of a path to the
;            file should be given
;       cscale=colour scale (0-to-1; normally uses 0.99)
;       
;SEXI
;Setting this switch allows you to input the data as [hh,mm,ss.s], [dd,mm,ss.s]
;_pm will then be interpreted as arcmins
;box size will be interpreted as being in arcsecs
;circle radius will be interpreted as being in arcsecs

;DSS
;Rather than inputting your own file on your system, you can get make makefinder
;obtain the desired image from the DSS database.  If you set DSS as a
;switch (e.g. /DSS) then the DSS is queried using the provided
;coordinates (with or without /sexi).  If you pass a string to DSS, it
;will query NED to obtain the coordinates of your object then obtain
;the corresponding DSS image (e.g. DSS='Akn564').  An error message is
;retured if no image is returned from the DSS image.
;
;Remember, if you are behind a firewall or use proxy settings, this
;could effect whether the DSS querying process can speak to the outside. 
;
;Example: makefinder, 180, 0, 0.02, 0.02, ps='Chart.ps', Title='Andromeda'
;       This would make a plot showing the 0.04x0.04 degrees area around the
;       point (180, 0).  The output would be sent to Chart.ps, and it would
;       have the title, Andromeda. 
;
;Written by James Mullaney, University of Durham
;22 July, 2008.
;Report bugs to j.r.mullaney@dur.ac.uk


IF n_elements(sexi) GT 0 THEN BEGIN
   ra = 15*ten(ra[0],ra[1],ra[2])
   dec = ten(dec[0],dec[1],dec[2])
   
   ra_pm = ra_pm/60.
   dec_pm = dec_pm/60.
ENDIF

IF n_elements(ps) GT 0 THEN window, xsi=600, ysi=600, /PIXMAP ELSE window, xsi=600, ysi=600

loadct,0

IF n_elements(ps) GT 0 THEN PSON, filename=ps, PAGE_SIZE=[21,20]

IF n_elements(file) GT 1 THEN BEGIN
   print, 'Filename must be a scalar string'
   RETURN
ENDIF

IF n_elements(file) EQ 1 THEN BEGIN
   fitsfile = file_search(file)
ENDIF ELSE BEGIN
   fitsfile = file_search('*.fi*', count = count)
   fitsfile = (count GT 1)?fitsfile[0]:fitsfile
ENDELSE

IF n_elements(ID) EQ 0 THEN BEGIN
   im = readfits(fitsfile, header)
ENDIF ELSE BEGIN
   IF size(ID, /TYPE) EQ 7 THEN BEGIN
      ;If ID is a string, then query NED to match the ID.
      querydss, ID, im, header, /NED
   ENDIF ELSE BEGIN
      ;If not, then use the RA and DEC to query the DSS
      querydss, [ra, dec], im, header
   ENDELSE
ENDELSE

adxy, header, ra, dec, x, y

;GETROT, cut_header, rot, cdelt
;HROT, cut_im, cut_header, rot_im, rot_header, rot,$
;      x, y, 2, MISSING=min(cut_im), /pivot

GETROT, header, rot, cdelt
HROT, im, header, rot_im, rot_header, rot,$
      x, y, 0, MISSING=min(im),/pivot

adxy, rot_header, ra-ra_pm, dec-dec_pm, xmax, ymin
adxy, rot_header, ra+ra_pm, dec+dec_pm, xmin, ymax

hextract, rot_im, rot_header, cut_im, cut_header, xmin, xmax, ymin, ymax

rot_im = cut_im
rot_header = cut_header

EXTAST, rot_header, astr

im_size = size(rot_im)

xyad, rot_header, 0, 0, ra_min, dec_min 
xyad, rot_header, im_size[1]-1, im_size[2]-1, ra_max, dec_max

plot, [-1000,1000],[-1000,1000], $
      XRANGE=[ra_min,ra_max], XSTYLE=1, $
      YRANGE=[dec_min,dec_max], YSTYLE=1, $
      POSITION=[0.1,0.1,0.9,0.9]

IF n_elements(cscale) EQ 0 THEN cscale=0.99

plotimage = bytscl(sigrange(alog10(rot_im),fraction=cscale))

IF n_elements(ps) GT 0 THEN plotimage = max(plotimage)-plotimage
tvimage, plotimage, $
         POSITION=[0.1,0.1,0.9,0.9]


IF n_elements(title) GT 0 THEN BEGIN
   plot, [-1000,1000],[-1000,1000], $
         XRANGE=[ra_min,ra_max], XSTYLE=1, $
         YRANGE=[dec_min,dec_max], YSTYLE=1, $
         POSITION=[0.1,0.1,0.9,0.9], $
         /NODATA, /NOERASE,$
         XTITLE='RA (Degrees; J2000)', YTITLE='Dec (Degrees; J2000)',$
         TITLE=title
ENDIF ELSE BEGIN
   plot, [-1000,1000],[-1000,1000], $
         XRANGE=[ra_min,ra_max], XSTYLE=1, $
         YRANGE=[dec_min,dec_max], YSTYLE=1, $
         POSITION=[0.1,0.1,0.9,0.9], $
         /NODATA, /NOERASE,$
         XTITLE='RA (Degrees; J2000)', YTITLE='Dec (Degrees; J2000)'
ENDELSE

IF n_elements(mark) GT 0 THEN BEGIN
   IF n_elements(ps) GT 0 THEN BEGIN
      oplot, [ra,ra], [dec,dec], psym=6, color=1
   ENDIF ELSE BEGIN
      oplot, [ra,ra], [dec,dec], psym=6, color=1
   ENDELSE
ENDIF

IF n_elements(circle) GT 0 THEN BEGIN
   IF n_elements(circle) EQ 1 THEN BEGIN
      circle = (n_elements(sexi) GT 0)?circle/3600.:circle
      tvcircle, circle, ra, dec, /DATA, color=255
   ENDIF ELSE IF n_elements(circle) EQ 3 THEN BEGIN
      circle = (n_elements(sexi) GT 0)?circle/3600.:circle
      tvcircle, circle[2], circle[0], circle[1], /DATA, color=255
   ENDIF ELSE BEGIN
      print, 'CIRCLE must either be an array of 3 values: [xc,yc,rad]'  
      print, 'OR a scaler representing radius in degrees'
      return
   ENDELSE
ENDIF

IF n_elements(box) GT 0 THEN BEGIN
   CASE n_elements(box) OF
      1:BEGIN
         ;Only the width of the SQUARE supplied:
         box = (n_elements(sexi) GT 0)?box/3600.:box
         tvbox, (n_elements(sexi) GT 0)?box/3600.:box, ra, dec, color=255, /DATA
      END
      2:BEGIN
         ;Width and length of RECTANGLE supplied:
         box[0:1] = (n_elements(sexi) GT 0)?box[0:1]/3600.:box[0:1]
         tvbox, [box[0],box[1]], ra, dec, color=255, /DATA
      END
      3:BEGIN
         ;Width and length AND angle of RECTANGLE supplied:
         box[0:1] = (n_elements(sexi) GT 0)?box[0:1]/3600.:box[0:1]
         tvbox, [box[0],box[1]], ra, dec, ANGLE=box[2], color=255, /DATA
         line_start = [ra, dec]
         line_fin = line_start+5*box[0]*[sin(box[2]*!PI/180.),cos(box[2]*!PI/180.)]
         plots, [line_start[0], line_fin[0]], [line_start[1],line_fin[1]]
      END
      4:BEGIN
         ;Position of box and width and length of RECTANGLE supplied:
         box[2:3] = (n_elements(sexi) GT 0)?box[0:3]/3600.:box[0:3]
         tvbox, [box[2],box[3]],box[0],box[1], color=255, /DATA
      END
      5:BEGIN
         ;Position of box and width and length AND angle of RECTANGLE supplied:
         box[0:3] = (n_elements(sexi) GT 0)?box[0:3]/3600.:box[0:3]
         tvbox, [box[2],box[3]],box[0],box[1],ANGLE=box[4], color=255, /DATA
         line_start = [box[0], box[1]]
         line_fin = line_start+5*box[2]*[sin(box[4]*!PI/180.),cos(box[4]*!PI/180.)]
      END
      ELSE: print, 'BOX must be upto a 5 element vector'
   ENDCASE
ENDIF

IF n_elements(secscale) GT 0 THEN BEGIN
   scale = secscale*(1./3600.)
   ra_start_scale = ra_min+0.5*(ra_max-ra_min)
   ra_end_scale = ra_min+0.5*(ra_max-ra_min) + scale
   
   dec_start_scale = dec_min+0.2*(dec_max-dec_min)
   dec_end_scale = dec_min+0.2*(dec_max-dec_min) + scale
   
   oplot, [ra_start_scale, ra_end_scale], $
          [dec_start_scale, dec_start_scale], $
          color=(n_elements(ps) GT 0)?1:255
   xyouts, (ra_start_scale + ra_end_scale)/2.,$
           dec_start_scale-0.03*(dec_max-dec_min), $
           strtrim(string(secscale,FORMAT='(F4.1)'),2)+'"', $
           alignment=0.5, color=(n_elements(ps) GT 0)?1:255
   
   oplot, [ra_end_scale, ra_end_scale], $
          [dec_start_scale, dec_end_scale], $
          color=(n_elements(ps) GT 0)?1:255
   xyouts, ra_end_scale - (ra_max-ra_min)*0.01,$
           (dec_start_scale+dec_end_scale)/2., $
           strtrim(string(secscale, FORMAT='(F4.1)'),2)+'"', $
           alignment=0.5, orientation=90, color=(n_elements(ps) GT 0)?1:255
ENDIF
   
plotsym, 2, 4
oplot, [ra_min+0.95*(ra_max-ra_min)], $
       [dec_min+0.9*(dec_max-dec_min)],$
       psym=8

xyouts, [ra_min+0.92*(ra_max-ra_min)], $
        [dec_min+0.87*(dec_max-dec_min)], 'E'

plotsym, 6, 4
oplot, [ra_min+0.95*(ra_max-ra_min)], $
       [dec_min+0.9*(dec_max-dec_min)],$
       psym=8

xyouts, [ra_min+0.96*(ra_max-ra_min)], $
        [dec_min+0.92*(dec_max-dec_min)], 'N'

IF n_elements(DSS) GT 0 THEN BEGIN
   xyouts, [ra_min+0.8*(ra_max-ra_min)], $
           [dec_min+0.05*(dec_max-dec_min)], 'DSS image'
ENDIF

IF n_elements(ps) GT 0 THEN PSOFF

END