#pragma TextEncoding = "UTF-8"
#pragma rtGlobals=3
#pragma version=1.03

#include <Readback ModifyStr>

static constant kCHECK_FOR_UPDATES=1  // set to 1 to check for updates on init
static constant kFREQUENCY=604800 // (seconds). updater won't check more frequently than this. weekly=604800.
static strconstant ksVERSION="1.03"
static strconstant ksPROJECT="Encircle"
static strconstant ksLOCATION="http://www.igorexchange.com/node/7895/release"
static strconstant ksSHORTNAME="Encircle"
static strconstant ksPREFSFILE="EncirclePrefs.bin"

// by tony.withers@uwo.ca

// Description
// draws a polgon around a cloud of data points on a plot.

// How to use
// Right click on a data point to invoke the trace menu and select 
// "Encircle points".
// If you reshape the graph window you can repeat this procedure to 
// update the shape of the polygon.
// Selecting "Remove Encircled Region" deletes the shaded region.
// Polygons are drawn in the user back layer and can be edited with the 
// drawing tools.

// Details
// Polygon surrounds all of the data in the wave.
// It makes most sense to apply this to x-y data.

// Version history
// 1.03 if Package Updater is installed, will check IgorExchange for new releases 
// 1.02 deal with offsets
// 1.01 bug fix: remove gratuitous print statement

Menu "TracePopup"
	"Encircle points", /Q, DrawCatmullRomBoundingSpline()
	"Remove encircled region", /Q, RemoveBoundingSpline()
end

function DrawCatmullRomBoundingSpline()
	
#if (exists("Updater#CheckForUpdates")==6)
	if(Updater#CheckForUpdates(kCHECK_FOR_UPDATES,kFREQUENCY,ksVERSION,ksPROJECT,ksLOCATION,ksSHORTNAME,ksPREFSFILE))
		return 0 // allow user to reload new version of this file before continuing
	endif
#endif
	
	// initialize (or reinitialize) the package data folder
	NewDataFolder /O root:Packages
	NewDataFolder /O root:Packages:Encircle
	
	// Create a data folder reference variable
	DFREF dfr = root:Packages:Encircle
	
	// figure out graph and trace names
	GetLastUserMenuInfo
	string s_info=traceinfo(S_graphName, S_traceName, 0)
	//print s_info
	if(numtype(NumberByKey("AXISZ",s_info))==0)
		doalert 1, "you've selected a contour trace. Continue?"
		if (v_flag==2)
			return 0
		endif
	endif	
	
	if(strlen(stringbykey("XWAVE",s_info))==0)
		duplicate /free tracenametowaveref(S_graphName, S_traceName), xdata
		xdata=x
	else
		wave xdata=XWaveRefFromTrace(S_graphName, S_traceName)
	endif	
	wave ydata=tracenametowaveref(S_graphName, S_traceName)
	
	// number of points to calculate between nodes
	variable nPoints=20	
	variable alpha=0.5
	// alpha=0 for standard (uniform) Catmull-Rom spline
	// alpha=0.5 for centripetal Catmull-Rom spline
	// alpha=1 for chordal Catmull-Rom spline
	
	if (numpnts(ydata)<3)
		doalert 0, "must have more than two data points"
		return 0
	endif
	
	Convexhull xdata,ydata
	wave W_XHull=W_XHull, W_YHull=W_YHull
	
	// deal with offsets
	W_XHull+=GetNumFromModifyStr(s_info,"offset","{",0)
	W_YHull+=GetNumFromModifyStr(s_info,"offset","{",1)
	
	variable vertices=numpnts(W_XHull)	
	
	make /free/n=(vertices+3) x_hull, y_hull
	x_hull=W_XHull[p-(vertices)*(p>(vertices-1))]
	y_hull=W_YHull[p-(vertices)*(p>(vertices-1))]	

	// find the centroid of the hull
	variable x0=0, y0=0, varA=0
	variable i, factor
	for(i=0;i<vertices; i+=1)
		factor=x_hull[i]*y_hull[i+1] - x_hull[i+1]*y_hull[i]
		x0+=(x_hull[i]+x_hull[i+1])*factor
		y0+=(y_hull[i]+y_hull[i+1])*factor
		varA+=factor
	endfor
	varA*=0.5
	x0/=6*varA
	y0/=6*varA
	
	string newXname, newYname
	newXname=cleanupname(S_traceName, 0)+"_xCR"
	newYname=cleanupname(S_traceName, 0)+"_yCR"

	// could end up overwriting a different trace's shaded region
	// if you use long or liberal wave names	
	make /o/n=0 dfr:$newXname /WAVE=CR_x, dfr:$newYname /WAVE=CR_y	 
	make /free /n=(nPoints) x_out, y_out
	make /free /n=4 x_in, y_in	
	
	for(i=1;i<=vertices;i+=1)
		x_in=x_hull[i-1+p]
		y_in=y_hull[i-1+p]
		CatmullRomSegment(x_in, y_in, x_out, y_out, alpha)
		concatenate /NP {x_out}, CR_x
		concatenate /NP {y_out}, CR_y
	endfor
	// wrap around all the way back to first point
	CR_x[numpnts(CR_x)]={CR_x[0]}
	CR_y[numpnts(CR_y)]={CR_y[0]}
	
	// add a center point at the start of polygon
	// so that we can scale outwards
	insertpoints 0, 2, CR_x, CR_y
	CR_x[0]={x0,nan}
	CR_y[0]={y0,nan}	
		
	// need to know the aspect ratio of the window...
	getwindow /Z $S_graphName psize
	variable hpixels=V_right-V_left
	variable vpixels=V_bottom-V_top
	
	string s_Xax=stringbykey("XAXIS",s_info)
	string s_Yax=stringbykey("YAXIS",s_info)
	getaxis /Q $s_Xax
	variable xaxrange=abs(V_max-V_min)
	getaxis /Q $s_Yax
	variable yaxrange=abs(V_max-V_min)
	variable xrange=(wavemax(x_hull)-wavemin(x_hull))
	variable yrange=(wavemax(y_hull)-wavemin(y_hull))
	
	// ... to figure out the aspect ratio of the bounding shape 
	hpixels=xrange*hpixels/xaxrange // horizontal spread of data in pixels
	vpixels=yrange*vpixels/yaxrange // vertical spread of data in pixels
	
	// expand a bit to inset symbols
	variable v_expand=0.2 // expansion factor
	variable hscale=1+v_expand*(vpixels/hpixels)/(1+vpixels/hpixels)
	variable vscale=1+v_expand/(1+vpixels/hpixels)
		
	// prepare to draw the polygon
	string objName=cleanupname(S_traceName, 0)+"loop"	
	DrawAction /W=$S_graphName/L=UserBack getgroup=$objName, delete
	
	SetDrawLayer /W=$S_graphName UserBack
	SetDrawEnv /W=$S_graphName gstart,gname=$objName	
	SetDrawEnv /W=$S_graphName xcoord=$s_Xax,ycoord=$s_Yax
	
	variable v_opacity=4000
	variable v_thickness=2
	
	string s_color=""
	Variable v1,v2
	v1= strsearch(s_info, "rgb(x)=", 0)
	if( v1 >= 0 )
		v1+=7 // characters in key
		v2=strsearch(s_info, ")", v1)	
		s_color=s_info[v1,v2]
		variable red, green, blue
		sscanf s_color, "(%d,%d,%d)", red, green, blue
		SetDrawEnv /W=$S_graphName linefgc= (red, green, blue),fillfgc= (red, green, blue, v_opacity)
	endif
	SetDrawEnv /W=$S_graphName fillpat= 1,linethick= v_thickness
	drawpoly /W=$S_graphName x0, y0, hscale, vscale, CR_x, CR_y
	SetDrawEnv /W=$S_graphName gstop,gname=$objName
end

function RemoveBoundingSpline()
	GetLastUserMenuInfo
	string objName=cleanupname(S_traceName, 0)+"loop"	
	
	DrawAction /L=UserBack /W=$S_graphName getgroup=$objName, delete
	DFREF dfr = root:Packages:Encircle
	string basename=cleanupname(S_traceName, 0)
	killwaves /Z dfr:$basename+"_xCR", dfr:$basename+"_yCR"	
end

// Create a spline that passes through a series of x-y points,
// looping back to the start point.
// Can be used with ConvexHull to make a loop around a cloud of points.
function CatmullRomSplineLoop(xdata, ydata [,segPoints, alpha])
	wave xdata, ydata	
	variable segPoints, alpha
	
	if( ParamIsDefault(segPoints) )
		segPoints=20
	endif
	if( ParamIsDefault(alpha) )
		alpha=0.5
	endif
	
	variable np=numpnts(ydata)	

	string newXname, newYname
	newXname=nameofwave(xdata)+"_CR"
	newYname=nameofwave(ydata)+"_CR"
	
	if(exists(newXname)||exists(newYname))
		doalert 1, "overwrite existing wave?"
		if (V_Flag==2)
			return 0
		endif
	endif
	
	make /o/n=0 $newXname /WAVE=CR_x, $newYname /WAVE=CR_y	
	make /free/n=(np+2) x_temp, y_temp
	x_temp=xdata[p-(np-1)*(p>(np-1))]
	y_temp=ydata[p-(np-1)*(p>(np-1))]	
	make  /free /n=(segPoints) x_out, y_out
	make /free /n=4 x_in, y_in	
	
	variable i
	for(i=1;i<np;i+=1)
		x_in=x_temp[i-1+p]
		y_in=y_temp[i-1+p]
		CatmullRomSegment(x_in, y_in, x_out, y_out, alpha)
		concatenate /NP {x_out}, CR_x
		concatenate /NP {y_out}, CR_y
	endfor
	// wrap around all the way back to first point
	CR_x[numpnts(CR_x)]={CR_x[0]}
	CR_y[numpnts(CR_y)]={CR_y[0]}
end

function CatmullRomSegment(w_x, w_y, xout, yout, alpha)
	wave w_x, w_y // x and y values of four control points
	wave xout, yout // output waves
	variable alpha
	
	// alpha=0 for standard (uniform) Catmull-Rom spline
	// alpha=0.5 for centripetal Catmull-Rom spline
	// alpha=1 for chordal Catmull-Rom spline
	
	variable nPoints=numpnts(xout)
	variable t0, t1, t2, t3

	t0=0
	t1=sqrt( (w_x[1]-w_x[0])^2 + (w_y[1]-w_y[0])^2 )^alpha + t0	
	t2=sqrt( (w_x[2]-w_x[1])^2 + (w_y[2]-w_y[1])^2 )^alpha + t1
	t3=sqrt( (w_x[3]-w_x[2])^2 + (w_y[3]-w_y[2])^2 )^alpha + t2
	
	make /free/n=(nPoints) w_t,A1x,A1y,A2x,A2y,A3x,A3y,B1x,B1y,B2x,B2y
	
	w_t=t1+p/(nPoints)*(t2-t1)
	A1x=(t1-w_t)/(t1-t0)*w_x[0] + (w_t-t0)/(t1-t0)*w_x[1]
	A1y=(t1-w_t)/(t1-t0)*w_y[0] + (w_t-t0)/(t1-t0)*w_y[1]
	A2x=(t2-w_t)/(t2-t1)*w_x[1] + (w_t-t1)/(t2-t1)*w_x[2]
	A2y=(t2-w_t)/(t2-t1)*w_y[1] + (w_t-t1)/(t2-t1)*w_y[2]	
	A3x=(t3-w_t)/(t3-t2)*w_x[2] + (w_t-t2)/(t3-t2)*w_x[3]
	A3y=(t3-w_t)/(t3-t2)*w_y[2] + (w_t-t2)/(t3-t2)*w_y[3]	
	B1x=(t2-w_t)/(t2-t0)*A1x + (w_t-t0)/(t2-t0)*A2x
	B1y=(t2-w_t)/(t2-t0)*A1y + (w_t-t0)/(t2-t0)*A2y
	B2x=(t3-w_t)/(t3-t1)*A2x + (w_t-t1)/(t3-t1)*A3x
	B2y=(t3-w_t)/(t3-t1)*A2y + (w_t-t1)/(t3-t1)*A3y
	
	xout=(t2-w_t)/(t2-t1)*B1x + (w_t-t1)/(t2-t1)*B2x
	yout=(t2-w_t)/(t2-t1)*B1y + (w_t-t1)/(t2-t1)*B2y
end

// this function isn't needed for drawing shaded regions.

// Create a spline that passes through a series of x-y points
// Catmull E. and Rom R. (1974) A class of local interpolating splines. 
// In Computer  Aided  Geometric Design, R.E. Barnhill and R.F. Reisenfeld, Eds. 
// Academic Press, New York, 1974, pp. 317–326.
function CatmullRomSpline(xdata, ydata [,segPoints, alpha])
	wave xdata, ydata	
	variable segPoints, alpha
	
	if( ParamIsDefault(segPoints) )
		segPoints=20
	endif
	if( ParamIsDefault(alpha) )
		alpha=0.5
	endif
	
	variable np=numpnts(ydata)	

	string newXname, newYname
	newXname=nameofwave(xdata)+"_CR"
	newYname=nameofwave(ydata)+"_CR"
	
	if(exists(newXname)||exists(newYname))
		doalert 1, "overwrite existing wave?"
		if (V_Flag==2)
			return 0
		endif
	endif
	
	make /o/n=0 $newXname /WAVE=CR_x, $newYname /WAVE=CR_y	
	make /free /n=(segPoints) x_out, y_out
	make /free /n=4 x_in, y_in	
	
	variable i
	for(i=0;i<(np-1);i+=1)
		if (i==0) // first segment
			// reflect through end points
			x_in[0]=xdata[0]-(xdata[1]-xdata[0])
			x_in[1,3]=xdata[p-1]
			y_in[0]=ydata[0]-(ydata[1]-ydata[0])
			y_in[1,3]=ydata[p-1]
		elseif(i==(np-2)) // last segment
			x_in[0,2]=xdata[i-1+p]
			x_in[3]=xdata[np-1]+(xdata[np-1]-xdata[np-2])
			y_in[0,2]=ydata[i-1+p]
			y_in[3]=ydata[np-1]+(ydata[np-1]-ydata[np-2])
		else
			x_in=xdata[i-1+p]
			y_in=ydata[i-1+p]
		endif
		CatmullRomSegment(x_in, y_in, x_out, y_out, alpha)
		concatenate /NP {x_out}, CR_x
		concatenate /NP {y_out}, CR_y
	endfor
	// insert end point
	CR_x[numpnts(CR_x)]={xdata[np-1]}
	CR_y[numpnts(CR_y)]={ydata[np-1]}
end


// --------------------  package updater code ----------------------

static structure PackagePrefs
	uint32	 version		// Preferences structure version number. 100 means 1.00.
	uint32 lastUpdate	// good for another 21 years?!
	uint32 reserved[100]	// Reserved for future use
endStructure

static function checkforupdates()	
		
	if(kCHECK_FOR_UPDATES==0)
		return 0
	endif
	
	variable chatty=0 // set to 1 for a lot of dialog 
	string s_currentVersion=ksVERSION
	string s_filename=ksSHORTNAME+".ipf"
	
	STRUCT PackagePrefs prefs
	LoadPackagePreferences ksSHORTNAME, ksPREFSFILE, 0, prefs
	if (V_flag!=0 || V_bytesRead==0)
		prefs.lastUpdate=0
		SavePackagePreferences ksSHORTNAME, ksPREFSFILE, 0, prefs		
	endif
	
	if ((datetime-prefs.lastUpdate)<kFREQUENCY)
		return 0 // silently give up
	endif
		
	variable v_major, v_minor
	sscanf s_currentVersion, "%d.%d", v_major, v_minor
	
	string cmd=""
	
	if(chatty)
		sprintf cmd, "Check for updates to %s?", ksPROJECT
		doalert 1, cmd
		if (v_flag==2)
			return 0
		endif
	endif
	printf "%s: Checking for updates\r", ksPROJECT
	
	prefs.lastUpdate=datetime // will be rounded to 32 bit integer
	SavePackagePreferences ksSHORTNAME, ksPREFSFILE, 0, prefs

	String webPageText = FetchURL(ksLOCATION)
	Variable error = GetRTError(1)	
	if (error != 0)
		sprintf cmd, "%s: Could not load %s", ksPROJECT, ksLOCATION
		print cmd
		if(chatty)
			doalert 0, cmd
		endif
		return 0
	endif
	
	variable v_start=-1, v_end=-1, v_success=1
	string s_releaseVersion="", s_url="", s_file="", s_currentIgor=""
	
	v_start= strsearch(webPageText, ksPROJECT+" IGOR", 0)
	v_success=v_success&&(v_start>-1)
	v_end= strsearch(webPageText, "</a>", v_start)
	v_start+=strlen(ksPROJECT)+6
	v_end-=1
	v_success=v_success&&(v_end>v_start)&&((v_end-v_start)<30)
	if(v_success)		
		s_releaseVersion=webPageText[v_start, v_end]
	else
		sprintf cmd, "%s: Could not find packages at %s.\r", ksPROJECT, ksLOCATION
		cmd+="Please check for new versions at http://www.igorexchange.com/"
		if (chatty)
			doalert 0, cmd
		endif
		print cmd
		return 0
	endif
		
	variable ReleaseIgorMajor, ReleaseIgorMinor, ReleaseMajor, ReleaseMinor, IgorCurrent, IgorCurrentMinor
	sscanf s_releaseVersion, "%d.%d.x-%d.%d", ReleaseIgorMajor, ReleaseIgorMinor, ReleaseMajor, ReleaseMinor
	
	s_currentIgor=StringByKey("IGORFILEVERSION", igorinfo(3))
	sscanf s_currentIgor, "%d.%d", IgorCurrent, IgorCurrentMinor
	if (IgorCurrent<ReleaseIgorMajor || (IgorCurrent==ReleaseIgorMajor &&  IgorCurrentMinor<ReleaseIgorMinor))
		sprintf cmd, "%s: Most recent version of %s requires Igor Pro %d.%d or later. Please update Igor!", ksSHORTNAME, ksPROJECT, ReleaseIgorMajor, ReleaseIgorMinor
		if (chatty)
			doalert 0, cmd
		endif
		print cmd
		return 0
	endif
	
	if((ReleaseMajor<v_major || (ReleaseMajor==v_major && ReleaseMinor<=v_minor)))
		sprintf cmd, "%s: Current version (%s) is the most recent version", ksSHORTNAME, s_currentVersion
		if (chatty)
			doalert 0, cmd
		endif
		print cmd
		return 0
	endif
	
	v_start=strsearch(webPageText, "http://www.igorexchange.com/project_download/files/"+ksSHORTNAME, v_end)
	v_end=strsearch(webPageText, ".ipf", v_start)
	if (v_start==-1 || v_end==-1 || ((v_end-v_start)>200))
		sprintf cmd, "%s: Could not resolve download location", ksSHORTNAME
		print cmd
		if (chatty)
			doalert 0, cmd
		endif
		return 0
	endif
	v_end+=3
	
	s_url=webPageText[v_start, v_end]
	
	if (chatty)
		sprintf cmd, "%s: Found version %d.%d. Download now?", ksSHORTNAME, ReleaseMajor, ReleaseMinor
		doalert 1, cmd
		if (v_flag==2)
			return 0
		endif
	endif
		
	webPageText=FetchURL(s_url)	
	killwindow /Z UpdaterTempNB
	NewNotebook /N=UpdaterTempNB/F=0/OPTS=3/V=0/ENCG=1 as "Updater Temp File" // create a hidden notebook
	notebook UpdaterTempNB text=replacestring("\r\n",webPageText,"\r")	// strip newlines from FetchURL output
		
	v_start= strsearch(s_releaseVersion, "-", 0)
	printf "%s: loaded %s%s from %s\r", ksPROJECT, ksSHORTNAME, s_releaseVersion[v_start+1,strlen(s_releaseVersion)-1], s_url
	
	// try to find current location of procedure file
	pathinfo IgorUserFiles
	s_path=FindFile(s_path, s_filename) // look in subfolders of user files folder for this procedure file
		
	if (strlen(s_path)) // s_path will be "" if file not found		
		cmd=""
		if (!chatty)
			sprintf cmd, "New update found for %s\r", ksPROJECT
		endif
		sprintf cmd, "%sOverwrite %s%s with new version?\r\r", cmd, s_path, s_filename
		doalert 1, cmd
		if(v_flag==2)
			printf "%s: update cancelled\r", ksPROJECT
			return 0
		endif
		newPath /O/Q/Z ipfPath, s_path
		SaveNotebook /O/S=6/P=ipfPath UpdaterTempNB as s_filename
	else
		sprintf cmd, "Select location of %s and click Save to overwrite", s_filename
		SaveNotebook /M=cmd/I/S=6/P=IgorUserFiles UpdaterTempNB as s_filename
	endif

	if (strlen(s_path))	
		printf "%s: saved release to %s\r", ksPROJECT, s_path
		sprintf cmd, "Allow Igor to reload %s to complete installation", s_filename
		doalert 0, cmd
	else
		printf "%s: update failed\r", ksPROJECT
	endif
	
	// clean up	
	killwindow /Z UpdaterTempNB	
	killpath /Z ipfPath
	return 1
end


// search recursively for s_file within folder described by s_path.
// don't expect this to be quick.
static function /T FindFile(s_path, s_file)
	String s_path, s_file
	
	variable folderLimit=5000 // quits after looking in this many folders.
	Variable folderIndex, fileIndex, i, startFolder, stopFolder, folderCount
	string s_folderList, fileName

	if (strlen(s_path)==0)
		// it would be nice to default to the current save location, but I couldn't find a way to do that
		s_path=SpecialDirPath("Documents", 0, 0, 0)
	endif
		
	make /o/T/n=1 w_folders
	w_folders={s_path}
	folderCount=1
	startFolder=0
	
	do
		stopFolder=numpnts(w_folders)
		for(folderIndex=startFolder;folderIndex<stopFolder; folderIndex+=1)
			newPath /O/Q/Z tempPath, w_folders[folderIndex]
			s_folderlist=IndexedDir(tempPath, -1, 0)
			
			for (i=0;i<itemsinList(s_folderlist);i+=1)
				w_folders[numpnts(w_folders)]={w_folders[folderIndex]+stringfromlist(i,s_folderlist)+":"}				
				folderCount+=1
			endfor
		endfor
		if(numpnts(w_folders)==startFolder || folderCount>folderLimit)
			break
		endif
		startFolder=stopFolder
	while(1)

	for(folderIndex=0;folderIndex<folderCount;folderIndex+=1)
		newPath /O/Q/Z tempPath, w_folders[folderIndex]
		fileindex=0
		do
			fileName=IndexedFile(tempPath, fileindex, ".ipf")
			if(stringmatch(s_file, fileName))
				killpath /Z tempPath
				return w_folders[folderIndex]
			endif
			fileindex+=1
		while(strlen(fileName))
	endfor
	killpath /Z tempPath
	return ""
end