Index: /trunk/anuga_core/source/anuga/config.py
===================================================================
--- /trunk/anuga_core/source/anuga/config.py	(revision 8807)
+++ /trunk/anuga_core/source/anuga/config.py	(revision 8808)
@@ -201,5 +201,5 @@
 optimised_gradient_limiter = True # Use hardwired gradient limiter
 
-points_file_block_line_size = 5000 # Number of lines read in from a points file
+points_file_block_line_size = 5e7 # Number of lines read in from a points file
                                   # when blocking
 
Index: /trunk/anuga_core/source/anuga/fit_interpolate/fit.py
===================================================================
--- /trunk/anuga_core/source/anuga/fit_interpolate/fit.py	(revision 8807)
+++ /trunk/anuga_core/source/anuga/fit_interpolate/fit.py	(revision 8808)
@@ -263,66 +263,6 @@
             self.Atz = Atz
         else:
-            fitsmooth.combine_partial_AtA_Atz(self.AtA, AtA,\
+            fitsmooth.combine_partial_AtA_Atz(self.AtA, AtA, \
                     self.Atz, Atz, zdim, self.mesh.number_of_nodes)
-
-    def _build_matrix_AtA_Atz_file(self, filename, attribute_name=None, max_read_lines=500,\
-                                   verbose=False):
-        """Build:
-        AtA  m x m  interpolation matrix, and,
-        Atz  m x a  interpolation matrix where,
-        m is the number of basis functions phi_k (one per vertex)
-        a is the number of data attributes
-
-        This algorithm uses a quad tree data structure for fast binning of
-        data points.
-
-        If Ata is None, the matrices AtA and Atz are created.
-
-        This function can be called again and again, with sub-sets of
-        the point coordinates.  Call fit to get the results.
-
-        Preconditions
-        z and points are numeric
-        Point_coordindates and mesh vertices have the same origin.
-
-        The number of attributes of the data points does not change
-        """
-        # PADARN NOTE: THe following block of code is translated from
-        # things that were being done in the Geospatial_data object
-        # which is no longer required if data from file in C.
-        #---------------------------------------------------------
-        # Default attribute name here seems to only have possibility
-        # of being 'None'.
-        G_data = Geospatial_data(filename,
-                                         max_read_lines=1,
-                                         load_file_now=True,
-                                         verbose=verbose)
-
-
-        if attribute_name is None:
-            if G_data.default_attribute_name is not None:
-                attribute_name = G_data.default_attribute_name
-            else:
-                attribute_name = G_data.attributes.keys()[0]
-                # above line takes the first one from keys
-
-        if verbose is True:
-            log.critical('Using attribute %s' % attribute_name)
-            log.critical('Available attributes: %s' % (G_data.attributes.keys()))
-
-        msg = 'Attribute name %s does not exist in data set' % attribute_name
-        assert G_data.attributes.has_key(attribute_name), msg
-        #---------------------------------------------------------
-
-        # MAKE SURE READING ABSOLUTE POINT COORDINATES
-        [AtA, Atz, npts] = fitsmooth.build_matrix_AtA_Atz_fileread(self.root.root, \
-               self.mesh.number_of_nodes, \
-               self.mesh.triangles, \
-               filename, \
-               attribute_name, \
-               max_read_lines)
-        self.point_count = npts
-        self.AtA = AtA
-        self.Atz = Atz
 
     def fit(self, point_coordinates_or_filename=None, z=None,
@@ -330,6 +270,5 @@
             point_origin=None,
             attribute_name=None,
-            max_read_lines=500,
-            use_blocking_option2=True):
+            max_read_lines=1e7):
         """Fit a smooth surface to given 1d array of data points z.
 
@@ -350,57 +289,40 @@
                 use_blocking_option2 = False
 
-        # NOTE PADARN 12/12/12: Currently trying to get all blocking to be done
-        # in C, but blocking option 1, which does everything using the python
-        # interface can be used in the meantime (much slower).
-
-        #-----------------------BLOCKING OPTION 1----------------------------
-        if use_blocking_option2 is False:
+        # NOTE PADARN 29/03/13: File reading from C has been removed. Now 
+        # the input is either a set of points, or a filename which is then
+        # handled by the Geospatial_data object
+
+        if verbose:
+            print 'Fit.fit: Initializing'
+
+        # Use blocking to load in the point info
+        if isinstance(point_coordinates_or_filename, basestring):
+            msg = "Don't set a point origin when reading from a file"
+            assert point_origin is None, msg
+            filename = point_coordinates_or_filename
+
+            G_data = Geospatial_data(filename,
+                                     max_read_lines=max_read_lines,
+                                     load_file_now=False,
+                                     verbose=verbose)
+
+            for i, geo_block in enumerate(G_data):
+
+               # Build the array
+                points = geo_block.get_data_points(absolute=True)
+                z = geo_block.get_attributes(attribute_name=attribute_name)
+
+                self._build_matrix_AtA_Atz(points, z, attribute_name, verbose)
+
+            point_coordinates = None
+
             if verbose:
-                print 'Fit.fit: Initializing'
-
-            # Use blocking to load in the point info
-            if isinstance(point_coordinates_or_filename, basestring):
-                msg = "Don't set a point origin when reading from a file"
-                assert point_origin is None, msg
-                filename = point_coordinates_or_filename
-
-                G_data = Geospatial_data(filename,
-                                         max_read_lines=max_read_lines,
-                                         load_file_now=False,
-                                         verbose=verbose)
-
-                for i, geo_block in enumerate(G_data):
-
-                   # Build the array
-                    points = geo_block.get_data_points(absolute=True)
-                    z = geo_block.get_attributes(attribute_name=attribute_name)
-
-                    self._build_matrix_AtA_Atz(points, z, attribute_name, verbose)
-
-                point_coordinates = None
-
-                if verbose:
-                    print ''
-            else:
-                point_coordinates = point_coordinates_or_filename
-
-        #-----------------------BLOCKING OPTION 2----------------------------
+                print ''
         else:
-            if verbose:
-                print 'Fit.fit: Initializing'
-
-            if isinstance(point_coordinates_or_filename, basestring):
-                msg = "Don't set a point origin when reading from a file"
-                assert point_origin is None, msg
-                filename = point_coordinates_or_filename
-
-                self._build_matrix_AtA_Atz_file(filename, attribute_name=attribute_name,\
-                                                verbose=verbose)
-
-                point_coordinates = None
-            else:
-                point_coordinates = point_coordinates_or_filename
-        #--------------------------------------------------------------------
-
+            point_coordinates = point_coordinates_or_filename
+
+
+        # This condition either means a filename was read or the function
+        # recieved a None as input
         if point_coordinates is None:
             if verbose:
@@ -410,5 +332,5 @@
             assert self.Atz is not None
 
-            # FIXME (DSG) - do  a message
+
         else:
             point_coordinates = ensure_absolute(point_coordinates,
Index: /trunk/anuga_core/source/anuga/fit_interpolate/fitsmooth.c
===================================================================
--- /trunk/anuga_core/source/anuga/fit_interpolate/fitsmooth.c	(revision 8807)
+++ /trunk/anuga_core/source/anuga/fit_interpolate/fitsmooth.c	(revision 8808)
@@ -11,5 +11,4 @@
 #include "sparse_dok.h" /*in utilities */
 #include "quad_tree.h"
-#include <netcdf.h>
 #include "omp.h"
 
@@ -24,97 +23,4 @@
     #define PYVERSION273
 #endif
-
-//------------------------ NET CDF READING ----------------------------
-
-// Note Padarn 06/12/12: I am removing the use of these functions from
-// the workflow of the new code. I think they should be reintroduced, but
-// to make sure all the unit tests work properly more checking of filetypes
-// etc needs to be done.
-// Note Padarn 12/12/12: These have been moved back into the workflow, 
-// but still there is no way of handling .csv or .txt files efficiently.
-
-
-// Reads a block of a netcdf file into the arrays 'point_coordiantes' and 'point_values'. The
-// block is specified by a start point (np_start) and and an end point (np_end)
-// NOTE PADARN: FIX - only reads elevation at the moment
-int _read_net_cdf_points_block(char * filename, char * attribute_name, int np_start,int np_end,double *point_coordinates,double *point_values){
-
-    //point_coordinates and _values are intitialised outside
-    int np=np_end-np_start+1;
-    int status;                       /* error status */
-    int ncid;                          /* netCDF ID */
-    int rh_id;                         /* variable ID */
-    size_t start_points[] = {np_start, 0}; /* start at first value */
-    size_t count_points[] = {np, 2};
-    size_t start[] = {np_start};
-    size_t count[] = {np};
-    
-    if ((status = nc_open(filename, NC_NOWRITE, &ncid))) ERR(status);
-
-    if ((status = nc_inq_varid(ncid, "points", &rh_id))) ERR(status);
-
-    /* read coordinates from netCDF variable */
-    if ((status = nc_get_vara_double(ncid, rh_id, start_points, count_points, point_coordinates))) ERR(status);
-
-    if ((status = nc_inq_varid(ncid, attribute_name, &rh_id))) ERR(status);
-
-    /* read values from netCDF variable */
-    if ((status = nc_get_vara_double(ncid, rh_id, start, count, point_values))) ERR(status);
-
-    if ((status = nc_close(ncid))) ERR(status);
-
-    return 0;
-
-}
-
-// Takes a netcdf file and gets the number of points stored in it (also the number of values for
-// attributes specified)
-// NOTE PADARN: FIX - uses 'points' to get points, this may be a problem.
-int _read_net_cdf_entries(char * filename,int * x){
-
-    int status;                       /* error status */
-    int ncid;
-    int p_id;
-    int ndims;
-    size_t length = {0};
-
-    if ((status = nc_open(filename, NC_NOWRITE, &ncid))) ERR(status);
-
-    if ((status = nc_inq_varid(ncid, "points", &p_id))) ERR(status);
-    
-    if((status = nc_inq_varndims(ncid, p_id, &ndims))) ERR(status);
-    
-    int * ndimsize = malloc(sizeof(int)*ndims);
-
-    if((status = nc_inq_vardimid(ncid, p_id, ndimsize))) ERR(status);
-
-    if((status = nc_inq_dimlen(ncid, ndimsize[0], &length))) ERR(status);
-
-    *x = (int)length;
-
-    if ((status = nc_close(ncid))) ERR(status);
-
-    return 0;
-
-}
-
-// Takes a netcdf file and gets the x and y corners
-int _read_net_cdf_corners(char * filename,double *xllcorner, double *yllcorner){
-
-    int status;                       /* error status */
-    int ncid;
-
-    if ((status = nc_open(filename, NC_NOWRITE, &ncid))) ERR(status);
-
-    if ((status = nc_get_att_double(ncid, NC_GLOBAL, "xllcorner", xllcorner))) ERR(status);
-    if ((status = nc_get_att_double(ncid, NC_GLOBAL, "yllcorner", yllcorner))) ERR(status);
-
-    if ((status = nc_close(ncid))) ERR(status);
-
-    return 0;
-
-}
-
-//--------------------------------------------------------------------------
 
 //-------------------------- QUANTITY FITTING ------------------------------
@@ -272,96 +178,4 @@
 // Builds the AtA and Atz interpolation matrix
 // and residual. Uses a quad_tree for fast access to the triangles of the mesh.
-// This function takes a filename and an attribute name, in netcdf format .pts, 
-// and reads the attribute value and point location from this file.
-int _build_matrix_AtA_Atz_fileread(int N, long * triangles,
-                      char * filename,
-                      char * attribute_name,
-                      int blocksize,
-                      sparse_dok * AtA,
-                      double * Atz,quad_tree * quadtree)
-              {
-
-
-    int np;
-    _read_net_cdf_entries(filename,&np);
-    //n=10;
-
-    double xllcorner, yllcorner;
-    _read_net_cdf_corners(filename,&xllcorner,&yllcorner);
-
-    //printf(" xllcorner %g yllcorner %g\n",xllcorner,yllcorner);
-
-
-    double * point_coordinates = malloc(2*sizeof(double)*blocksize);
-    double * point_values = malloc(sizeof(double)*blocksize);
-    int k;
-    int left;
-
-    int i;
-    for(i=0;i<N;i++){
-        Atz[i]=0;
-    } 
-
-    edge_key_t key;
-    int w;
-    int start = 0;
-    left = np;
-
-    
-    while(left > 0){
-        if(left-blocksize<=0){
-            blocksize = left;
-        }
-
-        // Read data into arrays
-        _read_net_cdf_points_block(filename,attribute_name,start,
-        start+blocksize-1,point_coordinates,point_values);
-        
-        start = start + blocksize;
-
-            #pragma omp parallel for private(k,i,key,w)
-            for(k=0;k<blocksize;k++){
-
-
-                double x = point_coordinates[2*k]+xllcorner;
-                double y = point_coordinates[2*k+1]+yllcorner;
-                triangle * T = search(quadtree,x,y);
-
-                if(T!=NULL){
-                    double * sigma = calculate_sigma(T,x,y);
-                    int js[3];
-                    for(i=0;i<3;i++){
-                        js[i]=triangles[3*(T->index)+i];
-                    }
-                    
-                    #pragma omp critical
-                    { 
-                    for(i=0;i<3;i++){
-
-                        Atz[js[i]] += sigma[i]*point_values[k];
-
-                        for(w=0;w<3;w++){
-                            
-                            key.i=js[i];
-                            key.j=js[w];
-
-                            add_dok_entry(AtA,key,sigma[i]*sigma[w]);
-                        }                        
-                    }
-                    }
-                    free(sigma);
-                    sigma=NULL;
-                }
-
-            } 
-        
-        left=left-blocksize;
-    }
-    
-    return 0;
-}
-
-// Builds the AtA and Atz interpolation matrix
-// and residual. Uses a quad_tree for fast access to the triangles of the mesh.
 // This function takes a list of point coordinates, and associated point values
 // (for any number of attributes).
@@ -391,5 +205,5 @@
 
 
-    //#pragma omp parallel for private(k,i,key,w)
+    #pragma omp parallel for private(k,i,key,w)
     for(k=0;k<npts;k++){
 
@@ -651,80 +465,4 @@
     
     #endif
-
-}
-
-// Read a data file (netcdf format) to build
-// AtA and Atz. Returns a pointer to the sparse dok matrix AtA wrapped in a capsule
-// object and a python list for the array Atz.
-PyObject *build_matrix_AtA_Atz_fileread(PyObject *self, PyObject *args) {
-
-
-    PyArrayObject *triangles;
-    char * file_name;
-    char * attribute_name;
-    int block_size;
-    int N; // Number of triangles
-    int err;
-    PyObject *tree;
-
-    // Convert Python arguments to C
-    if (!PyArg_ParseTuple(args, "OiOssi",&tree, &N,
-                                            &triangles,
-                                            &file_name,
-                                            &attribute_name,
-                                            &block_size
-                                            )) {
-      PyErr_SetString(PyExc_RuntimeError,
-              "fitsmooth.c: could not parse input");
-      return NULL;
-    }
-
-    
-
-    CHECK_C_CONTIG(triangles);
-
-    #ifdef PYVERSION273
-    quad_tree * quadtree = (quad_tree*) PyCapsule_GetPointer(tree,"quad tree");
-    #else
-    quad_tree * quadtree = (quad_tree*) PyCObject_AsVoidPtr(tree);
-    #endif
-
-    sparse_dok * dok_AtA; // Should be an input argument?
-    dok_AtA = make_dok();
-    double * Atz = malloc(sizeof(double)*N);
-
-    err = _build_matrix_AtA_Atz_fileread(N,(long*) triangles->data,
-                      file_name,
-                      attribute_name,
-                      block_size,
-                      dok_AtA,
-                      Atz,
-                      quadtree);
-
-    if (err != 0) {
-      PyErr_SetString(PyExc_RuntimeError,
-              "Unknown Error");
-      return NULL;
-    }
-
-    #ifdef PYVERSION273
-    PyObject * AtA_cap =  PyCapsule_New((void*) dok_AtA,
-                  "sparse dok",
-                  &delete_dok_cap);
-    #else
-    PyObject * AtA_cap =  PyCObject_FromVoidPtr((void*) dok_AtA,
-                  &delete_dok_cobj);
-    #endif
-
-    PyObject *Atz_ret = c_double_array_to_list(Atz,N); 
-
-    int npts;
-    err = _read_net_cdf_entries(file_name, &npts);
-
-    PyObject *lst = PyList_New(3);
-    PyList_SET_ITEM(lst, 0, AtA_cap);
-    PyList_SET_ITEM(lst, 1, Atz_ret);
-    PyList_SET_ITEM(lst, 2, PyInt_FromLong((long)npts));
-    return lst;
 
 }
@@ -1121,5 +859,4 @@
     {"build_quad_tree",build_quad_tree, METH_VARARGS, "Print out"},
     {"build_smoothing_matrix",build_smoothing_matrix, METH_VARARGS, "Print out"},
-    {"build_matrix_AtA_Atz_fileread",build_matrix_AtA_Atz_fileread, METH_VARARGS, "Print out"},
     {"build_matrix_AtA_Atz_points",build_matrix_AtA_Atz_points, METH_VARARGS, "Print out"},
     {"combine_partial_AtA_Atz",combine_partial_AtA_Atz, METH_VARARGS, "Print out"},
Index: /trunk/anuga_core/source/anuga/utilities/compile.py
===================================================================
--- /trunk/anuga_core/source/anuga/utilities/compile.py	(revision 8807)
+++ /trunk/anuga_core/source/anuga/utilities/compile.py	(revision 8808)
@@ -15,20 +15,5 @@
 import os, string, sys
 
-
-NETCDF_LIB_DIR = os.getenv('NETCDF_LIB_DIR', '')
-NETCDF_INCLUDE_DIR = os.getenv('NETCDF_INCLUDE_DIR', '')
-
-#print 'NETCDF_LIB_DIR: ',NETCDF_LIB_DIR
-#print 'NETCDF_INCLUDE_DIR: ',NETCDF_INCLUDE_DIR
-
 I_dirs = ''
-if NETCDF_INCLUDE_DIR != '' :
-    I_dirs = ' -I"%s" ' % NETCDF_INCLUDE_DIR
-
-netcdf_lib_dirs = ''
-if NETCDF_LIB_DIR != '' :
-    netcdf_lib_dirs = ' -L"%s" ' % NETCDF_LIB_DIR
-
-#print 'netcdf_lib_dirs: ',netcdf_lib_dirs
 
 #NumPy ------------------------------------
@@ -333,15 +318,15 @@
   # Make shared library (*.so or *.dll)
   if FN=="fitsmooth.c":
-    fitlibs = libs + netcdf_lib_dirs
+    fitlibs = libs 
     if sys.platform == 'win32':	  
       if fitlibs is "":
-        s = '%s -%s %s ../utilities/quad_tree.o ../utilities/sparse_dok.o ../utilities/sparse_csr.o -o %s.%s -lm  -fopenmp netcdf.dll' %(loader, sharedflag, object_files, root1, libext)
+        s = '%s -%s %s ../utilities/quad_tree.o ../utilities/sparse_dok.o ../utilities/sparse_csr.o -o %s.%s -lm  -fopenmp' %(loader, sharedflag, object_files, root1, libext)
       else:
-        s = '%s -%s %s ../utilities/quad_tree.o ../utilities/sparse_dok.o ../utilities/sparse_csr.o -o %s.%s %s -lm  -fopenmp netcdf.dll' %(loader, sharedflag, object_files, root1, libext, fitlibs)
+        s = '%s -%s %s ../utilities/quad_tree.o ../utilities/sparse_dok.o ../utilities/sparse_csr.o -o %s.%s %s -lm  -fopenmp' %(loader, sharedflag, object_files, root1, libext, fitlibs)
     else:    
       if fitlibs is "":
-        s = '%s -%s %s ../utilities/quad_tree.o ../utilities/sparse_dok.o ../utilities/sparse_csr.o -o %s.%s -lm  -fopenmp -lnetcdf' %(loader, sharedflag, object_files, root1, libext)
+        s = '%s -%s %s ../utilities/quad_tree.o ../utilities/sparse_dok.o ../utilities/sparse_csr.o -o %s.%s -lm  -fopenmp' %(loader, sharedflag, object_files, root1, libext)
       else:
-        s = '%s -%s %s ../utilities/quad_tree.o ../utilities/sparse_dok.o ../utilities/sparse_csr.o -o %s.%s %s -lm  -fopenmp -lnetcdf' %(loader, sharedflag, object_files, root1, libext, fitlibs)
+        s = '%s -%s %s ../utilities/quad_tree.o ../utilities/sparse_dok.o ../utilities/sparse_csr.o -o %s.%s %s -lm  -fopenmp' %(loader, sharedflag, object_files, root1, libext, fitlibs)
   elif FN=="quad_tree_ext.c":
     if libs is "":
