;+
;NAME:
;	mk_eis_map.pro
;PURPOSE:
;	Make a map structure from an EIS data window structure.
;	Corrects 40" slot observations to only include ~40" per raster step (and not the full data window width on the CCD)
;	Contains options to sum an observation's data cube over entire or part of any data dimension (wavelength, x, y)
;NOTES:
;	1. If ccd and whisker keywords are unset, summation will be over the wavelength bins.
;	2. If index keyword is unset, summation will be over the entire dimension.
;	3. Use the ccd keyword to help figure out what value of shift to use.
;
;INPUTS:
;	windata (struct) == EIS data window structure output by eis_getwindata.pro
;KEYWORDS:
;	ccd (/) == set to sum over multiple/all exposures (i.e. sum over raster positions in x-direction)
;	whisker (/) == set to sum over part/all of y-direction along slit (whisker plot)
;	index (int/intarr) == index or boundary of indices to include in the summation
;	shift (int) == manually set the x-shift of a 40" window from its left-hand-side
;OUTPUTS:
;	eis_map (struct) == contains the data array and the most common map tags
;
;EXAMPLES:
;IDL> eis_map = mk_eis_map(windata)
;IDL> eis_whisker_map = mk_eis_map(windata, /whisker, index=[300,511])
;
;CONTACT:
;	Gary Kilper (GSFC/ORAU; garyk1@gmail.com)
;	Dec 13, 2010: Initial public version
;	Feb 10, 2011: Corrected for negative values of windata.hdr.cdelt1
;-

function mk_eis_map, windata, ccd=ccd, whisker=whisker, index=index, shift=shift
	if ~tag_exist(windata, 'hdr') then begin
		print, '% MK_EIS_MAP: Error: Input must be an EIS data window from eis_getwindata.pro. '
		return, 0
	endif

	;set up data array and units for wavelength maps
	if keyword_set(ccd) then begin		;sum over time/exposures (x-direction of a raster)
		if size(index, /n_elements) EQ 1 then data = total(windata.int[*,index,*], 2) $
		else if size(index, /n_elements) EQ 2 then data = total(windata.int[*,index(0):index(1),*], 2) $
		else data = total(windata.int, 2)
		xcyc = float([windata.nl/2., windata.ycen])
		dxdy = [1., float(windata.scale(1))]
		xyunits = ['Wavelength Bin', 'arcsecs']
	endif else if keyword_set(whisker) then begin		;sum over slit (y-direction)
		if size(index, /n_elements) EQ 1 then data = windata.int[*,*,index] $		;when a matrix has last dimension with size of 1, it automatically reverts to a smaller dimension matrix
		else if size(index, /n_elements) EQ 2 then data = total(windata.int[*,*,index(0):index(1)], 3) $
		else data = total(windata.int, 3)
		xcyc = float([windata.nl/2., windata.nx/2.])
		dxdy = [1., 1.]
		xyunits = ['Wavelength Bin', 'Exposure']

	;set up data array and units for a normal spatial map
	endif else begin
		xcyc = float([windata.xcen, windata.ycen])
		dxdy = float(windata.scale)
		xyunits = ['arcsecs', 'arcsecs']
		slit_id = windata.hdr.slit_id
		if slit_id EQ '40"' or slit_id EQ '266"' then begin		;40" or 266" slot
			if size(windata.int, /n_dimensions) EQ 4 then begin		;sit-n-stare slot observation
				data = rotate(total(windata.int[*,*,*,0], 2), 5)		;use first observation since the second is usually saturated
			endif else begin
				d1 = abs(windata.hdr.cdelt1)
				if ~keyword_set(shift) and slit_id EQ '40"' then $
					shift = windata.nl - round(d1) - 2 $		;do not include first x-pixel of window on right-hand-side
				else shift = 0
				data = fltarr(round(d1*windata.nx), windata.ny)
				for n=0, windata.nx-1 do $
					data[round(n*d1):round((n+1)*d1)-1,*] = rotate(total(windata.int[1+shift:round((n+1)*d1)-round(n*d1)+shift,n,*], 2), 5)
			endelse
		endif else begin			;1" or 2" slit, raster or sit-n-stare	;sum over wavelength bins
			if size(index, /n_elements) EQ 1 then data = total(windata.int[index,*,*], 1) $
			else if size(index, /n_elements) EQ 2 then data = total(windata.int[index(0):index(1),*,*], 1) $
			else data = total(windata.int, 1)
			if windata.hdr.nraster EQ 1 then begin		;sit-n-stare
				xcyc(0) = windata.nx/2.
				dxdy(0) = 1.
				xyunits(0) = 'Exposure'
			endif
		endelse
	endelse

	;make and return map
	eis_map = make_map(data, xc=xcyc(0), yc=xcyc(1), dx=dxdy(0), dy=dxdy(1), xunits=xyunits(0), yunits=xyunits(1), $
			time=windata.hdr.date_obs, dur=float(anytim2tai(windata.hdr.date_end)-anytim2tai(windata.hdr.date_obs)), $
			id='EIS '+windata.line_id)
	ang = pb0r(windata.hdr.date_obs, /arcsec)
	add_prop, eis_map, b0=ang[1], rsun=ang[2]
	return, eis_map
end