;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;
;;	Reads an enviheader file.  
;;	    returns a struct :
;;		    ns:number of samples in the file
;;		    nl:number of lines in the file
;;		    nb:number of bands in the file
;;		    map: an envi map structure (includes utm coords and pixelsize)
;;		    desc: The description field in the .hdr
;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
function getFileInfo, name

	 envi_open_file, name, r_fid=id
	 
	 if id eq -1 then return, {errorstruct, name:name, ns:-1, nl:-1, nb:-1}

	 envi_file_query, id, nb=nb, nl=nl, ns=ns, h_map=maph, descrip=desc, $
	 	interleave=interleave, data_type=type
	 if maph eq -1 then return, {errorstruct, name:name, ns:-1, nl:-1, nb:-1}
	 
	 handle_value, maph, map
	 envi_file_mng, id=id, /remove

	 return, {imagestruct, name:name, ns:ns, nl:nl, nb:nb, map:map, desc:desc, $
	 	interleave:interleave, type:type}
end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;
;;	Sets envi header info (map and ns, nl, nb, type, interleave)
;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro setHdr, info, fname
	envi_setup_head, fname=fname, ns=info.ns, nl=info.nl, nb=info.nb, $
		interleave=info.interleave, data_type=info.type, $
		descrip=info.desc, map_info=info.map, /write
end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;
;;	return a byte mask ns x nl based on the roi specified by roi_id
;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
function buildmask,roi_id, ns, nl
	tot=ns*nl

	data=envi_get_roi(roi_id)
	
	;index = where(data ge 0 and data lt tot)
	;if index(0) eq -1 then return, -1
	
	;if n_elements(data) ne n_elements(index) then $
	;	data=data(index)

	;index = 0

	mask=bytarr(ns, nl)
	if n_elements(data) gt 10000000 then begin
		top = n_elements(data) / 100
		for i=0l, top - 1 do mask(data(i*100:((i+1)*100)-1)) = 1
		for i=0l, (n_elements(data) mod 100) -1 do mask(data(i+100*top)) = 1
	endif else $
		mask(data) = 1

	return, mask
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;
;;	return a byte mask ns x nl based on the roi specified by roi_id
;;		in addition, cut out the roi specified by cut_roi
;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
function buildfinalmask, roi_id, cut_roi, ns, nl
	tot=ns*nl

	data=envi_get_roi(roi_id)
	
	index = where(data ge 0 and data lt tot)
	if index(0) eq -1 then begin
		print, 'No valid data in ROI'
		return, -1
	endif
	
	if n_elements(data) ne n_elements(index) then $
		data=data(index)
	index = 0

	mask=bytarr(ns, nl)
	mask(data) = 1
	data=0
	
	hole=envi_get_roi(cut_roi)
	holdex = where(hole ge 0 and hole lt tot)
	if holdex(0) eq -1 then return, mask
	
	mask(hole(holdex)) = 0
	return, mask	
end


;;given a list of UTM coordinates and a ENVI map struct converts UTM to image coordinates
function imagespace, mapcoords, mapinfo	

;; need to check and convert vector coordinates to the current map
;; system...
  iproj = ENVI_PROJ_CREATE(/utm,zone=14)
  oproj = mapinfo.proj

  ENVI_CONVERT_PROJECTION_COORDINATES, transpose(mapcoords[0,*]), transpose(mapcoords[1,*]), $
    iproj, x, y, oproj

  mapcoords[0,*]=x
  mapcoords[1,*]=y
  
  mapcoords[0,*]=round((mapcoords[0,*] - mapinfo.mc[2])/mapinfo.ps[0]) - mapinfo.mc[0]
  mapcoords[1,*]=round((mapinfo.mc[3] - mapcoords[1,*])/mapinfo.ps[1]) - mapinfo.mc[1]
  return, long(mapcoords)
end


pro shptoMask, shpfile, imgfile, maskfile
	envistart
;; open and read basic parameters from the image file
	info=getFileInfo(imgfile)
	mapinfo = info.map

	tmproi = envi_create_roi(ns=info.ns, nl=info.nl)
	cutroi = envi_create_roi(ns=info.ns, nl=info.nl)
	
	myshape=OBJ_NEW('IDLffShape', shpfile)
	myshape->GetProperty, n_entities=n_ents
        
	for i=0, n_ents-1 do begin
		cur_ent = myshape->GetEntity(i)

		points = *(cur_ent.vertices)
		parts = *(cur_ent.parts)

		imgcoords = imagespace(points, mapinfo)
		myshape->DestroyEntity, cur_ent
		for j=0, n_elements(parts)-1 do begin
	
			if j lt n_elements(parts)-1 then begin
				xpts=transpose(points(0,parts(j):parts(j+1)-1))
				ypts=transpose(points(1,parts(j):parts(j+1)-1))
				
				if not isInternal(points(*,parts(j):parts(j+1)-1)) then begin
					envi_define_roi, tmproi, /polygon, xpts=xpts,ypts=ypts
				endif else $
					envi_define_roi, cutroi, /polygon, xpts=xpts,ypts=ypts

			endif else begin
				if not isInternal(points(*,parts(j):n_elements(points(0,*))-1)) then begin
					envi_define_roi, tmproi, /polygon, $
						xpts=transpose(points(0,parts(j):n_elements(points(0,*))-1)), $
						ypts=transpose(points(1,parts(j):n_elements(points(0,*))-1))
				endif else $
					envi_define_roi, cutroi, /polygon, $
						xpts=transpose(points(0,parts(j):n_elements(points(0,*))-1)), $
						ypts=transpose(points(1,parts(j):n_elements(points(0,*))-1))
			endelse
		
		endfor
	endfor
	mask=buildfinalmask(tmproi, cutroi, info.ns, info.nl)
	help, mask
	envi_delete_rois, tmproi
	envi_delete_rois, cutroi

	openw, oun, maskfile, /get
	writeu, oun, mask
	close, oun
	free_lun, oun
	info.type= 1
	info.nb  = 1
	info.desc='shpToMask Output from '+shpfile+' based on '+imgfile
	setHdr, info, maskFile
	
	myshape->Close
	
	obj_destroy, myshape
	
end