PRO eis_make_image, l1name, wvl, image, origin=origin, iexp=iexp, $
                    scale=scale, map=map, show=show, yip=yip, $
                    xoffset=xoffset, yoffset=yoffset, shift=shift, fix=fix, $
                    wpix=wpix

;+
; NAME:
;     EIS_MAKE_IMAGE
;
; PURPOSE:
;     Given a level-1 FITS file name and a wavelength, this routine
;     creates an image at the specified wavelength by summing 7 pixels
;     across the line profile.
;
; CATEGORY:
;     Hinode/EIS
;
; INPUTS:
;     L1NAME:  Can either be a level-1 EIS FITS file name, or a window
;              data structure produced by eis_getwindata.
;     WVL:     A wavelength (in Angstroms)
;
; OPTIONAL INPUTS:
;     IEXP:    If the data-set has multiple exposures per raster position,
;              then setting IEXP allows a particular exposure to be
;              selected. Otherwise the different exposures are
;              averaged to produce the final image.
;     XOFFSET: Used to correct the coordinates of the EIS map. E.g.,
;              if the EIS coordinate is 200, but the real solar
;              coordinate is 180, then xoffset=20.
;     YOFFSET: Used to correct the coordinates of the EIS map. E.g.,
;              if the EIS coordinate is 200, but the real solar
;              coordinate is 180, then yoffset=20.
;     WPIX:    Specifies the number of wavelength pixels to use when
;              creating the image. If the central pixel is P, then the
;              region P-WPIX to P+WPIX will be averaged. The default value
;              is 3.
;
; OUTPUTS:
;     IMAGE:   Either a 2D image array, or an IDL map structure (if
;              the keyword /MAP was set).
;
; OPTIONAL OUTPUTS:
;     ORIGIN:  A 2 element array giving the origin in heliocentric
;              coordinates of the image.
;     SCALE:   A 2 element array giving the X and Y scales for the
;              image. 
;     YIP:     The Y initial position (YIP) for the image, i.e., the
;              detector Y-pixel number that corresponds to the bottom
;              row of the raster.
;
; KEYWORD PARAMETERS:
;     MAP:     If set, then the output image is returned as an IDL
;              map, so it can be plotted with plot_map. The coordinates
;              are corrected with the values from eis_aia_offsets to
;              give better accuracy.
;     SHOW:    If set, then a plot is displayed showing the average
;              spectrum over the whole raster, with a thick vertical
;              line over-plotted showing the central wavelength that
;              has been selected by the user. Based on this, the user
;              may choose a slightly different wavelength to better
;              select the line center.
;     SHIFT:   The windata structure is processed with eis_shift_spec
;              before the image is created. This prevents distortions
;              to the image caused by the orbit variation and slit
;              tilt. 
;
; CALLS:
;     EIS_GETWINDATA(), EIS_SHIFT_SPEC(), EIS_AIA_OFFSETS(),
;     EIS_FIX_WINDATA()
;
; MODIFICATION HISTORY:
;     Ver.1, 29-Aug-2008, Peter Young
;     Ver.2, 3-Dec-2008, Peter Young
;       Corrected problem if the 7 pixel averaging region extends
;       outside of the wavelength window.
;     Ver.3, 5-Dec-2008, Peter Young
;       Added ORIGIN= keyword.
;     Ver.4, 23-Oct-2009, Peter Young
;       Allowed l1name to be a window data structure.
;     Ver.5, 28-Jun-2011, Peter Young
;       For data-sets with multiple exposures per raster position, the
;       routine now averages the exposures.
;     Ver.6, 2-Apr-2012, Peter Young
;       Added IEXP= keyword.
;     Ver.7, 20-Aug-2012, Peter Young
;       Now processes the windata structure with eis_shift_spec; added
;       /MAP keyword; added /SHOW keyword.
;     Ver.8, 6-Sep-2012, Peter Young
;       Added YIP= optional output.
;     Ver.9, 11-Feb-2013, Peter Young
;       Added XOFFSET= and YOFFSET= keywords.
;     Ver.10, 21-Aug-2013, Peter Young
;       The call to eis_shift_spec is only performed now if /shift is
;       set, as the call was taking a long time for some data-sets.
;     Ver.11, 9-Jul-2014, Peter Young
;       Modified pointing information in the output map (when /map is
;       set) to include the new EIS-AIA offset
;       information. Also the 2" X-offset for the period prior
;       to the grating position change is now applied.
;     Ver.12, 22-Oct-2014, Peter Young
;       I've changed the definition of origin to fix problem
;       when there are strange coordinate values.
;     Ver.13, 26-Jun-2017, Peter Young
;       I've added the tags 'wvl' and 'slit_ind' to the map
;       output structure.
;     Ver.14, 22-Jul-2022, Peter Young
;       Added optional input wpix=.
;-

IF n_params() LT 3 THEN BEGIN
  print,'Use:  IDL> eis_make_image, l1name, wvl, image [, origin=, iexp=, scale=, /map'
  print,'                             /show, yip=, xoffset=, yoffset=, /shift, wpix=]'
  return
ENDIF

IF n_elements(wpix) EQ 0 THEN wpix=3

IF datatype(l1name) EQ 'STR' THEN BEGIN
  wd=eis_getwindata(l1name,wvl,/refill)
  IF keyword_set(fix) THEN wd=eis_fix_windata(wd,/quiet)
  IF keyword_set(shift) THEN wd=eis_shift_spec(temporary(wd))
ENDIF ELSE BEGIN
  wd=l1name
  IF keyword_set(fix) THEN wd=eis_fix_windata(wd,/quiet)
  IF wd.wave_corr_set EQ 1 AND keyword_set(shift) THEN BEGIN
    wd=eis_shift_spec(temporary(wd))
  ENDIF 
  IF wvl GT max(wd.wvl) OR wvl LT min(wd.wvl) THEN BEGIN
    print,'%EIS_MAKE_IMAGE: The specified wavelength is not consistent with the window data structure.'
    print,'Returning...'
    return
  ENDIF 
ENDELSE

yip=wd.hdr.yws
nexp_prp=wd.hdr.nexp_prp

IF nexp_prp GT 1 AND n_elements(iexp) EQ 0 THEN BEGIN
  exp_times=reform(wd.exposure_time[0,*])
  print,'%EIS_MAKE_IMAGE: This data-set has multiple exposures per raster position. The exposure times are: '
  FOR i=0,nexp_prp-1 DO BEGIN
    print,format='(i5,f6.1,"s")',i,exp_times[i]
  ENDFOR
  print,'The exposures will be averaged together to produce the output image. Please use the keyword IEXP= to '
  print,'select an individual exposure.'
  print,''
ENDIF 

IF nexp_prp EQ 1 AND n_elements(iexp) NE 0 THEN BEGIN
  print,'%EIS_MAKE_IMAGE: nexp_prp=1, so IEXP keyword is ignored.'
ENDIF 

;
; PRY, 22-Oct-2014: I've changed definition of origin.
;
origin=[wd.solar_x[0],wd.solar_y[0]]
scale=wd.scale

int=wd.int
wvl_arr=wd.wvl

getmin=min(abs(wvl-wvl_arr),imin)

n=n_elements(wvl_arr)

;
; i0:i1 is the wavelength pixel range that is averaged.
;
i0=max([imin-wpix,0])
i1=min([imin+wpix,n-1])

missing=wd.missing

IF wd.hdr.nexp_prp GT 1 THEN BEGIN
  IF n_elements(iexp) NE 0 THEN BEGIN
    IF iexp GE wd.hdr.nexp_prp THEN BEGIN
      print,'%EIS_MAKE_IMAGE: iexp must take a value between 0 and '+trim(wd.hdr.nexp_prp-1)
      print,'Returning...'
      return
    ENDIF ELSE BEGIN
      IF keyword_set(show) THEN BEGIN
        spec=average(int[*,*,*,iexp],[2,3],missing=missing)
        plot,wd.wvl,spec,psym=10,/xsty
        oplot,wvl_arr[imin]*[1,1],[-1e10,1e10],th=2
      ENDIF 
      image=average(int[i0:i1,*,*,iexp],1,missing=missing)
    ENDELSE 
  ENDIF ELSE BEGIN 
    IF keyword_set(show) THEN BEGIN
      spec=average(int,[2,3,4],missing=missing)
      plot,wd.wvl,spec,psym=10,/xsty
      oplot,wvl_arr[imin]*[1,1],[-1e10,1e10],th=2
    ENDIF 
    image=average(int[i0:i1,*,*,*],[1,4],missing=missing)
  ENDELSE 
ENDIF ELSE BEGIN
  IF keyword_set(show) THEN BEGIN
    spec=average(int,[2,3],missing=missing)
    plot,wd.wvl,spec,psym=10,/xsty
    oplot,wvl_arr[imin]*[1,1],[-1e10,1e10],th=2
  ENDIF 
  image=average(int[i0:i1,*,*],1,missing=missing)
ENDELSE 

IF keyword_set(map) THEN BEGIN
  xcen=wd.xcen
  ycen=wd.ycen
 ;
 ; Get EIS offset
 ;
  xy=eis_aia_offsets(wd.hdr.date_obs)
  xcen=xcen+xy[0]
  ycen=ycen+xy[1]
 ;
 ; Check if X-offset prior to grating move needs to be applied
 ;
  date_obs_tai=anytim2tai(wd.hdr.date_obs)
  check_tai=anytim2tai('24-aug-2008')
  IF date_obs_tai LT check_tai THEN BEGIN
    IF wvl LT 230 THEN xcen=xcen-2.0
  ENDIF 
 ;
 ; Apply manual offsets (if any)
 ;
  IF n_elements(xoffset) NE 0 THEN xcen=xcen-xoffset ELSE xcen=xcen
  IF n_elements(yoffset) NE 0 THEN ycen=ycen-yoffset ELSE ycen=ycen
 ;
 ; Get slit index
 ;
  slit_ind=wd.hdr.slit_ind
 ;
  map=make_map(image, xc=xcen, yc=ycen, dx=wd.scale[0], dy=wd.scale[1], $
               time=wd.hdr.date_obs, id='Hinode/EIS', $
               missing=missing, wvl=wvl, slit_ind= slit_ind)
  image=map
ENDIF



END