Index: /issm/trunk-jpl/src/m/parameterization/setmask.m
===================================================================
--- /issm/trunk-jpl/src/m/parameterization/setmask.m	(revision 17431)
+++ /issm/trunk-jpl/src/m/parameterization/setmask.m	(revision 17432)
@@ -1,3 +1,3 @@
-function md=setmask(md,floatingicename,groundedicename)
+function md=setmask(md,floatingicename,groundedicename,varargin)
 %SETMASK - establish boundaries between grounded and floating ice.
 %
@@ -15,7 +15,18 @@
 
 %some checks on list of arguments
-if ((nargin~=3) | (nargout~=1)),
+if ((mod(nargin,2)==0) | (nargout~=1))
 	help mask
 	error('mask error message');
+end
+
+if(nargin>0)
+	if(varargin(1)='icedomain')
+		icedomainname=varargin(2);	
+	else
+		error('mask error message: wrong field specified. Only icedomain allowed for now.');
+	end
+	if ~exist(icedomainname),
+        error(['setmask error message: file ' icedomainname ' not found!']);
+	end
 end
 
@@ -42,5 +53,18 @@
 
 %level sets
-md.mask.ice_levelset=-1.*ones(md.mesh.numberofvertices,1);
 md.mask.groundedice_levelset=vertexongroundedice;
 md.mask.groundedice_levelset(find(vertexongroundedice==0.))=-1.;
+
+if(nargin>0)
+	if(varargin(1)='icedomain')
+		md.mask.ice_levelset = 1.*ones(md.mesh.numberofvertices,1);
+		%use contourtomesh to set ice values inside ice domain
+		[vertexinsideicedomain,elementinsideicedomain]=ContourToMesh(elements,x,y,icedomainname,'node',1);
+		pos=find(vertexinsideicedomain==1.);
+		md.mask.ice_levelset(pos) = -1.;
+	end
+else
+	md.mask.ice_levelset = -1.*ones(md.mesh.numberofvertices,1);
+end
+
+
Index: /issm/trunk-jpl/src/m/parameterization/setmask.py
===================================================================
--- /issm/trunk-jpl/src/m/parameterization/setmask.py	(revision 17431)
+++ /issm/trunk-jpl/src/m/parameterization/setmask.py	(revision 17432)
@@ -1,6 +1,8 @@
 import numpy
+import os
+from model import *
 from FlagElements import *
 
-def setmask(md, floatingicename, groundedicename):
+def setmask(md, floatingicename, groundedicename, *args):
 	"""
 	SETMASK - establish boundaries between grounded and floating ice.
@@ -18,4 +20,13 @@
 	      md=setmask(md,'Iceshelves.exp','Islands.exp');
 	"""
+	#some checks on list of arguments
+	if not isinstance(md,model):
+		raise TypeError("setmask error message")
+
+	if len(args)%2:
+		raise TypeError("odd number of arguments provided in setmask")
+
+	#process options
+	options=pairoptions(*args)
 
 	#Get assigned fields
@@ -42,7 +53,17 @@
 
 	#level sets
-	md.mask.ice_levelset         = -1.*numpy.ones((md.mesh.numberofvertices,1))
 	md.mask.groundedice_levelset = -1.*numpy.ones((md.mesh.numberofvertices,1))
 	md.mask.groundedice_levelset[md.mesh.elements[numpy.nonzero(elementongroundedice),:]-1]=1.
 
+	if(len(args)):
+		md.mask.ice_levelset = 1.*numpy.ones((md.mesh.numberofvertices,1))
+		icedomainfile = options.getfieldvalue('icedomain','none')
+		if not os.path.exists(icedomainfile):
+			raise IOError("setmask error message: ice domain file '%s' not found." % icedomainfile)
+		#use contourtomesh to set ice values inside ice domain
+		[vertexinsideicedomain,elementinsideicedomain]=ContourToMesh(elements,x,y,icedomainfile,'node',1)
+		md.mask.ice_levelset[numpy.nonzero(vertexinsideicedomain)[0]] = -1.
+	else:
+		md.mask.ice_levelset = -1.*numpy.ones((md.mesh.numberofvertices,1))
+
 	return md