Index: /issm/trunk-jpl/src/m/mesh/triangle.py
===================================================================
--- /issm/trunk-jpl/src/m/mesh/triangle.py	(revision 13519)
+++ /issm/trunk-jpl/src/m/mesh/triangle.py	(revision 13520)
@@ -1,56 +1,66 @@
-from numpy import *
-import TriMesh as tm
-import NodeConnectivity as nc
-import ElementConnectivity as ec
+import numpy
+from TriMesh import *
+from NodeConnectivity import *
+from ElementConnectivity import *
+from MatlabFuncs import *
 
-def triangle(md, domainname, resolution,riftname=''):
-	#TRIANGLE - create model mesh using the triangle package
-	#
-	#   This routine creates a model mesh using TriMesh and a domain outline, to within a certain resolution
-	#   where md is a @model object, domainname is the name of an Argus domain outline file, 
-	#   and resolution is a characteristic length for the mesh (same unit as the domain outline
-	#   unit). Riftname is an optional argument (Argus domain outline) describing rifts.
-	#
-	#   Usage:
-	#      md=triangle(md,domainname,resolution)
-	#   or md=triangle(md,domainname, resolution, riftname)
-	#
-	#   Examples:
-	#      md=triangle(md,'DomainOutline.exp',1000);
-	#      md=triangle(md,'DomainOutline.exp',1000,'Rifts.exp');
+def triangle(md,domainname,*args):
+	"""
+	TRIANGLE - create model mesh using the triangle package
 
+	   This routine creates a model mesh using TriMesh and a domain outline, to within a certain resolution
+	   where md is a @model object, domainname is the name of an Argus domain outline file, 
+	   and resolution is a characteristic length for the mesh (same unit as the domain outline
+	   unit). Riftname is an optional argument (Argus domain outline) describing rifts.
+
+	   Usage:
+	      md=triangle(md,domainname,resolution)
+	   or md=triangle(md,domainname, resolution, riftname)
+
+	   Examples:
+	      md=triangle(md,'DomainOutline.exp',1000);
+	      md=triangle(md,'DomainOutline.exp',1000,'Rifts.exp');
+	"""
 
 	#Figure out a characteristic area. Resolution is a node oriented concept (ex a 1000m  resolution node would 
 	#be made of 1000*1000 area squares). 
 
+	if len(args)==1:
+		resolution=args[0]
+		riftname=''
+	if len(args)==2:
+		riftname=args[0]
+		resolution=args[1]
+
 	#Check that mesh was not already run, and warn user: 
-	if md.mesh.numberofelements != 0.:
-		choice = input('This model already has a mesh. Are you sure you want to go ahead? (y/n)')
-		if choice != 'y':
+	if md.mesh.numberofelements:
+		choice = raw_input('This model already has a mesh. Are you sure you want to go ahead? (y/n)')
+		if not strcmp(choice,'y'):
 			print 'no meshing done ... exiting'
-			return []
-		
-	area = resolution**2.
+			return None
+
+	area = resolution**2
 
 	#Mesh using TriMesh
-	[md.mesh.elements,md.mesh.x,md.mesh.y,md.mesh.segments,md.mesh.segmentmarkers]=tm.TriMesh(domainname,riftname,area)
-
+	[md.mesh.elements,md.mesh.x,md.mesh.y,md.mesh.segments,md.mesh.segmentmarkers]=TriMesh(domainname,riftname,area)
 
 	#Fill in rest of fields:
-	md.mesh.numberofelements = len(md.mesh.elements)
-	md.mesh.numberofvertices = len(md.mesh.x)
-	md.mesh.z = zeros(md.mesh.numberofvertices)
-	md.mesh.vertexonboundary = zeros(md.mesh.numberofvertices)
+	md.mesh.numberofelements = numpy.size(md.mesh.elements,axis=0)
+	md.mesh.numberofvertices = numpy.size(md.mesh.x)
+	md.mesh.z = numpy.zeros(md.mesh.numberofvertices)
+	md.mesh.vertexonboundary = numpy.zeros(md.mesh.numberofvertices)
 	md.mesh.vertexonboundary[md.mesh.segments[:,0:2].astype(int)-1] = 1.
-	md.mesh.vertexonbed = ones(md.mesh.numberofvertices)
-	md.mesh.vertexonsurface = ones(md.mesh.numberofvertices)
-	md.mesh.elementonbed = ones(md.mesh.numberofelements)
-	md.mesh.elementonsurface = ones(md.mesh.numberofelements)
+	md.mesh.vertexonbed = numpy.ones(md.mesh.numberofvertices)
+	md.mesh.vertexonsurface = numpy.ones(md.mesh.numberofvertices)
+	md.mesh.elementonbed = numpy.ones(md.mesh.numberofelements)
+	md.mesh.elementonsurface = numpy.ones(md.mesh.numberofelements)
 
 	#Now, build the connectivity tables for this mesh.
-	[md.mesh.vertexconnectivity]= nc.NodeConnectivity(md.mesh.elements, md.mesh.numberofvertices)
-	[md.mesh.elementconnectivity] = ec.ElementConnectivity(md.mesh.elements, md.mesh.vertexconnectivity)
-	
+	[md.mesh.vertexconnectivity]= NodeConnectivity(md.mesh.elements, md.mesh.numberofvertices)
+	[md.mesh.elementconnectivity] = ElementConnectivity(md.mesh.elements, md.mesh.vertexconnectivity)
+
 	#type of model
 	md.mesh.dimension = 2.
+
 	return md
+