; NAME:
;       function readpiece.pro (v1.0)
;
; PURPOSE:
;		This function returns a 3-dimensional subarray of pixel values from a larger input cube.
; CALLING SEQUENCE:
;       x=readpiece,filename,section
;
; INPUTS:
;       filename    - name of FITS image
;		section		- lower-left corner and upper-right corner of subarray
;
; KEYWORD PARAMETERS:

; EXAMPLE
;
;		Return a 10X10 subarray cube for all planes in the image.
;
;		IDL> x=readpiece('myimage.fits',[21,30,500,509])
;
; REFERENCE:
;
;
; MODIFICATION HISTORY:
;       Written by:	Eddie Bergeron & Don Figer, IDTL, April 8, 2003
;
;		Modified by:
;		Ernie Morse, IDTL, April 10, 2003
;			Fixed bug in calculation of # of pixels to read in each call to fits_read
;-

function readpiece,filename,section

	;	section is in form [x1,x2,y1,y2]

	; retrieve header
	fits_read,filename,data,header,/header_only

	; get dimensions of data
	xsize=sxpar(header,'NAXIS1')
	ysize=sxpar(header,'NAXIS2')
	nreads=sxpar(header,'NAXIS3')

	; define output image cube
	im=intarr(section(1)-section(0)+1,section(3)-section(2)+1,nreads)

	; loop through all images in cube and load the new image cube
	for i=0,nreads-1 do begin
		print,i
		fits_read,filename,cc,first=(double(xsize)*ysize)*i,last=((double(xsize)*ysize)*(i+1))-1
		cim=reform(cc,xsize,ysize)
		im(*,*,i)=cim(section(0):section(1),section(2):section(3))
	endfor

	return,im

end