FLOWPATH.SML

  Download

More scripts: Display Toolbar

Syntax Highlighing:

comments, key words, predefined symbols, class members & methods, functions & classes
            
######################
#
# FLOWPATH.SML
#
# Demostration SML ToolScript
#
# Requires TNTmips Version 7.4 or later
# Get the latest version of this toolscript from our website www.microimages.com
#####################
# View ToolScript
#
# The following symbols are predefined
#    class VIEW View            {use to access the view the tool script is attached to}
#    class GROUP Group          {use to access the group being viewed if the script is run from a group view}
#    class LAYOUT Layout        {use to access the layout being viewed if the script is run from a layout view}
#    number ToolIsActive        Will be 0 if tool is inactive or 1 if tool is active
#
# To use this ToolScript add one or more DEM's to your Group Display
# The DEM you actually want to use in the watershed process should be the first layer in the group
# i.e. the layer at the bottom of the Group Controls window if you want to change the DEM that you were using
# switch from the ToolScript tool to some other tool and move that DEM to the bottom of the Group Controls
# and reactivate the Toolscript it will ask if you want to switch the input DEM
# you can run several successive flowpaths and buffer zones keeping the same DEM
array numeric seedx[10];
array numeric seedy[10];
class WATERSHED w;
numeric numpts;
class RVC_RASTER DEM;
class RVC_GEOREFERENCE rGeoref;
class TRANS2D_MAPGEN transparm;
class RECT3D extentsDEM;
numeric firstpass;
class GRE_LAYER_VECTOR VecBuf;
class GRE_LAYER_VECTOR VecFlow;
class GRE_LAYER_VECTOR BasinLayer;
class GRE_LAYER_VECTOR BoundaryLayer;
class GRE_LAYER_RASTER DEMLayer;
numeric haslayers;
class XmForm dlgform;
numeric setDefaultWhenClose;
class RVC_VECTOR VectIn;
class RVC_VECTOR VecBoundary;
array numeric xPoints[10],yPoints[10];
numeric xMax,yMax,xMin,yMin;
class PromptNum PromptDistance;
class RVC_VECTOR VECFLOW;
string tempfilename$;
numeric tempinode;
string tempobjname$;
class Color bocolor;
class Color fcolor;
class Color bucolor;
class Color bacolor;
vector Buffer2, VECFLOW2, Buffer, TempBuffer, BasinVector, TempFlowBuffer;
string userflowpathFilename$, userflowpathObjname$;
string userBasinFilename$, userBasinObjname$;
string demFilename$, demObjname$;
string flowpathFilename$, flowpathObjname$;
numeric userbasinInode, demInode, flowpathInode, i;
# function to create a status to attach popup messages to
# the advantage is that if the process is quick dialog won't pop in/out
# requires tntdisp 20001012 or later
proc StatusStart() {
	class StatusHandle status;
	class StatusContext context;
	status = StatusDialogCreate(dlgform,1);
	context = StatusContextCreate(status);
	}
proc StatusStop() {
	StatusContextDestroy(context);
	StatusDialogDestroy(status);
	}
# called when user presses save button on dialog
# if something wasn't calculated by choice output object is empty
proc DoSave() {
	if (haslayers) {
		numeric answer = PopupYesNo("Save Output Objects?");
		if (answer) {
			# Get output Filename
			string destFilename$ = GetOutputFileName(" ","Watershed output","rvc");
			CreateProjectFile(destFilename$,"Watershed output");
			numeric destParentInode = 0;
			numeric userflowpathInode = ObjectNumber(userflowpathFilename$,userflowpathObjname$,"VECTOR");
			CopyObject(userflowpathFilename$,userflowpathInode,destFilename$,destParentInode);
			CreateVector(Buffer2,destFilename$,"Buffer","Flow path buffer zone","Polygonal");
			Buffer2 = VectorToBufferZone(VectIn,"line",PromptDistance.value,"meters");
			userbasinInode = ObjectNumber(userBasinFilename$,userBasinObjname$,"VECTOR");
			CopyObject(userBasinFilename$,userbasinInode,destFilename$,destParentInode);
		}
	}
}
# function to remove layers from display
proc DoRemove() {
	if (haslayers) {
		View.DisableRedraw = 1;
		# save user's changes to colors;
		fcolor.red = VecFlow.Line.NormalStyle.Color.red;
		fcolor.green = VecFlow.Line.NormalStyle.Color.green;
		fcolor.blue = VecFlow.Line.NormalStyle.Color.blue;
		bucolor.red = VecBuf.Line.NormalStyle.Color.red;
		bucolor.green = VecBuf.Line.NormalStyle.Color.green;
		bucolor.blue = VecBuf.Line.NormalStyle.Color.blue;
		bacolor.red = BasinLayer.Line.NormalStyle.Color.red;
		bacolor.green = BasinLayer.Line.NormalStyle.Color.green;
		bacolor.blue = BasinLayer.Line.NormalStyle.Color.blue;
		bocolor.red = BoundaryLayer.Line.NormalStyle.Color.red;
		bocolor.green = BoundaryLayer.Line.NormalStyle.Color.green;
		bocolor.blue = BoundaryLayer.Line.NormalStyle.Color.blue;
		VecFlow = LayerDestroy(VecFlow);
		VecBuf = LayerDestroy(VecBuf);
		BasinLayer = LayerDestroy(BasinLayer);
		View.DisableRedraw = 0;
		haslayers = 0;
		}
}
# called when user presses remove button on dialog
proc cbDoRemove() {
	DoRemove();
	ViewRedrawIfNeeded(View);
}
# called when the user presses the set button on dialog
# sets the number of seedpoints to use
proc DoSet() {
	numpts = PopupNum("Enter number of seed points", numpts,1,10);
	i = 1;
	}
# Close the window, switching to default tool
# called when user presses close button and when clicks x in dialog
proc DoClose() {
	if (setDefaultWhenClose) {
		setDefaultWhenClose = false;
		View.SetDefaultTool();
	}
}
# Called the first time the tool is activated.
# If the tool implements a dialog it should be created (but not displayed) here.
func OnInitialize () {
# setup intial color for vector layers
class Color bocolor;
class Color fcolor;
class Color bucolor;
class Color bacolor;
	bocolor.red = 65535;		bocolor.green = 0;			bocolor.blue = 0;
	fcolor.red = 0;			fcolor.green = 0;				fcolor.blue = 65535;
	bucolor.red = 65535;		bucolor.green = 65535;		bucolor.blue = 0;
	bacolor.red = 0;			bacolor.green = 65535;		bacolor.blue = 0;
	#StatusStart();
	# get the raster input DEM
	if (Group.FirstLayer.Type == "Raster") {
		DispGetRasterFromLayer(DEM,Group.FirstLayer);
		DEMLayer = Group.FirstLayer;
		# get georeference and extents for the DEM
		DEM.GetDefaultGeoref(rGeoref);
		rGeoref.GetTransParm(transparm, 0, rGeoref.GetCalibModel() );
		extentsDEM = DEMLayer.Extents;		# extents in object coordinates
		extentsDEM = transparm.ConvertRectFwd(extentsDEM);		# convert to map coordinates
		}
	else {
		PopupString("First Layer must be a raster object for Watershed Toolscript");
		WaitForExit();
		}
	demFilename$ = GetObjectFileName(DEM);
	demInode = GetObjectNumber(DEM);
	demObjname$ = GetObjectName(demFilename$,demInode);
	# Initialize watershed object (assigns object handle)
	w = WatershedInit(demFilename$,demObjname$);
	# Fill all depressions in DEM and compute watersheds.
	# a depressionless version of the DEM is automatically created
	# as a temporary internal object associated with watershed handle w
	WatershedCompute(w,"FillAllDepressions,FlowPath");
	# get the Flow Paths for the entire raster
	WatershedGetObject(w,"VectorFlowPath",flowpathFilename$,flowpathObjname$);
	flowpathInode = ObjectNumber(flowpathFilename$,flowpathObjname$,"VECTOR");
	CreateTempVector(VECFLOW);
	tempfilename$ = GetObjectFileName(VECFLOW);
	tempinode = GetObjectNumber(VECFLOW);
	tempobjname$ = GetObjectName(tempfilename$,tempinode);
	tempinode = CopyObject(flowpathFilename$,flowpathInode,tempfilename$);
	tempfilename$ = GetObjectFileName(VECFLOW);
	tempobjname$ = GetObjectName(tempfilename$,tempinode);
	CloseVector(VECFLOW);
	OpenVector(VECFLOW2,tempfilename$,tempobjname$);
	firstpass = 1;
	haslayers = 0;
	numpts = 1;
	# setup the dialog
	dlgform = CreateFormDialog("FlowPath and Buffer Zone");
	WidgetAddCallback(dlgform.Shell.PopdownCallback,DoClose);
	class PushButtonItem btnItemSave;
	class PushButtonItem btnItemRemove;
	class PushButtonItem btnItemSet;
	class PushButtonItem btnItemClose;
	btnItemSave = CreatePushButtonItem("Save Output Layers...",DoSave);
	btnItemSave.IconName = "save";
	btnItemRemove = CreatePushButtonItem("Remove Output Layers",cbDoRemove);
	btnItemRemove.IconName = "remove_sel";
	btnItemSet = CreatePushButtonItem("Set Number of Seedpoints...",DoSet);
	btnItemSet.IconName = "apply";
	btnItemClose = CreatePushButtonItem("Close",DoClose);
	btnItemClose.IconName = "delete";
	class XmRowColumn btnrowaction;
	btnrowaction = CreateIconButtonRow(dlgform,btnItemSave,btnItemRemove,btnItemSet,btnItemClose);
	btnrowaction.TopWidget = dlgform;
	btnrowaction.RightWidget = dlgform;
	btnrowaction.LeftWidget = dlgform;
	class XmToggleButton btnSnap;
	class XmToggleButton btnFlow;
	class XmToggleButton btnBasin;
	class XmToggleButton btnBuffer;
	btnSnap = CreateToggleButton(dlgform,"Move Seed Point to FlowPath");
	btnFlow = CreateToggleButton(dlgform,"Compute Flow Path");
	btnBasin = CreateToggleButton(dlgform,"Compute Upstream Basin");
	btnBuffer = CreateToggleButton(dlgform,"Compute Buffer Zone");
	PromptDistance = CreatePromptNum(dlgform,"Buffer Distance",5,0,100);
	btnSnap.Set = 1;
	btnFlow.Set = 1;
	btnBasin.Set = 1;
	btnBuffer.Set = 1;
	btnSnap.TopWidget = btnrowaction;
	btnSnap.LeftWidget = dlgform;
	btnFlow.TopWidget = btnSnap;
	btnFlow.LeftWidget = dlgform;
	btnBasin.TopWidget = btnFlow;
	btnBasin.LeftWidget = dlgform;
	btnBuffer.TopWidget = btnBasin;
	btnBuffer.LeftWidget = dlgform;
	PromptDistance.TopWidget = btnBuffer;
	PromptDistance.LeftWidget = dlgform;
	class XmSeparator btnsep;
	btnsep = CreateHorizontalSeparator(dlgform);
	btnsep.TopWidget = PromptDistance;
	btnsep.TopOffset = 4;
	btnsep.LeftWidget = dlgform;
	btnsep.RightWidget = dlgform;
	#StatusStop();
	}  # end of OnInitialize
#compute and display flowpath,buffer and buffer zone if chosen by user
proc DoFlowPath() {
	#StatusStart();
	# compute vector flow paths originating at seed point location.
	# This step requires the previous computation  of the depressionless DEM
	# don't draw view until we are done
	View.DisableRedraw = 1;
	# only calculate what we have to
	if ((btnFlow.Set == 0) and (btnBuffer.Set == 0) and (btnBasin.Set == 0)) {
		return;
	}
	if ((btnFlow.Set == 1) or (btnBuffer.Set == 1) and (btnBasin.Set == 1)) {
		WatershedComputeElements(w,seedx,seedy,numpts,"FlowPath,Basin");
	}
	if ((btnBasin.Set == 1) and (btnBuffer.Set == 0) and (btnFlow.Set == 0)){
		WatershedComputeElements(w,seedx,seedy,numpts,"Basin");
	}
	if (btnBasin.Set == 0) {
		WatershedComputeElements(w,seedx,seedy,numpts,"FlowPath");
	}
	if  ((btnFlow.Set == 1) or (btnBuffer.Set == 1)) {
	WatershedGetObject(w,"VectorUserFlowPath",userflowpathFilename$,userflowpathObjname$);
	OpenVector(VectIn,userflowpathFilename$,userflowpathObjname$);
	}
	if (btnFlow.Set == 1) { # had to calculate it for buffer now add it if they want it
		VecFlow = GroupQuickAddVectorVar(Group,VectIn);
		# change the displayed color according to user (previous) choice
		VecFlow.Line.NormalStyle.Color.red = fcolor.red;
		VecFlow.Line.NormalStyle.Color.green = fcolor.green;
		VecFlow.Line.NormalStyle.Color.blue = fcolor.blue;
	}
	# Compute Buffer zone around flow path
	if (btnBuffer.Set == 1) {
		CreateTempVector(Buffer);
		CreateTempVector(TempBuffer);
		TempBuffer = VectorToBufferZone(VectIn,"line",PromptDistance.value,"meters");
		Buffer = VectorExtract(VecBoundary,TempBuffer,"InsideClip");
		VecBuf = GroupQuickAddVectorVar(Group,Buffer);
		VecBuf.Line.NormalStyle.Color.red = bucolor.red;
		VecBuf.Line.NormalStyle.Color.green = bucolor.green;
		VecBuf.Line.NormalStyle.Color.blue = bucolor.blue;
		}
	if (btnBasin.Set == 1) {
		WatershedGetObject(w,"VectorUserBasin",userBasinFilename$,userBasinObjname$);
		OpenVector(BasinVector,userBasinFilename$,userBasinObjname$);
		BasinLayer = GroupQuickAddVectorVar(Group,BasinVector);
		BasinLayer.Line.NormalStyle.Color.red = bacolor.red;
		BasinLayer.Line.NormalStyle.Color.green = bacolor.green;
		BasinLayer.Line.NormalStyle.Color.blue = bacolor.blue;
	}
	BoundaryLayer.Line.NormalStyle.Color.red = bocolor.red;
	BoundaryLayer.Line.NormalStyle.Color.green = bocolor.green;
	BoundaryLayer.Line.NormalStyle.Color.blue = bocolor.blue;
	View.DisableRedraw = 0;
	ViewRedrawIfNeeded(View);
	haslayers = 1;
	#StatusStop();
	}  # end of DoFlowPath
# Called when tool is to be destroyed, will not be called if tool was never activated.
# If the tool implements a dialog it should be destroyed here.
# don't destroy added layers here since they are already destroyed before this gets called
func OnDestroy () {
	# WatershedClose(w);
	DestroyWidget(dlgform);
	}  # end of OnDestroy
# Called when tool is activated.
# If the tool implements a dialog it should be "managed" (displayed) here.
func OnActivate () {
	#StatusStart();
	if (!firstpass) {
		answer = PopupYesNo("Keep Same Input DEM?");
		if (!answer) {
			CloseRaster(DEM);
			if (Group.FirstLayer.Type == "Raster") {
				DispGetRasterFromLayer(DEM,Group.FirstLayer);
				DEMLayer = Group.FirstLayer;
				}
		else {
			PopupString("First Layer must be a raster object for Watershed Toolscript");
			WaitForExit();
			}
		demFilename$ = GetObjectFileName(DEM);
		demInode = GetObjectNumber(DEM);
		demObjname$ = GetObjectName(demFilename$,demInode);
		# Initialize watershed object (assigns object handle)
		w = WatershedInit(demFilename$,demObjname$);
		# Fill all depressions in DEM and compute watersheds.
		# a depressionless version of the DEM is automatically created
		# as a temporary internal object associated with watershed handle w
		WatershedCompute(w,"FillAllDepressions,FlowPath");
		# get the Flow Paths for the entire raster
		WatershedGetObject(w,"VectorFlowPath",flowpathFilename$,flowpathObjname$);
		flowpathInode = ObjectNumber(flowpathFilename$,flowpathObjname$,"VECTOR");
		CreateTempVector(VECFLOW);
		tempfilename$ = GetObjectFileName(VECFLOW);
		tempinode = GetObjectNumber(VECFLOW);
		tempobjname$ = GetObjectName(tempfilename$,tempinode);
		tempinode = CopyObject(flowpathFilename$,flowpathInode,tempfilename$);
		tempfilename$ = GetObjectFileName(VECFLOW);
		tempobjname$ = GetObjectName(tempfilename$,tempinode);
		CloseVector(VECFLOW);
		OpenVector(VECFLOW2,tempfilename$,tempobjname$);
		}
	}
	# draw vector box around DEM so can show and clip buffer zone to it
	CreateTempVector(VecBoundary,"VectorToolkit", "", extentsDEM);
	class RVC_GEOREFERENCE boxGeoref;
	boxGeoref.SetCoordRefSys(rGeoref.GetCoordRefSys());
	boxGeoref.SetImplied();
	boxGeoref.Make(VecBoundary);
	VectorToolkitInit(VecBoundary);
	xPoints[1] = extentsDEM.x1;
	yPoints[1] = extentsDEM.y1;
	xPoints[2] = extentsDEM.x2;
	yPoints[2] = extentsDEM.y1;
	xPoints[3] = extentsDEM.x2;
	yPoints[3] = extentsDEM.y2;
	xPoints[4] = extentsDEM.x1;
	yPoints[4] = extentsDEM.y2;
	xPoints[5] = extentsDEM.x1;
	yPoints[5] = extentsDEM.y1;
	VectorAddLine(VecBoundary, 5, xPoints, yPoints);
	VectorValidate(VecBoundary);
	View.DisableRedraw = 1;
	BoundaryLayer = GroupQuickAddVectorVar(Group,VecBoundary);
	BoundaryLayer.IgnoreExtents = 1;
	BoundaryLayer.Line.NormalStyle.Color.red = bocolor.red;
	BoundaryLayer.Line.NormalStyle.Color.green = bocolor.green;
	BoundaryLayer.Line.NormalStyle.Color.blue = bocolor.blue;
	View.DisableRedraw = 0;
	ViewRedrawIfNeeded(View);
	# set i
	i = 1;
	dlgform.managed = 1;
	SetPopupDialogParent(dlgform);
	setDefaultWhenClose = true;
	#StatusStop();
	} # end of OnActivate
# Called when tool is deactivated (usually when switching to another tool).
# If the tool implements a dialog it should be "unmanaged" (hidden) here.
func OnDeactivate () {
	DoRemove();
	View.DisableRedraw = 1;
	BoundaryLayer = LayerDestroy(BoundaryLayer);
	View.DisableRedraw = 0;
	ViewRedrawIfNeeded(View);
	firstpass = 0;
	setDefaultWhenClose = false;
	dlgform.managed = 0;
	}  # end of OnDeactivate
# Called when user releases 'left' pointer/mouse button.
func OnLeftButtonPress() {
	DoRemove();
	class POINT2D pt;
	numeric xhold;
	numeric yhold;
	if (i <= numpts) {
	pt.x = PointerX;
	pt.y = PointerY;
	xhold = PointerX;
	yhold = PointerY;
	pt = TransPoint2D(pt,ViewGetTransLayerToScreen(View, DEMLayer, 1));
	seedx[i] = pt.x;
	seedy[i] = pt.y;
	# make sure seed points are within bounds of raster
	if (seedx[i] > (DEM.$Info.NumCols-1))
		seedx[i] = (DEM.$Info.NumCols-1);
	if (seedx[i] < 0)
		seedx[i] = 0;
	if (seedy[i] > (DEM.$Info.NumLins-1))
		seedy[i] = (DEM.$Info.NumLins-1);
	if (seedy[i] < 0)
		seedy[i] = 0;
	# if user wants to move seedpoint to flowpath
	if (btnSnap.Set) {
		class RVC_VECTOR tempflowvector;
		class RVC_VECTOR stub;
		numeric linenumber;
		numeric a;
		numeric b;
		numeric returndistance;
		numeric tempx;
		numeric tempy;
		array numeric tempseedx[1];
		array numeric tempseedy[1];
		CreateTempVector(stub,"VectorToolkit", "", extentsDEM);
		tempseedx[1] = seedx[i];
		tempseedy[1] = seedy[i];
		#StatusStart();
		WatershedComputeElements(w,tempseedx,tempseedy,1,"FlowPath");
		#StatusStop();
		WatershedGetObject(w,"VectorUserFlowPath", userflowpathFilename$, userflowpathObjname$)
		OpenVector(tempflowvector, userflowpathFilename$, userflowpathObjname$);
		CreateTempVector(TempFlowBuffer);
		TempFlowBuffer = VectorToBufferZone(tempflowvector, "line", 1, "meters");
		stub = VectorExtract(TempFlowBuffer, VECFLOW2, "InsideClip");
		# if this vector has no lines then the user defined flowpath runs straight off DEM
		# and doesn't connect with the flowpath so skip these calculations, thus leaving
		# the seedpoint where it was originally
		if (NumVectorLines(stub) > 0)
			{
			pt.x = xhold;
			pt.y = yhold;
			pt = TransPoint2D(pt,ViewGetTransLayerToScreen(View,DEMLayer,1));
			pt = TransPoint2D(pt,ViewGetTransLayerToView(View,DEMLayer));
			pt = TransPoint2D(pt,ViewGetTransMapToView(View,DEMLayer.projection,1));
			MapToObject(GetLastUsedGeorefObject(stub), pt.x, pt.y, stub, tempx, tempy);
			linenumber = FindClosestLine(stub,pt.x,pt.y,GetLastUsedGeorefObject(stub),9999999.9,returndistance);
			ClosestPointOnLine(stub,linenumber,tempx,tempy,a,b);
			ObjectToMap(stub,a,b,GetLastUsedGeorefObject(stub),tempx,tempy);
			MapToObject(GetLastUsedGeorefObject(DEM),tempx,tempy,DEM,a,b);
			seedx[i] = a;
			seedy[i] = b;
			}
		CloseVector(tempflowvector);
		CloseVector(TempFlowBuffer);
		CloseVector(stub);
		}
	i = i + 1;
	}
	if (i > numpts) {
		DoFlowPath();
	i = 1;
	}
}