Index: /anuga_core/source/anuga/abstract_2d_finite_volumes/general_mesh.py
===================================================================
--- /anuga_core/source/anuga/abstract_2d_finite_volumes/general_mesh.py	(revision 4778)
+++ /anuga_core/source/anuga/abstract_2d_finite_volumes/general_mesh.py	(revision 4779)
@@ -259,5 +259,28 @@
                 
         return V
-        
+    
+    def get_node(self, i,
+                 absolute=False):
+        """Return node coordinates.
+
+        The nodes are ordered in an Nx2 array where N is the number of nodes.
+        This is the same format they were provided in the constructor
+        i.e. without any duplication.
+
+        Boolean keyword argument absolute determines whether coordinates
+        are to be made absolute by taking georeference into account
+        Default is False as many parts of ANUGA expects relative coordinates.
+        (To see which, switch to default absolute=True and run tests).        
+        """
+
+       
+        V = self.nodes[i,:]
+        if absolute is True:
+            if not self.geo_reference.is_absolute():
+                V[0] += self.geo_reference.xllcorner
+                V[1] += self.geo_reference.yllcorner
+                
+        return V
+    
         
 
Index: /anuga_core/source/anuga/abstract_2d_finite_volumes/quantity.py
===================================================================
--- /anuga_core/source/anuga/abstract_2d_finite_volumes/quantity.py	(revision 4778)
+++ /anuga_core/source/anuga/abstract_2d_finite_volumes/quantity.py	(revision 4779)
@@ -726,6 +726,11 @@
 
         # Call fit_interpolate.fit function
-        args = (coordinates, triangles, points, values)
-        kwargs = {'data_origin': data_georef.get_origin(),
+        #args = (coordinates, triangles, points, values)
+        args = (points, )
+        kwargs = {'vertex_coordinates': coordinates,
+                  'triangles': triangles,
+                  'mesh': None,
+                  'point_attributes': values,
+                  'data_origin': data_georef.get_origin(),
                   'mesh_origin': mesh_georef.get_origin(),
                   'alpha': alpha,
@@ -776,8 +781,9 @@
             raise msg
 
-        coordinates = self.domain.get_nodes(absolute=True)
-        triangles = self.domain.triangles      #FIXME
-            
-        vertex_attributes = fit_to_mesh(coordinates, triangles, filename,
+        #coordinates = self.domain.get_nodes(absolute=True)
+        #triangles = self.domain.triangles      #FIXME
+        #vertex_attributes = fit_to_mesh(coordinates, triangles, filename, 
+        vertex_attributes = fit_to_mesh(filename,
+                                        mesh = self.domain, 
                                         alpha=alpha,
                                         attribute_name=attribute_name,
Index: /anuga_core/source/anuga/abstract_2d_finite_volumes/test_general_mesh.py
===================================================================
--- /anuga_core/source/anuga/abstract_2d_finite_volumes/test_general_mesh.py	(revision 4778)
+++ /anuga_core/source/anuga/abstract_2d_finite_volumes/test_general_mesh.py	(revision 4779)
@@ -9,4 +9,5 @@
 from Numeric import allclose, array, ones, Float
 from general_mesh import General_mesh
+from anuga.coordinate_transforms.geo_reference import Geo_reference
 
 
@@ -216,10 +217,44 @@
         assert unique_vertices == [0,2,4,5,6,7]
 
-
-        
+    def test_get_node(self):
+        """test_get_triangles_and_vertices_per_node -
+
+        Test that tuples of triangle, vertex can be extracted
+        from inverted triangles structure
+
+        """
+        
+        x0 = 314036.58727982
+        y0 = 6224951.2960092
+        geo = Geo_reference(56, x0, y0)
+        
+        a = [0.0, 0.0]
+        b = [0.0, 2.0]
+        c = [2.0, 0.0]
+        d = [0.0, 4.0]
+        e = [2.0, 2.0]
+        f = [4.0, 0.0]
+
+        nodes = array([a, b, c, d, e, f])
+
+        nodes_absolute = geo.get_absolute(nodes)
+        
+        #bac, bce, ecf, dbe, daf, dae
+        triangles = array([[1,0,2], [1,2,4], [4,2,5], [3,1,4]])
+
+        domain = General_mesh(nodes, triangles,
+                       geo_reference = geo)
+        node = domain.get_node(2)        
+        self.assertEqual(c, node)
+        
+        node = domain.get_node(2, absolute=True)     
+        self.assertEqual(nodes_absolute[2], node)
+        
+
 
 #-------------------------------------------------------------
 if __name__ == "__main__":
-    suite = unittest.makeSuite(Test_General_Mesh,'test')    
+    suite = unittest.makeSuite(Test_General_Mesh,'test') 
+    #suite = unittest.makeSuite(Test_General_Mesh,'test_get_node')    
     runner = unittest.TextTestRunner()
     runner.run(suite)
Index: /anuga_core/source/anuga/abstract_2d_finite_volumes/test_quantity.py
===================================================================
--- /anuga_core/source/anuga/abstract_2d_finite_volumes/test_quantity.py	(revision 4778)
+++ /anuga_core/source/anuga/abstract_2d_finite_volumes/test_quantity.py	(revision 4779)
@@ -541,8 +541,8 @@
 
         #Now try by setting the same values directly
-        vertex_attributes = fit_to_mesh(quantity.domain.get_nodes(),
+        vertex_attributes = fit_to_mesh(data_points,
+                                        quantity.domain.get_nodes(),
                                         quantity.domain.triangles, #FIXME
-                                        data_points,
-                                        z,
+                                        point_attributes=z,
                                         alpha = 0,
                                         verbose=False)
@@ -584,5 +584,6 @@
 
         #Reference
-        ref = fit_to_mesh(vertex_coordinates, triangles, data_points, z,
+        ref = fit_to_mesh(data_points, vertex_coordinates, triangles,
+                          point_attributes=z,
                           data_origin = data_georef.get_origin(),
                           mesh_origin = mesh_georef.get_origin(),
@@ -1048,4 +1049,6 @@
         x0 = 314036.58727982
         y0 = 6224951.2960092
+        #x0 = 0.0
+        #y0 = 0.0
 
         a = [0.0, 0.0]
Index: /anuga_core/source/anuga/fit_interpolate/fit.py
===================================================================
--- /anuga_core/source/anuga/fit_interpolate/fit.py	(revision 4778)
+++ /anuga_core/source/anuga/fit_interpolate/fit.py	(revision 4779)
@@ -27,6 +27,7 @@
 import types
 
-from Numeric import zeros, Float, ArrayType,take 
-
+from Numeric import zeros, Float, ArrayType,take, Int
+
+from anuga.abstract_2d_finite_volumes.neighbour_mesh import Mesh
 from anuga.caching import cache            
 from anuga.geospatial_data.geospatial_data import Geospatial_data, \
@@ -49,6 +50,7 @@
     
     def __init__(self,
-                 vertex_coordinates,
-                 triangles,
+                 vertex_coordinates=None,
+                 triangles=None,
+                 mesh=None,
                  mesh_origin=None,
                  alpha = None,
@@ -98,4 +100,5 @@
                  vertex_coordinates,
                  triangles,
+                 mesh,
                  mesh_origin,
                  verbose,
@@ -434,7 +437,8 @@
 ############################################################################
 
-def fit_to_mesh(vertex_coordinates,
-                triangles,
-                point_coordinates, # this can also be a points file name
+def fit_to_mesh(point_coordinates, # this can also be a points file name
+                vertex_coordinates=None,
+                triangles=None,
+                mesh=None,
                 point_attributes=None,
                 alpha=DEFAULT_ALPHA,
@@ -450,6 +454,9 @@
     """
     
-    args = (vertex_coordinates, triangles, point_coordinates, )
-    kwargs = {'point_attributes': point_attributes,
+    args = (point_coordinates, )
+    kwargs = {'vertex_coordinates': vertex_coordinates,
+              'triangles': triangles,
+              'mesh': mesh,
+              'point_attributes': point_attributes,
               'alpha': alpha,
               'verbose': verbose,
@@ -480,7 +487,8 @@
                      args, kwargs)
 
-def _fit_to_mesh(vertex_coordinates,
-                 triangles,
-                 point_coordinates, # this can also be a points file name
+def _fit_to_mesh(point_coordinates, # this can also be a points file name
+                 vertex_coordinates=None,
+                 triangles=None,
+                 mesh=None,
                  point_attributes=None,
                  alpha=DEFAULT_ALPHA,
@@ -513,8 +521,11 @@
           alpha: Smoothing parameter.
 
-          acceptable overshoot: controls the allowed factor by which fitted values
+          acceptable overshoot: NOT IMPLEMENTED
+          controls the allowed factor by which
+          fitted values
           may exceed the value of input data. The lower limit is defined
           as min(z) - acceptable_overshoot*delta z and upper limit
           as max(z) + acceptable_overshoot*delta z
+          
 
           mesh_origin: A geo_reference object or 3-tuples consisting of
@@ -531,8 +542,20 @@
     # Duncan and Ole think that this isn't worth caching.
     # Caching happens at the higher level anyway.
-    interp = Fit(vertex_coordinates,
-                 triangles,
+    
+    if mesh is None:
+        # Fixme (DSG) Throw errors if triangles or vertex_coordinates
+        # are None
+            
+        #Convert input to Numeric arrays
+        triangles = ensure_numeric(triangles, Int)
+        vertex_coordinates = ensure_absolute(vertex_coordinates,
+                                             geo_reference = mesh_origin)
+
+        if verbose: print 'FitInterpolate: Building mesh'        
+        mesh = Mesh(vertex_coordinates, triangles)
+        mesh.check_integrity()
+    
+    interp = Fit(mesh=mesh,
                  verbose=verbose,
-                 mesh_origin=mesh_origin,
                  alpha=alpha)
 
@@ -624,7 +647,8 @@
     if verbose: print "points file loaded"
     if verbose: print "fitting to mesh"
-    f = fit_to_mesh(vertex_coordinates,
+    f = fit_to_mesh(point_coordinates,
+                    vertex_coordinates,
                     triangles,
-                    point_coordinates,
+                    None,
                     point_attributes,
                     alpha = alpha,
Index: /anuga_core/source/anuga/fit_interpolate/general_fit_interpolate.py
===================================================================
--- /anuga_core/source/anuga/fit_interpolate/general_fit_interpolate.py	(revision 4778)
+++ /anuga_core/source/anuga/fit_interpolate/general_fit_interpolate.py	(revision 4779)
@@ -44,6 +44,7 @@
     
     def __init__(self,
-                 vertex_coordinates,
-                 triangles,
+                 vertex_coordinates=None,
+                 triangles=None,
+                 mesh=None,
                  mesh_origin=None,
                  verbose=False,
@@ -54,4 +55,7 @@
         function values at vertices to function values at data points
 
+        Pass in a mesh instance or vertex_coordinates and triangles
+        and optionally mesh_origin
+        
         Inputs:
 
@@ -64,4 +68,6 @@
           triangles: List of 3-tuples (or a Numeric array) of
               integers representing indices of all vertices in the mesh.
+
+        mesh: A mesh instance describing the mesh.
 
           mesh_origin: A geo_reference object or 3-tuples consisting of
@@ -79,15 +85,19 @@
         if max_vertices_per_cell == None:
             max_vertices_per_cell = MAX_VERTICES_PER_CELL
-        
-        #Convert input to Numeric arrays
-        triangles = ensure_numeric(triangles, Int)
-        vertex_coordinates = ensure_absolute(vertex_coordinates,
-                                             geo_reference = mesh_origin)
 
-        #Don't pass geo_reference to mesh.  It doesn't work. (Still??)
-        
-        if verbose: print 'FitInterpolate: Building mesh'        
-        self.mesh = Mesh(vertex_coordinates, triangles)
-        self.mesh.check_integrity()
+        if mesh is None:
+            # Fixme (DSG) Throw errors if triangles or vertex_coordinates
+            # are None
+            
+            #Convert input to Numeric arrays
+            triangles = ensure_numeric(triangles, Int)
+            vertex_coordinates = ensure_absolute(vertex_coordinates,
+                                                 geo_reference = mesh_origin)
+
+            if verbose: print 'FitInterpolate: Building mesh'        
+            self.mesh = Mesh(vertex_coordinates, triangles)
+            self.mesh.check_integrity()
+        else:
+            self.mesh = mesh
         
         if verbose: print 'FitInterpolate: Building quad tree'
Index: /anuga_core/source/anuga/fit_interpolate/interpolate.py
===================================================================
--- /anuga_core/source/anuga/fit_interpolate/interpolate.py	(revision 4778)
+++ /anuga_core/source/anuga/fit_interpolate/interpolate.py	(revision 4779)
@@ -84,9 +84,9 @@
         
         FitInterpolate.__init__(self,
-                                vertex_coordinates,
-                                triangles,
-                                mesh_origin,
-                                verbose,
-                                max_vertices_per_cell)
+                                vertex_coordinates=vertex_coordinates,
+                                triangles=triangles,
+                                mesh_origin=mesh_origin,
+                                verbose=verbose,
+                                max_vertices_per_cell=max_vertices_per_cell)
 
     # FIXME: What is a good start_blocking_len value?
Index: /anuga_core/source/anuga/fit_interpolate/test_fit.py
===================================================================
--- /anuga_core/source/anuga/fit_interpolate/test_fit.py	(revision 4778)
+++ /anuga_core/source/anuga/fit_interpolate/test_fit.py	(revision 4779)
@@ -270,5 +270,5 @@
         points = [a, b, c]
         elements = [[0,2,1]]
-        f = fit_to_mesh(points, elements, txt_file,
+        f = fit_to_mesh(txt_file, points, elements, 
                         alpha=0.0, max_read_lines=2)
         answer = linear_function(points)
@@ -312,5 +312,5 @@
         points = geo.get_data_points(absolute=True)
         atts = geo.get_attributes()
-        f = fit_to_mesh(vertices, triangles,points,atts,
+        f = fit_to_mesh(points,atts, vertices, triangles,
                                 alpha=0.0, max_read_lines=2,
                         use_cache=True, verbose=True)
@@ -352,5 +352,45 @@
         fileName_pts = tempfile.mktemp(".pts")
         geo.export_points_file(fileName_pts)
-        f = fit_to_mesh(vertices, triangles,fileName_pts,
+        f = fit_to_mesh(fileName_pts, vertices, triangles,
+                                alpha=0.0, max_read_lines=2)
+        answer = linear_function(vertices)
+        #print "f\n",f
+        #print "answer\n",answer
+        assert allclose(f, answer)
+        os.remove(fileName)
+        os.remove(fileName_pts)
+        
+    def test_fit_to_mesh_pts_passing_mesh_in(self):
+        a = [-1.0, 0.0]
+        b = [3.0, 4.0]
+        c = [4.0,1.0]
+        d = [-3.0, 2.0] #3
+        e = [-1.0,-2.0]
+        f = [1.0, -2.0] #5
+
+        vertices = [a, b, c, d,e,f]
+        triangles = [[0,1,3], [1,0,2], [0,4,5], [0,5,2]] #abd bac aef afc
+
+
+        fileName = tempfile.mktemp(".txt")
+        file = open(fileName,"w")
+        file.write(" x, y, elevation \n\
+-2.0, 2.0, 0.\n\
+-1.0, 1.0, 0.\n\
+0.0, 2.0 , 2.\n\
+1.0, 1.0 , 2.\n\
+ 2.0,  1.0 ,3. \n\
+ 0.0,  0.0 , 0.\n\
+ 1.0,  0.0 , 1.\n\
+ 0.0,  -1.0, -1.\n\
+ -0.2, -0.5, -0.7\n\
+ -0.9, -1.5, -2.4\n\
+ 0.5,  -1.9, -1.4\n\
+ 3.0,  1.0 , 4.\n")
+        file.close()
+        geo = Geospatial_data(fileName)
+        fileName_pts = tempfile.mktemp(".pts")
+        geo.export_points_file(fileName_pts)
+        f = fit_to_mesh(fileName_pts, vertices, triangles,
                                 alpha=0.0, max_read_lines=2)
         answer = linear_function(vertices)
@@ -391,5 +431,5 @@
         file.close()
         
-        f = fit_to_mesh(vertices, triangles,fileName,
+        f = fit_to_mesh(fileName, vertices, triangles,
                                 alpha=0.0, max_read_lines=2)
                         #use_cache=True, verbose=True)
@@ -432,5 +472,5 @@
         file.close()
         
-        f = fit_to_mesh(vertices, triangles,fileName,
+        f = fit_to_mesh(fileName, vertices, triangles,
                         alpha=0.0, 
                         attribute_name='elevation', max_read_lines=2)
@@ -686,5 +726,7 @@
         z = [z1, z2, z3]
 
-        f = fit_to_mesh(points, triangles, data_coords, z, alpha=0.0)
+        f = fit_to_mesh(data_coords, vertex_coordinates=points,
+                        triangles=triangles,
+                        point_attributes=z, alpha=0.0)
         answer = [[0, 0], [5., 10.], [5., 10.]]
 
@@ -720,5 +762,7 @@
         
         #Fit
-        zz = fit_to_mesh(vertex_coordinates, triangles, data_points, z,
+        zz = fit_to_mesh(data_points, vertex_coordinates=vertex_coordinates,
+                         triangles=triangles,
+                         point_attributes=z,
                          data_origin = data_geo.get_origin(),
                          mesh_origin = mesh_geo.get_origin(),
@@ -1052,5 +1096,5 @@
 if __name__ == "__main__":
     suite = unittest.makeSuite(Test_Fit,'test')
-    #suite = unittest.makeSuite(Test_Fit,'cache_test_fit_to_mesh_pts')
+    #suite = unittest.makeSuite(Test_Fit,'test_fit_to_mesh_pts_passing_mesh_in')
     #suite = unittest.makeSuite(Test_Fit,'')
     runner = unittest.TextTestRunner(verbosity=1)
Index: /anuga_core/source/anuga/fit_interpolate/test_interpolate.py
===================================================================
--- /anuga_core/source/anuga/fit_interpolate/test_interpolate.py	(revision 4778)
+++ /anuga_core/source/anuga/fit_interpolate/test_interpolate.py	(revision 4779)
@@ -14,5 +14,6 @@
 
 from Scientific.IO.NetCDF import NetCDFFile
-from Numeric import allclose, array, transpose, zeros, Float, sometrue, alltrue, take, where
+from Numeric import allclose, array, transpose, zeros, Float, sometrue, \
+     alltrue, take, where
 
 
Index: /anuga_core/source/anuga/utilities/quad.py
===================================================================
--- /anuga_core/source/anuga/utilities/quad.py	(revision 4778)
+++ /anuga_core/source/anuga/utilities/quad.py	(revision 4779)
@@ -142,5 +142,5 @@
 	if len(args) == 2:
 	    point_id = int(args[1])
-            x, y = self.mesh.get_nodes()[point_id]
+            x, y = self.mesh.get_node(point_id, absolute=True)
 
             #print point_id, x, y
@@ -183,7 +183,10 @@
                     self.store(point)
                 else:
-                    x = self.mesh.coordinates[point][0]
-                    y = self.mesh.coordinates[point][1]
-                    print "(" + str(x) + "," + str(y) + ")"
+                    # Have to take into account of georef.
+                    #x = self.mesh.coordinates[point][0]
+                    #y = self.mesh.coordinates[point][1]
+                    node = self.mesh.get_node(point, absolute=True)
+                    print "node", node
+                    print "(" + str(node[0]) + "," + str(node[1]) + ")"
                     raise 'point not in region: %s' %str(point)
 		
@@ -229,5 +232,6 @@
                     if not triangles.has_key(k):
                         # print 'k',k
-                        tri = self.mesh.get_vertex_coordinates(k)
+                        tri = self.mesh.get_vertex_coordinates(k,
+                                                               absolute=True)
                         n0 = self.mesh.get_normal(k, 0)
                         n1 = self.mesh.get_normal(k, 1)
@@ -441,5 +445,6 @@
     #print mesh.coordinates
 
-    nodes = mesh.get_nodes()
+    
+    nodes = mesh.get_nodes(absolute=True)
     xmin = min(nodes[:,0])
     xmax = max(nodes[:,0])
@@ -447,4 +452,6 @@
     ymax = max(nodes[:,1])
 
+    # Don't know why this didn't work
+    #xmin, xmax, ymin, ymax = mesh.get_extent(absolute=True)
     
     # Ensure boundary points are fully contained in region
