Index: issm/trunk-jpl/src/m/parameterization/setmask.py
===================================================================
--- issm/trunk-jpl/src/m/parameterization/setmask.py	(revision 13009)
+++ issm/trunk-jpl/src/m/parameterization/setmask.py	(revision 13449)
@@ -1,21 +1,23 @@
-from numpy import *
-import FlagElements as fe
+import numpy
+from FlagElements import *
 
 def setmask(md, floatingicename, groundedicename):
-	#SETMASK - establish boundaries between grounded and floating ice.
-	#
-	#   By default, ice is considered grounded. The contour floatingicename defines nodes 
-	#   for which ice is floating. The contour groundedicename defines nodes inside an floatingice, 
-	#   that are grounded (ie: ice rises, islands, etc ...)
-	#   All input files are in the Argus format (extension .exp).
-	#
-	#   Usage:
-	#      md=setmask(md,floatingicename,groundedicename)
-	#
-	#   Examples:
-	#      md=setmask(md,'all','');
-	#      md=setmask(md,'Iceshelves.exp','Islands.exp');
+	"""
+	SETMASK - establish boundaries between grounded and floating ice.
 
-	#%Get assigned fields
+	   By default, ice is considered grounded. The contour floatingicename defines nodes 
+	   for which ice is floating. The contour groundedicename defines nodes inside an floatingice, 
+	   that are grounded (ie: ice rises, islands, etc ...)
+	   All input files are in the Argus format (extension .exp).
+
+	   Usage:
+	      md=setmask(md,floatingicename,groundedicename)
+
+	   Examples:
+	      md=setmask(md,'all','');
+	      md=setmask(md,'Iceshelves.exp','Islands.exp');
+	"""
+
+	#Get assigned fields
 	x = md.mesh.x
 	y = md.mesh.y
@@ -23,33 +25,27 @@
 
 	#Assign elementonfloatingice, elementongroundedice, vertexongroundedice and vertexonfloatingice. Only change at your own peril! This is synchronized heavily with the GroundingLineMigration module. {{{
-	elementonfloatingice = fe.FlagElements(md, floatingicename)
-	elementongroundedice = fe.FlagElements(md, groundedicename) 
+	elementonfloatingice = FlagElements(md, floatingicename)
+	elementongroundedice = FlagElements(md, groundedicename) 
 
 	#Because groundedice nodes and elements can be included into an floatingice, we need to update. Remember, all the previous 
 	#arrays come from domain outlines that can intersect one another: 
 
-	elementonfloatingice = logical_and(elementonfloatingice,~elementongroundedice)
-	elementongroundedice = ~elementonfloatingice
+	elementonfloatingice = numpy.logical_and(elementonfloatingice,numpy.logical_not(elementongroundedice))
+	elementongroundedice = numpy.logical_not(elementonfloatingice)
 
 	#the order here is important. we choose vertexongroundedice as default on the grounding line.
-	vertexonfloatingice = zeros(md.mesh.numberofvertices,'bool')
-	vertexongroundedice = zeros(md.mesh.numberofvertices,'bool')
-
-	pos=argwhere(elementongroundedice==1)
-	pos=md.mesh.elements[pos,:]-1
-	if pos.size:
-		vertexongroundedice[pos]=True
-
-	pos=argwhere(~vertexongroundedice)
-	if pos.size:
-		vertexonfloatingice[pos]=True;
-	#%}}}
+	vertexonfloatingice = numpy.zeros(md.mesh.numberofvertices,'bool')
+	vertexongroundedice = numpy.zeros(md.mesh.numberofvertices,'bool')
+	vertexongroundedice[md.mesh.elements[numpy.nonzero(elementongroundedice),:].astype(int)-1]=True
+	vertexonfloatingice[numpy.nonzero(numpy.logical_not(vertexongroundedice))]=True
+	#}}}
 
 	#Return: 
-	md.mask.elementonfloatingice = double(elementonfloatingice)
-	md.mask.vertexonfloatingice = double(vertexonfloatingice)
-	md.mask.elementongroundedice = double(elementongroundedice)
-	md.mask.vertexongroundedice = double(vertexongroundedice)
-	md.mask.vertexonwater = zeros(md.mesh.numberofvertices)
-	md.mask.elementonwater = zeros(md.mesh.numberofelements)
+	md.mask.elementonfloatingice = elementonfloatingice.astype(float)
+	md.mask.vertexonfloatingice = vertexonfloatingice.astype(float)
+	md.mask.elementongroundedice = elementongroundedice.astype(float)
+	md.mask.vertexongroundedice = vertexongroundedice.astype(float)
+	md.mask.vertexonwater = numpy.zeros(md.mesh.numberofvertices)
+	md.mask.elementonwater = numpy.zeros(md.mesh.numberofelements)
+
 	return md
