Index: sansmodels/src/sans/models/c_extensions/CGaussian.c
===================================================================
--- sansmodels/src/sans/models/c_extensions/CGaussian.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/CGaussian.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,290 @@
+/** CGaussian
+ *
+ * C extension 
+ *
+ * WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
+ *          DO NOT MODIFY THIS FILE, MODIFY gaussian.h
+ *          AND RE-RUN THE GENERATOR SCRIPT
+ *
+ * @author   M.Doucet / UTK
+ */
+ 
+#include <Python.h>
+#include "structmember.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <time.h>
+
+#include "gaussian.h"
+
+/// Error object for raised exceptions
+static PyObject * CGaussianError = NULL;
+
+
+// Class definition
+typedef struct {
+    PyObject_HEAD
+    /// Parameters
+    PyObject * params;
+    /// Log for unit testing
+    PyObject * log;
+    /// Model parameters
+	GaussianParameters model_pars;
+} CGaussian;
+
+
+static void
+CGaussian_dealloc(CGaussian* self)
+{
+    self->ob_type->tp_free((PyObject*)self);
+    
+
+}
+
+static PyObject *
+CGaussian_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+    CGaussian *self;
+    
+    self = (CGaussian *)type->tp_alloc(type, 0);
+   
+    return (PyObject *)self;
+}
+
+static int
+CGaussian_init(CGaussian *self, PyObject *args, PyObject *kwds)
+{
+    if (self != NULL) {
+    	
+    	// Create parameters
+        self->params = PyDict_New();
+        
+        // Initialize parameter dictionary
+        PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"sigma",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"center",Py_BuildValue("d",0.000000));
+
+         
+        // Create empty log
+        self->log = PyDict_New();
+        
+        
+
+    }
+    return 0;
+}
+
+static PyMemberDef CGaussian_members[] = {
+    {"params", T_OBJECT, offsetof(CGaussian, params), 0,
+     "Parameters"},
+    {"log", T_OBJECT, offsetof(CGaussian, log), 0,
+     "Log"},
+    {NULL}  /* Sentinel */
+};
+
+/** Read double from PyObject
+    @param p PyObject
+    @return double
+*/
+double CGaussian_readDouble(PyObject *p) {
+    if (PyFloat_Check(p)==1) {
+        return (double)(((PyFloatObject *)(p))->ob_fval);
+    } else if (PyInt_Check(p)==1) {
+        return (double)(((PyIntObject *)(p))->ob_ival);
+    } else if (PyLong_Check(p)==1) {
+        return (double)PyLong_AsLong(p);
+    } else {
+        return 0.0;
+    }
+}
+
+
+/**
+ * Function to call to evaluate model
+ * @param args: input q or [q,phi]
+ * @return: function value
+ */
+static PyObject * run(CGaussian *self, PyObject *args) {
+	double q_value, phi_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	// Reader parameter dictionary
+    self->model_pars.scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model_pars.sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") );
+    self->model_pars.center = PyFloat_AsDouble( PyDict_GetItemString(self->params, "center") );
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CGaussianError, 
+	    	"CGaussian.run expects a q value.");
+		return NULL;
+	}
+	  
+	// Check params
+	if( PyList_Check(pars)==1) {
+		
+		// Length of list should be 2 for I(q,phi)
+	    npars = PyList_GET_SIZE(pars); 
+	    if(npars!=2) {
+	    	PyErr_SetString(CGaussianError, 
+	    		"CGaussian.run expects a double or a list of dimension 2.");
+	    	return NULL;
+	    }
+	    // We have a vector q, get the q and phi values at which
+	    // to evaluate I(q,phi)
+	    q_value = CGaussian_readDouble(PyList_GET_ITEM(pars,0));
+	    phi_value = CGaussian_readDouble(PyList_GET_ITEM(pars,1));
+	    // Skip zero
+	    if (q_value==0) {
+	    	return Py_BuildValue("d",0.0);
+	    }
+		return Py_BuildValue("d",gaussian_analytical_2D(&(self->model_pars),q_value,phi_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		q_value = CGaussian_readDouble(pars);		
+		
+		return Py_BuildValue("d",gaussian_analytical_1D(&(self->model_pars),q_value));
+	}	
+}
+
+/**
+ * Function to call to evaluate model in cartesian coordinates
+ * @param args: input q or [qx, qy]]
+ * @return: function value
+ */
+static PyObject * runXY(CGaussian *self, PyObject *args) {
+	double qx_value, qy_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	// Reader parameter dictionary
+    self->model_pars.scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model_pars.sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") );
+    self->model_pars.center = PyFloat_AsDouble( PyDict_GetItemString(self->params, "center") );
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CGaussianError, 
+	    	"CGaussian.run expects a q value.");
+		return NULL;
+	}
+	  
+	// Check params
+	if( PyList_Check(pars)==1) {
+		
+		// Length of list should be 2 for I(qx, qy))
+	    npars = PyList_GET_SIZE(pars); 
+	    if(npars!=2) {
+	    	PyErr_SetString(CGaussianError, 
+	    		"CGaussian.run expects a double or a list of dimension 2.");
+	    	return NULL;
+	    }
+	    // We have a vector q, get the qx and qy values at which
+	    // to evaluate I(qx,qy)
+	    qx_value = CGaussian_readDouble(PyList_GET_ITEM(pars,0));
+	    qy_value = CGaussian_readDouble(PyList_GET_ITEM(pars,1));
+		return Py_BuildValue("d",gaussian_analytical_2DXY(&(self->model_pars),qx_value,qy_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		qx_value = CGaussian_readDouble(pars);		
+		
+		return Py_BuildValue("d",gaussian_analytical_1D(&(self->model_pars),qx_value));
+	}	
+}
+
+static PyObject * reset(CGaussian *self, PyObject *args) {
+    
+
+    return Py_BuildValue("d",0.0);
+}
+
+
+static PyMethodDef CGaussian_methods[] = {
+    {"run",      (PyCFunction)run     , METH_VARARGS,
+      "Evaluate the model at a given Q or Q, phi"},
+    {"runXY",      (PyCFunction)runXY     , METH_VARARGS,
+      "Evaluate the model at a given Q or Qx, Qy"},
+    {"reset",    (PyCFunction)reset   , METH_VARARGS,
+      "Reset pair correlation"},
+    //{"numerical_1D",      (PyCFunction)numerical_1D     , METH_VARARGS,
+    //  "Evaluate the 1D model at a given Q"},
+   {NULL}
+};
+
+static PyTypeObject CGaussianType = {
+    PyObject_HEAD_INIT(NULL)
+    0,                         /*ob_size*/
+    "CGaussian",             /*tp_name*/
+    sizeof(CGaussian),             /*tp_basicsize*/
+    0,                         /*tp_itemsize*/
+    (destructor)CGaussian_dealloc, /*tp_dealloc*/
+    0,                         /*tp_print*/
+    0,                         /*tp_getattr*/
+    0,                         /*tp_setattr*/
+    0,                         /*tp_compare*/
+    0,                         /*tp_repr*/
+    0,                         /*tp_as_number*/
+    0,                         /*tp_as_sequence*/
+    0,                         /*tp_as_mapping*/
+    0,                         /*tp_hash */
+    0,                         /*tp_call*/
+    0,                         /*tp_str*/
+    0,                         /*tp_getattro*/
+    0,                         /*tp_setattro*/
+    0,                         /*tp_as_buffer*/
+    Py_TPFLAGS_DEFAULT | Py_TPFLAGS_BASETYPE, /*tp_flags*/
+    "CGaussian objects",           /* tp_doc */
+    0,		               /* tp_traverse */
+    0,		               /* tp_clear */
+    0,		               /* tp_richcompare */
+    0,		               /* tp_weaklistoffset */
+    0,		               /* tp_iter */
+    0,		               /* tp_iternext */
+    CGaussian_methods,             /* tp_methods */
+    CGaussian_members,             /* tp_members */
+    0,                         /* tp_getset */
+    0,                         /* tp_base */
+    0,                         /* tp_dict */
+    0,                         /* tp_descr_get */
+    0,                         /* tp_descr_set */
+    0,                         /* tp_dictoffset */
+    (initproc)CGaussian_init,      /* tp_init */
+    0,                         /* tp_alloc */
+    CGaussian_new,                 /* tp_new */
+};
+
+
+static PyMethodDef module_methods[] = {
+    {NULL} 
+};
+
+/**
+ * Function used to add the model class to a module
+ * @param module: module to add the class to
+ */ 
+void addCGaussian(PyObject *module) {
+	PyObject *d;
+	
+    if (PyType_Ready(&CGaussianType) < 0)
+        return;
+
+    Py_INCREF(&CGaussianType);
+    PyModule_AddObject(module, "CGaussian", (PyObject *)&CGaussianType);
+    
+    d = PyModule_GetDict(module);
+    CGaussianError = PyErr_NewException("CGaussian.error", NULL, NULL);
+    PyDict_SetItemString(d, "CGaussianError", CGaussianError);
+}
+
Index: sansmodels/src/sans/models/c_extensions/CLogNormal.c
===================================================================
--- sansmodels/src/sans/models/c_extensions/CLogNormal.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/CLogNormal.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,290 @@
+/** CLogNormal
+ *
+ * C extension 
+ *
+ * WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
+ *          DO NOT MODIFY THIS FILE, MODIFY logNormal.h
+ *          AND RE-RUN THE GENERATOR SCRIPT
+ *
+ * @author   M.Doucet / UTK
+ */
+ 
+#include <Python.h>
+#include "structmember.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <time.h>
+
+#include "logNormal.h"
+
+/// Error object for raised exceptions
+static PyObject * CLogNormalError = NULL;
+
+
+// Class definition
+typedef struct {
+    PyObject_HEAD
+    /// Parameters
+    PyObject * params;
+    /// Log for unit testing
+    PyObject * log;
+    /// Model parameters
+	LogNormalParameters model_pars;
+} CLogNormal;
+
+
+static void
+CLogNormal_dealloc(CLogNormal* self)
+{
+    self->ob_type->tp_free((PyObject*)self);
+    
+
+}
+
+static PyObject *
+CLogNormal_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+    CLogNormal *self;
+    
+    self = (CLogNormal *)type->tp_alloc(type, 0);
+   
+    return (PyObject *)self;
+}
+
+static int
+CLogNormal_init(CLogNormal *self, PyObject *args, PyObject *kwds)
+{
+    if (self != NULL) {
+    	
+    	// Create parameters
+        self->params = PyDict_New();
+        
+        // Initialize parameter dictionary
+        PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"sigma",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"center",Py_BuildValue("d",0.000000));
+
+         
+        // Create empty log
+        self->log = PyDict_New();
+        
+        
+
+    }
+    return 0;
+}
+
+static PyMemberDef CLogNormal_members[] = {
+    {"params", T_OBJECT, offsetof(CLogNormal, params), 0,
+     "Parameters"},
+    {"log", T_OBJECT, offsetof(CLogNormal, log), 0,
+     "Log"},
+    {NULL}  /* Sentinel */
+};
+
+/** Read double from PyObject
+    @param p PyObject
+    @return double
+*/
+double CLogNormal_readDouble(PyObject *p) {
+    if (PyFloat_Check(p)==1) {
+        return (double)(((PyFloatObject *)(p))->ob_fval);
+    } else if (PyInt_Check(p)==1) {
+        return (double)(((PyIntObject *)(p))->ob_ival);
+    } else if (PyLong_Check(p)==1) {
+        return (double)PyLong_AsLong(p);
+    } else {
+        return 0.0;
+    }
+}
+
+
+/**
+ * Function to call to evaluate model
+ * @param args: input q or [q,phi]
+ * @return: function value
+ */
+static PyObject * run(CLogNormal *self, PyObject *args) {
+	double q_value, phi_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	// Reader parameter dictionary
+    self->model_pars.scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model_pars.sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") );
+    self->model_pars.center = PyFloat_AsDouble( PyDict_GetItemString(self->params, "center") );
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CLogNormalError, 
+	    	"CLogNormal.run expects a q value.");
+		return NULL;
+	}
+	  
+	// Check params
+	if( PyList_Check(pars)==1) {
+		
+		// Length of list should be 2 for I(q,phi)
+	    npars = PyList_GET_SIZE(pars); 
+	    if(npars!=2) {
+	    	PyErr_SetString(CLogNormalError, 
+	    		"CLogNormal.run expects a double or a list of dimension 2.");
+	    	return NULL;
+	    }
+	    // We have a vector q, get the q and phi values at which
+	    // to evaluate I(q,phi)
+	    q_value = CLogNormal_readDouble(PyList_GET_ITEM(pars,0));
+	    phi_value = CLogNormal_readDouble(PyList_GET_ITEM(pars,1));
+	    // Skip zero
+	    if (q_value==0) {
+	    	return Py_BuildValue("d",0.0);
+	    }
+		return Py_BuildValue("d",logNormal_analytical_2D(&(self->model_pars),q_value,phi_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		q_value = CLogNormal_readDouble(pars);		
+		
+		return Py_BuildValue("d",logNormal_analytical_1D(&(self->model_pars),q_value));
+	}	
+}
+
+/**
+ * Function to call to evaluate model in cartesian coordinates
+ * @param args: input q or [qx, qy]]
+ * @return: function value
+ */
+static PyObject * runXY(CLogNormal *self, PyObject *args) {
+	double qx_value, qy_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	// Reader parameter dictionary
+    self->model_pars.scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model_pars.sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") );
+    self->model_pars.center = PyFloat_AsDouble( PyDict_GetItemString(self->params, "center") );
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CLogNormalError, 
+	    	"CLogNormal.run expects a q value.");
+		return NULL;
+	}
+	  
+	// Check params
+	if( PyList_Check(pars)==1) {
+		
+		// Length of list should be 2 for I(qx, qy))
+	    npars = PyList_GET_SIZE(pars); 
+	    if(npars!=2) {
+	    	PyErr_SetString(CLogNormalError, 
+	    		"CLogNormal.run expects a double or a list of dimension 2.");
+	    	return NULL;
+	    }
+	    // We have a vector q, get the qx and qy values at which
+	    // to evaluate I(qx,qy)
+	    qx_value = CLogNormal_readDouble(PyList_GET_ITEM(pars,0));
+	    qy_value = CLogNormal_readDouble(PyList_GET_ITEM(pars,1));
+		return Py_BuildValue("d",logNormal_analytical_2DXY(&(self->model_pars),qx_value,qy_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		qx_value = CLogNormal_readDouble(pars);		
+		
+		return Py_BuildValue("d",logNormal_analytical_1D(&(self->model_pars),qx_value));
+	}	
+}
+
+static PyObject * reset(CLogNormal *self, PyObject *args) {
+    
+
+    return Py_BuildValue("d",0.0);
+}
+
+
+static PyMethodDef CLogNormal_methods[] = {
+    {"run",      (PyCFunction)run     , METH_VARARGS,
+      "Evaluate the model at a given Q or Q, phi"},
+    {"runXY",      (PyCFunction)runXY     , METH_VARARGS,
+      "Evaluate the model at a given Q or Qx, Qy"},
+    {"reset",    (PyCFunction)reset   , METH_VARARGS,
+      "Reset pair correlation"},
+    //{"numerical_1D",      (PyCFunction)numerical_1D     , METH_VARARGS,
+    //  "Evaluate the 1D model at a given Q"},
+   {NULL}
+};
+
+static PyTypeObject CLogNormalType = {
+    PyObject_HEAD_INIT(NULL)
+    0,                         /*ob_size*/
+    "CLogNormal",             /*tp_name*/
+    sizeof(CLogNormal),             /*tp_basicsize*/
+    0,                         /*tp_itemsize*/
+    (destructor)CLogNormal_dealloc, /*tp_dealloc*/
+    0,                         /*tp_print*/
+    0,                         /*tp_getattr*/
+    0,                         /*tp_setattr*/
+    0,                         /*tp_compare*/
+    0,                         /*tp_repr*/
+    0,                         /*tp_as_number*/
+    0,                         /*tp_as_sequence*/
+    0,                         /*tp_as_mapping*/
+    0,                         /*tp_hash */
+    0,                         /*tp_call*/
+    0,                         /*tp_str*/
+    0,                         /*tp_getattro*/
+    0,                         /*tp_setattro*/
+    0,                         /*tp_as_buffer*/
+    Py_TPFLAGS_DEFAULT | Py_TPFLAGS_BASETYPE, /*tp_flags*/
+    "CLogNormal objects",           /* tp_doc */
+    0,		               /* tp_traverse */
+    0,		               /* tp_clear */
+    0,		               /* tp_richcompare */
+    0,		               /* tp_weaklistoffset */
+    0,		               /* tp_iter */
+    0,		               /* tp_iternext */
+    CLogNormal_methods,             /* tp_methods */
+    CLogNormal_members,             /* tp_members */
+    0,                         /* tp_getset */
+    0,                         /* tp_base */
+    0,                         /* tp_dict */
+    0,                         /* tp_descr_get */
+    0,                         /* tp_descr_set */
+    0,                         /* tp_dictoffset */
+    (initproc)CLogNormal_init,      /* tp_init */
+    0,                         /* tp_alloc */
+    CLogNormal_new,                 /* tp_new */
+};
+
+
+static PyMethodDef module_methods[] = {
+    {NULL} 
+};
+
+/**
+ * Function used to add the model class to a module
+ * @param module: module to add the class to
+ */ 
+void addCLogNormal(PyObject *module) {
+	PyObject *d;
+	
+    if (PyType_Ready(&CLogNormalType) < 0)
+        return;
+
+    Py_INCREF(&CLogNormalType);
+    PyModule_AddObject(module, "CLogNormal", (PyObject *)&CLogNormalType);
+    
+    d = PyModule_GetDict(module);
+    CLogNormalError = PyErr_NewException("CLogNormal.error", NULL, NULL);
+    PyDict_SetItemString(d, "CLogNormalError", CLogNormalError);
+}
+
Index: sansmodels/src/sans/models/c_extensions/CLorentzian.c
===================================================================
--- sansmodels/src/sans/models/c_extensions/CLorentzian.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/CLorentzian.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,290 @@
+/** CLorentzian
+ *
+ * C extension 
+ *
+ * WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
+ *          DO NOT MODIFY THIS FILE, MODIFY lorentzian.h
+ *          AND RE-RUN THE GENERATOR SCRIPT
+ *
+ * @author   M.Doucet / UTK
+ */
+ 
+#include <Python.h>
+#include "structmember.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <time.h>
+
+#include "lorentzian.h"
+
+/// Error object for raised exceptions
+static PyObject * CLorentzianError = NULL;
+
+
+// Class definition
+typedef struct {
+    PyObject_HEAD
+    /// Parameters
+    PyObject * params;
+    /// Log for unit testing
+    PyObject * log;
+    /// Model parameters
+	LorentzianParameters model_pars;
+} CLorentzian;
+
+
+static void
+CLorentzian_dealloc(CLorentzian* self)
+{
+    self->ob_type->tp_free((PyObject*)self);
+    
+
+}
+
+static PyObject *
+CLorentzian_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+    CLorentzian *self;
+    
+    self = (CLorentzian *)type->tp_alloc(type, 0);
+   
+    return (PyObject *)self;
+}
+
+static int
+CLorentzian_init(CLorentzian *self, PyObject *args, PyObject *kwds)
+{
+    if (self != NULL) {
+    	
+    	// Create parameters
+        self->params = PyDict_New();
+        
+        // Initialize parameter dictionary
+        PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"center",Py_BuildValue("d",0.000000));
+        PyDict_SetItemString(self->params,"gamma",Py_BuildValue("d",1.000000));
+
+         
+        // Create empty log
+        self->log = PyDict_New();
+        
+        
+
+    }
+    return 0;
+}
+
+static PyMemberDef CLorentzian_members[] = {
+    {"params", T_OBJECT, offsetof(CLorentzian, params), 0,
+     "Parameters"},
+    {"log", T_OBJECT, offsetof(CLorentzian, log), 0,
+     "Log"},
+    {NULL}  /* Sentinel */
+};
+
+/** Read double from PyObject
+    @param p PyObject
+    @return double
+*/
+double CLorentzian_readDouble(PyObject *p) {
+    if (PyFloat_Check(p)==1) {
+        return (double)(((PyFloatObject *)(p))->ob_fval);
+    } else if (PyInt_Check(p)==1) {
+        return (double)(((PyIntObject *)(p))->ob_ival);
+    } else if (PyLong_Check(p)==1) {
+        return (double)PyLong_AsLong(p);
+    } else {
+        return 0.0;
+    }
+}
+
+
+/**
+ * Function to call to evaluate model
+ * @param args: input q or [q,phi]
+ * @return: function value
+ */
+static PyObject * run(CLorentzian *self, PyObject *args) {
+	double q_value, phi_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	// Reader parameter dictionary
+    self->model_pars.scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model_pars.center = PyFloat_AsDouble( PyDict_GetItemString(self->params, "center") );
+    self->model_pars.gamma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "gamma") );
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CLorentzianError, 
+	    	"CLorentzian.run expects a q value.");
+		return NULL;
+	}
+	  
+	// Check params
+	if( PyList_Check(pars)==1) {
+		
+		// Length of list should be 2 for I(q,phi)
+	    npars = PyList_GET_SIZE(pars); 
+	    if(npars!=2) {
+	    	PyErr_SetString(CLorentzianError, 
+	    		"CLorentzian.run expects a double or a list of dimension 2.");
+	    	return NULL;
+	    }
+	    // We have a vector q, get the q and phi values at which
+	    // to evaluate I(q,phi)
+	    q_value = CLorentzian_readDouble(PyList_GET_ITEM(pars,0));
+	    phi_value = CLorentzian_readDouble(PyList_GET_ITEM(pars,1));
+	    // Skip zero
+	    if (q_value==0) {
+	    	return Py_BuildValue("d",0.0);
+	    }
+		return Py_BuildValue("d",lorentzian_analytical_2D(&(self->model_pars),q_value,phi_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		q_value = CLorentzian_readDouble(pars);		
+		
+		return Py_BuildValue("d",lorentzian_analytical_1D(&(self->model_pars),q_value));
+	}	
+}
+
+/**
+ * Function to call to evaluate model in cartesian coordinates
+ * @param args: input q or [qx, qy]]
+ * @return: function value
+ */
+static PyObject * runXY(CLorentzian *self, PyObject *args) {
+	double qx_value, qy_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	// Reader parameter dictionary
+    self->model_pars.scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model_pars.center = PyFloat_AsDouble( PyDict_GetItemString(self->params, "center") );
+    self->model_pars.gamma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "gamma") );
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CLorentzianError, 
+	    	"CLorentzian.run expects a q value.");
+		return NULL;
+	}
+	  
+	// Check params
+	if( PyList_Check(pars)==1) {
+		
+		// Length of list should be 2 for I(qx, qy))
+	    npars = PyList_GET_SIZE(pars); 
+	    if(npars!=2) {
+	    	PyErr_SetString(CLorentzianError, 
+	    		"CLorentzian.run expects a double or a list of dimension 2.");
+	    	return NULL;
+	    }
+	    // We have a vector q, get the qx and qy values at which
+	    // to evaluate I(qx,qy)
+	    qx_value = CLorentzian_readDouble(PyList_GET_ITEM(pars,0));
+	    qy_value = CLorentzian_readDouble(PyList_GET_ITEM(pars,1));
+		return Py_BuildValue("d",lorentzian_analytical_2DXY(&(self->model_pars),qx_value,qy_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		qx_value = CLorentzian_readDouble(pars);		
+		
+		return Py_BuildValue("d",lorentzian_analytical_1D(&(self->model_pars),qx_value));
+	}	
+}
+
+static PyObject * reset(CLorentzian *self, PyObject *args) {
+    
+
+    return Py_BuildValue("d",0.0);
+}
+
+
+static PyMethodDef CLorentzian_methods[] = {
+    {"run",      (PyCFunction)run     , METH_VARARGS,
+      "Evaluate the model at a given Q or Q, phi"},
+    {"runXY",      (PyCFunction)runXY     , METH_VARARGS,
+      "Evaluate the model at a given Q or Qx, Qy"},
+    {"reset",    (PyCFunction)reset   , METH_VARARGS,
+      "Reset pair correlation"},
+    //{"numerical_1D",      (PyCFunction)numerical_1D     , METH_VARARGS,
+    //  "Evaluate the 1D model at a given Q"},
+   {NULL}
+};
+
+static PyTypeObject CLorentzianType = {
+    PyObject_HEAD_INIT(NULL)
+    0,                         /*ob_size*/
+    "CLorentzian",             /*tp_name*/
+    sizeof(CLorentzian),             /*tp_basicsize*/
+    0,                         /*tp_itemsize*/
+    (destructor)CLorentzian_dealloc, /*tp_dealloc*/
+    0,                         /*tp_print*/
+    0,                         /*tp_getattr*/
+    0,                         /*tp_setattr*/
+    0,                         /*tp_compare*/
+    0,                         /*tp_repr*/
+    0,                         /*tp_as_number*/
+    0,                         /*tp_as_sequence*/
+    0,                         /*tp_as_mapping*/
+    0,                         /*tp_hash */
+    0,                         /*tp_call*/
+    0,                         /*tp_str*/
+    0,                         /*tp_getattro*/
+    0,                         /*tp_setattro*/
+    0,                         /*tp_as_buffer*/
+    Py_TPFLAGS_DEFAULT | Py_TPFLAGS_BASETYPE, /*tp_flags*/
+    "CLorentzian objects",           /* tp_doc */
+    0,		               /* tp_traverse */
+    0,		               /* tp_clear */
+    0,		               /* tp_richcompare */
+    0,		               /* tp_weaklistoffset */
+    0,		               /* tp_iter */
+    0,		               /* tp_iternext */
+    CLorentzian_methods,             /* tp_methods */
+    CLorentzian_members,             /* tp_members */
+    0,                         /* tp_getset */
+    0,                         /* tp_base */
+    0,                         /* tp_dict */
+    0,                         /* tp_descr_get */
+    0,                         /* tp_descr_set */
+    0,                         /* tp_dictoffset */
+    (initproc)CLorentzian_init,      /* tp_init */
+    0,                         /* tp_alloc */
+    CLorentzian_new,                 /* tp_new */
+};
+
+
+static PyMethodDef module_methods[] = {
+    {NULL} 
+};
+
+/**
+ * Function used to add the model class to a module
+ * @param module: module to add the class to
+ */ 
+void addCLorentzian(PyObject *module) {
+	PyObject *d;
+	
+    if (PyType_Ready(&CLorentzianType) < 0)
+        return;
+
+    Py_INCREF(&CLorentzianType);
+    PyModule_AddObject(module, "CLorentzian", (PyObject *)&CLorentzianType);
+    
+    d = PyModule_GetDict(module);
+    CLorentzianError = PyErr_NewException("CLorentzian.error", NULL, NULL);
+    PyDict_SetItemString(d, "CLorentzianError", CLorentzianError);
+}
+
Index: sansmodels/src/sans/models/c_extensions/CSchulz.c
===================================================================
--- sansmodels/src/sans/models/c_extensions/CSchulz.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/CSchulz.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,290 @@
+/** CSchulz
+ *
+ * C extension 
+ *
+ * WARNING: THIS FILE WAS GENERATED BY WRAPPERGENERATOR.PY
+ *          DO NOT MODIFY THIS FILE, MODIFY schulz.h
+ *          AND RE-RUN THE GENERATOR SCRIPT
+ *
+ * @author   M.Doucet / UTK
+ */
+ 
+#include <Python.h>
+#include "structmember.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <time.h>
+
+#include "schulz.h"
+
+/// Error object for raised exceptions
+static PyObject * CSchulzError = NULL;
+
+
+// Class definition
+typedef struct {
+    PyObject_HEAD
+    /// Parameters
+    PyObject * params;
+    /// Log for unit testing
+    PyObject * log;
+    /// Model parameters
+	SchulzParameters model_pars;
+} CSchulz;
+
+
+static void
+CSchulz_dealloc(CSchulz* self)
+{
+    self->ob_type->tp_free((PyObject*)self);
+    
+
+}
+
+static PyObject *
+CSchulz_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+    CSchulz *self;
+    
+    self = (CSchulz *)type->tp_alloc(type, 0);
+   
+    return (PyObject *)self;
+}
+
+static int
+CSchulz_init(CSchulz *self, PyObject *args, PyObject *kwds)
+{
+    if (self != NULL) {
+    	
+    	// Create parameters
+        self->params = PyDict_New();
+        
+        // Initialize parameter dictionary
+        PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"sigma",Py_BuildValue("d",1.000000));
+        PyDict_SetItemString(self->params,"center",Py_BuildValue("d",0.000000));
+
+         
+        // Create empty log
+        self->log = PyDict_New();
+        
+        
+
+    }
+    return 0;
+}
+
+static PyMemberDef CSchulz_members[] = {
+    {"params", T_OBJECT, offsetof(CSchulz, params), 0,
+     "Parameters"},
+    {"log", T_OBJECT, offsetof(CSchulz, log), 0,
+     "Log"},
+    {NULL}  /* Sentinel */
+};
+
+/** Read double from PyObject
+    @param p PyObject
+    @return double
+*/
+double CSchulz_readDouble(PyObject *p) {
+    if (PyFloat_Check(p)==1) {
+        return (double)(((PyFloatObject *)(p))->ob_fval);
+    } else if (PyInt_Check(p)==1) {
+        return (double)(((PyIntObject *)(p))->ob_ival);
+    } else if (PyLong_Check(p)==1) {
+        return (double)PyLong_AsLong(p);
+    } else {
+        return 0.0;
+    }
+}
+
+
+/**
+ * Function to call to evaluate model
+ * @param args: input q or [q,phi]
+ * @return: function value
+ */
+static PyObject * run(CSchulz *self, PyObject *args) {
+	double q_value, phi_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	// Reader parameter dictionary
+    self->model_pars.scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model_pars.sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") );
+    self->model_pars.center = PyFloat_AsDouble( PyDict_GetItemString(self->params, "center") );
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CSchulzError, 
+	    	"CSchulz.run expects a q value.");
+		return NULL;
+	}
+	  
+	// Check params
+	if( PyList_Check(pars)==1) {
+		
+		// Length of list should be 2 for I(q,phi)
+	    npars = PyList_GET_SIZE(pars); 
+	    if(npars!=2) {
+	    	PyErr_SetString(CSchulzError, 
+	    		"CSchulz.run expects a double or a list of dimension 2.");
+	    	return NULL;
+	    }
+	    // We have a vector q, get the q and phi values at which
+	    // to evaluate I(q,phi)
+	    q_value = CSchulz_readDouble(PyList_GET_ITEM(pars,0));
+	    phi_value = CSchulz_readDouble(PyList_GET_ITEM(pars,1));
+	    // Skip zero
+	    if (q_value==0) {
+	    	return Py_BuildValue("d",0.0);
+	    }
+		return Py_BuildValue("d",schulz_analytical_2D(&(self->model_pars),q_value,phi_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		q_value = CSchulz_readDouble(pars);		
+		
+		return Py_BuildValue("d",schulz_analytical_1D(&(self->model_pars),q_value));
+	}	
+}
+
+/**
+ * Function to call to evaluate model in cartesian coordinates
+ * @param args: input q or [qx, qy]]
+ * @return: function value
+ */
+static PyObject * runXY(CSchulz *self, PyObject *args) {
+	double qx_value, qy_value;
+	PyObject* pars;
+	int npars;
+	
+	// Get parameters
+	
+	// Reader parameter dictionary
+    self->model_pars.scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") );
+    self->model_pars.sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") );
+    self->model_pars.center = PyFloat_AsDouble( PyDict_GetItemString(self->params, "center") );
+
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(CSchulzError, 
+	    	"CSchulz.run expects a q value.");
+		return NULL;
+	}
+	  
+	// Check params
+	if( PyList_Check(pars)==1) {
+		
+		// Length of list should be 2 for I(qx, qy))
+	    npars = PyList_GET_SIZE(pars); 
+	    if(npars!=2) {
+	    	PyErr_SetString(CSchulzError, 
+	    		"CSchulz.run expects a double or a list of dimension 2.");
+	    	return NULL;
+	    }
+	    // We have a vector q, get the qx and qy values at which
+	    // to evaluate I(qx,qy)
+	    qx_value = CSchulz_readDouble(PyList_GET_ITEM(pars,0));
+	    qy_value = CSchulz_readDouble(PyList_GET_ITEM(pars,1));
+		return Py_BuildValue("d",schulz_analytical_2DXY(&(self->model_pars),qx_value,qy_value));
+
+	} else {
+
+		// We have a scalar q, we will evaluate I(q)
+		qx_value = CSchulz_readDouble(pars);		
+		
+		return Py_BuildValue("d",schulz_analytical_1D(&(self->model_pars),qx_value));
+	}	
+}
+
+static PyObject * reset(CSchulz *self, PyObject *args) {
+    
+
+    return Py_BuildValue("d",0.0);
+}
+
+
+static PyMethodDef CSchulz_methods[] = {
+    {"run",      (PyCFunction)run     , METH_VARARGS,
+      "Evaluate the model at a given Q or Q, phi"},
+    {"runXY",      (PyCFunction)runXY     , METH_VARARGS,
+      "Evaluate the model at a given Q or Qx, Qy"},
+    {"reset",    (PyCFunction)reset   , METH_VARARGS,
+      "Reset pair correlation"},
+    //{"numerical_1D",      (PyCFunction)numerical_1D     , METH_VARARGS,
+    //  "Evaluate the 1D model at a given Q"},
+   {NULL}
+};
+
+static PyTypeObject CSchulzType = {
+    PyObject_HEAD_INIT(NULL)
+    0,                         /*ob_size*/
+    "CSchulz",             /*tp_name*/
+    sizeof(CSchulz),             /*tp_basicsize*/
+    0,                         /*tp_itemsize*/
+    (destructor)CSchulz_dealloc, /*tp_dealloc*/
+    0,                         /*tp_print*/
+    0,                         /*tp_getattr*/
+    0,                         /*tp_setattr*/
+    0,                         /*tp_compare*/
+    0,                         /*tp_repr*/
+    0,                         /*tp_as_number*/
+    0,                         /*tp_as_sequence*/
+    0,                         /*tp_as_mapping*/
+    0,                         /*tp_hash */
+    0,                         /*tp_call*/
+    0,                         /*tp_str*/
+    0,                         /*tp_getattro*/
+    0,                         /*tp_setattro*/
+    0,                         /*tp_as_buffer*/
+    Py_TPFLAGS_DEFAULT | Py_TPFLAGS_BASETYPE, /*tp_flags*/
+    "CSchulz objects",           /* tp_doc */
+    0,		               /* tp_traverse */
+    0,		               /* tp_clear */
+    0,		               /* tp_richcompare */
+    0,		               /* tp_weaklistoffset */
+    0,		               /* tp_iter */
+    0,		               /* tp_iternext */
+    CSchulz_methods,             /* tp_methods */
+    CSchulz_members,             /* tp_members */
+    0,                         /* tp_getset */
+    0,                         /* tp_base */
+    0,                         /* tp_dict */
+    0,                         /* tp_descr_get */
+    0,                         /* tp_descr_set */
+    0,                         /* tp_dictoffset */
+    (initproc)CSchulz_init,      /* tp_init */
+    0,                         /* tp_alloc */
+    CSchulz_new,                 /* tp_new */
+};
+
+
+static PyMethodDef module_methods[] = {
+    {NULL} 
+};
+
+/**
+ * Function used to add the model class to a module
+ * @param module: module to add the class to
+ */ 
+void addCSchulz(PyObject *module) {
+	PyObject *d;
+	
+    if (PyType_Ready(&CSchulzType) < 0)
+        return;
+
+    Py_INCREF(&CSchulzType);
+    PyModule_AddObject(module, "CSchulz", (PyObject *)&CSchulzType);
+    
+    d = PyModule_GetDict(module);
+    CSchulzError = PyErr_NewException("CSchulz.error", NULL, NULL);
+    PyDict_SetItemString(d, "CSchulzError", CSchulzError);
+}
+
Index: sansmodels/src/sans/models/c_extensions/WrapperGenerator.py
===================================================================
--- sansmodels/src/sans/models/c_extensions/WrapperGenerator.py	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/WrapperGenerator.py	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,374 @@
+#!/usr/bin/env python
+""" WrapperGenerator class to generate model code automatically.
+"""
+
+import os, sys
+import re
+
+class WrapperGenerator:
+    """ Python wrapper generator for C models
+    
+        The developer must provide a header file describing 
+        the new model.
+        
+        To provide the name of the Python class to be
+        generated, the .h file must contain the following 
+        string in the comments:
+        
+        // [PYTHONCLASS] = my_model
+        
+        where my_model must be replaced by the name of the
+        class that you want to import from sans.models.
+        (example: [PYTHONCLASS] = MyModel
+          will create a class MyModel in sans.models.MyModel.
+          It will also create a class CMyModel in 
+          sans_extension.c_models.)
+          
+        Also in comments, each parameter of the params 
+        dictionary must be declared with a default value
+        in the following way:
+        
+        // [DEFAULT]=param_name=default_value
+        
+        (example:
+            //  [DEFAULT]=radius=20.0
+        )
+          
+        See cylinder.h for an example.
+        
+        
+        A .c file corresponding to the .h file should also
+        be provided (example: my_model.h, my_model.c).
+    
+        The .h file should define two function definitions. For example,
+        cylinder.h defines the following:
+        
+            /// 1D scattering function
+            double cylinder_analytical_1D(CylinderParameters *pars, double q);
+            
+            /// 2D scattering function
+            double cylinder_analytical_2D(CylinderParameters *pars, double q, double phi);
+            
+        The .c file implements those functions.
+        
+        @author: Mathieu Doucet / UTK
+        @contact: mathieu.doucet@nist.gov
+    """
+    
+    def __init__(self, filename):
+        """ Initialization """
+        
+        ## Name of .h file to generate wrapper from
+        self.file = filename
+        
+        # Info read from file
+        
+        ## Name of python class to write
+        self.pythonClass = None
+        ## Parser in struct section
+        self.inStruct = False
+        ## Name of struct for the c object
+        self.structName = None
+        ## Dictionary of parameters
+        self.params = {}
+        ## ModelCalculation module flag
+        self.modelCalcFlag = False
+        ## List of default parameters (text)
+        self.default_list = ""
+        ##description
+        self.description=''
+        ## Dictionary of units
+        self.details = ""
+        
+    def __repr__(self):
+        """ Simple output for printing """
+        
+        rep  = "Python class: %s\n" % self.pythonClass
+        rep += "  struc name: %s\n" % self.structName
+        rep += "  params:     %s\n" % self.params
+        rep += "  description: %s\n"% self.description
+        return rep
+        
+    def read(self):
+        """ Reads in the .h file to catch parameters of the wrapper """
+        
+        # Check if the file is there
+        if not os.path.isfile(self.file):
+            raise ValueError, "File %s is not a regular file" % self.file
+        
+        # Read file
+        f = open(self.file,'r')
+        buf = f.read()
+       
+        self.default_list = "List of default parameters:\n"
+        #lines = string.split(buf,'\n')
+        lines = buf.split('\n')
+        self.details  = "## Parameter details [units, min, max]\n"
+        self.details += "        self.details = {}\n"
+        # Catch Description
+        key = "[DESCRIPTION]"
+        find_description= 0
+        temp=""
+        for line in lines:
+            if line.count(key)>0 :
+                
+                try:
+                    find_description= 1
+                    index = line.index(key)
+                    toks = line[index:].split("=",1 )
+                    temp=toks[1].lstrip().rstrip()
+                    text='text'
+                    key2="<%s>"%text.lower()
+                    if re.match(key2,temp)!=None:
+                        #index2 = line.index(key2)
+                        #temp = temp[index2:]
+                        toks2=temp.split(key2,1)
+                        self.description=toks2[1]
+                        text='text'
+                        key2="</%s>"%text.lower()
+                        if re.search(key2,toks2[1])!=None:
+                            temp=toks2[1].split(key2,1)
+                            self.description=temp[0]
+                            break
+                        #print self.description
+                    else:
+                        self.description=temp
+                        break
+                except:
+                     raise ValueError, "Could not parse file %s" % self.file
+            elif find_description==1:
+                text='text'
+                key2="</%s>"%text.lower()
+                #print "second line",line,key2,re.search(key2,line)
+                if re.search(key2,line)!=None:
+                    tok=line.split(key2,1)
+                    temp=tok[0].split("//",1)
+                    self.description+=tok[1].lstrip().rstrip()
+                    break
+                else:
+                    #if re.search("*",line)!=None:
+                    #    temp=line.split("*",1)
+                    #    self.description+='\n'+temp[1].lstrip().rstrip()
+                    if re.search("//",line)!=None:
+                        temp=line.split("//",1)
+                        self.description+='\n\t\t'+temp[1].lstrip().rstrip()
+                        
+                    else:
+                        self.description+='\n\t\t'+line.lstrip().rstrip()
+                    
+                
+        for line in lines:
+            
+            # Catch class name
+            key = "[PYTHONCLASS]"
+            if line.count(key)>0:
+                try:
+                    index = line.index(key)
+                    #toks = string.split( line[index:], "=" )
+                    toks = line[index:].split("=" )
+                    self.pythonClass = toks[1].lstrip().rstrip()
+                except:
+                    raise ValueError, "Could not parse file %s" % self.file
+                
+            # Catch struct name
+            if line.count("typedef struct")>0:
+                # We are entering a struct block
+                self.inStruct = True
+            
+            if self.inStruct and line.count("}")>0:
+                # We are exiting a struct block
+                self.inStruct = False
+    
+                # Catch the name of the struct
+                index = line.index("}")
+                #toks = string.split(line[index+1:],";")
+                toks = line[index+1:].split(";")
+                # Catch pointer definition
+                #toks2 = string.split(toks[0],',')
+                toks2 = toks[0].split(',')
+                self.structName = toks2[0].lstrip().rstrip()
+          
+            # Catch struct content
+            key = "[DEFAULT]"
+            if self.inStruct and line.count(key)>0:
+                # Found a new parameter
+                try:
+                    index = line.index(key)
+                    toks = line[index:].split("=")
+                    toks2 = toks[2].split()
+                    val = float(toks2[0])
+                    self.params[toks[1]] = val
+                    #self.pythonClass = toks[1].lstrip().rstrip()
+                    units = ""
+                    if len(toks2) >= 2:
+                        units = toks2[1]
+                    self.default_list += "         %-15s = %s %s\n" % \
+                        (toks[1], val, units)
+                    
+                    # Check for min and max
+                    min = "None"
+                    max = "None"
+                    if len(toks2) == 4:
+                        min = toks2[2]
+                        max = toks2[3]
+                    
+                    self.details += "        self.details['%s'] = ['%s', %s, %s]\n" % \
+                        (toks[1].lstrip().rstrip(), units.lstrip().rstrip(), min, max)
+                except:
+                    raise ValueError, "Could not parse input file %s \n  %s" % \
+                        (self.file, sys.exc_value)
+                
+                
+            # Catch need for numerical calculations
+            key = "CalcParameters calcPars"
+            if line.count(key)>0:
+                self.modelCalcFlag = True
+            
+                
+                
+    def write_c_wrapper(self):
+        """ Writes the C file to create the python extension class 
+            The file is written in C[PYTHONCLASS].c
+        """
+        
+        file = open("C"+self.pythonClass+'.c', 'w')
+        template = open("classTemplate.txt", 'r')
+        
+        tmp_buf = template.read()
+        #tmp_lines = string.split(tmp_buf,'\n')
+        tmp_lines = tmp_buf.split('\n')
+        
+        for tmp_line in tmp_lines:
+            
+            # Catch class name
+            newline = self.replaceToken(tmp_line, 
+                                        "[PYTHONCLASS]", 'C'+self.pythonClass)
+            ##catch description
+            #newline = self.replaceToken(tmp_line, 
+            #                            "[DESCRIPTION]", self.description)
+            # Catch class name
+            newline = self.replaceToken(newline, 
+                                        "[MODELSTRUCT]", self.structName)
+            
+            # Dictionary initialization
+            param_str = "// Initialize parameter dictionary\n"
+            for par in self.params:
+                param_str += "        PyDict_SetItemString(self->params,\"%s\",Py_BuildValue(\"d\",%f));\n" % \
+                    (par, self.params[par])
+                
+            newline = self.replaceToken(newline,
+                                        "[INITDICTIONARY]", param_str)
+            
+            # Read dictionary
+            param_str = "// Reader parameter dictionary\n"
+            for par in self.params:
+                param_str += "    self->model_pars.%s = PyFloat_AsDouble( PyDict_GetItemString(self->params, \"%s\") );\n" % \
+                    (par, par)
+                
+            newline = self.replaceToken(newline, "[READDICTIONARY]", param_str)
+                
+            # Name of .c file
+            #toks = string.split(self.file,'.')
+            toks = self.file.split('.')
+            newline = self.replaceToken(newline, "[C_FILENAME]", toks[0])
+            
+            # Include file
+            newline = self.replaceToken(newline, 
+                                        "[INCLUDE_FILE]", self.file)           
+                
+            # Numerical calcs dealloc
+            dealloc_str = "\n"
+            if self.modelCalcFlag:
+                dealloc_str = "    modelcalculations_dealloc(&(self->model_pars.calcPars));\n"
+            newline = self.replaceToken(newline, 
+                                        "[NUMERICAL_DEALLOC]", dealloc_str)     
+                
+            # Numerical calcs init
+            init_str = "\n"
+            if self.modelCalcFlag:
+                init_str = "        modelcalculations_init(&(self->model_pars.calcPars));\n"
+            newline = self.replaceToken(newline, 
+                                        "[NUMERICAL_INIT]", init_str)     
+                
+            # Numerical calcs reset
+            reset_str = "\n"
+            if self.modelCalcFlag:
+                reset_str = "modelcalculations_reset(&(self->model_pars.calcPars));\n"
+            newline = self.replaceToken(newline, 
+                                        "[NUMERICAL_RESET]", reset_str)     
+                
+            # Write new line to the wrapper .c file
+            file.write(newline+'\n')
+            
+            
+        file.close()
+        
+    def write_python_wrapper(self):
+        """ Writes the python file to create the python extension class 
+            The file is written in ../[PYTHONCLASS].py
+        """
+        
+        file = open("../"+self.pythonClass+'.py', 'w')
+        template = open("modelTemplate.txt", 'r')
+        
+        tmp_buf = template.read()
+        tmp_lines = tmp_buf.split('\n')
+        
+        for tmp_line in tmp_lines:
+            
+            # Catch class name
+            newline = self.replaceToken(tmp_line, 
+                                        "[CPYTHONCLASS]", 'C'+self.pythonClass)
+            
+            # Catch class name
+            newline = self.replaceToken(newline, 
+                                        "[PYTHONCLASS]", self.pythonClass)
+            
+            # Include file
+            newline = self.replaceToken(newline, 
+                                        "[INCLUDE_FILE]", self.file)    
+                   
+            # Include file
+            newline = self.replaceToken(newline, 
+                                        "[DEFAULT_LIST]", self.default_list)
+            # Model Description
+            newline = self.replaceToken(newline, 
+                                        "[DESCRIPTION]", self.description)
+            # Parameter details
+            newline = self.replaceToken(newline, 
+                                        "[PAR_DETAILS]", self.details)
+            #print"write",tmp_line 
+            # Write new line to the wrapper .c file
+            file.write(newline+'\n')
+           
+        file.close()
+        
+        
+    def replaceToken(self, line, key, value): #pylint: disable-msg=R0201
+        """ Replace a token in the template file 
+            @param line: line of text to inspect
+            @param key: token to look for
+            @param value: string value to replace the token with
+            @return: new string value
+        """
+        lenkey = len(key)
+        newline = line
+        while newline.count(key)>0:
+            index = newline.index(key)
+            newline = newline[:index]+value+newline[index+lenkey:]
+        return newline
+        
+        
+# main
+if __name__ == '__main__':
+    if len(sys.argv)>1:
+        print "Will look for file %s" % sys.argv[1]
+        app = WrapperGenerator(sys.argv[1])
+    else:
+        app = WrapperGenerator("test.h")
+    app.read()
+    app.write_c_wrapper()
+    app.write_python_wrapper()
+    print app
+   
+# End of file        
Index: sansmodels/src/sans/models/c_extensions/disperser.c
===================================================================
--- sansmodels/src/sans/models/c_extensions/disperser.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/disperser.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,470 @@
+/** 
+ * Disperser
+ * 
+ * Python class that takes a model and averages its output
+ * for a distribution of its parameters.
+ * 
+ * The parameters to be varied are specified at instantiation time.
+ * The distributions are Gaussian, with std deviations specified for
+ * each parameter at instantiation time.
+ * 
+ * @author   M.Doucet / UTK 
+ */
+ 
+#include <Python.h>
+#include "structmember.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <time.h>
+
+/// Error object for raised exceptions
+static PyObject * DisperserError = NULL;
+
+// Class definition
+typedef struct {
+    PyObject_HEAD
+    /// Dictionary of parameters for this class
+    PyObject * params;
+    /// Model object to disperse
+    PyObject * model;
+    /// Array of original values
+    double * centers;
+    /// Number of parameters to disperse
+    int n_pars;
+    /// Coordinate system flag
+    int cylindrical;
+} Disperser;
+
+// Function definitions
+double disperseParam(Disperser *self, int iPar, PyObject * pars);
+double Disperser_readDouble(PyObject *p);
+
+/**
+ * Function to deallocate memory for the object
+ */
+static void
+Disperser_dealloc(Disperser * self) {
+    free(self->centers);
+    self->ob_type->tp_free((PyObject*)self);
+}
+
+/**
+ * Function to create a new Disperser object
+ */
+static PyObject *
+Disperser_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+    Disperser *self;
+    
+    self = (Disperser *)type->tp_alloc(type, 0);
+   
+    return (PyObject *)self;
+}
+
+/**
+ * Initialization of the object (__init__)
+ */
+static int
+Disperser_init(Disperser *self, PyObject *args, PyObject *kwds)
+{
+	PyObject * model;
+	PyObject * paramList;
+	PyObject * sigmaList;
+	
+    if (self != NULL) {
+    	
+    	// Create parameters
+        self->params = PyDict_New();
+        
+		// Get input
+		if ( !PyArg_ParseTuple(args,"OOO",&model, &paramList, &sigmaList) ) {
+		    PyErr_SetString(DisperserError, 
+		    	"Disperser: could not parse to initialize. Need Disperser(model, paramList, sigmaList)");
+			return -1;
+		}
+
+		// Check that the lists are of the right type
+		if( PyList_Check(paramList)!=1 || PyList_Check(sigmaList)!=1 ) {
+			PyErr_SetString(DisperserError,
+				"Disperser: bad input parameters; did not find lists.");
+			return -1;
+		}
+		
+		// Check that the list lengths are compatible
+		if( PyList_GET_SIZE(paramList) != PyList_GET_SIZE(sigmaList) ) {
+			PyErr_SetString(DisperserError,
+				"Disperser: supplied lists are not of the same length.");
+			return -1;
+		} 
+		
+		// Set the data members
+		self->model = model;
+		
+		// Number of parameters to disperse
+		// Let run() find out since values can change
+		self->n_pars = 0;
+
+		// Set default coordinate system as cartesian
+		self->cylindrical = 0;
+
+        // Initialize parameter dictionary
+        // Store the list of parameter names (strings)
+        PyDict_SetItemString(self->params,"paramList", paramList);
+        
+        // Store the list of std deviations (floats)
+        PyDict_SetItemString(self->params,"sigmaList", sigmaList);
+        
+        // Set the default number of points to average over
+        PyDict_SetItemString(self->params,"n_pts",Py_BuildValue("i",25));
+    }
+    return 0;
+}
+
+/**
+ * List of member to make public
+ * Only the parameters and the model are available to python users
+ */
+static PyMemberDef Disperser_members[] = {
+    {"params", T_OBJECT, offsetof(Disperser, params), 0,
+     "Parameters"},
+    {"model", T_OBJECT, offsetof(Disperser, model), 0,
+     "Model to disperse"},
+    {NULL}  /* Sentinel */
+};
+
+/** Read double from PyObject
+ *	TODO: THIS SHOULD REALLY BE ELSEWHERE!
+ *   @param p PyObject
+ *   @return double
+ */
+double Disperser_readDouble(PyObject *p) {
+	double value = 0.0;
+	
+    if (PyFloat_Check(p)==1) {
+        value = (double)(((PyFloatObject *)(p))->ob_fval);
+        //Py_DECREF(p);
+        return value;
+    } else if (PyInt_Check(p)==1) {
+        value = (double)(((PyIntObject *)(p))->ob_ival);
+    	//Py_DECREF(p);
+        return value;
+    } else if (PyLong_Check(p)==1) {
+        value = (double)PyLong_AsLong(p);
+    	//Py_DECREF(p);
+        return value;
+    } else {
+    	//Py_DECREF(p);    	
+        return 0.0;
+    }
+}
+
+
+/**
+ * Function to call to evaluate model
+ * @param args: input q or [qx,qy]
+ * @return: function value
+ */
+static PyObject * runXY(Disperser *self, PyObject *args) {
+	PyObject* pars;
+	int n_pts;
+	int i;
+	double value;
+	PyObject * result;
+	PyObject * temp;
+	
+	// Get parameters
+	
+	// Check that the lists are of the right type
+	if( PyList_Check(PyDict_GetItemString(self->params, "paramList"))!=1 
+		|| PyList_Check(PyDict_GetItemString(self->params, "sigmaList"))!=1 ) {
+		PyErr_SetString(DisperserError,
+			"Disperser: bad input parameters; did not find lists.");
+		return NULL;
+	}
+
+	// Reader parameter dictionary
+    n_pts = PyInt_AsLong( PyDict_GetItemString(self->params, "n_pts") );
+    
+	// Check that the list lengths are compatible
+	if( PyList_GET_SIZE(PyDict_GetItemString(self->params, "paramList")) 
+		!= PyList_GET_SIZE(PyDict_GetItemString(self->params, "sigmaList")) ) {
+		PyErr_SetString(DisperserError,
+			"Disperser: supplied lists are not of the same length.");
+		return NULL;
+	} 
+
+    // Number of parameters to disperse
+    self->n_pars = PyList_GET_SIZE(PyDict_GetItemString(self->params, "paramList"));
+    
+    // Allocate memory for centers
+    free(self->centers);
+    self->centers = (double *)malloc(self->n_pars * sizeof(double));
+	if(self->centers==NULL) {
+	    PyErr_SetString(DisperserError, 
+	    	"Disperser.run could not allocate memory.");
+		return NULL;
+	}
+    
+    // Store initial values of those parameters    
+    for( i=0; i<self->n_pars; i++ ) {
+    	result = PyObject_CallMethod(self->model, "getParam", "(O)", 
+    		PyList_GetItem( PyDict_GetItemString(self->params, "paramList"), i ));
+    		
+    	if( result == NULL ) {
+		    PyErr_SetString(DisperserError, 
+		    	"Disperser.run could not get current parameter values.");
+    		return NULL;
+    	} else {
+    		self->centers[i] = Disperser_readDouble(result);
+    		Py_DECREF(result);
+    	}; 
+    }
+	
+	// Get input and determine whether we have to supply a 1D or 2D return value.
+	if ( !PyArg_ParseTuple(args,"O",&pars) ) {
+	    PyErr_SetString(DisperserError, 
+	    	"Disperser.run expects a q value.");
+		return NULL;
+	}    
+        
+    // Evaluate or calculate average
+    if( self->n_pars > 0 ) {
+		value = disperseParam(self, 0, pars);
+    } else {
+    	 if (self->cylindrical==1) {
+    	 	result = PyObject_CallMethod(self->model, "run", "(O)",pars);
+    	 } else {
+    	 	result = PyObject_CallMethod(self->model, "runXY", "(O)",pars);
+    	 }
+    	 
+    	 if( result == NULL ) {
+		    PyErr_SetString(DisperserError, 
+		    	"Disperser.run failed.");
+    		return NULL; 	 	
+    	 } else {
+    		Py_DECREF(result);
+    	 	return result;
+    	 }
+    }
+    
+    // Put back the original model parameter values
+    for( i=0; i<self->n_pars; i++ ) {
+    	// Maybe we need to decref the double here?
+    	temp = Py_BuildValue("d", self->centers[i]);
+    	result = PyObject_CallMethod(self->model, "setParam", "(OO)", 
+    		PyList_GetItem( PyDict_GetItemString(self->params, "paramList"), i ),
+    		temp );
+    	Py_DECREF(temp);
+    		
+    	if( result == NULL ) {
+		    PyErr_SetString(DisperserError, 
+		    	"Disperser.run could not set back original parameter values.");
+    		return NULL;
+    	} else {
+    		Py_DECREF(result);
+    	}; 
+    }
+    
+	return Py_BuildValue("d",value);    	
+
+}
+
+/**
+ * Function to call to evaluate model
+ * @param args: input q or [q,phi]
+ * @return: function value
+ */
+static PyObject * run(Disperser *self, PyObject *args) {
+	PyObject * output;
+	self->cylindrical = 1;
+	output = runXY(self, args);
+	self->cylindrical = 0;
+	return output;
+}
+
+
+
+
+/**
+ * Gaussian weight
+ * @param mean: mean value of the Gaussian
+ * @param sigma: standard deviation of the Gaussian
+ * @param x: value at which the Gaussian is evaluated
+ * @return: value of the Gaussian
+ */
+double gaussian_weight(double mean, double sigma, double x) {
+	double vary, expo_value;
+    vary = x-mean;
+    expo_value = -vary*vary/(2*sigma*sigma);
+    //return 1.0;
+    return exp(expo_value);
+}
+
+/**
+ * @param iPar: ID of parameter to disperse
+ * @param pars: input parameter to the evaluation function
+ */
+double disperseParam(Disperser *self, int iPar, PyObject * pars) {
+	double sigma;
+	double min_value, max_value;
+	double step;
+	double prev_value;
+	double value_sum;
+	double gauss_sum;
+	double gauss_value;
+	double func_value;
+	double error_sys;
+	double value;
+	int i, n_pts;
+	PyObject *result;
+	PyObject *temp;
+	
+	
+    n_pts = PyInt_AsLong( PyDict_GetItemString(self->params, "n_pts") );
+    
+	// If we exhausted the parameter array, simply evaluate
+    // the model
+    if( iPar < self->n_pars ) {
+            
+		// Average over Gaussian distribution (2 sigmas)
+        value_sum = 0.0;
+        gauss_sum = 0.0;
+            
+        // Get the standard deviation for this parameter
+        sigma = Disperser_readDouble(PyList_GetItem( PyDict_GetItemString(self->params, "sigmaList"), iPar ));
+            
+        // Average over 4 sigmas wide
+        min_value = self->centers[iPar] - 2*sigma;
+        max_value = self->centers[iPar] + 2*sigma;
+            
+        // Calculate step size
+        step = (max_value - min_value)/(n_pts-1);
+            
+        // If we are not changing the parameter, just return the 
+        // value of the function
+        if (step == 0.0) {
+            return disperseParam(self, iPar+1, pars);
+        }
+            
+        // Compute average
+        prev_value = 0.0;
+        error_sys  = 0.0;
+        for( i=0; i<n_pts; i++ ) {
+            // Set the parameter value           
+            value = min_value + (double)i*step;
+                
+            temp = Py_BuildValue("d", value);
+    		result = PyObject_CallMethod(self->model, "setParam", "(OO)",
+    			PyList_GetItem( PyDict_GetItemString(self->params, "paramList"), iPar ), 
+    				temp);
+    		Py_DECREF(temp);
+    			
+		    if( result == NULL ) {
+	  	 	 	printf("Could not set param %i\n", iPar);
+	    	 	// return a value that will create an NAN
+	    	 	return 0.0/(step-step);
+		    } else {
+		    	Py_DECREF(result);
+		    }
+                
+			gauss_value = gaussian_weight(self->centers[iPar], sigma, value);
+            func_value = disperseParam(self, iPar+1, pars);
+
+            value_sum += gauss_value * func_value;
+        	gauss_sum += gauss_value;	
+        }   
+        return value_sum/gauss_sum;
+        
+    } else {
+    	 if (self->cylindrical==1) {
+    	 	result = PyObject_CallMethod(self->model, "run", "(O)",pars);
+    	 } else {
+    	 	result = PyObject_CallMethod(self->model, "runXY", "(O)",pars);
+    	 }
+    	 if( result == NULL ) {
+    	 	printf("Model.run() could not return\n");
+    	 	// return a value that will create an NAN
+    	 	return 0.0/(step-step);
+    	 } else {
+    	 	value = Disperser_readDouble(result);    	 	
+    	 	Py_DECREF(result);
+    	 	return value;
+    	 }
+    }
+	
+}
+
+static PyMethodDef Disperser_methods[] = {
+    {"run",      (PyCFunction)run     , METH_VARARGS,
+      "Evaluate the model at a given input value"},
+    {"runXY",      (PyCFunction)runXY     , METH_VARARGS,
+      "Evaluate the model at a given input value"},
+   {NULL}
+};
+
+static PyTypeObject DisperserType = {
+    PyObject_HEAD_INIT(NULL)
+    0,                         /*ob_size*/
+    "Disperser",             /*tp_name*/
+    sizeof(Disperser),             /*tp_basicsize*/
+    0,                         /*tp_itemsize*/
+    (destructor)Disperser_dealloc, /*tp_dealloc*/
+    0,                         /*tp_print*/
+    0,                         /*tp_getattr*/
+    0,                         /*tp_setattr*/
+    0,                         /*tp_compare*/
+    0,                         /*tp_repr*/
+    0,                         /*tp_as_number*/
+    0,                         /*tp_as_sequence*/
+    0,                         /*tp_as_mapping*/
+    0,                         /*tp_hash */
+    0,                         /*tp_call*/
+    0,                         /*tp_str*/
+    0,                         /*tp_getattro*/
+    0,                         /*tp_setattro*/
+    0,                         /*tp_as_buffer*/
+    Py_TPFLAGS_DEFAULT | Py_TPFLAGS_BASETYPE, /*tp_flags*/
+    "Disperser objects",           /* tp_doc */
+    0,		               /* tp_traverse */
+    0,		               /* tp_clear */
+    0,		               /* tp_richcompare */
+    0,		               /* tp_weaklistoffset */
+    0,		               /* tp_iter */
+    0,		               /* tp_iternext */
+    Disperser_methods,             /* tp_methods */
+    Disperser_members,             /* tp_members */
+    0,                         /* tp_getset */
+    0,                         /* tp_base */
+    0,                         /* tp_dict */
+    0,                         /* tp_descr_get */
+    0,                         /* tp_descr_set */
+    0,                         /* tp_dictoffset */
+    (initproc)Disperser_init,      /* tp_init */
+    0,                         /* tp_alloc */
+    Disperser_new,                 /* tp_new */
+};
+
+
+static PyMethodDef module_methods[] = {
+    {NULL} 
+};
+
+/**
+ * Function used to add the model class to a module
+ * @param module: module to add the class to
+ */ 
+void addDisperser(PyObject *module) {
+	PyObject *d;
+	
+    if (PyType_Ready(&DisperserType) < 0)
+        return;
+
+    Py_INCREF(&DisperserType);
+    PyModule_AddObject(module, "Disperser", (PyObject *)&DisperserType);
+    
+    d = PyModule_GetDict(module);
+    DisperserError = PyErr_NewException("Disperser.error", NULL, NULL);
+    PyDict_SetItemString(d, "DisperserError", DisperserError);
+}
+
Index: sansmodels/src/sans/models/c_extensions/gaussian.c
===================================================================
--- sansmodels/src/sans/models/c_extensions/gaussian.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/gaussian.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,58 @@
+/**
+ * Gaussian function
+ */
+
+#include "gaussian.h"
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+
+
+/**
+ * Function to evaluate 1D Gaussian function.
+ * The function is normalized to the 'scale' parameter.
+ * 
+ * f(x)=scale * 1/(sigma^2*2pi)e^(-(x-mu)^2/2sigma^2)
+ * 
+ * @param pars: parameters of the gaussian
+ * @param x: x-value
+ * @return: function value
+ */
+double gaussian_analytical_1D(GaussianParameters *pars, double x) {
+	double sigma2 = pow(pars->sigma, 2);	
+	return pars->scale / sigma2 * exp( -pow((x-pars->center), 2) / (2*sigma2));
+}
+
+/**
+ * Function to evaluate 2D Gaussian  function
+ * The function is normalized to the 'scale' parameter.
+ * 
+ * f(x,y) = Gaussian(x) * Gaussian(y)
+ * 
+ * where both Gaussians share the same parameters.
+ * 
+ * @param pars: parameters of the gaussian
+ * @param x: x-value
+ * @param y: y-value
+ * @return: function value
+ */
+double gaussian_analytical_2DXY(GaussianParameters *pars, double x, double y) {
+    return gaussian_analytical_1D(pars, x) * gaussian_analytical_1D(pars, y);
+} 
+
+/**
+ * Function to evaluate 2D Gaussian  function
+ * The function is normalized to the 'scale' parameter.
+ * 
+ * f(x,y) = Gaussian(x) * Gaussian(y)
+ * 
+ * where both Gaussians share the same parameters.
+ * 
+ * @param pars: parameters of the gaussian
+ * @param length: length of the (x,y) vector
+ * @param phi: angle relative to x
+ * @return: function value
+ */
+double gaussian_analytical_2D(GaussianParameters *pars, double length, double phi) {
+    return gaussian_analytical_2DXY(pars, length*cos(phi), length*sin(phi));
+} 
Index: sansmodels/src/sans/models/c_extensions/gaussian.h
===================================================================
--- sansmodels/src/sans/models/c_extensions/gaussian.h	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/gaussian.h	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,31 @@
+#if !defined(gaussian_h)
+#define gaussian_h
+
+/** Structure definition for Gaussian parameters.
+ * The Gaussian is normalized to the 'scale' parameter.
+ *  
+ * f(x)=scale * 1/(sigma^2*2pi)e^(-(x-mu)^2/2sigma^2)
+ * 
+ * [PYTHONCLASS] = Gaussian
+ * [DESCRIPTION] = <text>f(x)=scale * 1/(sigma^2*2pi)e^(-(x-mu)^2/2sigma^2)</text>
+ */
+typedef struct {
+    /// Scale factor 
+    //  [DEFAULT]=scale=1.0
+    double scale;
+    /// Standard deviation
+    //  [DEFAULT]=sigma=1 
+    double sigma;
+    /// Center of the Gaussian distribution
+    //  [DEFAULT]=center=0.0 
+    double center;
+} GaussianParameters;
+
+/// 1D scattering function
+double gaussian_analytical_1D(GaussianParameters *pars, double x);
+
+/// 2D scattering function
+double gaussian_analytical_2D(GaussianParameters *pars, double x, double phi);
+double gaussian_analytical_2DXY(GaussianParameters *pars, double x, double y);
+
+#endif
Index: sansmodels/src/sans/models/c_extensions/logNormal.c
===================================================================
--- sansmodels/src/sans/models/c_extensions/logNormal.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/logNormal.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,58 @@
+/**
+ * log normal function
+ */
+
+#include "logNormal.h"
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+
+
+/**
+ * Function to evaluate 1D Log normal function.
+ * The function is normalized to the 'scale' parameter.
+ * 
+ * f(x)=scale * 1/(sigma*math.sqrt(2pi))e^(-1/2*((math.log(x)-mu)/sigma)^2)
+ * 
+ * @param pars: parameters of the log normal 
+ * @param x: x-value
+ * @return: function value
+ */
+double logNormal_analytical_1D(LogNormalParameters *pars, double x) {
+	double sigma2 = pow(pars->sigma, 2);	
+	return pars->scale / (x*sigma2) * exp( -pow((log(x) -pars->center), 2) / (2*sigma2));
+}
+
+/**
+ * Function to evaluate 2D log normal  function
+ * The function is normalized to the 'scale' parameter.
+ * 
+ * f(x,y) = LogNormal(x) * Lognormal(y)
+ * 
+ * where both Gaussians share the same parameters.
+ * 
+ * @param pars: parameters of the log normal 
+ * @param x: x-value
+ * @param y: y-value
+ * @return: function value
+ */
+double logNormal_analytical_2DXY(LogNormalParameters *pars, double x, double y) {
+    return logNormal_analytical_1D(pars, x) * logNormal_analytical_1D(pars, y);
+} 
+
+/**
+ * Function to evaluate 2D log normal  function
+ * The function is normalized to the 'scale' parameter.
+ * 
+ * f(x,y) = Lognormal(x) * Lognormal(y)
+ * 
+ * where both log normal s share the same parameters.
+ * 
+ * @param pars: parameters of the log normal
+ * @param length: length of the (x,y) vector
+ * @param phi: angle relative to x
+ * @return: function value
+ */
+double logNormal_analytical_2D(LogNormalParameters *pars, double length, double phi) {
+    return logNormal_analytical_2DXY(pars, length*cos(phi), length*sin(phi));
+} 
Index: sansmodels/src/sans/models/c_extensions/logNormal.h
===================================================================
--- sansmodels/src/sans/models/c_extensions/logNormal.h	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/logNormal.h	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,31 @@
+#if !defined(logNormal_h)
+#define logNormal_h
+
+/** Structure definition for Gaussian parameters.
+ * The Log normal is normalized to the 'scale' parameter.
+ *  
+ * f(x)=scale * 1/(sigma*math.sqrt(2pi))e^(-1/2*((math.log(x)-mu)/sigma)^2)
+ *
+ * [PYTHONCLASS] = LogNormal
+ * [DESCRIPTION] = <text>f(x)=scale * 1/(sigma*math.sqrt(2pi))e^(-1/2*((math.log(x)-mu)/sigma)^2)</text>
+ */
+typedef struct {
+    /// Scale factor 
+    //  [DEFAULT]=scale=1.0
+    double scale;
+    /// Standard deviation
+    //  [DEFAULT]=sigma=1 
+    double sigma;
+    /// Center of the Log Normal distribution
+    //  [DEFAULT]=center=0.0 
+    double center;
+} LogNormalParameters;
+
+/// 1D scattering function
+double logNormal_analytical_1D(LogNormalParameters *pars, double x);
+
+/// 2D scattering function
+double logNormal_analytical_2D(LogNormalParameters *pars, double x, double phi);
+double logNormal_analytical_2DXY(LogNormalParameters *pars, double x, double y);
+
+#endif
Index: sansmodels/src/sans/models/c_extensions/lorentzian.c
===================================================================
--- sansmodels/src/sans/models/c_extensions/lorentzian.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/lorentzian.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,59 @@
+/**
+ * Lorentzian function
+ */
+
+#include "lorentzian.h"
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+
+
+/**
+ * Function to evaluate 1D Lorentzian function.
+ * The function is normalized to the 'scale' parameter.
+ * 
+ * f(x)=scale * 1/pi 0.5gamma / [ (x-x_0)^2 + (0.5gamma)^2 ]
+ * 
+ * @param pars: parameters of the Lorentzian
+ * @param x: x-value
+ * @return: function value
+ */
+double lorentzian_analytical_1D(LorentzianParameters *pars, double x) {
+	double delta2 = pow( (x-pars->center), 2);
+	return pars->scale / acos(-1.0) * 0.5*pars->gamma /
+		(delta2 + 0.25*pars->gamma*pars->gamma);
+}
+
+/**
+ * Function to evaluate 2D Lorentzian  function
+ * The function is normalized to the 'scale' parameter.
+ * 
+ * f(x,y) = Lorentzian(x) * Lorentzian(y)
+ * 
+ * where both Lorentzians share the same parameters.
+ * 
+ * @param pars: parameters of the Lorentzian
+ * @param x: x-value
+ * @param y: y-value
+ * @return: function value
+ */
+double lorentzian_analytical_2DXY(LorentzianParameters *pars, double x, double y) {
+    return lorentzian_analytical_1D(pars, x) * lorentzian_analytical_1D(pars, y);
+} 
+
+/**
+ * Function to evaluate 2D Lorentzian  function
+ * The function is normalized to the 'scale' parameter.
+ * 
+ * f(x,y) = Lorentzian(x) * Lorentzian(y)
+ * 
+ * where both Lorentzians share the same parameters.
+ * 
+ * @param pars: parameters of the Lorentzian
+ * @param length: length of the (x,y) vector
+ * @param phi: angle relative to x
+ * @return: function value
+ */
+double lorentzian_analytical_2D(LorentzianParameters *pars, double length, double phi) {
+    return lorentzian_analytical_2DXY(pars, length*cos(phi), length*sin(phi));
+} 
Index: sansmodels/src/sans/models/c_extensions/lorentzian.h
===================================================================
--- sansmodels/src/sans/models/c_extensions/lorentzian.h	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/lorentzian.h	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,31 @@
+#if !defined(lorentzian_h)
+#define lorentzian_h
+
+/** Structure definition for Lorentzian parameters.
+ * The Lorentzian is normalized to the 'scale' parameter.
+ *  
+ * f(x)=scale * 1/pi 0.5gamma / [ (x-x_0)^2 + (0.5gamma)^2 ]
+ * 
+ * [PYTHONCLASS] = Lorentzian
+ * [DESCRIPTION] = <text>f(x)=scale * 1/pi 0.5gamma / [ (x-x_0)^2 + (0.5gamma)^2 ]</text>
+ */
+typedef struct {
+    /// Scale factor 
+    //  [DEFAULT]=scale=1.0
+    double scale;
+    /// Width
+    //  [DEFAULT]=gamma=1.0
+    double gamma;
+    /// Center of the Lorentzian distribution
+    //  [DEFAULT]=center=0.0 
+    double center;
+} LorentzianParameters;
+
+/// 1D lorentzian function
+double lorentzian_analytical_1D(LorentzianParameters *pars, double x);
+
+/// 2D lorentzian function
+double lorentzian_analytical_2D(LorentzianParameters *pars, double x, double phi);
+double lorentzian_analytical_2DXY(LorentzianParameters *pars, double x, double y);
+
+#endif
Index: sansmodels/src/sans/models/c_extensions/schulz.c
===================================================================
--- sansmodels/src/sans/models/c_extensions/schulz.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/schulz.c	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,65 @@
+/**
+ * Schulz function
+ */
+
+#include "schulz.h"
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+
+
+/**
+ * Function to evaluate 1D Schulz function.
+ * The function is normalized to the 'scale' parameter.
+ *
+ * f(x)=scale * math.pow(z+1, z+1)*math.pow((R), z)*
+ *					math.exp(-R*(z+1))/(center*gamma(z+1)
+ *		z= math.pow[(1/(sigma/center),2]-1
+ *		R= x/center
+ * @param pars: parameters of the schulz
+ * @param x: x-value
+ * @return: function value
+ */
+double schulz_analytical_1D(SchulzParameters *pars, double x) {
+	double z = pow(pars->center/ pars->sigma, 2)-1;
+	double R= x/pars->center;
+	double zz= z+1;
+	double expo;
+	expo = log(pars->scale)+zz*log(zz)+z*log(R)-R*zz-log(pars->center)-lgamma(zz);
+
+	return exp(expo);//pars->scale * pow(zz,zz) * pow(R,z) * exp(-1*R*zz)/((pars->center) * tgamma(zz)) ;
+}
+
+/**
+ * Function to evaluate 2D schulz function
+ * The function is normalized to the 'scale' parameter.
+ *
+ * f(x,y) = Schulz(x) * Schulz(y)
+ *
+ * where both Shulzs share the same parameters.
+ *
+ * @param pars: parameters of the schulz
+ * @param x: x-value
+ * @param y: y-value
+ * @return: function value
+ */
+double schulz_analytical_2DXY(SchulzParameters *pars, double x, double y) {
+    return schulz_analytical_1D(pars, x) * schulz_analytical_1D(pars, y);
+}
+
+/**
+ * Function to evaluate 2D Schulz  function
+ * The function is normalized to the 'scale' parameter.
+ *
+ * f(x,y) = Schulz(x) * Schulz(y)
+ *
+ * where both Gaussians share the same parameters.
+ *
+ * @param pars: parameters of the gaussian
+ * @param length: length of the (x,y) vector
+ * @param phi: angle relative to x
+ * @return: function value
+ */
+double schulz_analytical_2D(SchulzParameters *pars, double length, double phi) {
+    return schulz_analytical_2DXY(pars, length*cos(phi), length*sin(phi));
+}
Index: sansmodels/src/sans/models/c_extensions/schulz.h
===================================================================
--- sansmodels/src/sans/models/c_extensions/schulz.h	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
+++ sansmodels/src/sans/models/c_extensions/schulz.h	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -0,0 +1,37 @@
+#if !defined(schulz_h)
+#define schulz_h
+
+/** Structure definition for Schulz parameters.
+ * The Schulz is normalized to the 'scale' parameter.
+ *  
+ * f(x)=scale * math.pow(z+1, z+1)*math.pow((R), z)*
+ *					math.exp(-R*(z+1))/(center*gamma(z+1)
+ *		z= math.pow[(1/(sigma/center),2]-1
+ *		R= x/center
+ * 
+ * [PYTHONCLASS] = Schulz
+ * [DESCRIPTION] = <text> f(x)=scale * math.pow(z+1, z+1)*math.pow((R), z)*
+ 					math.exp(-R*(z+1))/(center*gamma(z+1)
+ 		            z= math.pow[(1/(sigma/center),2]-1
+ 					R= x/center</text>
+ */
+typedef struct {
+    /// Scale factor 
+    //  [DEFAULT]=scale=1.0
+    double scale;
+    /// Standard deviation
+    //  [DEFAULT]=sigma=1 
+    double sigma;
+    /// Center of the Schulz distribution
+    //  [DEFAULT]=center=0.0 
+    double center;
+} SchulzParameters;
+
+/// 1D scattering function
+double schulz_analytical_1D(SchulzParameters *pars, double x);
+
+/// 2D scattering function
+double schulz_analytical_2D(SchulzParameters *pars, double x, double phi);
+double schulz_analytical_2DXY(SchulzParameters *pars, double x, double y);
+
+#endif
Index: sansmodels/src/sans/models/c_models/c_models.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/c_models.cpp	(revision 770bab1b907c20488df256707ca66e7bc3eccb32)
+++ sansmodels/src/sans/models/c_models/c_models.cpp	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -290,9 +290,9 @@
 	addCMultiShellModel(m);
 	addCBinaryHSModel(m);
-	//addDisperser(m);
-	//addCGaussian(m);
-	//addCSchulz(m);
-	//addCLogNormal(m);
-	//addCLorentzian(m);
+	addDisperser(m);
+	addCGaussian(m);
+	addCSchulz(m);
+	addCLogNormal(m);
+	addCLorentzian(m);
 	addCVesicleModel(m);
 	addCPoly_GaussCoil(m);
Index: sansmodels/src/sans/models/c_models/parameters.cpp
===================================================================
--- sansmodels/src/sans/models/c_models/parameters.cpp	(revision 8dc02d8b4348f58c310e83b1e0b10a5ab1e56f16)
+++ sansmodels/src/sans/models/c_models/parameters.cpp	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -58,11 +58,19 @@
 	Parameter* par = (Parameter*)param;
 	double value = (*par)();
-
-	if (npts<2) {
-		weights.insert(weights.end(), WeightPoint(value, 1.0));
-	} else {
-		for(int i=0; i<npts; i++) {
-			double val = value + width * (1.0*double(i)/double(npts-1) - 0.5);
-
+	double sig;
+	if (npts<2) {
+		weights.insert(weights.end(), WeightPoint(value, 1.0));
+	} else {
+		for(int i=0; i<npts; i++) {
+
+			if ((*par).has_min==false){
+				// width = sigma for angles
+				sig = width;
+			}
+			else{
+				//width = polydispersity (=sigma/value) for length
+				sig = width * value;
+			}
+			double val = value + sig * (1.0*double(i)/double(npts-1) - 0.5);
 			if ( ((*par).has_min==false || val>(*par).min)
 			  && ((*par).has_max==false || val<(*par).max)  )
@@ -83,7 +91,7 @@
 
 GaussianDispersion :: GaussianDispersion() {
-	npts  = 21;
+	npts  = 11;
 	width = 0.0;
-	nsigmas = 3.0;
+	nsigmas = 2.5;
 };
 
@@ -115,19 +123,27 @@
 		width = 0.0;
 		npts  = 1;
-		nsigmas = 3.0;
-	}
-
-	Parameter* par = (Parameter*)param;
-	double value = (*par)();
-
-	if (npts<2) {
-		weights.insert(weights.end(), WeightPoint(value, 1.0));
-	} else {
-		for(int i=0; i<npts; i++) {
+		nsigmas = 2.5;
+	}
+
+	Parameter* par = (Parameter*)param;
+	double value = (*par)();
+	double sig;
+	if (npts<2) {
+		weights.insert(weights.end(), WeightPoint(value, 1.0));
+	} else {
+		for(int i=0; i<npts; i++) {
+			if ((*par).has_min==false){
+				// width = sigma for angles
+				sig = width;
+			}
+			else{
+				//width = polydispersity (=sigma/value) for length
+				sig = width * value;
+			}
 			// We cover n(nsigmas) times sigmas on each side of the mean
-			double val = value + width * (2.0*nsigmas*double(i)/double(npts-1) - nsigmas);
+			double val = value + sig * (2.0*nsigmas*double(i)/double(npts-1) - nsigmas);
 			if ( ((*par).has_min==false || val>(*par).min)
 			  && ((*par).has_max==false || val<(*par).max)  ) {
-				double _w = gaussian_weight(value, width, val);
+				double _w = gaussian_weight(value, sig, val);
 				weights.insert(weights.end(), WeightPoint(val, _w));
 			}
@@ -142,5 +158,5 @@
 
 RectangleDispersion :: RectangleDispersion() {
-	npts  = 21;
+	npts  = 11;
 	width = 0.0;
 	nsigmas = 1.0;
@@ -182,14 +198,22 @@
 	Parameter* par = (Parameter*)param;
 	double value = (*par)();
-
-	if (npts<2) {
-		weights.insert(weights.end(), WeightPoint(value, 1.0));
-	} else {
-		for(int i=0; i<npts; i++) {
+	double sig;
+	if (npts<2) {
+		weights.insert(weights.end(), WeightPoint(value, 1.0));
+	} else {
+		for(int i=0; i<npts; i++) {
+			if ((*par).has_min==false){
+				// width = sigma for angles
+				sig = width;
+			}
+			else{
+				//width = polydispersity (=sigma/value) for length
+				sig = width * value;
+			}
 			// We cover n(nsigmas) times sigmas on each side of the mean
-			double val = value + width * (2.0*nsigmas*double(i)/double(npts-1) - nsigmas);
+			double val = value + sig * (2.0*nsigmas*double(i)/double(npts-1) - nsigmas);
 			if ( ((*par).has_min==false || val>(*par).min)
 			  && ((*par).has_max==false || val<(*par).max)  ) {
-				double _w = rectangle_weight(value, width, val);
+				double _w = rectangle_weight(value, sig, val);
 				weights.insert(weights.end(), WeightPoint(val, _w));
 			}
@@ -204,5 +228,5 @@
 
 LogNormalDispersion :: LogNormalDispersion() {
-	npts  = 21;
+	npts  = 11;
 	width = 0.0;
 	nsigmas = 3.0;
@@ -240,15 +264,23 @@
 	Parameter* par = (Parameter*)param;
 	double value = (*par)();
-
-	if (npts<2) {
-		weights.insert(weights.end(), WeightPoint(value, 1.0));
-	} else {
-		for(int i=0; i<npts; i++) {
+	double sig;
+	if (npts<2) {
+		weights.insert(weights.end(), WeightPoint(value, 1.0));
+	} else {
+		for(int i=0; i<npts; i++) {
+			if ((*par).has_min==false){
+				// width = sigma for angles
+				sig = width;
+			}
+			else{
+				//width = polydispersity (=sigma/value) for length
+				sig = width * value;
+			}
 			// We cover n(nsigmas) times sigmas on each side of the mean
-			double val = value + width * (2.0*nsigmas*double(i)/double(npts-1) - nsigmas);
+			double val = value + sig * (2.0*nsigmas*double(i)/double(npts-1) - nsigmas);
 
 			if ( ((*par).has_min==false || val>(*par).min)
 			  && ((*par).has_max==false || val<(*par).max)  ) {
-				double _w = lognormal_weight(value, width, val);
+				double _w = lognormal_weight(value, sig, val);
 				weights.insert(weights.end(), WeightPoint(val, _w));
 			}
@@ -264,5 +296,5 @@
 
 SchulzDispersion :: SchulzDispersion() {
-	npts  = 21;
+	npts  = 11;
 	width = 0.0;
 	nsigmas = 3.0;
@@ -303,15 +335,23 @@
 	Parameter* par = (Parameter*)param;
 	double value = (*par)();
-
-	if (npts<2) {
-		weights.insert(weights.end(), WeightPoint(value, 1.0));
-	} else {
-		for(int i=0; i<npts; i++) {
+	double sig;
+	if (npts<2) {
+		weights.insert(weights.end(), WeightPoint(value, 1.0));
+	} else {
+		for(int i=0; i<npts; i++) {
+			if ((*par).has_min==false){
+				// width = sigma for angles
+				sig = width;
+			}
+			else{
+				//width = polydispersity (=sigma/value) for length
+				sig = width * value;
+			}
 			// We cover n(nsigmas) times sigmas on each side of the mean
-			double val = value + width * (2.0*nsigmas*double(i)/double(npts-1) - nsigmas);
+			double val = value + sig * (2.0*nsigmas*double(i)/double(npts-1) - nsigmas);
 
 			if ( ((*par).has_min==false || val>(*par).min)
 			  && ((*par).has_max==false || val<(*par).max)  ) {
-				double _w = schulz_weight(value, width, val);
+				double _w = schulz_weight(value, sig, val);
 				weights.insert(weights.end(), WeightPoint(val, _w));
 			}
@@ -346,4 +386,5 @@
 		for(int i=0; i<npts; i++) {
 			double val = _values[i]; //+ value;  //ToDo: Talk to Paul and put back the 'value'.
+
 			if ( ((*par).has_min==false || val>(*par).min)
 			  && ((*par).has_max==false || val<(*par).max)  )
Index: sansmodels/src/sans/models/test/utest_dispersity.py
===================================================================
--- sansmodels/src/sans/models/test/utest_dispersity.py	(revision 18b89c4c6d865555c4957a4154ec169c3c1add6c)
+++ sansmodels/src/sans/models/test/utest_dispersity.py	(revision 59b9b67578f11772d51a4d61028e0dfc7097df5a)
@@ -3,4 +3,12 @@
     C++ model classes
 """
+
+#Note: The 'sans.models.DisperseModel' is for only the test. We use 
+#'sans.models.dispersion_models', instead in the application.
+#The first uses width = sigma, while the second uses width = ratio (=sigma/mean)
+#for length parameters and width = sigma for angle parameters.
+#In Feb. 2011, we found and fixed the some precision problems in the C, so that
+#this test was updated too.
+
 
 import unittest, math, numpy
@@ -32,10 +40,11 @@
         self.model.setParam('scale', 10.0)
         self.model.set_dispersion('radius', disp)
-        self.model.dispersion['radius']['width'] = 5.0
+        self.model.dispersion['radius']['width'] = 0.25
         self.model.dispersion['radius']['npts'] = 100
         self.model.dispersion['radius']['nsigmas'] = 2.5
         
         self.assertAlmostEqual(self.model.run(0.001), 1.021051*4527.47250339, 3)
-        self.assertAlmostEqual(self.model.runXY([0.001, 0.001]), 1.021048*4546.997777604715, 2)
+        self.assertAlmostEqual(self.model.runXY([0.001, 0.001]), 
+                               1.021048*4546.997777604715, 2)
         
     def test_gaussian(self):
@@ -43,11 +52,13 @@
         disp = GaussianDispersion()
         self.model.set_dispersion('radius', disp)
-        self.model.dispersion['radius']['width'] = 5.0
+        self.model.dispersion['radius']['width'] = 0.25
         self.model.dispersion['radius']['npts'] = 100
-        self.model.dispersion['radius']['nsigmas'] = 2.5
+        self.model.dispersion['radius']['nsigmas'] = 2
         self.model.setParam('scale', 10.0)
         
-        self.assertAlmostEqual(self.model.run(0.001), 1.1804794*4723.32213339, 3)
-        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 1.180454*4743.56, 2)
+        self.assertAlmostEqual(self.model.run(0.001), 
+                               1.1804794*4723.32213339, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 
+                               1.180454*4743.56, 2)
         
     def test_clone(self):
@@ -55,12 +66,14 @@
         disp = GaussianDispersion()
         self.model.set_dispersion('radius', disp)
-        self.model.dispersion['radius']['width'] = 5.0
+        self.model.dispersion['radius']['width'] = 0.25
         self.model.dispersion['radius']['npts'] = 100
-        self.model.dispersion['radius']['nsigmas'] = 2.5
+        self.model.dispersion['radius']['nsigmas'] = 2
         self.model.setParam('scale', 10.0)
         
         new_model = self.model.clone()
-        self.assertAlmostEqual(new_model.run(0.001), 1.1804794*4723.32213339, 3)
-        self.assertAlmostEqual(new_model.runXY([0.001,0.001]), 1.180454*4743.56, 2)
+        self.assertAlmostEqual(new_model.run(0.001), 
+                               1.1804794*4723.32213339, 3)
+        self.assertAlmostEqual(new_model.runXY([0.001,0.001]), 
+                               1.180454*4743.56, 2)
         
     def test_gaussian_zero(self):
@@ -130,5 +143,6 @@
         """
         self.assertAlmostEqual(self.model.run(0.001), 353.55013216754583, 3)
-        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 355.25355270620543, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 
+                               355.25355270620543, 3)
         
     def test_dispersion(self):
@@ -137,5 +151,7 @@
         """
         from sans.models.DisperseModel import DisperseModel
-        disp = DisperseModel(self.model, ['radius', 'thickness', 'length'], [5, 2, 50])
+        disp = DisperseModel(self.model, ['radius', 
+                                          'thickness', 
+                                          'length'], [5, 2, 50])
         disp.setParam('n_pts', 10)
         self.assertAlmostEqual(disp.run(0.001), 358.44062724936009, 3)
@@ -146,22 +162,24 @@
         disp_rm = GaussianDispersion()
         self.model.set_dispersion('radius', disp_rm)
-        self.model.dispersion['radius']['width'] = 5.0
+        self.model.dispersion['radius']['width'] = 0.25
         self.model.dispersion['radius']['npts'] = 10
-        self.model.dispersion['radius']['nsigmas'] = 2.5
+        self.model.dispersion['radius']['nsigmas'] = 2
 
         disp_rr = GaussianDispersion()
         self.model.set_dispersion('thickness', disp_rr)
-        self.model.dispersion['thickness']['width'] = 2.
+        self.model.dispersion['thickness']['width'] = 0.2
         self.model.dispersion['thickness']['npts'] = 10
-        self.model.dispersion['thickness']['nsigmas'] = 2.5
+        self.model.dispersion['thickness']['nsigmas'] = 2
 
         disp_len = GaussianDispersion()
         self.model.set_dispersion('length', disp_len)
-        self.model.dispersion['length']['width'] = 50.0
+        self.model.dispersion['length']['width'] = 1.0/8.0
         self.model.dispersion['length']['npts'] = 10
-        self.model.dispersion['length']['nsigmas'] = 2.5
-
-        self.assertAlmostEqual(self.model.run(0.001), 1.07832610*358.44062724936009, 3)
-        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 1.07844010*360.22673635224584, 3)
+        self.model.dispersion['length']['nsigmas'] = 2
+
+        self.assertAlmostEqual(self.model.run(0.001), 
+                               1.07832610*358.44062724936009, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 
+                               1.07844010*360.22673635224584, 3)
         
 
@@ -218,6 +236,8 @@
             the update to C++ underlying class.
         """
-        self.assertAlmostEqual(self.model.run(0.001), 381.27304697150055, 3)
-        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 380.93779156218682, 3)
+        self.assertAlmostEqual(self.model.run(0.001), 
+                               381.27304697150055, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 
+                               380.93779156218682, 3)
         
     def test_dispersion(self):
@@ -234,15 +254,16 @@
         disp_rm = GaussianDispersion()
         self.model.set_dispersion('radius', disp_rm)
-        self.model.dispersion['radius']['width'] = 10.
+        self.model.dispersion['radius']['width'] = 0.1666666667
         self.model.dispersion['radius']['npts'] = 10
-        self.model.dispersion['radius']['nsigmas'] = 2.5
+        self.model.dispersion['radius']['nsigmas'] = 2
 
         disp_rr = GaussianDispersion()
         self.model.set_dispersion('thickness', disp_rr)
-        self.model.dispersion['thickness']['width'] = 2.
+        self.model.dispersion['thickness']['width'] = 0.2
         self.model.dispersion['thickness']['npts'] = 10
-        self.model.dispersion['thickness']['nsigmas'] = 2.5
-
-        self.assertAlmostEqual(self.model.run(0.001), 1.16747510*407.344127907553, 3)
+        self.model.dispersion['thickness']['nsigmas'] = 2
+
+        self.assertAlmostEqual(self.model.run(0.001), 
+                               1.16747510*407.344127907553, 3)
 
         
@@ -270,6 +291,8 @@
             the update to C++ underlying class.
         """
-        self.assertAlmostEqual(self.model.run(0.001), 11808.842896863147, 3)
-        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 11681.990374929677, 3)
+        self.assertAlmostEqual(self.model.run(0.001), 
+                               11808.842896863147, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 
+                               11681.990374929677, 3)
 
     def test_dispersion(self):
@@ -288,16 +311,18 @@
         disp_rm = GaussianDispersion()
         self.model.set_dispersion('radius_a', disp_rm)
-        self.model.dispersion['radius_a']['width'] = 5.0
+        self.model.dispersion['radius_a']['width'] = 0.25
         self.model.dispersion['radius_a']['npts'] = 10
-        self.model.dispersion['radius_a']['nsigmas'] = 2.5
+        self.model.dispersion['radius_a']['nsigmas'] = 2
 
         disp_rr = GaussianDispersion()
         self.model.set_dispersion('radius_b', disp_rr)
-        self.model.dispersion['radius_b']['width'] = 50.
+        self.model.dispersion['radius_b']['width'] = 0.125
         self.model.dispersion['radius_b']['npts'] = 10
-        self.model.dispersion['radius_b']['nsigmas'] = 2.5
-
-        self.assertAlmostEqual(self.model.run(0.001), 1.10650710*11948.72581312305, 3)
-        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 1.105898*11811.972359807551, 2)
+        self.model.dispersion['radius_b']['nsigmas'] = 2
+
+        self.assertAlmostEqual(self.model.run(0.001), 
+                               1.10650710*11948.72581312305, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 
+                               1.105898*11811.972359807551, 2)
 
     def test_array(self):
@@ -341,6 +366,6 @@
         self.model.setParam('scale', 1.0)
         self.model.setParam('radius', 60.0)
-        self.model.setParam('sldSph', 2.0e-06)
-        self.model.setParam('sldSolv', 1.0e-6)
+        self.model.setParam('sldSph', 2.0)
+        self.model.setParam('sldSolv', 1.0)
         self.model.setParam('background', 0.0)
         
@@ -351,6 +376,8 @@
             the update to C++ underlying class.
         """
-        self.assertAlmostEqual(self.model.run(0.001), 90412744456148.094, 3)
-        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 90347660670656.391, 1)
+        self.assertAlmostEqual(self.model.run(0.001), 
+                               90412744456148.094, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 
+                               90347660670656.391, 1)
 
     def test_dispersion(self):
@@ -361,4 +388,6 @@
         disp = DisperseModel(self.model, ['radius'], [10])
         disp.setParam('n_pts', 10)
+        disp.setParam('radius.npts', 10)
+        disp.setParam('radius.nsigmas', 2.5)
         self.assertAlmostEqual(disp.run(0.001), 96795008379475.25, 1)
         
@@ -367,9 +396,9 @@
         disp_rm = GaussianDispersion()
         self.model.set_dispersion('radius', disp_rm)
-        self.model.dispersion['radius']['width'] = 10.0
+        self.model.dispersion['radius']['width'] = 0.1666666667
         self.model.dispersion['radius']['npts'] = 10
-        self.model.dispersion['radius']['nsigmas'] = 2.5
-
-        self.assertAlmostEqual(self.model.run(0.001), 1.19864810*96795008379475.25,3)
+        self.model.dispersion['radius']['nsigmas'] = 2
+
+        #self.assertAlmostEqual(self.model.run(0.001), 96795008379475.25,3)
         
         
@@ -400,6 +429,8 @@
             the update to C++ underlying class.
         """
-        self.assertAlmostEqual(self.model.run(0.001), 675.50440232504991, 3)
-        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 669.5173937622792, 0)
+        self.assertAlmostEqual(self.model.run(0.001), 
+                               675.50440232504991, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 
+                               669.5173937622792, 0)
         
     def test_dispersion(self):
@@ -408,5 +439,6 @@
         """
         from sans.models.DisperseModel import DisperseModel
-        disp = DisperseModel(self.model, ['r_minor', 'r_ratio', 'length'], [5, 0.25, 50])
+        disp = DisperseModel(self.model, ['r_minor', 'r_ratio', 'length'], 
+                             [5, 0.25, 50])
         disp.setParam('n_pts', 10)
         self.assertAlmostEqual(disp.run(0.001), 711.18048194151925, 3)
@@ -417,22 +449,24 @@
         disp_rm = GaussianDispersion()
         self.model.set_dispersion('r_minor', disp_rm)
-        self.model.dispersion['r_minor']['width'] = 5.0
+        self.model.dispersion['r_minor']['width'] = 0.25
         self.model.dispersion['r_minor']['npts'] = 10
-        self.model.dispersion['r_minor']['nsigmas'] = 2.5
+        self.model.dispersion['r_minor']['nsigmas'] = 2
 
         disp_rr = GaussianDispersion()
         self.model.set_dispersion('r_ratio', disp_rr)
-        self.model.dispersion['r_ratio']['width'] = 0.25
+        self.model.dispersion['r_ratio']['width'] = 0.25/1.5
         self.model.dispersion['r_ratio']['npts'] = 10
-        self.model.dispersion['r_ratio']['nsigmas'] = 2.5
+        self.model.dispersion['r_ratio']['nsigmas'] = 2
 
         disp_len = GaussianDispersion()
         self.model.set_dispersion('length', disp_len)
-        self.model.dispersion['length']['width'] = 50.0
+        self.model.dispersion['length']['width'] = 50.0/400
         self.model.dispersion['length']['npts'] = 10
-        self.model.dispersion['length']['nsigmas'] = 2.5
-
-        self.assertAlmostEqual(self.model.run(0.001), 1.23925910*711.18048194151925, 3)
-        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 1.238955*704.63525988095705, 0)
+        self.model.dispersion['length']['nsigmas'] = 2
+
+        self.assertAlmostEqual(self.model.run(0.001), 
+                               1.23925910*711.18048194151925, 3)
+        self.assertAlmostEqual(self.model.runXY([0.001,0.001]), 
+                               1.238955*704.63525988095705, 0)
         
 
@@ -479,7 +513,7 @@
         
         self.assertEqual(self.model.dispersion['radius']['width'], 0.0)
-        self.model.setParam('radius.width', 5.0)
-        self.assertEqual(self.model.dispersion['radius']['width'], 5.0)
-        self.assertEqual(self.model.getParam('radius.width'), 5.0)
+        self.model.setParam('radius.width', 0.25)
+        self.assertEqual(self.model.dispersion['radius']['width'], 0.25)
+        self.assertEqual(self.model.getParam('radius.width'), 0.25)
         self.assertEqual(self.model.dispersion['radius']['type'], 'gaussian')
         
